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     // !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
 351     // the particle indexes are 0 based but fParticles is a 1 based array
 352     // -1 is the trivial code (when it does not exist)
 353     // !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
 354     TParticle *iparticle = (TParticle*)fParticles.At(i);
 355     Int_t kf = iparticle->GetPdgCode();
 356     Bool_t hasMother = (iparticle->GetFirstMother() >= 0);
 357     Bool_t hasDaughter = (iparticle->GetNDaughters() > 0);
 358
 359     if(hasDaughter) {
 360       // This particle has decayed
 361       // It will not be tracked
 362       // Add it only once with coordinates not
 363       // smeared with primary vertex position
 364       Float_t p[3] = {p[0] = iparticle->Px(),
 365                       p[1] = iparticle->Py(),
 366                       p[2] = iparticle->Pz()};
 367       mass = TDatabasePDG::Instance()->GetParticle(kf)->Mass();
 368       energy = sqrt(mass*mass + p[0]*p[0] + p[1]*p[1] + p[2]*p[2]);
 369       v[0] = iparticle->Vx();
 370       v[1] = iparticle->Vy();
 371       v[2] = iparticle->Vz();
 372       Float_t time = iparticle->T();
 373
 374       Int_t imo = -1;
 375       if(hasMother) {
 376         imo = iparticle->GetFirstMother(); //index of mother particle in fParticles
 377       } // if has mother
 378       Bool_t trackFlag = kFALSE;   // tFlag is kFALSE --> do not track the particle
 379
 380       PushTrack(trackFlag, (imo>=0 ? idsOnStack[imo+1] : imo), kf,
 381                 p[0], p[1], p[2], energy,
 382                 v[0], v[1], v[2], time,
 383                 polar[0], polar[1], polar[2],
 384                 (hasMother ? kPDecay : kPNoProcess), nt);
 385       idsOnStack[i] = nt;
 386       fNprimaries++;
 387       KeepTrack(nt);
 388     }
 389     else {
 390       // This is a final state particle
 391       // It will be tracked
 392       // Add it TWICE to the stack !!!
 393       // First time with event-wide coordinates (for femtoscopy) -
 394       //   this one will not be tracked
 395       // Second time with event-wide c0ordinates and vertex smearing
 396       //   this one will be tracked
 397       Float_t p[3] = {p[0] = iparticle->Px(),
 398                       p[1] = iparticle->Py(),
 399                       p[2] = iparticle->Pz()};
 400       mass = TDatabasePDG::Instance()->GetParticle(kf)->Mass();
 401       energy = sqrt(mass*mass + p[0]*p[0] + p[1]*p[1] + p[2]*p[2]);
 402       v[0] = iparticle->Vx();
 403       v[1] = iparticle->Vy();
 404       v[2] = iparticle->Vz();
 405
 406       Int_t type    = iparticle->GetStatusCode(); // 1-from jet / 0-from hydro
 407       Int_t coeffT=1;
 408       if(type==1) coeffT=-1; //to separate particles from jets
 409
 410       Int_t imo = -1;
 411
 412       if(hasMother) {
 413         imo = iparticle->GetFirstMother();
 414       } // if has mother
 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       fNprimaries++;
 439       KeepTrack(nt);
 440     }
 441   }
 442
 443   SetHighWaterMark(fNprimaries);
 444
 445   TArrayF eventVertex;
 446   eventVertex.Set(3);
 447   eventVertex[0] = origin0[0];
 448   eventVertex[1] = origin0[1];
 449   eventVertex[2] = origin0[2];
 450
 451   // Builds the event header, to be called after each event
 452   AliGenEventHeader* header = new AliGenHijingEventHeader("UHKM");
 453
 454   ((AliGenHijingEventHeader*) header)->SetNProduced(fNprimaries);
 455   ((AliGenHijingEventHeader*) header)->SetPrimaryVertex(eventVertex);
 456   ((AliGenHijingEventHeader*) header)->SetImpactParameter(0.0);
 457   ((AliGenHijingEventHeader*) header)->SetTotalEnergy(0.0);
 458   ((AliGenHijingEventHeader*) header)->SetHardScatters(0);
 459   ((AliGenHijingEventHeader*) header)->SetParticipants(0, 0);
 460   ((AliGenHijingEventHeader*) header)->SetCollisions(0, 0, 0, 0);
 461   ((AliGenHijingEventHeader*) header)->SetSpectators(0, 0, 0, 0);
 462   ((AliGenHijingEventHeader*) header)->SetReactionPlaneAngle(0);//evrot);
 463
 464   header->SetPrimaryVertex(fVertex);
 465   AddHeader(header);
 466   fCollisionGeometry = (AliGenHijingEventHeader*)  header;
 467
 468   delete idsOnStack;
 469
 470 }
 471
 472 void AliGenUHKM::Copy(TObject &) const
 473 {
 474   Fatal("Copy","Not implemented!\n");
 475 }
 476
 477 void AliGenUHKM::SetAllParameters() {
 478   // Forward all input parameters to the TGenerator::TUHKMgen object
 479
 480   fUHKMgen->SetEcms(fHydjetParams.fSqrtS);
 481   fUHKMgen->SetBmin(fHydjetParams.fBmin);
 482   fUHKMgen->SetBmax(fHydjetParams.fBmax);
 483   fUHKMgen->SetAw(fHydjetParams.fAw);
 484   fUHKMgen->SetSeed(fHydjetParams.fSeed);
 485
 486   fUHKMgen->SetChFrzTemperature(fHydjetParams.fT);
 487   fUHKMgen->SetMuB(fHydjetParams.fMuB);
 488   fUHKMgen->SetMuS(fHydjetParams.fMuS);
 489   fUHKMgen->SetMuQ(fHydjetParams.fMuI3);
 490   fUHKMgen->SetTauB(fHydjetParams.fTau);
 491   fUHKMgen->SetThFrzTemperature(fHydjetParams.fThFO);
 492   fUHKMgen->SetMuPionThermal(fHydjetParams.fMu_th_pip);
 493
 494   fUHKMgen->SetSigmaTau(fHydjetParams.fSigmaTau);
 495   fUHKMgen->SetRmaxB(fHydjetParams.fR);
 496   fUHKMgen->SetYlMax(fHydjetParams.fYlmax);
 497   fUHKMgen->SetEtaRMax(fHydjetParams.fUmax);
 498   fUHKMgen->SetMomAsymmPar(fHydjetParams.fDelta);
 499   fUHKMgen->SetCoordAsymmPar(fHydjetParams.fEpsilon);
 500
 501   fUHKMgen->SetGammaS(fHydjetParams.fCorrS);
 502   fUHKMgen->SetEtaType(fHydjetParams.fEtaType);
 503   fUHKMgen->SetFlagWeakDecay(fHydjetParams.fWeakDecay);
 504
 505   //PYQUEN parameters
 506
 507   fUHKMgen->SetPyquenNhsel(fHydjetParams.fNhsel);
 508   fUHKMgen->SetPyquenShad(fHydjetParams.fIshad);
 509   fUHKMgen->SetPyquenPtmin(fHydjetParams.fPtmin);
 510   fUHKMgen->SetPyquenT0(fHydjetParams.fT0);
 511   fUHKMgen->SetPyquenTau0(fHydjetParams.fTau0);
 512   fUHKMgen->SetPyquenNf(fHydjetParams.fNf);
 513   fUHKMgen->SetPyquenIenglu(fHydjetParams.fIenglu);
 514   fUHKMgen->SetPyquenIanglu(fHydjetParams.fIanglu);
 515
 516   fUHKMgen->SetPDGParticleFile(fParticleFilename);
 517   fUHKMgen->SetPDGDecayFile(fDecayFilename);
 518   //  fUHKMgen->SetUseCharmParticles(fUseCharmParticles);
 519   //  fUHKMgen->SetMinimumWidth(fMinWidth);
 520   //  fUHKMgen->SetMaximumWidth(fMaxWidth);
 521   //  fUHKMgen->SetMinimumMass(fMinMass);
 522   //  fUHKMgen->SetMaximumMass(fMaxMass);
 523
 524   cout << "AliGenUHKM::Init() no. stable flagged particles = " << fStableFlagged << endl;
 525   for(Int_t i=0; i<fStableFlagged; i++) {
 526     cout << "AliGenUHKM::Init() flag no. " << i
 527          << " PDG = " << fStableFlagPDG[i]
 528          << " flag = " << fStableFlagStatus[i] << endl;
 529     fUHKMgen->SetPDGParticleStable(fStableFlagPDG[i], fStableFlagStatus[i]);
 530   }
 531
 532  cout<<" Print all parameters "<<endl;
 533  cout<<" SqrtS = "<<fHydjetParams.fSqrtS<<endl;
 534  cout<<" Bmin  = "<< fHydjetParams.fBmin<<endl;
 535  cout<<" Bmax= "<<fHydjetParams.fBmax<<endl;
 536  cout<<" Aw= "<<fHydjetParams.fAw<<endl;
 537  cout<<" Seed= "<<fHydjetParams.fSeed<<endl;
 538
 539  cout<<" ---Stat-model parameters----------- "<<endl;
 540
 541  cout<<" ChFrzTemperature= "<<fHydjetParams.fT<<endl;
 542  cout<<" MuB= "<<fHydjetParams.fMuB<<endl;
 543  cout<<" MuS= "<<fHydjetParams.fMuS<<endl;
 544  cout<<" MuQ= "<<fHydjetParams.fMuI3<<endl;
 545  cout<<" TauB= "<<fHydjetParams.fTau<<endl;
 546  cout<<" ThFrzTemperature= "<<fHydjetParams.fThFO<<endl;
 547  cout<<" MuPionThermal= "<<fHydjetParams.fMu_th_pip<<endl;
 548
 549 cout<<"-----Volume parameters -------------- "<<endl;
 550
 551  cout<<" SigmaTau= "<<fHydjetParams.fSigmaTau<<endl;
 552  cout<<" RmaxB= "<<fHydjetParams.fR<<endl;
 553  cout<<" YlMax= "<<fHydjetParams.fYlmax<<endl;
 554  cout<<" EtaRMax= "<<fHydjetParams.fUmax<<endl;
 555  cout<<" MomAsymmPar= "<<fHydjetParams.fDelta<<endl;
 556  cout<<" CoordAsymmPar= "<<fHydjetParams.fEpsilon<<endl;
 557
 558 cout<<" --------Flags------ "<<endl;
 559
 560  cout<<" GammaS= "<<fHydjetParams.fCorrS<<endl;
 561  cout<<" EtaType= "<<fHydjetParams.fEtaType<<endl;
 562  cout<<" FlagWeakDecay= "<<fHydjetParams.fWeakDecay<<endl;
 563
 564   cout<<"----PYQUEN parameters---"<<endl;
 565
 566   cout<<" Nhsel= "<<fHydjetParams.fNhsel<<endl;
 567   cout<<" Shad= "<<fHydjetParams.fIshad<<endl;
 568   cout<<" Ptmin= "<<fHydjetParams.fPtmin<<endl;
 569   cout<<" T0= "<<fHydjetParams.fT0<<endl;
 570   cout<<" Tau0= "<<fHydjetParams.fTau0<<endl;
 571   cout<<" Nf= "<<fHydjetParams.fNf<<endl;
 572   cout<<" Ienglu= "<<fHydjetParams.fIenglu<<endl;
 573   cout<<" Ianglu= "<<fHydjetParams.fIanglu<<endl;
 574
 575   //  cout<<"----PDG table parameters---"<<endl;
 576
 577   //  cout<<" UseCharmParticles= "<<fUseCharmParticles<<endl;
 578   //  cout<<" MinimumWidth= "<<fMinWidth<<endl;
 579   //  cout<<" MaximumWidth= "<<fMaxWidth<<endl;
 580   //  cout<<" MinimumMass= "<<fMinMass<<endl;
 581   //  cout<<" MaximumMass= "<<fMaxMass<<endl;
 582
 583   //  cout << "AliGenUHKM::SetAllParameters() OUT" << endl;
 584 }
 585
 586 // add the additional PDG codes from UHKM(SHARE table) to ROOT's table
 587 void AliGenUHKM::CheckPDGTable() {
 588   // Add temporarely all particle definitions from HYDJET++ which miss in the ROOT's PDG tables
 589   // to the TDatabasePDG table.
 590
 591   DatabasePDG *uhkmPDG = fUHKMgen->PDGInfo();              // UHKM's PDG table
 592   TParticlePDG *rootTestParticle;
 593   ParticlePDG *uhkmTestParticle;
 594
 595   // loop over all particles in the SHARE table
 596   for(Int_t i=0; i<uhkmPDG->GetNParticles(); i++) {
 597     // get a particle specie
 598     uhkmTestParticle = uhkmPDG->GetPDGParticleByIndex(i);
 599     // check if this code exists in ROOT's table
 600     rootTestParticle = TDatabasePDG::Instance()->GetParticle(uhkmTestParticle->GetPDG());
 601     if(!rootTestParticle) {    // if not then add it to the ROOT's PDG database
 602
 603       TDatabasePDG::Instance()->AddParticle(uhkmTestParticle->GetName(), uhkmTestParticle->GetName(),
 604                                             uhkmTestParticle->GetMass(), uhkmTestParticle->GetElectricCharge(),
 605                                             (uhkmTestParticle->GetWidth()<1e-10 ? kTRUE : kFALSE),
 606                                             uhkmTestParticle->GetWidth(),
 607                                             (Int_t(uhkmTestParticle->GetBaryonNumber())==0 ? "meson" : "baryon"),
 608                                             uhkmTestParticle->GetPDG());
 609       if(uhkmTestParticle->GetWidth()<1e-10)
 610         cout << uhkmTestParticle->GetPDG() << "its mass "
 611              << TDatabasePDG::Instance()->GetParticle(uhkmTestParticle->GetPDG())->Mass()
 612              << TDatabasePDG::Instance()->GetParticle(uhkmTestParticle->GetPDG())->Width() << endl;
 613     }
 614   }  // end for
 615 }