If not needed do not create the histograms in the neutral meson selection task
authorgconesab <gconesab@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 15 Jul 2010 14:43:18 +0000 (14:43 +0000)
committergconesab <gconesab@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 15 Jul 2010 14:43:18 +0000 (14:43 +0000)
PWG4/PartCorrBase/AliNeutralMesonSelection.cxx
PWG4/PartCorrBase/AliNeutralMesonSelection.h

index ab866b5d8457a4ecd5abfd91811a96e9ec7ff2d4..d6f03f67c326f8f32aa95e9f9751331e35fb6ce6 100755 (executable)
@@ -134,57 +134,59 @@ TList *  AliNeutralMesonSelection::GetCreateOutputObjects()
   
   TList * outputContainer = new TList() ; 
   outputContainer->SetName("MesonDecayHistos") ; 
-       if(fKeepNeutralMesonHistos) outputContainer->SetOwner(kFALSE);
-       
-  fhAnglePairNoCut  = new TH2F
-    ("AnglePairNoCut",
-     "Angle between all #gamma pair vs E_{#pi^{0}}",fHistoNEBins,fHistoEMin,fHistoEMax,fHistoNAngleBins,fHistoAngleMin,fHistoAngleMax); 
-  fhAnglePairNoCut->SetYTitle("Angle (rad)");
-  fhAnglePairNoCut->SetXTitle("E_{ #pi^{0}} (GeV)");
-  
-  fhAnglePairOpeningAngleCut  = new TH2F
-    ("AnglePairOpeningAngleCut",
-     "Angle between all #gamma pair (opening angle + azimuth cut) vs E_{#pi^{0}}"
-     ,fHistoNEBins,fHistoEMin,fHistoEMax,fHistoNAngleBins,fHistoAngleMin,fHistoAngleMax); 
-  fhAnglePairOpeningAngleCut->SetYTitle("Angle (rad)");
-  fhAnglePairOpeningAngleCut->SetXTitle("E_{ #pi^{0}} (GeV)");
-  
-  fhAnglePairAllCut  = new TH2F
-    ("AnglePairAllCut",
-     "Angle between all #gamma pair (opening angle + inv mass cut+azimuth) vs E_{#pi^{0}}"
-     ,fHistoNEBins,fHistoEMin,fHistoEMax,fHistoNAngleBins,fHistoAngleMin,fHistoAngleMax); 
-  fhAnglePairAllCut->SetYTitle("Angle (rad)");
-  fhAnglePairAllCut->SetXTitle("E_{ #pi^{0}} (GeV)");    
-  
-  //
-  fhInvMassPairNoCut  = new TH2F
-    ("InvMassPairNoCut","Invariant Mass of all #gamma pair vs E_{#pi^{0}}",
-     fHistoNPtBins,fHistoPtMin,fHistoPtMax,fHistoNIMBins,fHistoIMMin,fHistoIMMax); 
-  fhInvMassPairNoCut->SetYTitle("Invariant Mass (GeV/c^{2})");
-  fhInvMassPairNoCut->SetXTitle("E_{ #pi^{0}} (GeV)");
-  
-  fhInvMassPairOpeningAngleCut  = new TH2F
-    ("InvMassPairOpeningAngleCut",
-     "Invariant Mass of #gamma pair (angle cut) vs E_{#pi^{0}}",
-     fHistoNPtBins,fHistoPtMin,fHistoPtMax,fHistoNIMBins,fHistoIMMin,fHistoIMMax); 
-  fhInvMassPairOpeningAngleCut->SetYTitle("Invariant Mass (GeV/c^{2})");
-  fhInvMassPairOpeningAngleCut->SetXTitle(" E_{#pi^{0}}(GeV)");
-  
-  fhInvMassPairAllCut  = new TH2F
-    ("InvMassPairAllCut",
-     "Invariant Mass of #gamma pair (opening angle+invmass cut) vs E_{#pi^{0}}",
-     fHistoNPtBins,fHistoPtMin,fHistoPtMax,fHistoNIMBins,fHistoIMMin,fHistoIMMax); 
-  fhInvMassPairAllCut->SetYTitle("Invariant Mass (GeV/c^{2})");
-  fhInvMassPairAllCut->SetXTitle("E_{#pi^{0}}(GeV)");
-  
-  outputContainer->Add(fhAnglePairNoCut) ; 
-  outputContainer->Add(fhAnglePairOpeningAngleCut) ;
-  outputContainer->Add(fhAnglePairAllCut) ; 
-  
-  outputContainer->Add(fhInvMassPairNoCut) ; 
-  outputContainer->Add(fhInvMassPairOpeningAngleCut) ; 
-  outputContainer->Add(fhInvMassPairAllCut) ; 
-  
+  if(fKeepNeutralMesonHistos){
+         
+         outputContainer->SetOwner(kFALSE);
+         
+         fhAnglePairNoCut  = new TH2F
+         ("AnglePairNoCut",
+          "Angle between all #gamma pair vs E_{#pi^{0}}",fHistoNEBins,fHistoEMin,fHistoEMax,fHistoNAngleBins,fHistoAngleMin,fHistoAngleMax); 
+         fhAnglePairNoCut->SetYTitle("Angle (rad)");
+         fhAnglePairNoCut->SetXTitle("E_{ #pi^{0}} (GeV)");
+         
+         fhAnglePairOpeningAngleCut  = new TH2F
+         ("AnglePairOpeningAngleCut",
+          "Angle between all #gamma pair (opening angle + azimuth cut) vs E_{#pi^{0}}"
+          ,fHistoNEBins,fHistoEMin,fHistoEMax,fHistoNAngleBins,fHistoAngleMin,fHistoAngleMax); 
+         fhAnglePairOpeningAngleCut->SetYTitle("Angle (rad)");
+         fhAnglePairOpeningAngleCut->SetXTitle("E_{ #pi^{0}} (GeV)");
+         
+         fhAnglePairAllCut  = new TH2F
+         ("AnglePairAllCut",
+          "Angle between all #gamma pair (opening angle + inv mass cut+azimuth) vs E_{#pi^{0}}"
+          ,fHistoNEBins,fHistoEMin,fHistoEMax,fHistoNAngleBins,fHistoAngleMin,fHistoAngleMax); 
+         fhAnglePairAllCut->SetYTitle("Angle (rad)");
+         fhAnglePairAllCut->SetXTitle("E_{ #pi^{0}} (GeV)");    
+         
+         //
+         fhInvMassPairNoCut  = new TH2F
+         ("InvMassPairNoCut","Invariant Mass of all #gamma pair vs E_{#pi^{0}}",
+          fHistoNPtBins,fHistoPtMin,fHistoPtMax,fHistoNIMBins,fHistoIMMin,fHistoIMMax); 
+         fhInvMassPairNoCut->SetYTitle("Invariant Mass (GeV/c^{2})");
+         fhInvMassPairNoCut->SetXTitle("E_{ #pi^{0}} (GeV)");
+         
+         fhInvMassPairOpeningAngleCut  = new TH2F
+         ("InvMassPairOpeningAngleCut",
+          "Invariant Mass of #gamma pair (angle cut) vs E_{#pi^{0}}",
+          fHistoNPtBins,fHistoPtMin,fHistoPtMax,fHistoNIMBins,fHistoIMMin,fHistoIMMax); 
+         fhInvMassPairOpeningAngleCut->SetYTitle("Invariant Mass (GeV/c^{2})");
+         fhInvMassPairOpeningAngleCut->SetXTitle(" E_{#pi^{0}}(GeV)");
+         
+         fhInvMassPairAllCut  = new TH2F
+         ("InvMassPairAllCut",
+          "Invariant Mass of #gamma pair (opening angle+invmass cut) vs E_{#pi^{0}}",
+          fHistoNPtBins,fHistoPtMin,fHistoPtMax,fHistoNIMBins,fHistoIMMin,fHistoIMMax); 
+         fhInvMassPairAllCut->SetYTitle("Invariant Mass (GeV/c^{2})");
+         fhInvMassPairAllCut->SetXTitle("E_{#pi^{0}}(GeV)");
+         
+         outputContainer->Add(fhAnglePairNoCut) ; 
+         outputContainer->Add(fhAnglePairOpeningAngleCut) ;
+         outputContainer->Add(fhAnglePairAllCut) ; 
+         
+         outputContainer->Add(fhInvMassPairNoCut) ; 
+         outputContainer->Add(fhInvMassPairOpeningAngleCut) ; 
+         outputContainer->Add(fhInvMassPairAllCut) ; 
+  }
   return outputContainer;
 }
 
@@ -261,19 +263,24 @@ Bool_t  AliNeutralMesonSelection::SelectPair(TLorentzVector gammai, TLorentzVect
   Double_t e       = (gammai+gammaj).E();
   
   //Fill histograms with no cuts applied.
-  fhAnglePairNoCut->Fill(e,angle);
-  fhInvMassPairNoCut->Fill(e,invmass);
-    
+  if(fKeepNeutralMesonHistos){
+         fhAnglePairNoCut  ->Fill(e,angle);
+         fhInvMassPairNoCut->Fill(e,invmass);
+  }
   //Cut on the aperture of the pair
   if(IsAngleInWindow(angle,e)){
-    fhAnglePairOpeningAngleCut     ->Fill(e,angle);
-    fhInvMassPairOpeningAngleCut->Fill(e,invmass);
+       if(fKeepNeutralMesonHistos){
+               fhAnglePairOpeningAngleCut  ->Fill(e,angle);
+               fhInvMassPairOpeningAngleCut->Fill(e,invmass);
+       }
     //AliDebug(2,Form("Angle cut: pt %f, phi %f",pt,phi));
     
     //Cut on the invariant mass of the pair
     if((invmass>fInvMassMinCut) && (invmass<fInvMassMaxCut)){ 
-      fhInvMassPairAllCut  ->Fill(e,invmass);
-      fhAnglePairAllCut       ->Fill(e,angle);
+         if(fKeepNeutralMesonHistos){
+                 fhInvMassPairAllCut  ->Fill(e,invmass);
+                 fhAnglePairAllCut    ->Fill(e,angle);
+         }
       goodpair = kTRUE;
       //AliDebug(2,Form("IM cut: pt %f, phi %f",pt,phi));
     }//(invmass>0.125) && (invmass<0.145)
index 0c85577d2c62f899bbfbdaca5a86c75647bbacbd..ff5a96a1805ad405a54c1274e72e1f5a7ea25916 100755 (executable)
@@ -109,26 +109,26 @@ class AliNeutralMesonSelection : public TObject {
   Bool_t   fKeepNeutralMesonHistos ; // Keep neutral meson selection histograms
   
   //Histograms
-  TH2F * fhAnglePairNoCut  ;  //Aperture angle of decay photons, no cuts
-  TH2F * fhAnglePairOpeningAngleCut   ;  //Aperture angle of decay photons, cut on opening angle
-  TH2F * fhAnglePairAllCut   ;  //Aperture angle of decay photons, all cuts
-  TH2F * fhInvMassPairNoCut    ;  //Invariant mass of decay photons, no cuts
-  TH2F * fhInvMassPairOpeningAngleCut  ;  //Invariant mass of decay photons, cut on opening angle
-  TH2F * fhInvMassPairAllCut   ;  //Invariant mass of decay photons, all cuts
+  TH2F * fhAnglePairNoCut  ;              //! Aperture angle of decay photons, no cuts
+  TH2F * fhAnglePairOpeningAngleCut   ;   //! Aperture angle of decay photons, cut on opening angle
+  TH2F * fhAnglePairAllCut   ;            //! Aperture angle of decay photons, all cuts
+  TH2F * fhInvMassPairNoCut    ;          //! Invariant mass of decay photons, no cuts
+  TH2F * fhInvMassPairOpeningAngleCut  ;  //Invariant mass of decay photons, cut on opening angle
+  TH2F * fhInvMassPairAllCut   ;          //! Invariant mass of decay photons, all cuts
   
   //Histograms binning and range    
-  Int_t   fHistoNEBins ;  //Number of bins in pi0 E axis
-  Float_t fHistoEMax ;    //Maximum value of pi0 E histogram range
-  Float_t fHistoEMin ;    //Minimum value of pi0 E histogram range
-  Int_t   fHistoNPtBins ;  //Number of bins in Pt trigger axis
-  Float_t fHistoPtMax ;    //Maximum value of Pt trigger histogram range
-  Float_t fHistoPtMin ;    //Minimum value of Pt trigger histogram range               
-  Int_t   fHistoNAngleBins ; //Number of bins in angle axis
-  Float_t fHistoAngleMax ;//Maximum value of angle histogram range
-  Float_t fHistoAngleMin ;//Minimum value of angle histogram range
-  Int_t   fHistoNIMBins ; //Number of bins in Invariant Mass axis
-  Float_t fHistoIMMax ;   //Maximum value of Invariant Mass histogram range
-  Float_t fHistoIMMin ;   //Minimum value of Invariant Mass histogram range
+  Int_t   fHistoNEBins ;     // Number of bins in pi0 E axis
+  Float_t fHistoEMax ;       // Maximum value of pi0 E histogram range
+  Float_t fHistoEMin ;       // Minimum value of pi0 E histogram range
+  Int_t   fHistoNPtBins ;    // Number of bins in Pt trigger axis
+  Float_t fHistoPtMax ;      // Maximum value of Pt trigger histogram range
+  Float_t fHistoPtMin ;      // Minimum value of Pt trigger histogram range            
+  Int_t   fHistoNAngleBins ; // Number of bins in angle axis
+  Float_t fHistoAngleMax ;   // Maximum value of angle histogram range
+  Float_t fHistoAngleMin ;   // Minimum value of angle histogram range
+  Int_t   fHistoNIMBins ;    // Number of bins in Invariant Mass axis
+  Float_t fHistoIMMax ;      // Maximum value of Invariant Mass histogram range
+  Float_t fHistoIMMin ;      // Minimum value of Invariant Mass histogram range
   
   ClassDef(AliNeutralMesonSelection,3)