Fast simulation of pi/K -> mu X (H. Woehri, A. de Falco)
authormorsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 30 Apr 2007 09:39:21 +0000 (09:39 +0000)
committermorsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 30 Apr 2007 09:39:21 +0000 (09:39 +0000)
EVGEN/AliGenMUONCocktailpp.cxx
EVGEN/AliGenMUONCocktailpp.h

index b3aef2c..e96d1e4 100644 (file)
@@ -27,6 +27,7 @@
 #include <TObjArray.h>
 #include <TParticle.h>
 #include <TF1.h>
+#include <TVirtualMC.h>
 
 #include "AliGenCocktailEventHeader.h"
 
@@ -37,6 +38,7 @@
 #include "AliMC.h"
 #include "AliRun.h"
 #include "AliStack.h"
+#include "AliDecayer.h"
 #include "AliLog.h"
 #include "AliGenPythia.h"
 
@@ -46,11 +48,15 @@ ClassImp(AliGenMUONCocktailpp)
 //________________________________________________________________________
 AliGenMUONCocktailpp::AliGenMUONCocktailpp()
     :AliGenCocktail(),
+     fDecayer(0),
+     fDecayModeResonance(kAll),
+     fDecayModePythia(kAll),
      fTotalRate(0),
      fMuonMultiplicity(0),
      fMuonPtCut(0.),
      fMuonThetaMinCut(0.),
      fMuonThetaMaxCut(0.),
+     fMuonOriginCut(-999.),
      fNSucceded(0),
      fNGenerated(0)
 {
@@ -97,6 +103,8 @@ void AliGenMUONCocktailpp::Init()
     Double_t sigmaupsilonP = 0.502e-6;  
     Double_t sigmaupsilonPP = 0.228e-6;
     
+    AliInfo(Form("the parametrised resonances uses the decay mode %d",fDecayModeResonance));
+
 // Generation using CDF scaled distribution for pp@14TeV (from D. Stocco)
 //----------------------------------------------------------------------
 // Jpsi generator
@@ -105,7 +113,9 @@ void AliGenMUONCocktailpp::Init()
     genjpsi->SetPtRange(0.,100.);
     genjpsi->SetYRange(-8.,8.);
     genjpsi->SetPhiRange(0.,360.);
-    genjpsi->SetForceDecay(kAll);
+    genjpsi->SetForceDecay(fDecayModeResonance);
+    if (!gMC) genjpsi->SetDecayer(fDecayer);
+
     genjpsi->Init();  // generation in 4pi
     ratiojpsi = sigmajpsi / sigmaReaction * genjpsi->GetRelativeArea(ptMin,ptMax,yMin,yMax,phiMin,phiMax); // get weight
     AliInfo(Form("jpsi prod. cross-section in pp(14 TeV) %5.3g b",sigmajpsi));
@@ -125,7 +135,9 @@ void AliGenMUONCocktailpp::Init()
     genpsiP->SetPtRange(0.,100.);
     genpsiP->SetYRange(-8.,8.);
     genpsiP->SetPhiRange(0.,360.);
-    genpsiP->SetForceDecay(kAll);
+    genpsiP->SetForceDecay(fDecayModeResonance);
+    if (!gMC) genpsiP->SetDecayer(fDecayer);
+
     genpsiP->Init();  // generation in 4pi
     ratiopsiP = sigmapsiP / sigmaReaction * genpsiP->GetRelativeArea(ptMin,ptMax,yMin,yMax,phiMin,phiMax);
     AliInfo(Form("psiP prod. cross-section in pp(14 TeV) %5.3g b",sigmapsiP));
@@ -145,7 +157,8 @@ void AliGenMUONCocktailpp::Init()
     genupsilon->SetPtRange(0.,100.); 
     genupsilon->SetYRange(-8.,8.);
     genupsilon->SetPhiRange(0.,360.);
-    genupsilon->SetForceDecay(kAll);
+    genupsilon->SetForceDecay(fDecayModeResonance);
+    if (!gMC) genupsilon->SetDecayer(fDecayer);
     genupsilon->Init();  // generation in 4pi
     ratioupsilon = sigmaupsilon / sigmaReaction * genupsilon->GetRelativeArea(ptMin,ptMax,yMin,yMax,phiMin,phiMax);
     AliInfo(Form("Upsilon 1S prod. cross-section in pp(14 TeV) %5.3g b",sigmaupsilon));
@@ -165,7 +178,8 @@ void AliGenMUONCocktailpp::Init()
     genupsilonP->SetPtRange(0.,100.);
     genupsilonP->SetYRange(-8.,8.);
     genupsilonP->SetPhiRange(0.,360.);
-    genupsilonP->SetForceDecay(kAll);
+    genupsilonP->SetForceDecay(fDecayModeResonance);
+    if (!gMC) genupsilonP->SetDecayer(fDecayer);  
     genupsilonP->Init();  // generation in 4pi
     ratioupsilonP = sigmaupsilonP / sigmaReaction * genupsilonP->GetRelativeArea(ptMin,ptMax,yMin,yMax,phiMin,phiMax);
     AliInfo(Form("Upsilon 2S prod. cross-section in pp(14 TeV) %5.3g b",sigmaupsilonP));
@@ -185,7 +199,8 @@ void AliGenMUONCocktailpp::Init()
     genupsilonPP->SetPtRange(0.,100.); 
     genupsilonPP->SetYRange(-8.,8.);
     genupsilonPP->SetPhiRange(0.,360.);
-    genupsilonPP->SetForceDecay(kAll);
+    genupsilonPP->SetForceDecay(fDecayModeResonance);
+    if (!gMC) genupsilonPP->SetDecayer(fDecayer);  
     genupsilonPP->Init();  // generation in 4pi
     ratioupsilonPP = sigmaupsilonPP / sigmaReaction * genupsilonPP->GetRelativeArea(ptMin,ptMax,yMin,yMax,phiMin,phiMax);
     AliInfo(Form("Upsilon 3S prod. cross-section in pp(14 TeV) %5.3g b",sigmaupsilonPP));
@@ -197,14 +212,14 @@ void AliGenMUONCocktailpp::Init()
     genupsilonPP->Init(); // generation in selected kinematical range
     AddGenerator(genupsilonPP,"UpsilonPP", ratioupsilonPP); // Adding Generator
     fTotalRate+=ratioupsilonPP;
-    
 //------------------------------------------------------------------
 // Pythia generator
     AliGenPythia *pythia = new AliGenPythia(1);
     pythia->SetProcess(kPyMbMSEL1);
     pythia->SetStrucFunc(kCTEQ5L);
     pythia->SetEnergyCMS(14000.);
-    pythia->SetForceDecay(kAll);
+    AliInfo(Form("\n\npythia uses the decay mode %d",fDecayModePythia));
+    pythia->SetForceDecay(fDecayModePythia);
     pythia->SetPtRange(0.,100.);
     pythia->SetYRange(-8.,8.);
     pythia->SetPhiRange(0.,360.);
@@ -213,6 +228,13 @@ void AliGenMUONCocktailpp::Init()
     AddGenerator(pythia,"Pythia",1);
     fTotalRate+=1.;
 
+    TIter next(fEntries);
+    AliGenCocktailEntry *entry;
+    if (fStack) {
+       while((entry = (AliGenCocktailEntry*)next())) {
+           entry->Generator()->SetStack(fStack);
+       }
+    }
 }
 
 //_________________________________________________________________________
@@ -269,23 +291,25 @@ void AliGenMUONCocktailpp::Generate()
        Int_t iPart;
        fNGenerated++;
        Int_t numberOfMuons=0;  
-       
-       for(iPart=0; iPart<partArray->GetEntriesFast(); iPart++){      
+
+       Int_t maxPart = partArray->GetEntriesFast();
+       for(iPart=0; iPart<maxPart; iPart++){      
            
-           if ( (TMath::Abs(gAlice->GetMCApp()->Particle(iPart)->GetPdgCode())==13)  &&
-                (gAlice->GetMCApp()->Particle(iPart)->Theta()*180./TMath::Pi()>fMuonThetaMinCut) &&
-                (gAlice->GetMCApp()->Particle(iPart)->Theta()*180./TMath::Pi()<fMuonThetaMaxCut) &&
-                (gAlice->GetMCApp()->Particle(iPart)->Pt()>fMuonPtCut)                             ) { 
-               numberOfMuons++;
+         TParticle *part = gAlice->GetMCApp()->Particle(iPart);
+         if ( TMath::Abs(part->GetPdgCode()) == 13 ){  
+           if((part->Vz() > fMuonOriginCut) && //take only the muons that decayed before the abs + 1 int. length in C abs
+              (part->Theta()*180./TMath::Pi()>fMuonThetaMinCut) &&
+              (part->Theta()*180./TMath::Pi()<fMuonThetaMaxCut) &&
+              (part->Pt()>fMuonPtCut)) {
+             numberOfMuons++;
            }
-       }
-       
+         }
+       }       
        if (numberOfMuons >= fMuonMultiplicity) primordialTrigger = kTRUE;
     }
     fNSucceded++;
 
-    AliInfo(Form("Generated Events are %d and Succeeded Events are %d",fNGenerated,fNSucceded));
+//     AliInfo(Form("Generated Events are %d and Succeeded Events are %d",fNGenerated,fNSucceded));
     AliDebug(5,Form("Generated Events are %d and Succeeded Events are %d",fNGenerated,fNSucceded));
 }
 
index 02c269a..a46290e 100644 (file)
@@ -4,6 +4,7 @@
  * See cxx source for full Copyright notice                               */
 
 #include "AliGenCocktail.h"
+#include "AliDecayer.h"
 
 class AliGenCocktailEntry;
 
@@ -21,27 +22,37 @@ class AliGenMUONCocktailpp : public AliGenCocktail
     Float_t GetMuonPtCut()         const {return fMuonPtCut;}
     Float_t GetMuonThetaMin()      const {return fMuonThetaMinCut;}
     Float_t GetMuonThetaMax()      const {return fMuonThetaMaxCut;}        
+    Float_t GetMuonOriginCut()     const {return fMuonOriginCut;}          
+    Float_t GetDecayModeResonance()const {return fDecayModeResonance;}
+    Float_t GetDecayModePythia()   const {return fDecayModePythia;}
     
-    void    SetMuonMultiplicity(Int_t MuonMultiplicity) { fMuonMultiplicity= MuonMultiplicity;}
+    void    SetMuonMultiplicity(Int_t MuonMultiplicity) { fMuonMultiplicity = MuonMultiplicity;}
     void    SetMuonPtCut(Float_t PtCut) { fMuonPtCut = PtCut;}
+    void    SetMuonOriginCut(Float_t originCut) { fMuonOriginCut = originCut;}
     void    SetMuonThetaRange(Float_t ThetaMin, Float_t ThetaMax){
        fMuonThetaMinCut=ThetaMin;
        fMuonThetaMaxCut=ThetaMax; }    
-
+    void    SetDecayer(AliDecayer* decayer){fDecayer = decayer;}
+    void    SetDecayModeResonance(Decay_t decay){ fDecayModeResonance = decay;}
+    void    SetDecayModePythia(Decay_t decay){ fDecayModePythia = decay;}
  protected:
 
     //
  private:
     AliGenMUONCocktailpp(const AliGenMUONCocktailpp &cocktail); 
     AliGenMUONCocktailpp & operator=(const AliGenMUONCocktailpp &cocktail); 
-
+    AliDecayer* fDecayer;
+    Decay_t fDecayModeResonance; //decay mode in which resonances are forced to decay, default: kAll
+    Decay_t fDecayModePythia; //decay mode in which particles in Pythia are forced to decay, default: kAll
     Float_t fTotalRate;// Total rate of the full cocktail processes
     Int_t   fMuonMultiplicity; // Muon multiplicity for the primordial trigger
     Float_t fMuonPtCut;// Transverse momentum cut for muons
     Float_t fMuonThetaMinCut;// Minimum theta cut for muons
     Float_t fMuonThetaMaxCut; // Maximum theta cut for muons
+    Float_t fMuonOriginCut; //use only muons whose "part->Vz()" value is larger than fMuonOrigin
     Int_t   fNSucceded;// Number of Succes in the (di)-muon generation in the acceptance
     Int_t   fNGenerated;// Number of generated cocktails
+    
 
     ClassDef(AliGenMUONCocktailpp,1)  //  cocktail for physics in the Alice
 };