]> git.uio.no Git - u/mrichter/AliRoot.git/commitdiff
reject events with no tracks, add track multiplicity counter to the MCReader, add...
authorgconesab <gconesab@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 12 Nov 2010 09:57:25 +0000 (09:57 +0000)
committergconesab <gconesab@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 12 Nov 2010 09:57:25 +0000 (09:57 +0000)
PWG4/PartCorrBase/AliAnaPartCorrMaker.cxx
PWG4/PartCorrBase/AliAnaPartCorrMaker.h
PWG4/PartCorrBase/AliCaloTrackMCReader.cxx
PWG4/PartCorrBase/AliCaloTrackReader.cxx

index ddaa58a858cfcbe512fdee6e2bada55d28de6022..1d0da840d43164c293cc3f4bd64b0fc6e5b6c0a3 100755 (executable)
@@ -47,7 +47,7 @@ TObject(),
 fOutputContainer(new TList ), fAnalysisContainer(new TList ),
 fMakeHisto(kFALSE), fMakeAOD(kFALSE), fAnaDebug(0), 
 fReader(0), fCaloUtils(0), 
-fCuts(new TList), fhNEvents(0x0)
+fCuts(new TList), fhNEvents(0x0), fhTrackMult(0x0)
 {
   //Default Ctor
   if(fAnaDebug > 1 ) printf("*** Analysis Maker  Constructor *** \n");
@@ -64,7 +64,7 @@ fMakeHisto(maker.fMakeHisto), fMakeAOD(maker.fMakeAOD),
 fAnaDebug(maker.fAnaDebug),
 fReader(),//new AliCaloTrackReader(*maker.fReader)), 
 fCaloUtils(),//(new AliCalorimeterUtils(*maker.fCaloUtils)),
-fCuts(new TList()), fhNEvents(maker.fhNEvents)
+fCuts(new TList()),fhNEvents(maker.fhNEvents) ,fhTrackMult(maker.fhTrackMult)
 {
   // cpy ctor
        
@@ -195,9 +195,14 @@ TList *AliAnaPartCorrMaker::GetOutputContainer()
       }// Analysis with histograms as output on
     }//Loop on analysis defined
   }//Analysis list available
+  
+  //General event histograms
+  
   fhNEvents        = new TH1I("hNEvents", "Number of analyzed events"   , 1 , 0 , 1  ) ;
   fOutputContainer->Add(fhNEvents);
-       
+  fhTrackMult      = new TH1I("hTrackMult", "Number of tracks per events"   , 2000 , 0 , 2000  ) ;
+  fOutputContainer->Add(fhTrackMult);
+  
   return fOutputContainer;
   
 }
@@ -325,7 +330,8 @@ void AliAnaPartCorrMaker::ProcessEvent(const Int_t iEntry, const char * currentF
   }
        
   fhNEvents->Fill(0); //Event analyzed
-       
+  fhTrackMult->Fill(fReader->GetTrackMultiplicity()); 
+
   fReader->ResetLists();
   
   //printf(">>>>>>>>>> AFTER >>>>>>>>>>>\n");
index 07f085bc229a48f67cc9a7574a761658576e9d9c..461d8ac8765884223cbf301b5ba503c31a718a70 100755 (executable)
@@ -90,9 +90,10 @@ class AliAnaPartCorrMaker : public TObject {
 
   TList * fCuts ;                 //! List with analysis cuts
 
-  TH1I  * fhNEvents;           //! Number of events counter histogram
-       
-  ClassDef(AliAnaPartCorrMaker,7)
+  TH1I  * fhNEvents;         //! Number of events counter histogram
+  TH1I  * fhTrackMult;       //! Number of tracks per event histogram
+
+  ClassDef(AliAnaPartCorrMaker,8)
 } ;
  
 
index d44df993bc1b1be6af9c7c8aef184f45b1f777a1..aa59ae822bc2ad61aeaf0f31eb3f5ebb4aec056d 100755 (executable)
@@ -213,10 +213,10 @@ Bool_t AliCaloTrackMCReader::FillInputEvent(const Int_t iEntry, const char * cur
   
   fEventNumber = iEntry;
   fCurrentFileName = TString(currentFileName);
-       
+  fTrackMult = 0;
   //In case of analysis of events with jets, skip those with jet pt > 5 pt hard        
   if(fComparePtHardAndJetPt && GetStack()) {
-       if(!ComparePtHardAndJetPt()) return kFALSE ;
+    if(!ComparePtHardAndJetPt()) return kFALSE ;
   }
        
   //Fill Vertex array
@@ -242,66 +242,69 @@ Bool_t AliCaloTrackMCReader::FillInputEvent(const Int_t iEntry, const char * cur
       particle->Momentum(momentum);
       //---------- Charged particles ----------------------
       if(charge != 0){
-                if(fFillCTS && (momentum.Pt() > fCTSPtMin)){
-         //Particles in CTS acceptance
-               
-         if(fCheckFidCut && !fFiducialCut->IsInFiducialCut(momentum,"CTS")) continue;
-                
-    if(TMath::Abs(pdg) == 11 && GetStack()->Particle(particle->GetFirstMother())->GetPdgCode()==22) continue ;
-       
-         if(fDebug > 3 && momentum.Pt() > 0.2)
-           printf("AliCaloTrackMCReader::FillInputEvent() - CTS : Selected tracks E %3.2f, pt %3.2f, phi %3.2f, eta %3.2f\n",
-                  momentum.E(),momentum.Pt(),momentum.Phi()*TMath::RadToDeg(),momentum.Eta());
-         
-         x[0] = particle->Vx(); x[1] = particle->Vy(); x[2] = particle->Vz();
-         p[0] = particle->Px(); p[1] = particle->Py(); p[2] = particle->Pz();
-         //Create object and write it to file
-         AliAODTrack *aodTrack = new AliAODTrack(0, iParticle, p, kTRUE, x, kFALSE,NULL, 0, 0, 
-                       NULL,
-                       0x0,//primary,
-                       kFALSE, // No fit performed
-                       kFALSE, // No fit performed
-                       AliAODTrack::kPrimary, 
-                       0);
-         SetTrackChargeAndPID(pdg, aodTrack);
-         fAODCTS->Add(aodTrack);//reference the selected object to the list
-       }
-       //Keep some charged particles in calorimeters lists
-       if((fFillPHOS || fFillEMCAL) && KeepChargedParticles(pdg)) FillCalorimeters(iParticle, particle, momentum);
-       
+        
+        if(TMath::Abs(momentum.Eta())< fTrackMultEtaCut) fTrackMult++;
+        
+        if(fFillCTS && (momentum.Pt() > fCTSPtMin)){
+          //Particles in CTS acceptance
+          
+          if(fCheckFidCut && !fFiducialCut->IsInFiducialCut(momentum,"CTS")) continue;
+          
+          if(TMath::Abs(pdg) == 11 && GetStack()->Particle(particle->GetFirstMother())->GetPdgCode()==22) continue ;
+          
+          if(fDebug > 3 && momentum.Pt() > 0.2)
+            printf("AliCaloTrackMCReader::FillInputEvent() - CTS : Selected tracks E %3.2f, pt %3.2f, phi %3.2f, eta %3.2f\n",
+                   momentum.E(),momentum.Pt(),momentum.Phi()*TMath::RadToDeg(),momentum.Eta());
+          
+          x[0] = particle->Vx(); x[1] = particle->Vy(); x[2] = particle->Vz();
+          p[0] = particle->Px(); p[1] = particle->Py(); p[2] = particle->Pz();
+          //Create object and write it to file
+          AliAODTrack *aodTrack = new AliAODTrack(0, iParticle, p, kTRUE, x, kFALSE,NULL, 0, 0, 
+                                                  NULL,
+                                                  0x0,//primary,
+                                                  kFALSE, // No fit performed
+                                                  kFALSE, // No fit performed
+                                                  AliAODTrack::kPrimary, 
+                                                  0);
+          SetTrackChargeAndPID(pdg, aodTrack);
+          fAODCTS->Add(aodTrack);//reference the selected object to the list
+        }
+        //Keep some charged particles in calorimeters lists
+        if((fFillPHOS || fFillEMCAL) && KeepChargedParticles(pdg)) FillCalorimeters(iParticle, particle, momentum);
+        
       }//Charged
       
       //-------------Neutral particles ----------------------
       else if(charge == 0 && (fFillPHOS || fFillEMCAL)){
-       //Skip neutrinos or other neutral particles
-       //if(SkipNeutralParticles(pdg) || particle->IsPrimary()) continue ; // protection added (MG)
-       if(SkipNeutralParticles(pdg)) continue ;
-       //Fill particle/calocluster arrays
-       if(!fDecayPi0) {
-         FillCalorimeters(iParticle, particle, momentum);
-       }
-       else {
-         //Sometimes pi0 are stable for the generator, if needed decay it by hand
-         if(pdg == 111 ){
-           if(momentum.Pt() >  fPHOSPtMin || momentum.Pt() >  fEMCALPtMin){
-             TLorentzVector lvGamma1, lvGamma2 ;
-             //Double_t angle = 0;
-             
-             //Decay
-             MakePi0Decay(momentum,lvGamma1,lvGamma2);//,angle);
-             
-             //Check if Pi0 is in the acceptance of the calorimeters, if aperture angle is small, keep it
-             TParticle * pPhoton1 = new TParticle(22,1,iParticle,0,0,0,lvGamma1.Px(),lvGamma1.Py(),
-                                                  lvGamma1.Pz(),lvGamma1.E(),0,0,0,0);   
-             TParticle * pPhoton2 = new TParticle(22,1,iParticle,0,0,0,lvGamma2.Px(),lvGamma2.Py(),
-                                                  lvGamma2.Pz(),lvGamma2.E(),0,0,0,0);
-             //Fill particle/calocluster arrays
-             FillCalorimeters(iParticle,pPhoton1,lvGamma1);
-             FillCalorimeters(iParticle,pPhoton2,lvGamma2);
-           }//pt cut
-         }//pi0
-         else FillCalorimeters(iParticle,particle, momentum); //Add the rest
-       }
+        //Skip neutrinos or other neutral particles
+        //if(SkipNeutralParticles(pdg) || particle->IsPrimary()) continue ; // protection added (MG)
+        if(SkipNeutralParticles(pdg)) continue ;
+        //Fill particle/calocluster arrays
+        if(!fDecayPi0) {
+          FillCalorimeters(iParticle, particle, momentum);
+        }
+        else {
+          //Sometimes pi0 are stable for the generator, if needed decay it by hand
+          if(pdg == 111 ){
+            if(momentum.Pt() >  fPHOSPtMin || momentum.Pt() >  fEMCALPtMin){
+              TLorentzVector lvGamma1, lvGamma2 ;
+              //Double_t angle = 0;
+              
+              //Decay
+              MakePi0Decay(momentum,lvGamma1,lvGamma2);//,angle);
+              
+              //Check if Pi0 is in the acceptance of the calorimeters, if aperture angle is small, keep it
+              TParticle * pPhoton1 = new TParticle(22,1,iParticle,0,0,0,lvGamma1.Px(),lvGamma1.Py(),
+                                                   lvGamma1.Pz(),lvGamma1.E(),0,0,0,0);   
+              TParticle * pPhoton2 = new TParticle(22,1,iParticle,0,0,0,lvGamma2.Px(),lvGamma2.Py(),
+                                                   lvGamma2.Pz(),lvGamma2.E(),0,0,0,0);
+              //Fill particle/calocluster arrays
+              FillCalorimeters(iParticle,pPhoton1,lvGamma1);
+              FillCalorimeters(iParticle,pPhoton2,lvGamma2);
+            }//pt cut
+          }//pi0
+          else FillCalorimeters(iParticle,particle, momentum); //Add the rest
+        }
       }//neutral particles
     } //particle with correct status
   }//particle loop
@@ -309,7 +312,7 @@ Bool_t AliCaloTrackMCReader::FillInputEvent(const Int_t iEntry, const char * cur
   fIndex2ndPhoton = -1; //In case of overlapping studies, reset for each event 
        
   return kTRUE;        
-
+  
 }
 
 //________________________________________________________________
index c2a474946eb120548bea3308002db43c222fd55e..965c66c8e22c63eca198d88af5384cf0bad82e58 100755 (executable)
@@ -63,7 +63,7 @@ ClassImp(AliCaloTrackReader)
 //    fSecondInputFileName(""),fSecondInputFirstEvent(0), 
 //    fAODCTSNormalInputEntries(0), fAODEMCALNormalInputEntries(0), 
 //    fAODPHOSNormalInputEntries(0), 
-    fTrackStatus(0),   fESDtrackCuts(0), fTrackMult(0), fTrackMultEtaCut(0.9),
+    fTrackStatus(0),   fESDtrackCuts(0), fTrackMult(0), fTrackMultEtaCut(0.8),
     fReadStack(kFALSE), fReadAODMCParticles(kFALSE), 
     fDeltaAODFileName("deltaAODPartCorr.root"),fFiredTriggerClassName(""),
     fAnaLED(kFALSE),fTaskName(""),fCaloUtils(0x0), 
@@ -143,34 +143,34 @@ AliCaloTrackReader::~AliCaloTrackReader() {
 
 //_________________________________________________________________________
 Bool_t AliCaloTrackReader::ComparePtHardAndJetPt(){
-       // Check the event, if the requested ptHard is much larger than the jet pT, then there is a problem.
-       // Only for PYTHIA.
-       if(!fReadStack) return kTRUE; //Information not filtered to AOD
-       
-       if(!strcmp(GetGenEventHeader()->ClassName(), "AliGenPythiaEventHeader")){
-               TParticle * jet =  0;
-               AliGenPythiaEventHeader* pygeh= (AliGenPythiaEventHeader*) GetGenEventHeader();
-               Int_t nTriggerJets =  pygeh->NTriggerJets();
-               Float_t ptHard = pygeh->GetPtHard();
-               
-               //if(fDebug > 1) printf("AliMCAnalysisUtils::PythiaEventHeader: Njets: %d, pT Hard %f\n",nTriggerJets, ptHard);
+  // Check the event, if the requested ptHard is much larger than the jet pT, then there is a problem.
+  // Only for PYTHIA.
+  if(!fReadStack) return kTRUE; //Information not filtered to AOD
+  
+  if(!strcmp(GetGenEventHeader()->ClassName(), "AliGenPythiaEventHeader")){
+    TParticle * jet =  0;
+    AliGenPythiaEventHeader* pygeh= (AliGenPythiaEventHeader*) GetGenEventHeader();
+    Int_t nTriggerJets =  pygeh->NTriggerJets();
+    Float_t ptHard = pygeh->GetPtHard();
+    
+    //if(fDebug > 1) printf("AliMCAnalysisUtils::PythiaEventHeader: Njets: %d, pT Hard %f\n",nTriggerJets, ptHard);
     Float_t tmpjet[]={0,0,0,0};
-               for(Int_t ijet = 0; ijet< nTriggerJets; ijet++){
-                       pygeh->TriggerJet(ijet, tmpjet);
-                       jet = new TParticle(94, 21, -1, -1, -1, -1, tmpjet[0],tmpjet[1],tmpjet[2],tmpjet[3], 0,0,0,0);
-                       //Compare jet pT and pt Hard
-                       //if(fDebug > 1) printf("AliMCAnalysisUtils:: %d pycell jet pT %f\n",ijet, jet->Pt());
-                       if(jet->Pt() > fPtHardAndJetPtFactor * ptHard) {
-                               printf("AliMCAnalysisUtils::PythiaEventHeader: Njets: %d, pT Hard %2.2f, pycell jet pT %2.2f, rejection factor %1.1f\n",
-                                          nTriggerJets, ptHard, jet->Pt(), fPtHardAndJetPtFactor);
-                               return kFALSE;
-                       }
-               }
+    for(Int_t ijet = 0; ijet< nTriggerJets; ijet++){
+      pygeh->TriggerJet(ijet, tmpjet);
+      jet = new TParticle(94, 21, -1, -1, -1, -1, tmpjet[0],tmpjet[1],tmpjet[2],tmpjet[3], 0,0,0,0);
+      //Compare jet pT and pt Hard
+      //if(fDebug > 1) printf("AliMCAnalysisUtils:: %d pycell jet pT %f\n",ijet, jet->Pt());
+      if(jet->Pt() > fPtHardAndJetPtFactor * ptHard) {
+       printf("AliMCAnalysisUtils::PythiaEventHeader: Njets: %d, pT Hard %2.2f, pycell jet pT %2.2f, rejection factor %1.1f\n",
+              nTriggerJets, ptHard, jet->Pt(), fPtHardAndJetPtFactor);
+       return kFALSE;
+      }
+    }
     if(jet) delete jet; 
-       }
-         
-       return kTRUE ;
-       
+  }
+  
+  return kTRUE ;
+  
 }
 
 //____________________________________________________________________________
@@ -233,24 +233,24 @@ TClonesArray* AliCaloTrackReader::GetAODMCParticles(Int_t input) const {
 
 //____________________________________________________________________________
 AliAODMCHeader* AliCaloTrackReader::GetAODMCHeader(Int_t input) const {
-       //Return MC header in AOD. Do it for the corresponding input event.
+  //Return MC header in AOD. Do it for the corresponding input event.
   AliAODMCHeader *mch = NULL;
-       if(fDataType == kAOD){
-               //Normal input AOD
-               if(input == 0) {
+  if(fDataType == kAOD){
+    //Normal input AOD
+    if(input == 0) {
       mch = (AliAODMCHeader*)((AliAODEvent*)fInputEvent)->FindListObject("mcheader");
     }
-//             //Second input AOD
-//             else if(input == 1){ 
-//       mch = (AliAODMCHeader*) fSecondInputAODEvent->FindListObject("mcheader");
-//  }
-               else {
-                       printf("AliCaloTrackReader::GetAODMCHeader() - wrong AOD input index, %d\n",input);
-               }
-       }
-       else {
-               printf("AliCaloTrackReader::GetAODMCHeader() - Input are not AODs\n");
-       }
+    //         //Second input AOD
+    //         else if(input == 1){ 
+    //       mch = (AliAODMCHeader*) fSecondInputAODEvent->FindListObject("mcheader");
+    //  }
+    else {
+      printf("AliCaloTrackReader::GetAODMCHeader() - wrong AOD input index, %d\n",input);
+    }
+  }
+  else {
+    printf("AliCaloTrackReader::GetAODMCHeader() - Input are not AODs\n");
+  }
   
   return mch;
 }
@@ -258,17 +258,17 @@ AliAODMCHeader* AliCaloTrackReader::GetAODMCHeader(Int_t input) const {
 //_______________________________________________________________
 void AliCaloTrackReader::Init()
 {
-       //Init reader. Method to be called in AliAnaPartCorrMaker
-       
-       //Get the file with second input events if the filename is given
-       //Get the tree and connect the AODEvent. Only with AODs
-
-       if(fReadStack && fReadAODMCParticles){
-               printf("AliCaloTrackReader::Init() - Cannot access stack and mcparticles at the same time, change them \n");
-               fReadStack = kFALSE;
-               fReadAODMCParticles = kFALSE;
-       }
-       
+  //Init reader. Method to be called in AliAnaPartCorrMaker
+  
+  //Get the file with second input events if the filename is given
+  //Get the tree and connect the AODEvent. Only with AODs
+  
+  if(fReadStack && fReadAODMCParticles){
+    printf("AliCaloTrackReader::Init() - Cannot access stack and mcparticles at the same time, change them \n");
+    fReadStack = kFALSE;
+    fReadAODMCParticles = kFALSE;
+  }
+  
 //     if(fSecondInputFileName!=""){
 //             if(fDataType == kAOD){
 //                     TFile * input2   = new TFile(fSecondInputFileName,"read");
@@ -440,8 +440,11 @@ Bool_t AliCaloTrackReader::FillInputEvent(const Int_t iEntry, const char * curre
   if(fFillPHOSCells)  
     FillInputPHOSCells();
        
-  if(fFillCTS)   
+  if(fFillCTS){   
     FillInputCTS();
+    //Accept events with at least one track
+    if(fTrackMult == 0) return kFALSE;
+  }
   if(fFillEMCAL) 
     FillInputEMCAL();
   if(fFillPHOS)