]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - THijing/AliGenHijing.cxx
Draw charged tracks as straight lines when magnetic field is near zero.
[u/mrichter/AliRoot.git] / THijing / AliGenHijing.cxx
index c0e6c152273b7791711117c8c171068fee74aeda..41ccc49babaf07dee2bc87c450b9aa46a7738f5d 100644 (file)
 // Generator using HIJING as an external generator
 // The main HIJING options are accessable for the user through this interface.
 // Uses the THijing implementation of TGenerator.
+// Author:
+// Andreas Morsch    (andreas.morsch@cern.ch)
 //
-// andreas.morsch@cern.ch
 
-#include <TArrayI.h>
 #include <TGraph.h>
 #include <THijing.h>
 #include <TLorentzVector.h>
@@ -73,10 +73,11 @@ AliGenHijing::AliGenHijing(Int_t npart)
     fDnDb       =  0;
     fPtMinJet   = -2.5;        
     fRadiation  =  3;
-    fEventVertex.Set(3);
-//
+    //
     SetSimpleJets();
     SetNoGammas();
+    SetRandomPz();
+    SwitchOffHeavyQuarks(kFALSE);
 //
     fParticles = new TClonesArray("TParticle",10000);    
 //
@@ -86,7 +87,8 @@ AliGenHijing::AliGenHijing(Int_t npart)
 
 }
 
-AliGenHijing::AliGenHijing(const AliGenHijing & Hijing)
+AliGenHijing::AliGenHijing(const AliGenHijing & hijing):
+    AliGenMC(hijing)
 {
 // copy constructor
 }
@@ -111,7 +113,7 @@ void AliGenHijing::Init()
                      fAProjectile, fZProjectile, fATarget, fZTarget, 
                      fMinImpactParam, fMaxImpactParam));
 
-    fHijing=(THijing*) fgMCEvGen;
+    fHijing=(THijing*) fMCEvGen;
     fHijing->SetIHPR2(2,  fRadiation);
     fHijing->SetIHPR2(3,  fTrigger);
     fHijing->SetIHPR2(6,  fShadowing);
@@ -151,6 +153,15 @@ void AliGenHijing::Init()
        fHijing->SetHIPR1(11, 2.5);
     }
     
+//
+// Heavy quarks
+//    
+    if (fNoHeavyQuarks) {
+       fHijing->SetIHPR2(49, 1);
+    } else {
+       fHijing->SetIHPR2(49, 0);
+    }
+    
     
     AliGenMC::Init();
     
@@ -170,7 +181,7 @@ void AliGenHijing::Generate()
   Float_t polar[3]    =   {0,0,0};
   Float_t origin[3]   =   {0,0,0};
   Float_t origin0[3]  =   {0,0,0};
-  Float_t p[3], random[6];
+  Float_t p[3];
   Float_t tof;
 
 //  converts from mm/c to s
@@ -178,7 +189,7 @@ void AliGenHijing::Generate()
 //
   Int_t nt  = 0;
   Int_t jev = 0;
-  Int_t j, kf, ks, imo;
+  Int_t j, kf, ks, ksp, imo;
   kf = 0;
     
 
@@ -186,28 +197,20 @@ void AliGenHijing::Generate()
   fTrials = 0;
   for (j = 0;j < 3; j++) origin0[j] = fOrigin[j];
   if(fVertexSmear == kPerEvent) {
-      Float_t dv[3];
-      dv[2] = 1.e10;
-      while(TMath::Abs(dv[2]) > fCutVertexZ*fOsigma[2]) {
-         Rndm(random,6);
-         for (j=0; j < 3; j++) {
-             dv[j] = fOsigma[j]*TMath::Cos(2*random[2*j]*TMath::Pi())*
-                 TMath::Sqrt(-2*TMath::Log(random[2*j+1]));
-         }
-      }
-      for (j=0; j < 3; j++) origin0[j] += dv[j];
-  } else if (fVertexSmear == kPerTrack) {
-//         fHijing->SetMSTP(151,0);
-      for (j = 0; j < 3; j++) {
-//           fHijing->SetPARP(151+j, fOsigma[j]*10.);
-      }
-  }
+      Vertex();
+      for (j=0; j < 3; j++) origin0[j] = fVertex[j];
+  } 
+
+
+  Float_t sign = (fRandomPz && (Rndm() < 0.5))? -1. : 1.;
   while(1)
   {
 //    Generate one event
 // --------------------------------------------------------------------------
-      fSpecn   = 0;  
-      fSpecp   = 0;
+      fProjectileSpecn    = 0;  
+      fProjectileSpecp    = 0;
+      fTargetSpecn        = 0;  
+      fTargetSpecp        = 0;
 // --------------------------------------------------------------------------
       fHijing->GenerateEvent();
       fTrials++;
@@ -234,9 +237,9 @@ void AliGenHijing::Generate()
 //      Get event vertex
 //
       TParticle *  iparticle = (TParticle *) fParticles->At(0);
-      fEventVertex[0] = origin0[0];
-      fEventVertex[1] = origin0[1];    
-      fEventVertex[2] = origin0[2];
+      fVertex[0] = origin0[0];
+      fVertex[1] = origin0[1]; 
+      fVertex[2] = origin0[2];
       
 //
 //      First select parent particles
@@ -280,19 +283,24 @@ void AliGenHijing::Generate()
          Bool_t  selected             =  kTRUE;
          kf        = iparticle->GetPdgCode();
          ks        = iparticle->GetStatusCode();
+         ksp       = iparticle->GetUniqueID();
          
 // --------------------------------------------------------------------------
 // Count spectator neutrons and protons
-         if(ks == 0 || ks == 1 || ks == 10 || ks == 11){
-             if(kf == kNeutron) fSpecn += 1;
-             if(kf == kProton)  fSpecp += 1;
+         if(ksp == 0 || ksp == 1){
+             if(kf == kNeutron) fProjectileSpecn += 1;
+             if(kf == kProton)  fProjectileSpecp += 1;
+         }
+         else if(ksp == 10 || ksp == 11){
+             if(kf == kNeutron) fTargetSpecn += 1;
+             if(kf == kProton)  fTargetSpecp += 1;
          }
 // --------------------------------------------------------------------------
 //         
          if (!fSelectAll) {
              selected = KinematicSelection(iparticle,0)&&SelectFlavor(kf);
-             if (!fSpectators && selected) selected = (ks != 0 && ks != 1 && ks != 10
-                                                       && ks != 11);
+             if (!fSpectators && selected) selected = (ksp != 0 && ksp != 1 && ksp != 10
+                                                       && ksp != 11);
          }
 //
 // Put particle on the stack if selected
@@ -302,34 +310,42 @@ void AliGenHijing::Generate()
              pSelected[i] = 1;
          } // selected
       } // particle loop final state
+
 //
-// Write particles to stack
+//    Time of the interactions
+      Float_t tInt = 0.;
+      if (fPileUpTimeWindow > 0.) tInt = fPileUpTimeWindow * (2. * gRandom->Rndm() - 1.);
+
 //
+// Write particles to stack
+
       for (i = 0; i<np; i++) {
          TParticle *  iparticle = (TParticle *) fParticles->At(i);
          Bool_t  hasMother   = (iparticle->GetFirstMother()     >=0);
          Bool_t  hasDaughter = (iparticle->GetFirstDaughter()   >=0);
-
          if (pSelected[i]) {
              kf   = iparticle->GetPdgCode();
              ks   = iparticle->GetStatusCode();
              p[0] = iparticle->Px();
              p[1] = iparticle->Py();
-             p[2] = iparticle->Pz();
+             p[2] = iparticle->Pz() * sign;
              origin[0] = origin0[0]+iparticle->Vx()/10;
              origin[1] = origin0[1]+iparticle->Vy()/10;
              origin[2] = origin0[2]+iparticle->Vz()/10;
-             tof = kconv*iparticle->T();
+             tof = kconv * iparticle->T() + sign * origin0[2] / 3.e10;
+             if (fPileUpTimeWindow > 0.) tof += tInt;
+             
              imo = -1;
              TParticle* mother = 0;
              if (hasMother) {
                  imo = iparticle->GetFirstMother();
                  mother = (TParticle *) fParticles->At(imo);
-                 imo = (mother->GetPdgCode() != 92) ? imo = newPos[imo] : -1;
+                 imo = (mother->GetPdgCode() != 92) ? newPos[imo] : -1;
              } // if has mother   
              Bool_t tFlag = (fTrackIt && !hasDaughter);
-             SetTrack(tFlag,imo,kf,p,origin,polar,
-                      tof,kPNoProcess,nt, 1., ks);
+             PushTrack(tFlag,imo,kf,p,origin,polar,tof,kPNoProcess,nt, 1., ks);
+
+             
              KeepTrack(nt);
              newPos[i] = nt;
          } // if selected
@@ -474,7 +490,7 @@ Bool_t AliGenHijing::SelectFlavor(Int_t pid)
     return res;
 }
 
-Bool_t AliGenHijing::Stable(TParticle*  particle)
+Bool_t AliGenHijing::Stable(TParticle*  particle) const
 {
 // Return true for a stable particle
 //
@@ -502,7 +518,11 @@ void AliGenHijing::MakeHeader()
                                                       fHijing->GetN01(),
                                                       fHijing->GetN10(),
                                                       fHijing->GetN11());
-    ((AliGenHijingEventHeader*) header)->SetSpectators(fSpecn, fSpecp);
+    ((AliGenHijingEventHeader*) header)->SetSpectators(fProjectileSpecn, fProjectileSpecp,
+                                                      fTargetSpecn,fTargetSpecp);
+    ((AliGenHijingEventHeader*) header)->SetReactionPlaneAngle(fHijing->GetHINT1(20));
+    
+
 
 // 4-momentum vectors of the triggered jets.
 //
@@ -530,11 +550,22 @@ void AliGenHijing::MakeHeader()
 // Bookkeeping for kinematic bias
     ((AliGenHijingEventHeader*) header)->SetTrials(fTrials);
 // Event Vertex
-    header->SetPrimaryVertex(fEventVertex);
-    gAlice->SetGenEventHeader(header);   
+    header->SetPrimaryVertex(fVertex);
+    AddHeader(header);
     fCollisionGeometry = (AliGenHijingEventHeader*)  header;
 }
 
+void AliGenHijing::AddHeader(AliGenEventHeader* header)
+{
+    // Passes header either to the container or to gAlice
+    if (fContainer) {
+       fContainer->AddHeader(header);
+    } else {
+       gAlice->SetGenEventHeader(header);      
+    }
+}
+
+
 Bool_t AliGenHijing::CheckTrigger()
 {
 // Check the kinematic trigger condition
@@ -575,8 +606,8 @@ Bool_t AliGenHijing::CheckTrigger()
        for (Int_t i = 0; i < np; i++) {
            TParticle* part = (TParticle*) fParticles->At(i);
            Int_t kf = part->GetPdgCode();
-           Int_t ks = part->GetStatusCode();
-           if (kf == 22 && ks == 40) {
+           Int_t ksp = part->GetUniqueID();
+           if (kf == 22 && ksp == 40) {
                Float_t phi = part->Phi();
                Float_t eta = part->Eta();
                if  (eta < fEtaMaxJet && 
@@ -593,11 +624,14 @@ Bool_t AliGenHijing::CheckTrigger()
 }
 
 
-
+void AliGenHijing::Copy(TObject &) const
+{
+  Fatal("Copy","Not implemented!\n");
+}
 
 AliGenHijing& AliGenHijing::operator=(const  AliGenHijing& rhs)
 {
-// Assignment operator
-    return *this;
+    rhs.Copy(*this); 
+    return (*this);
 }