File indexing completed on 2025-09-17 09:08: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 setKinematicsCM();
0198 bool setKinematics(double eCMIn) override;
0199 bool setKinematics(double eAIn, double eBIn) override;
0200 bool setKinematics(double, double, double, double, double, double) override;
0201 bool setKinematics(Vec4, Vec4) override;
0202 bool setKinematics();
0203
0204
0205 bool setBeamIDs(int idAIn, int idBIn = 0) override;
0206
0207
0208 void unifyFrames();
0209
0210
0211 void banner() const;
0212
0213
0214 const SubCollisionSet& subCollisions() const {
0215 return subColls; }
0216
0217
0218 const SubCollisionModel& subCollisionModel() const {
0219 return *collPtr.get(); }
0220
0221 SubCollisionModel* subCollPtr() {
0222 return collPtr.get();
0223 }
0224
0225
0226 const ImpactParameterGenerator impactParameterGenerator() const {
0227 return *bGenPtr.get(); }
0228
0229
0230 const Nucleus& projectile() const {
0231 return proj; }
0232
0233
0234 const Nucleus& target() const {
0235 return targ; }
0236
0237
0238 const NucleusModel& projectileModel() const {
0239 return *projPtr.get(); }
0240
0241
0242 const NucleusModel& targetModel() const {
0243 return *targPtr.get(); }
0244
0245
0246 const SigmaTotal sigmaNN() const {
0247 return sigTotNN; }
0248
0249 protected:
0250
0251 virtual void onInitInfoPtr() override {
0252 registerSubObject(sigTotNN); }
0253
0254
0255 void setBeamKinematics(int idA, int idB);
0256
0257
0258
0259 bool init(PythiaObject sel, string name, int n = 0);
0260
0261
0262 EventInfo mkEventInfo(Pythia &, Info &, const SubCollision * coll = 0);
0263
0264
0265 EventInfo getSignal(const SubCollision& coll);
0266 EventInfo getND() { return getMBIAS(0, 101); }
0267 EventInfo getND(const SubCollision& coll) { return getMBIAS(&coll, 101); }
0268 EventInfo getEl(const SubCollision& coll) { return getMBIAS(&coll, 102); }
0269 EventInfo getSDP(const SubCollision& coll) { return getMBIAS(&coll, 103); }
0270 EventInfo getSDT(const SubCollision& coll) { return getMBIAS(&coll, 104); }
0271 EventInfo getDD(const SubCollision& coll) { return getMBIAS(&coll, 105); }
0272 EventInfo getCD(const SubCollision& coll) { return getMBIAS(&coll, 106); }
0273 EventInfo getSDabsP(const SubCollision& coll)
0274 { return getSASD(&coll, 103); }
0275 EventInfo getSDabsT(const SubCollision& coll)
0276 { return getSASD(&coll, 104); }
0277 EventInfo getMBIAS(const SubCollision * coll, int procid);
0278 EventInfo getSASD(const SubCollision * coll, int procid);
0279
0280 bool genAbs(SubCollisionSet& subCollsIn, list<EventInfo>& subEventsIn);
0281 void addSASD(const SubCollisionSet& subCollsIn);
0282 bool addDD(const SubCollisionSet& subCollsIn, list<EventInfo>& subEventsIn);
0283 bool addSD(const SubCollisionSet& subCollsIn, list<EventInfo>& subEventsIn);
0284 void addSDsecond(const SubCollisionSet& subCollsIn);
0285 bool addCD(const SubCollisionSet& subCollsIn, list<EventInfo>& subEventsIn);
0286 void addCDsecond(const SubCollisionSet& subCollsIn);
0287 bool addEL(const SubCollisionSet& subCollsIn, list<EventInfo>& subEventsIn);
0288 void addELsecond(const SubCollisionSet& subCollsIn);
0289
0290 void resetEvent();
0291 bool buildEvent(list<EventInfo>& subEventsIn);
0292
0293 bool setupFullCollision(EventInfo& ei, const SubCollision& coll,
0294 Nucleon::Status projStatus, Nucleon::Status targStatus);
0295 bool isRemnant(const EventInfo& ei, int i, int past = 1 ) const {
0296 int statNow = ei.event[i].status()*past;
0297 if ( statNow == 63 ) return true;
0298 if ( statNow > 70 && statNow < 80 )
0299 return isRemnant(ei, ei.event[i].mother1(), -1);
0300 return false;
0301 }
0302 bool fixIsoSpin(EventInfo& ei);
0303 EventInfo& shiftEvent(EventInfo& ei);
0304 static int getBeam(Event& ev, int i);
0305
0306
0307 bool nextSASD(int proc);
0308
0309
0310
0311 bool addNucleonExcitation(EventInfo& orig, EventInfo& add,
0312 bool colConnect = false);
0313
0314
0315
0316 vector<int> findRecoilers(const Event& e, bool tside, int beam, int end,
0317 const Vec4& pdiff, const Vec4& pbeam);
0318
0319
0320 void addSubEvent(Event& evnt, Event& sub);
0321 static void addJunctions(Event& evnt, Event& sub, int coloff);
0322
0323
0324
0325 bool addNucleusRemnants();
0326
0327 public:
0328
0329
0330
0331 static bool
0332 getTransforms(Vec4 p1, Vec4 p2, const Vec4& p1p,
0333 pair<RotBstMatrix,RotBstMatrix>& R12);
0334 static double mT2(const Vec4& p) { return p.pPos()*p.pNeg(); }
0335 static double mT(const Vec4& p) { return sqrt(max(mT2(p), 0.0)); }
0336
0337 private:
0338
0339
0340 struct ProcessSelectorHook: public UserHooks {
0341
0342 ProcessSelectorHook(): proc(0), b(-1.0) {}
0343
0344
0345 virtual bool canVetoProcessLevel() {
0346 return true;
0347 }
0348
0349
0350 virtual bool doVetoProcessLevel(Event&) {
0351 return proc > 0 && infoPtr->code() != proc;
0352 }
0353
0354
0355 virtual bool canSetImpactParameter() const {
0356 return b >= 0.0;
0357 }
0358
0359
0360 virtual double doSetImpactParameter() {
0361 return b;
0362 }
0363
0364
0365 int proc;
0366
0367
0368 double b;
0369
0370 };
0371
0372
0373 struct HoldProcess {
0374
0375
0376 HoldProcess(shared_ptr<ProcessSelectorHook> hook, int proc,
0377 double b = -1.0) : saveHook(hook), saveProc(hook->proc), saveB(hook->b) {
0378 hook->proc = proc;
0379 hook->b = b;
0380 }
0381
0382
0383 ~HoldProcess() {
0384 if ( saveHook ) {
0385 saveHook->proc = saveProc;
0386 saveHook->b = saveB;
0387 }
0388 }
0389
0390
0391 shared_ptr<ProcessSelectorHook> saveHook;
0392
0393
0394 int saveProc;
0395
0396
0397 double saveB;
0398
0399 };
0400
0401
0402 shared_ptr<ProcessSelectorHook> selectMB;
0403
0404
0405 shared_ptr<ProcessSelectorHook> selectSASD;
0406
0407 private:
0408
0409 static const int MAXTRY = 999;
0410 static const int MAXEVSAVE = 999;
0411
0412
0413
0414 bool hasSignal;
0415
0416
0417 bool doHadronLevel;
0418
0419
0420 bool doSDTest;
0421
0422
0423 bool glauberOnly;
0424
0425
0426 SubCollisionSet subColls;
0427
0428
0429
0430 shared_ptr<SubCollisionModel> collPtr;
0431
0432
0433 shared_ptr<ImpactParameterGenerator> bGenPtr;
0434
0435
0436 Nucleus proj;
0437 Nucleus targ;
0438
0439
0440
0441 shared_ptr<NucleusModel> projPtr;
0442 shared_ptr<NucleusModel> targPtr;
0443
0444
0445 bool doVarECM;
0446
0447
0448
0449 int recoilerMode;
0450
0451
0452 int bMode;
0453
0454
0455 bool doAbort;
0456
0457 };
0458
0459
0460
0461 }
0462
0463 #endif