File indexing completed on 2025-01-18 10:06:24
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012 #ifndef Pythia8_HeavyIons_H
0013 #define Pythia8_HeavyIons_H
0014
0015 #include "Pythia8/HIInfo.h"
0016 #include "Pythia8/PhysicsBase.h"
0017
0018 namespace Pythia8 {
0019
0020
0021 class Pythia;
0022
0023
0024
0025
0026
0027
0028
0029
0030
0031 class HeavyIons : public PhysicsBase {
0032
0033 public:
0034
0035
0036
0037
0038 HeavyIons(Pythia& mainPythiaIn)
0039 : mainPythiaPtr(&mainPythiaIn), HIHooksPtr(0),
0040 pythia(vector<Pythia*>(1, &mainPythiaIn)) {
0041 }
0042
0043
0044 virtual ~HeavyIons() {}
0045
0046
0047
0048
0049
0050
0051
0052 virtual bool init() = 0;
0053
0054
0055
0056
0057
0058 virtual bool next() = 0;
0059
0060
0061
0062 static void addSpecialSettings(Settings& settings);
0063
0064
0065
0066 static bool isHeavyIon(Settings& settings);
0067
0068
0069 bool setHIUserHooksPtr(HIUserHooksPtr userHooksPtrIn) {
0070 HIHooksPtr = userHooksPtrIn; return true;
0071 }
0072
0073
0074 virtual bool setKinematics(double ) {
0075 loggerPtr->ERROR_MSG("method not implemented for this heavy ion model");
0076 return false; }
0077 virtual bool setKinematics(double , double ) {
0078 loggerPtr->ERROR_MSG("method not implemented for this heavy ion model");
0079 return false; }
0080 virtual bool setKinematics(double, double, double, double, double, double) {
0081 loggerPtr->ERROR_MSG("method not implemented for this heavy ion model");
0082 return false; }
0083 virtual bool setKinematics(Vec4, Vec4) {
0084 loggerPtr->ERROR_MSG("method not implemented for this heavy ion model");
0085 return false; }
0086
0087
0088 virtual bool setBeamIDs(int , int = 0) {
0089 loggerPtr->ERROR_MSG("method not implemented for this heavy ion model");
0090 return false; }
0091
0092
0093
0094 HIInfo hiInfo;
0095
0096
0097 virtual void stat();
0098
0099 protected:
0100
0101
0102
0103 void updateInfo();
0104
0105
0106
0107
0108
0109
0110 void clearProcessLevel(Pythia& pyt);
0111
0112
0113 static void setupSpecials(Settings& settings, string match);
0114
0115
0116
0117 static void setupSpecials(Pythia& p, string match);
0118
0119
0120 int idProj, idTarg;
0121
0122
0123
0124 Pythia * mainPythiaPtr;
0125
0126
0127
0128 SigmaTotal sigTotNN;
0129
0130
0131
0132 HIUserHooksPtr HIHooksPtr;
0133
0134
0135
0136 vector<Pythia *> pythia;
0137
0138
0139 vector<string> pythiaNames;
0140
0141
0142
0143 vector<Info*> info;
0144
0145
0146
0147 struct InfoGrabber : public UserHooks {
0148
0149
0150 Info * getInfo() {
0151 return infoPtr;
0152 }
0153
0154 };
0155
0156 };
0157
0158
0159
0160
0161
0162 class Angantyr : public HeavyIons {
0163
0164 public:
0165
0166
0167 enum PythiaObject {
0168 HADRON = 0,
0169 MBIAS = 1,
0170 SASD = 2,
0171 SIGPP = 3,
0172 SIGPN = 4,
0173 SIGNP = 5,
0174 SIGNN = 6,
0175 ALL = 7
0176 };
0177
0178 public:
0179
0180
0181
0182
0183 Angantyr(Pythia& mainPythiaIn);
0184
0185 virtual ~Angantyr();
0186
0187
0188 virtual bool init() override;
0189
0190
0191 virtual bool next() override;
0192
0193
0194 bool setUserHooksPtr(PythiaObject sel, UserHooksPtr userHooksPtrIn);
0195
0196
0197 bool setKinematics(double eCMIn) override;
0198 bool setKinematics(double eAIn, double eBIn) override;
0199 bool setKinematics(double, double, double, double, double, double) override;
0200 bool setKinematics(Vec4, Vec4) override;
0201 bool setKinematics();
0202
0203
0204 bool setBeamIDs(int idAIn, int idBIn = 0) override;
0205
0206
0207 void unifyFrames();
0208
0209
0210 void banner() const;
0211
0212
0213 const SubCollisionSet& subCollisions() const {
0214 return subColls; }
0215
0216
0217 const SubCollisionModel& subCollisionModel() const {
0218 return *collPtr.get(); }
0219
0220
0221 const ImpactParameterGenerator impactParameterGenerator() const {
0222 return *bGenPtr.get(); }
0223
0224
0225 const Nucleus& projectile() const {
0226 return proj; }
0227
0228
0229 const Nucleus& target() const {
0230 return targ; }
0231
0232
0233 const NucleusModel& projectileModel() const {
0234 return *projPtr.get(); }
0235
0236
0237 const NucleusModel& targetModel() const {
0238 return *targPtr.get(); }
0239
0240
0241 const SigmaTotal sigmaNN() const {
0242 return sigTotNN; }
0243
0244 protected:
0245
0246 virtual void onInitInfoPtr() override {
0247 registerSubObject(sigTotNN); }
0248
0249
0250 void setBeamKinematics(int idA, int idB);
0251
0252
0253
0254 bool init(PythiaObject sel, string name, int n = 0);
0255
0256
0257 EventInfo mkEventInfo(Pythia &, Info &, const SubCollision * coll = 0);
0258
0259
0260 EventInfo getSignal(const SubCollision& coll);
0261 EventInfo getND() { return getMBIAS(0, 101); }
0262 EventInfo getND(const SubCollision& coll) { return getMBIAS(&coll, 101); }
0263 EventInfo getEl(const SubCollision& coll) { return getMBIAS(&coll, 102); }
0264 EventInfo getSDP(const SubCollision& coll) { return getMBIAS(&coll, 103); }
0265 EventInfo getSDT(const SubCollision& coll) { return getMBIAS(&coll, 104); }
0266 EventInfo getDD(const SubCollision& coll) { return getMBIAS(&coll, 105); }
0267 EventInfo getCD(const SubCollision& coll) { return getMBIAS(&coll, 106); }
0268 EventInfo getSDabsP(const SubCollision& coll)
0269 { return getSASD(&coll, 103); }
0270 EventInfo getSDabsT(const SubCollision& coll)
0271 { return getSASD(&coll, 104); }
0272 EventInfo getMBIAS(const SubCollision * coll, int procid);
0273 EventInfo getSASD(const SubCollision * coll, int procid);
0274
0275 bool genAbs(SubCollisionSet& subCollsIn, list<EventInfo>& subEventsIn);
0276 void addSASD(const SubCollisionSet& subCollsIn);
0277 bool addDD(const SubCollisionSet& subCollsIn, list<EventInfo>& subEventsIn);
0278 bool addSD(const SubCollisionSet& subCollsIn, list<EventInfo>& subEventsIn);
0279 void addSDsecond(const SubCollisionSet& subCollsIn);
0280 bool addCD(const SubCollisionSet& subCollsIn, list<EventInfo>& subEventsIn);
0281 void addCDsecond(const SubCollisionSet& subCollsIn);
0282 bool addEL(const SubCollisionSet& subCollsIn, list<EventInfo>& subEventsIn);
0283 void addELsecond(const SubCollisionSet& subCollsIn);
0284
0285 bool buildEvent(list<EventInfo>& subEventsIn);
0286
0287 bool setupFullCollision(EventInfo& ei, const SubCollision& coll,
0288 Nucleon::Status projStatus, Nucleon::Status targStatus);
0289 bool isRemnant(const EventInfo& ei, int i, int past = 1 ) const {
0290 int statNow = ei.event[i].status()*past;
0291 if ( statNow == 63 ) return true;
0292 if ( statNow > 70 && statNow < 80 )
0293 return isRemnant(ei, ei.event[i].mother1(), -1);
0294 return false;
0295 }
0296 bool fixIsoSpin(EventInfo& ei);
0297 EventInfo& shiftEvent(EventInfo& ei);
0298 static int getBeam(Event& ev, int i);
0299
0300
0301 bool nextSASD(int proc);
0302
0303
0304
0305 bool addNucleonExcitation(EventInfo& orig, EventInfo& add,
0306 bool colConnect = false);
0307
0308
0309
0310 vector<int> findRecoilers(const Event& e, bool tside, int beam, int end,
0311 const Vec4& pdiff, const Vec4& pbeam);
0312
0313
0314 void addSubEvent(Event& evnt, Event& sub);
0315 static void addJunctions(Event& evnt, Event& sub, int coloff);
0316
0317
0318
0319 bool addNucleusRemnants();
0320
0321 public:
0322
0323
0324
0325 static bool
0326 getTransforms(Vec4 p1, Vec4 p2, const Vec4& p1p,
0327 pair<RotBstMatrix,RotBstMatrix>& R12);
0328 static double mT2(const Vec4& p) { return p.pPos()*p.pNeg(); }
0329 static double mT(const Vec4& p) { return sqrt(max(mT2(p), 0.0)); }
0330
0331 private:
0332
0333
0334 struct ProcessSelectorHook: public UserHooks {
0335
0336 ProcessSelectorHook(): proc(0), b(-1.0) {}
0337
0338
0339 virtual bool canVetoProcessLevel() {
0340 return true;
0341 }
0342
0343
0344 virtual bool doVetoProcessLevel(Event&) {
0345 return proc > 0 && infoPtr->code() != proc;
0346 }
0347
0348
0349 virtual bool canSetImpactParameter() const {
0350 return b >= 0.0;
0351 }
0352
0353
0354 virtual double doSetImpactParameter() {
0355 return b;
0356 }
0357
0358
0359 int proc;
0360
0361
0362 double b;
0363
0364 };
0365
0366
0367 struct HoldProcess {
0368
0369
0370 HoldProcess(shared_ptr<ProcessSelectorHook> hook, int proc,
0371 double b = -1.0) : saveHook(hook), saveProc(hook->proc), saveB(hook->b) {
0372 hook->proc = proc;
0373 hook->b = b;
0374 }
0375
0376
0377 ~HoldProcess() {
0378 if ( saveHook ) {
0379 saveHook->proc = saveProc;
0380 saveHook->b = saveB;
0381 }
0382 }
0383
0384
0385 shared_ptr<ProcessSelectorHook> saveHook;
0386
0387
0388 int saveProc;
0389
0390
0391 double saveB;
0392
0393 };
0394
0395
0396 shared_ptr<ProcessSelectorHook> selectMB;
0397
0398
0399 shared_ptr<ProcessSelectorHook> selectSASD;
0400
0401 private:
0402
0403 static const int MAXTRY = 999;
0404 static const int MAXEVSAVE = 999;
0405
0406
0407
0408 bool hasSignal;
0409
0410
0411 bool doHadronLevel;
0412
0413
0414 bool doSDTest;
0415
0416
0417 bool glauberOnly;
0418
0419
0420 SubCollisionSet subColls;
0421
0422
0423
0424 shared_ptr<SubCollisionModel> collPtr;
0425
0426
0427 shared_ptr<ImpactParameterGenerator> bGenPtr;
0428
0429
0430 Nucleus proj;
0431 Nucleus targ;
0432
0433
0434
0435 shared_ptr<NucleusModel> projPtr;
0436 shared_ptr<NucleusModel> targPtr;
0437
0438
0439 bool doVarECM;
0440
0441
0442
0443 int recoilerMode;
0444
0445
0446 int bMode;
0447
0448
0449 bool doAbort;
0450
0451 };
0452
0453
0454
0455 }
0456
0457 #endif