Updates (N. Bastid)
authormorsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 31 Jul 2007 06:10:32 +0000 (06:10 +0000)
committermorsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 31 Jul 2007 06:10:32 +0000 (06:10 +0000)
EVGEN/AliGenCorrHF.cxx
EVGEN/AliGenMUONCocktailpp.cxx

index 5dfdcee..41e7c3c 100644 (file)
 // and quark fragmentation functions.
 // Is a generalisation of AliGenParam class for correlated pairs of hadrons.
 // In this version quark pairs and fragmentation functions are obtained from
-// Pythia6.124 using 100K events generated with kCharmppMNRwmi & kBeautyppMNRwmi 
+// Pythia6.124 using 100K events generated with kCharmppMNRwmi&kBeautyppMNRwmi 
 // in pp collisions at 14 TeV.
 // Decays are performed by Pythia. Used AliRoot version: v4-04-Release
 // Author: S. Grigoryan, LPC Clermont-Fd & YerPhI, Smbat.Grigoryan@cern.ch
-//
+// July 07: added quarks in the stack (B. Vulpescu)
 //-------------------------------------------------------------------------
 // How it works (for the given flavor):
 //
@@ -93,6 +93,7 @@
 #include <TRandom.h>
 #include <TTree.h>
 #include <TVirtualMC.h>
+#include <TVector3.h>
 
 #include "AliGenCorrHF.h"
 #include "AliLog.h"
@@ -206,7 +207,7 @@ void AliGenCorrHF::Init()
     }
 
     ComputeIntegral(fFile);
-
+    
     fParentWeight = 1./fNpart;   // fNpart is number of HF-hadron pairs
 
 // particle decay related initialization
@@ -239,14 +240,17 @@ void AliGenCorrHF::Generate()
   Float_t phi, theta;         // Phi and theta spherical angles of the parent particle momentum
   Float_t p[3], pc[3], och[3];// Momentum, polarisation and origin of the children particles from lujet
 
+  Float_t pq1[3], pq2[3];     // Momenta of the two quarks
 
-  Double_t dphi=0, ptq[2], yq[2], pth[2], plh[2], ph[2], phih[2];
+  Double_t dphi=0, ptq[2], yq[2], pth[2], plh[2], ph[2], phih[2], phiq[2];
   Int_t i, j, ipair, ihadron[2];
   for (i=0; i<2; i++) { 
-    ptq[i] =0; 
-    yq[i]  =0; 
-    pth[i] =0; 
-    plh[i] =0; 
+    ptq[i]     =0; 
+    yq[i]      =0; 
+    pth[i]     =0; 
+    plh[i]     =0;
+    phih[i]    =0; 
+    phiq[i]    =0;
     ihadron[i] =0; 
   }
 
@@ -266,6 +270,8 @@ void AliGenCorrHF::Generate()
   Float_t wgtp, wgtch, random[6];
   Int_t ipap = 0;
   Int_t nt   = 0;
+  Int_t ntq1 = 0;
+  Int_t ntq2 = 0;
 
 // Generating fNpart HF-hadron pairs per event
   while(ipap<fNpart) {
@@ -274,6 +280,35 @@ void AliGenCorrHF::Generate()
 
       GetQuarkPair(fFile, fgIntegral, yq[0], yq[1], ptq[0], ptq[1], dphi);
 
+// Add the quarks in the stack
+
+      phiq[0] = Rndm()*k2PI;
+      phiq[1] = phiq[0] + dphi*kDegrad; 
+      TVector3 qvect1 = TVector3();
+      TVector3 qvect2 = TVector3();
+      qvect1.SetPtEtaPhi(ptq[0],yq[0],phiq[0]);
+      qvect2.SetPtEtaPhi(ptq[1],yq[1],phiq[1]);
+      pq1[0] = qvect1.Px();
+      pq1[1] = qvect1.Py();
+      pq1[2] = qvect1.Pz();
+      pq2[0] = qvect2.Px();
+      pq2[1] = qvect2.Py();
+      pq2[2] = qvect2.Pz();
+      
+      wgtp  = fParentWeight;
+
+      PushTrack(0, -1, +fQuark, pq1, origin0, polar, 0, 
+               kPPrimary, nt, wgtp);
+      KeepTrack(nt);
+      ntq1 = nt;
+
+      PushTrack(0, -1, -fQuark, pq2, origin0, polar, 0, 
+               kPPrimary, nt, wgtp);
+      KeepTrack(nt);
+      ntq2 = nt;
+
+      nt = 2;
+      
       GetHadronPair(fFile, fQuark, yq[0], yq[1], ptq[0], ptq[1], ihadron[0], ihadron[1], plh[0], plh[1], pth[0], pth[1]);
 
 // Here we assume that  |phi_H1 - phi_H2| = |phi_Q1 - phi_Q2| = dphi
@@ -430,7 +465,7 @@ void AliGenCorrHF::Generate()
                      p1[0] = p[0]; p1[1] = p[1]; p1[2] = p[2];
                  } else {
                      ipap++;
-                     PushTrack(0, -1, ihadron[0], p1, origin0, polar, 0, 
+                     PushTrack(0, ntq1, ihadron[0], p1, origin0, polar, 0, 
                                kPPrimary, nt, wgtp);
                      KeepTrack(nt);
                      for (i = 1; i < np1; i++) {
@@ -444,7 +479,7 @@ void AliGenCorrHF::Generate()
                          KeepTrack(nt);
                          }
                      }
-                     PushTrack(0, -1, iPart, p, origin0, polar, 0, kPPrimary, nt, wgtp);
+                     PushTrack(0, ntq2, iPart, p, origin0, polar, 0, kPPrimary, nt, wgtp);
                      KeepTrack(nt);
                  } 
                  pParent[0] = nt;
@@ -502,9 +537,9 @@ void AliGenCorrHF::Generate()
            } else {
                ipap++;
                gAlice->GetMCApp()->
-               PushTrack(fTrackIt,-1,ihadron[0],p1,origin0,polar,0,kPPrimary,nt,wgtp);
+               PushTrack(fTrackIt,ntq1,ihadron[0],p1,origin0,polar,0,kPPrimary,nt,wgtp);
                gAlice->GetMCApp()->
-               PushTrack(fTrackIt,-1,iPart,p,origin0,polar,0,kPPrimary,nt,wgtp);
+               PushTrack(fTrackIt,ntq2,iPart,p,origin0,polar,0,kPPrimary,nt,wgtp);
            }
          }
          if (nhadron == 0) break;
index e96d1e4..b8f9553 100644 (file)
@@ -22,7 +22,8 @@
 // The free parameeters are :
 //      pp reaction cross-section
 //      production cross-sections in pp collisions 
-// 
+// July 07:added heavy quark production from AliGenCorrHF and heavy quark 
+//         production switched off in Pythia 
 
 #include <TObjArray.h>
 #include <TParticle.h>
@@ -41,7 +42,7 @@
 #include "AliDecayer.h"
 #include "AliLog.h"
 #include "AliGenPythia.h"
-
+#include "AliGenCorrHF.h"
 
 ClassImp(AliGenMUONCocktailpp)  
   
@@ -93,7 +94,9 @@ void AliGenMUONCocktailpp::Init()
     Double_t ratioupsilon;
     Double_t ratioupsilonP;
     Double_t ratioupsilonPP;
-    
+    Double_t ratioccbar;
+    Double_t ratiobbbar;
+
 // Cross sections in barns (from PPR Vol. II p: 552) pp - 14 TeV and 
 // corrected from feed down of higher resonances 
 
@@ -102,6 +105,8 @@ void AliGenMUONCocktailpp::Init()
     Double_t sigmaupsilon = 0.989e-6;  
     Double_t sigmaupsilonP = 0.502e-6;  
     Double_t sigmaupsilonPP = 0.228e-6;
+    Double_t sigmaccbar = 11.2e-3;
+    Double_t sigmabbbar = 0.51e-3;
     
     AliInfo(Form("the parametrised resonances uses the decay mode %d",fDecayModeResonance));
 
@@ -212,7 +217,30 @@ void AliGenMUONCocktailpp::Init()
     genupsilonPP->Init(); // generation in selected kinematical range
     AddGenerator(genupsilonPP,"UpsilonPP", ratioupsilonPP); // Adding Generator
     fTotalRate+=ratioupsilonPP;
+
+//------------------------------------------------------------------
+// Generator of charm
+    
+    AliGenCorrHF *gencharm = new AliGenCorrHF(1, 4);  
+          gencharm->SetMomentumRange(0,9999);
+          gencharm->SetForceDecay(kAll);
+         ratioccbar = sigmaccbar/sigmaReaction;
+          gencharm->Init();
+          AddGenerator(gencharm,"CorrHFCharm",ratioccbar);
+         fTotalRate+=ratioccbar;
+
 //------------------------------------------------------------------
+// Generator of beauty
+
+         AliGenCorrHF *genbeauty = new AliGenCorrHF(1, 5);  
+          genbeauty->SetMomentumRange(0,9999);
+          genbeauty->SetForceDecay(kAll);
+         ratiobbbar = sigmabbbar/sigmaReaction;
+         genbeauty->Init();
+          AddGenerator(genbeauty,"CorrHFBeauty",ratiobbbar); 
+         fTotalRate+=ratiobbbar;
+
+//--------------------------t----------------------------------------
 // Pythia generator
     AliGenPythia *pythia = new AliGenPythia(1);
     pythia->SetProcess(kPyMbMSEL1);
@@ -224,6 +252,7 @@ void AliGenMUONCocktailpp::Init()
     pythia->SetYRange(-8.,8.);
     pythia->SetPhiRange(0.,360.);
     pythia->SetPtHard(2.76,-1.0);
+    pythia->SwitchHFOff();
     pythia->Init(); 
     AddGenerator(pythia,"Pythia",1);
     fTotalRate+=1.;
@@ -290,9 +319,7 @@ void AliGenMUONCocktailpp::Generate()
 // in the muon spectrometer acceptance
        Int_t iPart;
        fNGenerated++;
-       Int_t numberOfMuons=0;  
-
-       Int_t maxPart = partArray->GetEntriesFast();
+       Int_t numberOfMuons=0;Int_t maxPart = partArray->GetEntriesFast();
        for(iPart=0; iPart<maxPart; iPart++){      
            
          TParticle *part = gAlice->GetMCApp()->Particle(iPart);