37 , fBoostType(
"CmsFixed")
38 , fUseLorentzBoost(kFALSE)
41 , fSpectatorsON(kTRUE)
42 , fInputFileVersion(1)
44 cout <<
"-I MpdDCMSMMGenerator: Opening input file " << fFileName.Data() << endl;
46 fInputFile = gzopen(fFileName.Data(),
"rb");
48 fInputFile = fopen(fFileName.Data(),
"r");
51 Fatal(
"MpdDCMSMMGenerator",
"Cannot open input file.");
58 Int_t header_size = 1;
60 for (Int_t
i = 0;
i < header_size;
i++) {
62 gzgets(fInputFile, read, 200);
64 fgets(read, 200, fInputFile);
66 cout <<
"-I MpdDCMSMMGenerator:" << read;
70 if (str0.find(
"Results of DCM-SMM") != string::npos) {
71 fInputFileVersion = 1;
73 }
else if (str0.find(
"Results of QGSM") != string::npos) {
74 fInputFileVersion = 0;
77 Fatal(
"MpdDCMSMMGenerator",
"Wrong input file format.");
81 A1 = stoi(str0.substr(str0.find(
"A1=") + 3, 3));
82 Z1 = stoi(str0.substr(str0.find(
"Z1=") + 3, 3));
83 A2 = stoi(str0.substr(str0.find(
"A2=") + 3, 3));
84 Z2 = stoi(str0.substr(str0.find(
"Z2=") + 3, 3));
88 T0 = stof(str0.substr(str0.find(
"T0=") + 3, 8));
89 sqS = stof(str0.substr(str0.find(
"sqrt(s)=") + 8, 8));
93 Double_t p = pow(e * e - mp * mp, 0.5);
94 Double_t mCMS = pow(2. * mp * (e + mp), 0.5);
95 cout <<
"mCMS=" << mCMS << endl;
97 if (fInputFileVersion == 0) {
98 fUseLorentzBoost = kTRUE;
99 if (fBoostType ==
"None")
100 fUseLorentzBoost = kFALSE;
101 else if (fBoostType ==
"CmsFixed") {
102 fBoostGamma = (e + mp) / mCMS;
103 fBoostBeta = pow(1. - 1. / (fBoostGamma * fBoostGamma), 0.5);
104 }
else if (fBoostType ==
"CmsFixedInverted") {
105 fBoostGamma = -(e + mp) / mCMS;
106 fBoostBeta = -pow(1. - 1. / (fBoostGamma * fBoostGamma), 0.5);
107 }
else if (fBoostType ==
"FixedFixedInverted") {
108 fBoostGamma = -e / mp;
114 cout <<
"-I MpdDCMSMMGenerator: A1=" << A1 <<
" Z1=" << Z1 <<
" A2=" << A2 <<
" Z2=" << Z2 <<
" T0=" << T0
115 <<
" sqS=" << sqS << endl;
143 cout <<
"-E MpdDCMSMMGenerator: Input file not open! " << endl;
149 cout <<
"-E- MpdDCMSMMGenerator::ReadEvent: "
150 <<
"No PrimaryGenerator!" << endl;
156 gzgets(fInputFile, read, 128);
158 fgets(read, 128, fInputFile);
161 Float_t b, bimpX, bimpY;
165 if (fInputFileVersion == 0)
166 sscanf(read,
"%d %f %f %f", &evnr, &b, &bimpX, &bimpY);
168 sscanf(read,
"%d %d %f %f %f", &evnr, &npp, &b, &bimpX, &bimpY);
171 if (gzeof(fInputFile)) {
173 if (feof(fInputFile)) {
175 cout <<
"-I MpdDCMSMMGenerator : End of input file reached." << endl;
176 const Bool_t ZeroSizeEvents = kFALSE;
190 Float_t phi =
atan2(bimpY, bimpX);
193 FairMCEventHeader*
event = primGen->GetEvent();
194 if (event && (!event->IsSet())) {
195 event->SetEventID(evnr);
197 event->MarkSet(kTRUE);
201 TDatabasePDG* dbPDG = TDatabasePDG::Instance();
204 for (Int_t ibeam = 0; ibeam < 3; ibeam++) {
207 gzgets(fInputFile, read, 128);
209 fgets(read, 128, fInputFile);
211 sscanf(read,
"%d", &np);
216 for (Int_t
i = 0;
i < np;
i++) {
218 gzgets(fInputFile, read, 128);
220 fgets(read, 128, fInputFile);
222 Int_t iN, iQ, iS = 0;
223 Float_t xxx = 0., mass;
224 Int_t pdgID_in = 0, iL;
226 if (fInputFileVersion == 0) {
228 sscanf(read,
"%d %d %f %f %f %f", &iN, &iQ, &xxx, &px, &py, &pz);
230 sscanf(read,
"%d %d %d %f %f %f %f", &iN, &iQ, &iS, &px, &py, &pz, &mass);
232 Float_t var10, var11;
233 Int_t res = sscanf(read,
"%d%d%d%d%d%f%f%f%f%f%f", &iQ, &iL, &iS, &iN, &pdgID_in, &px, &py, &pz, &pLabZ,
244 Fatal(__func__,
": data format mismatch, event %d\n", evnr);
249 Int_t pdgID = pdgID_in;
250 if (fInputFileVersion == 0) {
252 Double_t massFactor = 1.;
254 if (massFactor != 1.) {
268 if (pdgID == 310 || pdgID == 130 || TMath::Abs(pdgID) == 311) {
269 Double_t kSL = gRandom->Uniform(0., 1.);
277 TParticlePDG* particle = pdgID ? dbPDG->GetParticle(pdgID) : 0;
279 if (fInputFileVersion == 0) {
280 if (fUseLorentzBoost) {
282 Double_t
m = particle->Mass();
283 Double_t e = pow(px * px + py * py + pz * pz +
m *
m, 0.5);
284 Double_t pzF = fBoostGamma * (pz + fBoostBeta * e);
294 cout <<
"pz=" << pz <<
" N=" << iN <<
" Q=" << iQ <<
" pdg=" << pdgID <<
"\n";
295 primGen->AddTrack(pdgID, px, py, pz, 0., 0., 0.);
298 if (fInputFileVersion == 1 || ibeam == 2 || fSpectatorsON) {
300 cout <<
"-I MpdDCMSMMGenerator : unknown particle N=" << iN <<
" Q=" << iQ;
301 if (fInputFileVersion == 1)
302 cout <<
" PDG_in=" << pdgID_in <<
" ibeam=" << ibeam << endl;
317 for (Int_t ii = 0; ii < count; ii++) {
320 cout <<
"-E MpdDCMSMMGenerator: Input file not open! " << endl;
325 gzgets(fInputFile, read, 128);
327 fgets(read, 128, fInputFile);
331 for (Int_t ibeam = 0; ibeam < 3; ibeam++) {
334 gzgets(fInputFile, read, 128);
336 fgets(read, 128, fInputFile);
338 sscanf(read,
"%d", &np);
339 for (Int_t
i = 0;
i < np;
i++) {
341 gzgets(fInputFile, read, 128);
343 fgets(read, 128, fInputFile);
360 Int_t aA = abs(A), aS = abs(S), aZ = abs(Z);
365 const Int_t n000 = 7;
367 Float_t M000L[n000] = {0.0, 0.13400, 0.54700, 0.7680, 0.78200, 0.95700, 1.019390};
368 Float_t M000U[n000] = {0.0, 0.13510, 0.54810, 0.7711, 0.78310, 0.95810, 1.019610};
369 Int_t C000[n000] = {22, 111, 211, 113, 223, 331, 333};
370 for (Int_t
i = 0;
i < n000;
i++) {
371 if (mass < M000L[
i] || mass > M000U[
i])
377 if (k == 0 && mass > 0.4979 && mass < 0.4981)
381 if (mass > 0.1389 && mass < 0.1401)
383 else if (mass > 0.134 && mass < 0.136)
385 else if (mass == 0.776)
387 else if (mass > 0.0004 && mass < 0.0011)
389 else if (mass > 0.1055 && mass < 0.1065)
394 if (mass > 0.491 && mass < 0.4981) {
395 Double_t KSL = gRandom->Uniform(0., 1.);
400 }
else if (mass == 0.8922)
404 if (mass > 0.491 && mass < 0.4951)
406 else if (mass > 0.8815 && mass > 0.8895)
410 }
else if (aA == 1) {
413 if ((mass > 0.9389 && mass < 0.940) || mass < -1.)
415 else if (mass == 1.232)
417 }
else if (aZ == 1) {
418 if ((mass > 0.9379 && mass < 0.9401) || mass < -1.)
420 else if (mass == 1.232) {
423 else if (A * Z == -1)
426 }
else if (aZ == 2) {
430 }
else if (aS == 1) {
432 if (mass > 1.1149 && mass < 1.1161)
434 else if (mass > 1.189 && mass < 1.1921)
436 else if (mass == 1.382)
438 }
else if (aZ == 1) {
439 if (mass > 1.1889 && mass < 1.1921)
441 else if (mass > 1.196 && mass < 1.1971)
443 else if (mass == 1.3823)
445 else if (mass == 1.3875)
448 }
else if (aS == 2) {
450 if (mass > 1.3139 && mass < 1.3161)
452 else if (mass == 1.5318)
454 }
else if (aZ == 1) {
455 if (mass > 1.3139 && mass < 1.3221)
457 else if (mass == 1.535)
460 }
else if (aS == 3) {
462 if (mass > 1.6 && mass < 1.801)
469 if (Z == 1 && (mass < -1 || (mass > 1.875 && mass < 1.877)))
471 }
else if (S == -1) {
474 massFactor = 1.197436 / mass;
476 else if (Z == 0 && mass > 2.0545 && mass < 2.0555)
479 massFactor = 1.11568 / mass;
481 else if (Z == 1 && mass > 2.0535 && mass < 2.0545)
484 massFactor = 1.11568 / mass;
489 massFactor = 1.18937 / mass;
491 }
else if (S == -2) {
492 if (Z == -1 && mass > 2.2545 && mass < 2.2615) {
494 massFactor = 1.32132 / mass;
498 if (mass > 2.2305 && mass < 2.2315) {
500 massFactor = 1.11568 / mass;
502 else if (mass > 2.2525 && mass < 2.2555)
505 massFactor = 1.19255 / mass;
507 else if (mass > 2.2595 && mass < 2.2605)
510 massFactor = 1.197436 / mass;
512 }
else if (Z == 1 && mass > 2.2525 && mass < 2.2535) {
514 massFactor = 1.3149 / mass;
516 }
else if (S == -3) {
517 if (Z == -1 && mass > 2.4305 && mass < 2.4375) {
519 massFactor = 1.197436 / mass;
521 else if (Z == 0 && mass > 2.4301 && mass < 2.4315)
524 massFactor = 1.3149 / mass;
526 }
else if (S == -4 && Z == -1 && mass > 2.6295 && mass < 2.6365) {
528 massFactor = 1.32132 / mass;
533 if (Z == 1 && mass > 2.8085 && mass < 2.8095)
535 else if (Z == 2 && mass > 2.8085 && mass < 2.8095)
537 }
else if (S == -1) {
538 if (Z == 1 && mass > 2.9915 && mass < 2.9925)
543 if (Z == 2 && mass > 3.7275 && mass < 3.7285)
546 }
else if (S == -1) {
547 if (Z == 1 && mass > 3.9245 && mass < 3.9255)
549 if (Z == 2 && mass > 3.9245 && mass < 3.9255)
551 }
else if (S == -2) {
559 if (mass > -10 && k == 0) {
562 cout <<
"A=" << A <<
" S=" << S <<
" Z=" << Z <<
" mass= " << mass <<
" Pdg=" << k << endl;