Xrootd verbosity switched on for AliAnalysisManager debug level >2
[u/mrichter/AliRoot.git] / TUHKMgen / AliGenUHKM.cxx
1 /////////////////////////////////////////////////////////////////////////////
2 // Generator using UHKM 3.0 as an external generator.                      //
3 // ( only the HYDJET++ part is implemented for a moment)                   //
4 // temporary link:                                                         //
5 // http://lav01.sinp.msu.ru/~igor/hydjet++/hydjet++.txt                    //
6 // The main UHKM options are accessable  through this interface.           //
7 // Uses the TUHKMgen implementation of TGenerator.                         //
8 // Author of the first implementation: Sergey Zaporozhets                  //
9 // (zaporozh@sunhe.jinr.ru)                                                //
10 // Futhers modifications were made by                                      //
11 // Ionut Cristian Arsene (i.c.arsene@fys.uio.no)                           //
12 // & Malinina Liudmila(malinina@lav01.sinp.msu.ru) using as an example     //
13 //  AliGenTherminator.cxx created by Adam Kisiel                           //
14 //                                                                         //
15 ////////////////////////////////////////////////////////////////////////////
16
17 #include "TUHKMgen.h"
18 #ifndef DATABASE_PDG
19 #include "DatabasePDG.h"
20 #endif
21 #ifndef PARTICLE_PDG
22 #include "ParticlePDG.h"
23 #endif
24 #include <TLorentzVector.h>
25 #include <TPDGCode.h>
26 #include <TParticle.h>
27 #include <TClonesArray.h>
28 #include <TMCProcess.h>
29 #include "TDatabasePDG.h"
30 #include "TSystem.h"
31
32 #include "AliGenUHKM.h"
33 #include "AliRun.h"
34 #include "AliConst.h"
35 #include "AliDecayer.h"
36 #include "AliGenEventHeader.h"
37 #include "AliGenHijingEventHeader.h"
38 #include "AliLog.h"
39
40 #include <iostream>
41 #include <string>
42 using namespace std;
43
44
45 ClassImp(AliGenUHKM)
46
47 //_______________________________________
48 AliGenUHKM::AliGenUHKM()
49   :AliGenMC(),
50    fTrials(0),
51    fUHKMgen(0),
52    fHydjetParams(),
53    fStableFlagged(0)
54 {
55   // Default constructor setting up default reasonable parameter values
56   // for central Pb+Pb collisions at 5.5TeV 
57   
58   // LHC
59   fHydjetParams.fSqrtS=5500; //LHC
60   fHydjetParams.fAw=207;//Pb-Pb
61   fHydjetParams.fBmin=0.;
62   fHydjetParams.fBmax=0.5; //0-5% centrality
63   fHydjetParams.fT = 0.170;
64   fHydjetParams.fMuB = 0.0;
65   fHydjetParams.fMuS = 0.0;
66   fHydjetParams.fMuI3 = 0.0;
67   fHydjetParams.fThFO = 0.130;
68   fHydjetParams.fMu_th_pip = 0.0;
69   fHydjetParams.fSeed=0;
70   fHydjetParams.fTau=10.;
71   fHydjetParams.fSigmaTau=3.;
72   fHydjetParams.fR=11.;
73   fHydjetParams.fYlmax=4.0;
74   fHydjetParams.fUmax=1.1;
75   fHydjetParams.fDelta=0.;
76   fHydjetParams.fEpsilon=0.;
77   fHydjetParams.fWeakDecay=0; //>=0 on ,-1 off
78   fHydjetParams.fEtaType=1;//gaus
79   fHydjetParams.fCorrS=1.;
80   fHydjetParams.fNhsel=2;
81   fHydjetParams.fIshad=1;
82   fHydjetParams.fPtmin=7.0;
83   fHydjetParams.fT0=0.8;
84   fHydjetParams.fTau0=0.1;
85   fHydjetParams.fNf=0;
86   fHydjetParams.fIenglu=0;
87   fHydjetParams.fIanglu=0;
88
89
90 /* RHIC
91   fHydjetParams.fSqrtS=200; //RHIC
92   fHydjetParams.fAw=197;//Au-Au
93   fHydjetParams.fBmin=0.;
94   fHydjetParams.fBmax=0.5; //0-5% centrality
95   fHydjetParams.fT = 0.165;
96   fHydjetParams.fMuB = 0.0285;
97   fHydjetParams.fMuS = 0.007;
98   fHydjetParams.fMuI3 = -0.001;
99   fHydjetParams.fThFO = 0.100;
100   fHydjetParams.fMu_th_pip = 0.053;
101   fHydjetParams.fSeed=0;
102   fHydjetParams.fTau=8.;
103   fHydjetParams.fSigmaTau=2.;
104   fHydjetParams.fR=10.;
105   fHydjetParams.fYlmax=3.3;
106   fHydjetParams.fUmax=1.1;
107   fHydjetParams.fDelta=0.;
108   fHydjetParams.fEpsilon=0.;
109   fHydjetParams.fWeakDecay=0; //>=0 on ,-1 off
110   fHydjetParams.fEtaType=1;//gaus
111   fHydjetParams.fCorrS=1.;
112   fHydjetParams.fNhsel=2;
113   fHydjetParams.fIshad=0;
114   fHydjetParams.fPtmin=3.4;
115   fHydjetParams.fT0=0.3;
116   fHydjetParams.fTau0=0.4;
117   fHydjetParams.fNf=2;
118   fHydjetParams.fIenglu=0;
119   fHydjetParams.fIanglu=0;
120 */
121   strcpy(fParticleFilename, Form("%s/TUHKMgen/UHKM/particles.data", gSystem->Getenv("ALICE_ROOT")));
122   strcpy(fDecayFilename, Form("%s/TUHKMgen/UHKM/tabledecay.txt", gSystem->Getenv("ALICE_ROOT")));
123   for(Int_t i=0; i<500; i++) {
124     fStableFlagPDG[i] = 0;
125     fStableFlagStatus[i] = kFALSE;
126   }
127   fStableFlagged = 0;
128
129 }
130
131 //_______________________________________
132 AliGenUHKM::AliGenUHKM(Int_t npart)
133   :AliGenMC(npart),
134    fTrials(0),
135    fUHKMgen(0),
136    fHydjetParams(),
137    fStableFlagged(0)
138 {
139   // Constructor specifying the size of the particle table
140   // and setting up default reasonable parameter values
141   // for central Pb+Pb collisions at 5.5TeV
142
143   fName = "UHKM";
144   fTitle= "Particle Generator using UHKM 3.0";
145
146   fNprimaries = 0;   
147
148   //LHC
149   fHydjetParams.fSqrtS=5500; //LHC
150   fHydjetParams.fAw=207;//Pb-Pb
151   fHydjetParams.fBmin=0.;
152   fHydjetParams.fBmax=0.5; //0-5% centrality
153   fHydjetParams.fT = 0.170;
154   fHydjetParams.fMuB = 0.0;
155   fHydjetParams.fMuS = 0.0;
156   fHydjetParams.fMuI3 = 0.0;
157   fHydjetParams.fThFO = 0.130;
158   fHydjetParams.fMu_th_pip = 0.0;
159   fHydjetParams.fSeed=0;
160   fHydjetParams.fTau=10.;
161   fHydjetParams.fSigmaTau=3.;
162   fHydjetParams.fR=11.;
163   fHydjetParams.fYlmax=4.0;
164   fHydjetParams.fUmax=1.1;
165   fHydjetParams.fDelta=0.;
166   fHydjetParams.fEpsilon=0.;
167   fHydjetParams.fWeakDecay=0; //>=0 on ,-1 off
168   fHydjetParams.fEtaType=1;//gaus
169   fHydjetParams.fCorrS=1.;
170   fHydjetParams.fNhsel=2;
171   fHydjetParams.fIshad=1;
172   fHydjetParams.fPtmin=7.0;
173   fHydjetParams.fT0=0.8;
174   fHydjetParams.fTau0=0.1;
175   fHydjetParams.fNf=0;
176   fHydjetParams.fIenglu=0;
177   fHydjetParams.fIanglu=0;
178
179 /*RHIC
180   fHydjetParams.fSqrtS=200; //RHIC
181   fHydjetParams.fAw=197;//Au-Au
182   fHydjetParams.fBmin=0.;
183   fHydjetParams.fBmax=0.5; //0-5% centrality
184   fHydjetParams.fT = 0.165;
185   fHydjetParams.fMuB = 0.0285;
186   fHydjetParams.fMuS = 0.007;
187   fHydjetParams.fMuI3 = -0.001;
188   fHydjetParams.fThFO = 0.100;
189   fHydjetParams.fMu_th_pip = 0.053;
190   fHydjetParams.fSeed=0;
191   fHydjetParams.fTau=8.;
192   fHydjetParams.fSigmaTau=2.;
193   fHydjetParams.fR=10.;
194   fHydjetParams.fYlmax=3.3;
195   fHydjetParams.fUmax=1.1;
196   fHydjetParams.fDelta=0.;
197   fHydjetParams.fEpsilon=0.;
198   fHydjetParams.fWeakDecay=0;//>=0 on ,-1 off
199   fHydjetParams.fEtaType=1;//gaus
200   fHydjetParams.fCorrS=1.;
201   fHydjetParams.fNhsel=2;
202   fHydjetParams.fIshad=1;
203   fHydjetParams.fPtmin=3.4;
204   fHydjetParams.fT0=0.3;
205   fHydjetParams.fTau0=0.4;
206   fHydjetParams.fNf=2;
207   fHydjetParams.fIenglu=0;
208   fHydjetParams.fIanglu=0;
209 */
210
211   strcpy(fParticleFilename, Form("%s/TUHKMgen/UHKM/particles.data", gSystem->Getenv("ALICE_ROOT")));
212   strcpy(fDecayFilename, Form("%s/TUHKMgen/UHKM/tabledecay.txt", gSystem->Getenv("ALICE_ROOT")));
213   for(Int_t i=0; i<500; i++) {
214     fStableFlagPDG[i] = 0;
215     fStableFlagStatus[i] = kFALSE;
216   }
217   fStableFlagged = 0;  
218
219 }
220
221 //__________________________________________
222 AliGenUHKM::~AliGenUHKM()
223 {
224   // Destructor, do nothing
225   //  delete fParticles;
226 }
227
228 void AliGenUHKM::SetAllParametersRHIC() 
229 {
230   // Set reasonable default parameters for 0-5% central Au+Au collisions
231   // at 200 GeV at RHIC
232   SetEcms(200.0);                // RHIC top energy
233   SetAw(197);                    // Au+Au
234   SetBmin(0.0);                  // 0%
235   SetBmax(0.5);                  // 5%
236   SetChFrzTemperature(0.165);    // T_ch = 165 MeV
237   SetMuB(0.0285);                // mu_B = 28.5 MeV
238   SetMuS(0.007);                 // mu_S = 7 MeV
239   SetMuQ(-0.001);                // mu_Q = -1 MeV
240   SetThFrzTemperature(0.100);    // T_th = 100 MeV
241   SetMuPionThermal(0.053);       // mu_th_pion = 53 MeV 
242   SetSeed(0);                    // use UNIX time
243   SetTauB(8.0);                  // tau = 8 fm/c
244   SetSigmaTau(2.0);              // sigma_tau = 2 fm/c
245   SetRmaxB(10.0);                // fR = 10 fm
246   SetYlMax(3.3);                 // fYmax = 3.3
247   SetEtaRMax(1.1);               // Umax = 1.1
248   SetMomAsymmPar(0.0);           // delta = 0.0
249   SetCoordAsymmPar(0.0);         // epsilon = 0.0
250   //  SetFlagWeakDecay(0);           // weak decay on (<0 off !!!)
251   SetEtaType(1);                 // gaus distributed with fYmax dispersion (0 means boost invariant)
252   SetGammaS(1.0);                // gammaS = 1.0 (no strangeness canonical suppresion)
253   SetPyquenNhsel(2);             // hydro on, jets on, jet quenching on
254   SetPyquenShad(1);              // shadowing on (0 off)
255   SetPyquenPtmin(3.4);           // ptmin = 3.4 GeV/c
256   SetPyquenT0(0.3);              // T0 = 300 MeV
257   SetPyquenTau0(0.4);            // tau0 = 0.4 fm/c
258   SetPyquenNf(2);                // 2 flavours
259   SetPyquenIenglu(0);            // radiative and collisional energy loss
260   SetPyquenIanglu(0);            // small gluon angular distribution
261 }
262
263 void AliGenUHKM::SetAllParametersLHC()
264 {
265   // Set reasonable default parameters for 0-5% central Pb+Pb collisions
266   // at 5.5 TeV at LHC
267   SetEcms(5500.0);               // LHC
268   SetAw(207);                    // Pb+Pb
269   SetBmin(0.0);                  // 0%
270   SetBmax(0.5);                  // 5%
271   SetChFrzTemperature(0.170);    // T_ch = 170 MeV
272   SetMuB(0.0);                   // mu_B = 0 MeV
273   SetMuS(0.0);                   // mu_S = 0 MeV
274   SetMuQ(0.0);                   // mu_Q = 0 MeV
275   SetThFrzTemperature(0.130);    // T_th = 130 MeV
276   SetMuPionThermal(0.0);         // mu_th_pion = 0 MeV
277   SetSeed(0);                    // use UNIX time
278   SetTauB(10.0);                 // tau = 10 fm/c
279   SetSigmaTau(3.0);              // sigma_tau = 3 fm/c
280   SetRmaxB(11.0);                // fR = 11 fm
281   SetYlMax(4.0);                 // fYmax = 4.0
282   SetEtaRMax(1.1);               // Umax = 1.1
283   SetMomAsymmPar(0.0);           // delta = 0.0
284   SetCoordAsymmPar(0.0);         // epsilon = 0.0
285   //  SetFlagWeakDecay(0);           // weak decay on (<0 off !!!)
286   SetEtaType(1);                 // gaus distributed with fYmax dispersion (0 means boost invariant)
287   SetGammaS(1.0);                // gammaS = 1.0 (no strangeness canonical suppresion)
288   SetPyquenNhsel(2);             // hydro on, jets on, jet quenching on
289   SetPyquenShad(1);              // shadowing on (0 off)
290   SetPyquenPtmin(7.0);           // ptmin = 7.0 GeV/c
291   SetPyquenT0(0.8);              // T0 = 800 MeV
292   SetPyquenTau0(0.1);            // tau0 = 0.4 fm/c
293   SetPyquenNf(0);                // 0 flavours
294   SetPyquenIenglu(0);            // radiative and collisional energy loss
295   SetPyquenIanglu(0);            // small gluon angular distribution
296 }
297
298 //_________________________________________
299 void AliGenUHKM::Init()
300 {
301   // Initialization of the TGenerator::TUHKMgen interface object. 
302   // Model input parameters are transmited to the TUHKMgen object which forwards them
303   // further to the model.
304   // HYDJET++ is initialized (average multiplicities are calculated, particle species definitions and decay
305   // channels are loaded, etc.)
306
307   SetMC(new TUHKMgen());
308   fUHKMgen = (TUHKMgen*) fMCEvGen;
309   SetAllParameters();
310   
311   AliGenMC::Init();
312   fUHKMgen->Initialize();
313   CheckPDGTable();    
314   fUHKMgen->Print();  
315 }
316
317 //________________________________________
318 void AliGenUHKM::Generate()
319 {
320   // Generate one HYDJET++ event, get the output and push particles further 
321   // to AliRoot's stack 
322
323   Float_t polar[3] = {0,0,0};
324   Float_t origin[3]   = {0,0,0};
325   Float_t origin0[3]  = {0,0,0};
326   Float_t v[3];
327   Float_t mass, energy;
328
329   Vertex();
330   for(Int_t j=0; j<3; j++) origin0[j] = fVertex[j];
331
332   fTrials = 0;
333
334   Int_t nt  = 0;
335
336   fUHKMgen->GenerateEvent();
337   fTrials++;
338
339   fUHKMgen->ImportParticles(&fParticles,"All");
340
341   Int_t np = fParticles.GetEntriesFast();
342
343
344   Int_t* idsOnStack = new Int_t[np];
345   Int_t* newPos     = new Int_t[np];
346   for(Int_t i=0; i<np; i++) newPos[i] = i;
347
348   //_________ Loop for particle selection
349   for(Int_t i=0; i<np; i++) {
350     TParticle *iparticle = (TParticle*)fParticles.At(i);
351     Int_t kf = iparticle->GetPdgCode();
352     Bool_t hasMother = (iparticle->GetFirstMother() >= 0);
353     Bool_t hasDaughter = (iparticle->GetNDaughters() > 0);
354     
355     if(hasDaughter) {
356       // This particle has decayed
357       // It will not be tracked
358       // Add it only once with coordinates not
359       // smeared with primary vertex position
360       Float_t p[3] = {p[0] = iparticle->Px(),
361                       p[1] = iparticle->Py(),
362                       p[2] = iparticle->Pz()};
363       mass = TDatabasePDG::Instance()->GetParticle(kf)->Mass();
364       energy = sqrt(mass*mass + p[0]*p[0] + p[1]*p[1] + p[2]*p[2]);
365       v[0] = iparticle->Vx();
366       v[1] = iparticle->Vy();
367       v[2] = iparticle->Vz();
368       Float_t time = iparticle->T();
369
370       Int_t imo = -1;
371       if(hasMother) {
372         imo = iparticle->GetFirstMother(); //index of mother particle in fParticles
373       } // if has mother
374       Bool_t trackFlag = kFALSE;   // tFlag is kFALSE --> do not track the particle
375       
376       PushTrack(trackFlag, (imo>=0 ? idsOnStack[imo] : imo), kf,
377                 p[0], p[1], p[2], energy,
378                 v[0], v[1], v[2], time,
379                 polar[0], polar[1], polar[2],
380                 (hasMother ? kPDecay : kPNoProcess), nt);
381       idsOnStack[i] = nt;
382       
383       fNprimaries++;
384       KeepTrack(nt);
385     }
386     else {
387       // This is a final state particle
388       // It will be tracked
389       // Add it TWICE to the stack !!!
390       // First time with event-wide coordinates (for femtoscopy) -
391       //   this one will not be tracked
392       // Second time with event-wide c0ordinates and vertex smearing
393       //   this one will be tracked
394       
395       Float_t p[3] = {p[0] = iparticle->Px(),
396                       p[1] = iparticle->Py(),
397                       p[2] = iparticle->Pz()};
398       mass = TDatabasePDG::Instance()->GetParticle(kf)->Mass();
399       energy = sqrt(mass*mass + p[0]*p[0] + p[1]*p[1] + p[2]*p[2]);
400       v[0] = iparticle->Vx();
401       v[1] = iparticle->Vy();
402       v[2] = iparticle->Vz();
403
404       Int_t type    = iparticle->GetStatusCode(); // 1-from jet / 0-from hydro 
405       Int_t coeffT=1;
406       if(type==1) coeffT=-1; //to separate particles from jets
407       
408
409       Int_t imo = -1;
410       
411       if(hasMother) {
412         imo = iparticle->GetFirstMother();
413       } // if has mother
414       
415       Bool_t trackFlag = kFALSE;  // tFlag = kFALSE --> do not track this one, its for femtoscopy
416       PushTrack(trackFlag, (imo>=0 ? idsOnStack[imo] : imo), kf,
417                 p[0], p[1], p[2], energy,
418                 v[0], v[1], v[2], (iparticle->T())*coeffT,
419                 polar[0], polar[1], polar[2],
420                 hasMother ? kPDecay:kPNoProcess, nt);
421       
422       idsOnStack[i] = nt;
423       fNprimaries++;
424       KeepTrack(nt);
425
426       origin[0] = origin0[0]+v[0];
427       origin[1] = origin0[1]+v[1];
428       origin[2] = origin0[2]+v[2];
429       imo = nt;
430       
431       trackFlag = kTRUE;    // tFlag = kTRUE --> track this one
432
433       PushTrack(trackFlag, imo, kf,
434                 p[0], p[1], p[2], energy,
435                 origin[0], origin[1], origin[2], iparticle->T(),
436                 polar[0], polar[1], polar[2],
437                 hasMother ? kPDecay:kPNoProcess, nt);
438       
439       fNprimaries++;
440       KeepTrack(nt);
441     }
442   }
443
444   SetHighWaterMark(fNprimaries);
445
446   TArrayF eventVertex;
447   eventVertex.Set(3);
448   eventVertex[0] = origin0[0];
449   eventVertex[1] = origin0[1];
450   eventVertex[2] = origin0[2];
451
452   // Builds the event header, to be called after each event
453   AliGenEventHeader* header = new AliGenHijingEventHeader("UHKM");
454
455   ((AliGenHijingEventHeader*) header)->SetNProduced(fNprimaries);
456   ((AliGenHijingEventHeader*) header)->SetPrimaryVertex(eventVertex);
457   ((AliGenHijingEventHeader*) header)->SetImpactParameter(0.0);
458   ((AliGenHijingEventHeader*) header)->SetTotalEnergy(0.0);
459   ((AliGenHijingEventHeader*) header)->SetHardScatters(0);
460   ((AliGenHijingEventHeader*) header)->SetParticipants(0, 0);
461   ((AliGenHijingEventHeader*) header)->SetCollisions(0, 0, 0, 0);
462   ((AliGenHijingEventHeader*) header)->SetSpectators(0, 0, 0, 0);
463   ((AliGenHijingEventHeader*) header)->SetReactionPlaneAngle(0);//evrot);
464
465   header->SetPrimaryVertex(fVertex);
466   AddHeader(header);
467   fCollisionGeometry = (AliGenHijingEventHeader*)  header;
468
469   delete idsOnStack;
470
471 }
472
473 void AliGenUHKM::Copy(TObject &) const
474 {
475   Fatal("Copy","Not implemented!\n");
476 }
477
478 void AliGenUHKM::SetAllParameters() {
479   // Forward all input parameters to the TGenerator::TUHKMgen object
480
481   fUHKMgen->SetEcms(fHydjetParams.fSqrtS);
482   fUHKMgen->SetBmin(fHydjetParams.fBmin);
483   fUHKMgen->SetBmax(fHydjetParams.fBmax);
484   fUHKMgen->SetAw(fHydjetParams.fAw);
485   fUHKMgen->SetSeed(fHydjetParams.fSeed);
486
487   fUHKMgen->SetChFrzTemperature(fHydjetParams.fT);
488   fUHKMgen->SetMuB(fHydjetParams.fMuB);
489   fUHKMgen->SetMuS(fHydjetParams.fMuS);
490   fUHKMgen->SetMuQ(fHydjetParams.fMuI3);
491   fUHKMgen->SetTauB(fHydjetParams.fTau);
492   fUHKMgen->SetThFrzTemperature(fHydjetParams.fThFO);
493   fUHKMgen->SetMuPionThermal(fHydjetParams.fMu_th_pip);
494
495   fUHKMgen->SetSigmaTau(fHydjetParams.fSigmaTau);
496   fUHKMgen->SetRmaxB(fHydjetParams.fR);
497   fUHKMgen->SetYlMax(fHydjetParams.fYlmax);
498   fUHKMgen->SetEtaRMax(fHydjetParams.fUmax);
499   fUHKMgen->SetMomAsymmPar(fHydjetParams.fDelta);
500   fUHKMgen->SetCoordAsymmPar(fHydjetParams.fEpsilon);
501
502   fUHKMgen->SetGammaS(fHydjetParams.fCorrS);
503   fUHKMgen->SetEtaType(fHydjetParams.fEtaType);
504   fUHKMgen->SetFlagWeakDecay(fHydjetParams.fWeakDecay);
505
506   //PYQUEN parameters
507
508   fUHKMgen->SetPyquenNhsel(fHydjetParams.fNhsel);
509   fUHKMgen->SetPyquenShad(fHydjetParams.fIshad);
510   fUHKMgen->SetPyquenPtmin(fHydjetParams.fPtmin);
511   fUHKMgen->SetPyquenT0(fHydjetParams.fT0);
512   fUHKMgen->SetPyquenTau0(fHydjetParams.fTau0);
513   fUHKMgen->SetPyquenNf(fHydjetParams.fNf);
514   fUHKMgen->SetPyquenIenglu(fHydjetParams.fIenglu);
515   fUHKMgen->SetPyquenIanglu(fHydjetParams.fIanglu);
516
517   fUHKMgen->SetPDGParticleFile(fParticleFilename);
518   fUHKMgen->SetPDGDecayFile(fDecayFilename);
519   //  fUHKMgen->SetUseCharmParticles(fUseCharmParticles);
520   //  fUHKMgen->SetMinimumWidth(fMinWidth);
521   //  fUHKMgen->SetMaximumWidth(fMaxWidth);
522   //  fUHKMgen->SetMinimumMass(fMinMass);
523   //  fUHKMgen->SetMaximumMass(fMaxMass);
524
525   //  cout << "AliGenUHKM::Init() no. stable flagged particles = " << fStableFlagged << endl;
526   for(Int_t i=0; i<fStableFlagged; i++) {
527     //    cout << "AliGenUHKM::Init() flag no. " << i
528     //         << " PDG = " << fStableFlagPDG[i]
529     //         << " flag = " << fStableFlagStatus[i] << endl;
530     fUHKMgen->SetPDGParticleStable(fStableFlagPDG[i], fStableFlagStatus[i]);
531   }
532
533    cout<<" Print all parameters "<<endl;
534    cout<<" SqrtS = "<<fHydjetParams.fSqrtS<<endl;
535    cout<<" Bmin  = "<< fHydjetParams.fBmin<<endl;
536  cout<<" Bmax= "<<fHydjetParams.fBmax<<endl;
537  cout<<" Aw= "<<fHydjetParams.fAw<<endl;
538  cout<<" Seed= "<<fHydjetParams.fSeed<<endl;
539
540  cout<<" ---Stat-model parameters----------- "<<endl;
541
542  cout<<" ChFrzTemperature= "<<fHydjetParams.fT<<endl;
543  cout<<" MuB= "<<fHydjetParams.fMuB<<endl;
544  cout<<" MuS= "<<fHydjetParams.fMuS<<endl;
545  cout<<" MuQ= "<<fHydjetParams.fMuI3<<endl;
546  cout<<" TauB= "<<fHydjetParams.fTau<<endl;
547  cout<<" ThFrzTemperature= "<<fHydjetParams.fThFO<<endl;
548  cout<<" MuPionThermal= "<<fHydjetParams.fMu_th_pip<<endl;
549
550 cout<<"-----Volume parameters -------------- "<<endl;
551
552  cout<<" SigmaTau= "<<fHydjetParams.fSigmaTau<<endl;
553  cout<<" RmaxB= "<<fHydjetParams.fR<<endl;
554  cout<<" YlMax= "<<fHydjetParams.fYlmax<<endl;
555  cout<<" EtaRMax= "<<fHydjetParams.fUmax<<endl;
556  cout<<" MomAsymmPar= "<<fHydjetParams.fDelta<<endl;
557  cout<<" CoordAsymmPar= "<<fHydjetParams.fEpsilon<<endl;
558
559 cout<<" --------Flags------ "<<endl;
560
561  cout<<" GammaS= "<<fHydjetParams.fCorrS<<endl;
562  cout<<" EtaType= "<<fHydjetParams.fEtaType<<endl;
563  cout<<" FlagWeakDecay= "<<fHydjetParams.fWeakDecay<<endl;
564
565   cout<<"----PYQUEN parameters---"<<endl;
566
567   cout<<" Nhsel= "<<fHydjetParams.fNhsel<<endl;
568   cout<<" Shad= "<<fHydjetParams.fIshad<<endl;
569   cout<<" Ptmin= "<<fHydjetParams.fPtmin<<endl;
570   cout<<" T0= "<<fHydjetParams.fT0<<endl;
571   cout<<" Tau0= "<<fHydjetParams.fTau0<<endl;
572   cout<<" Nf= "<<fHydjetParams.fNf<<endl;
573   cout<<" Ienglu= "<<fHydjetParams.fIenglu<<endl;
574   cout<<" Ianglu= "<<fHydjetParams.fIanglu<<endl;
575
576   //  cout<<"----PDG table parameters---"<<endl;
577   
578   //  cout<<" UseCharmParticles= "<<fUseCharmParticles<<endl;
579   //  cout<<" MinimumWidth= "<<fMinWidth<<endl;
580   //  cout<<" MaximumWidth= "<<fMaxWidth<<endl;
581   //  cout<<" MinimumMass= "<<fMinMass<<endl;
582   //  cout<<" MaximumMass= "<<fMaxMass<<endl;
583
584   //  cout << "AliGenUHKM::SetAllParameters() OUT" << endl;
585 }
586
587 // add the additional PDG codes from UHKM(SHARE table) to ROOT's table
588 void AliGenUHKM::CheckPDGTable() {
589   // Add temporarely all particle definitions from HYDJET++ which miss in the ROOT's PDG tables
590   // to the TDatabasePDG table.
591
592   DatabasePDG *uhkmPDG = fUHKMgen->PDGInfo();              // UHKM's PDG table
593   TParticlePDG *rootTestParticle;
594   ParticlePDG *uhkmTestParticle;
595   
596   // loop over all particles in the SHARE table
597   for(Int_t i=0; i<uhkmPDG->GetNParticles(); i++) {
598     // get a particle specie
599     uhkmTestParticle = uhkmPDG->GetPDGParticleByIndex(i);
600     // check if this code exists in ROOT's table
601     rootTestParticle = TDatabasePDG::Instance()->GetParticle(uhkmTestParticle->GetPDG());
602     if(!rootTestParticle) {    // if not then add it to the ROOT's PDG database
603       
604       TDatabasePDG::Instance()->AddParticle(uhkmTestParticle->GetName(), uhkmTestParticle->GetName(), 
605                                             uhkmTestParticle->GetMass(), uhkmTestParticle->GetElectricCharge(),
606                                             (uhkmTestParticle->GetWidth()<1e-10 ? kTRUE : kFALSE),
607                                             uhkmTestParticle->GetWidth(), 
608                                             (Int_t(uhkmTestParticle->GetBaryonNumber())==0 ? "meson" : "baryon"),
609                                             uhkmTestParticle->GetPDG());    
610       if(uhkmTestParticle->GetWidth()<1e-10) 
611         cout << uhkmTestParticle->GetPDG() << " with mass "
612              << TDatabasePDG::Instance()->GetParticle(uhkmTestParticle->GetPDG())->Mass()
613              << TDatabasePDG::Instance()->GetParticle(uhkmTestParticle->GetPDG())->Width() << endl;  
614     }
615   }  // end for
616 }