[u/mrichter/AliRoot.git] / PWG4 / AliNeutralMesonSelection.cxx

 /**************************************************************************
 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
 *                                                                        *
 * Author: The ALICE Off-line Project.                                    *
 * Contributors are mentioned in the code where appropriate.              *
 *                                                                        *
 * Permission to use, copy, modify and distribute this software and its   *
 * documentation strictly for non-commercial purposes is hereby granted   *
 * without fee, provided that the above copyright notice appears in all   *
 * copies and that both the copyright notice and this permission notice   *
 * appear in the supporting documentation. The authors make no claims     *
 * about the suitability of this software for any purpose. It is          *
 * provided "as is" without express or implied warranty.                  *
 **************************************************************************/
/* $Id$ */

/* History of cvs commits:
 *
 * $Log$
 * Revision 1.2  2007/08/17 12:40:04  schutz
 * New analysis classes by Gustavo Conesa
 *
 * Revision 1.1.2.1  2007/07/26 10:32:09  schutz
 * new analysis classes in the the new analysis framework
 *
 *
 */

//_________________________________________________________________________
// Class that contains methods to select candidate pairs to neutral meson 
// 2 main selections, invariant mass around pi0 (also any other mass),
// apperture angle to distinguish from combinatorial.
// There is a 3rd cut based on the gamma correlation on phi or pt.
//-- Author: Gustavo Conesa (INFN-LNF)

// --- ROOT system ---
#include <TParticle.h>
#include <TLorentzVector.h>
#include <TH2.h>
#include <TList.h>
#include <TArrayD.h>
 
//---- AliRoot system ----
#include "AliNeutralMesonSelection.h" 
#include "AliLog.h"

ClassImp(AliNeutralMesonSelection)
  
  
//____________________________________________________________________________
  AliNeutralMesonSelection::AliNeutralMesonSelection() : 
    TObject(), fSelect(0), fM(0),
    fInvMassMaxCut(0.), fInvMassMinCut(0.),
    fAngleMaxParam(),  fMinPt(0),
    fDeltaPhiMaxCut(0.), fDeltaPhiMinCut(0.),
    fRatioMaxCut(0), fRatioMinCut(0),  fKeepNeutralMesonHistos(0),
    fhAnglePairNoCut(0),  fhAnglePairCorrelationCut(0),
    fhAnglePairOpeningAngleCut(0), fhAnglePairAllCut(0), 
    fhInvMassPairNoCut(0),   fhInvMassPairCorrelationCut(0), 
    fhInvMassPairOpeningAngleCut(0), fhInvMassPairAllCut(0) 
{
  //Default Ctor
  
  //Initialize parameters
  
  // kGammaHadron and kGammaJet 
  fAngleMaxParam.Set(4) ;
  fAngleMaxParam.Reset(0.);
  
  //Initialize parameters
  InitParameters();
}

//____________________________________________________________________________
AliNeutralMesonSelection::AliNeutralMesonSelection(const AliNeutralMesonSelection & g) :   
  TObject(),  
  fSelect(g.fSelect), fM(g.fM),
  fInvMassMaxCut(g.fInvMassMaxCut), fInvMassMinCut(g.fInvMassMinCut),
  fAngleMaxParam(g.fAngleMaxParam), fMinPt(g.fMinPt),
  fDeltaPhiMaxCut(g.fDeltaPhiMaxCut), fDeltaPhiMinCut(g.fDeltaPhiMinCut), 
  fRatioMaxCut(g.fRatioMaxCut), fRatioMinCut(g.fRatioMinCut), 
  fKeepNeutralMesonHistos(g.fKeepNeutralMesonHistos),
  fhAnglePairNoCut(g. fhAnglePairNoCut), 
  fhAnglePairCorrelationCut(g. fhAnglePairCorrelationCut), 
  fhAnglePairOpeningAngleCut(g. fhAnglePairOpeningAngleCut), 
  fhAnglePairAllCut(g. fhAnglePairAllCut), 
  fhInvMassPairNoCut(g.fhInvMassPairNoCut),  
  fhInvMassPairCorrelationCut(g.fhInvMassPairCorrelationCut), 
  fhInvMassPairOpeningAngleCut(g.fhInvMassPairOpeningAngleCut), 
  fhInvMassPairAllCut(g.fhInvMassPairAllCut)
{
  // cpy ctor
}

//_________________________________________________________________________
AliNeutralMesonSelection & AliNeutralMesonSelection::operator = (const AliNeutralMesonSelection & source)
{
  // assignment operator
  
  if(this == &source)return *this;
  ((TObject *)this)->operator=(source);

  fSelect = source.fSelect ;
  fM = source.fM ;
  fInvMassMaxCut = source.fInvMassMaxCut ; 
  fInvMassMinCut = source.fInvMassMinCut ;
  fAngleMaxParam = source.fAngleMaxParam ;
  fMinPt = source.fMinPt ;
  fDeltaPhiMaxCut = source.fDeltaPhiMaxCut ; 
  fDeltaPhiMinCut = source.fDeltaPhiMinCut ; 
  fRatioMaxCut = source.fRatioMaxCut ; 
  fRatioMinCut = source.fRatioMinCut ;
  fKeepNeutralMesonHistos = source.fKeepNeutralMesonHistos;
 
  fhAnglePairNoCut = source. fhAnglePairNoCut ; 
  fhAnglePairCorrelationCut = source. fhAnglePairCorrelationCut ; 
  fhAnglePairOpeningAngleCut = source. fhAnglePairOpeningAngleCut ; 
  fhAnglePairAllCut = source. fhAnglePairAllCut ; 
  fhInvMassPairNoCut = source.fhInvMassPairNoCut ; 
  fhInvMassPairCorrelationCut = source.fhInvMassPairCorrelationCut ; 
  fhInvMassPairOpeningAngleCut = source.fhInvMassPairOpeningAngleCut ; 
  fhInvMassPairAllCut = source.fhInvMassPairAllCut ; 
  
  return *this;
  
}

//____________________________________________________________________________
AliNeutralMesonSelection::~AliNeutralMesonSelection() 
{
  // Remove all pointers
    
  delete fhAnglePairNoCut  ; 
  delete fhAnglePairCorrelationCut  ; 
  delete fhAnglePairOpeningAngleCut   ; 
  delete fhAnglePairAllCut   ;  
  delete fhInvMassPairNoCut    ; 
  delete fhInvMassPairCorrelationCut    ; 
  delete fhInvMassPairOpeningAngleCut  ; 
  delete fhInvMassPairAllCut   ;    
  
}


//________________________________________________________________________
TList *  AliNeutralMesonSelection::GetCreateOutputObjects()
{  

  // Create histograms to be saved in output file and 
  // store them in outputContainer
  TList * outputContainer = new TList() ; 
  outputContainer->SetName("MesonDecayHistos") ; 
  
  fhAnglePairNoCut  = new TH2F
    ("AnglePairNoCut",
     "Angle between all #gamma pair vs E_{#pi^{0}}",200,0,50,200,0,0.2); 
  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}}"
     ,200,0,50,200,0,0.2); 
  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}}"
     ,200,0,50,200,0,0.2); 
  fhAnglePairAllCut->SetYTitle("Angle (rad)");
  fhAnglePairAllCut->SetXTitle("E_{ #pi^{0}} (GeV)");    
  
  //
  fhInvMassPairNoCut  = new TH2F
    ("InvMassPairNoCut","Invariant Mass of all #gamma pair vs E_{#pi^{0}}",
     120,0,120,360,0,0.5); 
  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}}",
     120,0,120,360,0,0.5); 
  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}}",
     120,0,120,360,0,0.5); 
  fhInvMassPairAllCut->SetYTitle("Invariant Mass (GeV/c^{2})");
  fhInvMassPairAllCut->SetXTitle("E_{#pi^{0}}(GeV)");

  fhAnglePairCorrelationCut  = new TH2F
    ("AnglePairCorrelationCut",
     "Angle between correlated #gamma pair vs E_{#pi^{0}}",200,0,50,200,0,0.2); 
  fhAnglePairCorrelationCut->SetYTitle("Angle (rad)");
  fhAnglePairCorrelationCut->SetXTitle("E_{ #pi^{0}} (GeV)");

  fhInvMassPairCorrelationCut  = new TH2F
    ("InvMassPairCorrelationCut","Invariant Mass of correlated #gamma pair vs E_{#pi^{0}}",
     120,0,120,360,0,0.5); 
  fhInvMassPairCorrelationCut->SetYTitle("Invariant Mass (GeV/c^{2})");
  fhInvMassPairCorrelationCut->SetXTitle("E_{ #pi^{0}} (GeV)");
  
  outputContainer->Add(fhAnglePairNoCut) ; 
  outputContainer->Add(fhAnglePairOpeningAngleCut) ;
  outputContainer->Add(fhAnglePairAllCut) ; 
  
  outputContainer->Add(fhInvMassPairNoCut) ; 
  outputContainer->Add(fhInvMassPairOpeningAngleCut) ; 
  outputContainer->Add(fhInvMassPairAllCut) ; 

  outputContainer->Add(fhAnglePairCorrelationCut) ; 
  outputContainer->Add(fhInvMassPairCorrelationCut) ; 
  
  return outputContainer;
}

 //____________________________________________________________________________
void AliNeutralMesonSelection::InitParameters()
{
 
  //Initialize the parameters of the analysis.
  fKeepNeutralMesonHistos = kTRUE ;

  //-------------kHadron, kJetLeadCone-----------------
  fAngleMaxParam.Set(4) ;
  fAngleMaxParam.AddAt(0.4,0);//={0.4,-0.25,0.025,-2e-4};
  fAngleMaxParam.AddAt(-0.25,1) ;
  fAngleMaxParam.AddAt(0.025,2) ;
  fAngleMaxParam.AddAt(-2e-4,3) ;

  fInvMassMaxCut  = 0.16 ;
  fInvMassMinCut  = 0.11 ;

  fM = 0.1349766;//neutralMeson mass

  fMinPt = 0.   ;
  fDeltaPhiMaxCut      = 4.5;
  fDeltaPhiMinCut      = 1.5 ;
  fRatioMaxCut    = 1.0 ;
  fRatioMinCut    = 0.1 ; 
}

//__________________________________________________________________________-
Bool_t AliNeutralMesonSelection::IsAngleInWindow(const Float_t angle,const Float_t e) const {
  //Check if the opening angle of the candidate pairs is inside 
  //our selection windowd

  Bool_t result = kFALSE;
  Double_t max =  fAngleMaxParam.At(0)*TMath::Exp(fAngleMaxParam.At(1)*e)
    +fAngleMaxParam.At(2)+fAngleMaxParam.At(3)*e;
  Double_t arg = (e*e-2*fM*fM)/(e*e);
  Double_t min = 100. ;
  if(arg>0.)
    min = TMath::ACos(arg);

  if((angle<max)&&(angle>=min))
    result = kTRUE;
 
  return result;
}

//____________________________________________________________________________
Bool_t  AliNeutralMesonSelection::CutPtPhi(Double_t ptg, Double_t phig, Double_t pt, Double_t phi)  const
{ 
  //Select pair if delta
  Bool_t cut = kFALSE ;
 
  if(fSelect == kNoSelectPhiPt) cut = kTRUE ;
  else if((phig-phi) > fDeltaPhiMinCut && ((phig-phi) < fDeltaPhiMaxCut)){
    //Cut on pt
    if((fSelect == kSelectPhiPtRatio && ptg > 0. && pt/ptg  > fRatioMinCut &&  pt/ptg  < fRatioMaxCut) ||
	(fSelect == kSelectPhiMinPt && pt > fMinPt)  )  cut = kTRUE ;
  }
  else cut = kFALSE ;
  
  return cut ;
  
}

//____________________________________________________________________________
Bool_t  AliNeutralMesonSelection::SelectPair(TParticle * pGamma, TLorentzVector gammai, TLorentzVector gammaj)  
{  
  
  //Search for the neutral pion within selection cuts
  
  Double_t ptg = pGamma->Pt();
  Double_t phig = pGamma->Phi() ;
  Bool_t goodpair = kFALSE ;
  
  Double_t pt  = (gammai+gammaj).Pt();
  Double_t phi = (gammai+gammaj).Phi();
  if(phi < 0)
    phi+=TMath::TwoPi();
  Double_t invmass = (gammai+gammaj).M();
  Double_t angle   = gammaj.Angle(gammai.Vect());
  Double_t e       = (gammai+gammaj).E();
  
  //Fill histograms with no cuts applied.
  fhAnglePairNoCut->Fill(e,angle);
  fhInvMassPairNoCut->Fill(e,invmass);
  
  //Cut on phig-phi meson
  if(CutPtPhi(ptg, phig, pt, phi)){
    
    fhAnglePairCorrelationCut     ->Fill(e,angle);
    fhInvMassPairCorrelationCut->Fill(e,invmass);
    
    //Cut on the aperture of the pair
    if(IsAngleInWindow(angle,e)){
      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);
	goodpair = kTRUE;
	AliDebug(2,Form("IM cut: pt %f, phi %f",pt,phi));
      }//(invmass>0.125) && (invmass<0.145)
    }//Opening angle cut
  } // cut on pt and phi
  
  
  return goodpair; 
  
}

//__________________________________________________________________
void AliNeutralMesonSelection::Print(const Option_t * opt) const
{

  //Print some relevant parameters set for the analysis
  if(! opt)
    return;
  
  Info("Print", "%s %s", GetName(), GetTitle() ) ;

  printf("mass : %f  \n", fM );
  printf("Invariant mass limits : %f < m < %f \n", fInvMassMinCut , fInvMassMinCut );
  printf("Angle selection param: \n");
  printf("p0 :     %f", fAngleMaxParam.At(0));
  printf("p1 :     %f", fAngleMaxParam.At(1));
  printf("p2 :     %f", fAngleMaxParam.At(2));
  printf("p3 :     %f", fAngleMaxParam.At(3));

  printf("pT meson       >    %f\n", fMinPt) ; 
  printf("Phi gamma-meson      <     %f\n", fDeltaPhiMaxCut) ; 
  printf("Phi gamma-meson      >     %f\n", fDeltaPhiMinCut) ;
  printf("pT meson / pT Gamma             <     %f\n", fRatioMaxCut) ; 
  printf("pT meson / pT Gamma             >     %f\n", fRatioMinCut) ;
  printf("Keep Neutral Meson Histos = %d\n",fKeepNeutralMesonHistos);

}
Commit	Line	Data
bdcfac30	1	/**************************************************************************
	2	* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
	3	* *
	4	* Author: The ALICE Off-line Project. *
	5	* Contributors are mentioned in the code where appropriate. *
	6	* *
	7	* Permission to use, copy, modify and distribute this software and its *
	8	* documentation strictly for non-commercial purposes is hereby granted *
	9	* without fee, provided that the above copyright notice appears in all *
	10	* copies and that both the copyright notice and this permission notice *
	11	* appear in the supporting documentation. The authors make no claims *
	12	* about the suitability of this software for any purpose. It is *
	13	* provided "as is" without express or implied warranty. *
	14	**************************************************************************/
	15	/* $Id$ */
	16
	17	/* History of cvs commits:
	18	*
	19	* $Log$
3bb2c538	20	* Revision 1.2 2007/08/17 12:40:04 schutz
	21	* New analysis classes by Gustavo Conesa
	22	*
bdcfac30	23	* Revision 1.1.2.1 2007/07/26 10:32:09 schutz
	24	* new analysis classes in the the new analysis framework
	25	*
	26	*
	27	*/
	28
	29	//_________________________________________________________________________
	30	// Class that contains methods to select candidate pairs to neutral meson
3bb2c538	31	// 2 main selections, invariant mass around pi0 (also any other mass),
	32	// apperture angle to distinguish from combinatorial.
	33	// There is a 3rd cut based on the gamma correlation on phi or pt.
bdcfac30	34	//-- Author: Gustavo Conesa (INFN-LNF)
	35
	36	// --- ROOT system ---
	37	#include <TParticle.h>
	38	#include <TLorentzVector.h>
	39	#include <TH2.h>
	40	#include <TList.h>
3bb2c538	41	#include <TArrayD.h>
3bb2c538	42
bdcfac30	43	//---- AliRoot system ----
bdcfac30	44	#include "AliNeutralMesonSelection.h"
bdcfac30	45	#include "AliLog.h"
	46
	47	ClassImp(AliNeutralMesonSelection)
	48
	49
	50	//____________________________________________________________________________
	51	AliNeutralMesonSelection::AliNeutralMesonSelection() :
	52	TObject(), fSelect(0), fM(0),
	53	fInvMassMaxCut(0.), fInvMassMinCut(0.),
	54	fAngleMaxParam(), fMinPt(0),
	55	fDeltaPhiMaxCut(0.), fDeltaPhiMinCut(0.),
	56	fRatioMaxCut(0), fRatioMinCut(0), fKeepNeutralMesonHistos(0),
	57	fhAnglePairNoCut(0), fhAnglePairCorrelationCut(0),
	58	fhAnglePairOpeningAngleCut(0), fhAnglePairAllCut(0),
	59	fhInvMassPairNoCut(0), fhInvMassPairCorrelationCut(0),
	60	fhInvMassPairOpeningAngleCut(0), fhInvMassPairAllCut(0)
	61	{
	62	//Default Ctor
	63
	64	//Initialize parameters
	65
	66	// kGammaHadron and kGammaJet
	67	fAngleMaxParam.Set(4) ;
	68	fAngleMaxParam.Reset(0.);
	69
	70	//Initialize parameters
	71	InitParameters();
	72	}
	73
	74	//____________________________________________________________________________
	75	AliNeutralMesonSelection::AliNeutralMesonSelection(const AliNeutralMesonSelection & g) :
	76	TObject(),
	77	fSelect(g.fSelect), fM(g.fM),
	78	fInvMassMaxCut(g.fInvMassMaxCut), fInvMassMinCut(g.fInvMassMinCut),
	79	fAngleMaxParam(g.fAngleMaxParam), fMinPt(g.fMinPt),
	80	fDeltaPhiMaxCut(g.fDeltaPhiMaxCut), fDeltaPhiMinCut(g.fDeltaPhiMinCut),
	81	fRatioMaxCut(g.fRatioMaxCut), fRatioMinCut(g.fRatioMinCut),
	82	fKeepNeutralMesonHistos(g.fKeepNeutralMesonHistos),
	83	fhAnglePairNoCut(g. fhAnglePairNoCut),
	84	fhAnglePairCorrelationCut(g. fhAnglePairCorrelationCut),
	85	fhAnglePairOpeningAngleCut(g. fhAnglePairOpeningAngleCut),
	86	fhAnglePairAllCut(g. fhAnglePairAllCut),
	87	fhInvMassPairNoCut(g.fhInvMassPairNoCut),
	88	fhInvMassPairCorrelationCut(g.fhInvMassPairCorrelationCut),
	89	fhInvMassPairOpeningAngleCut(g.fhInvMassPairOpeningAngleCut),
	90	fhInvMassPairAllCut(g.fhInvMassPairAllCut)
	91	{
	92	// cpy ctor
	93	}
	94
	95	//_________________________________________________________________________
	96	AliNeutralMesonSelection & AliNeutralMesonSelection::operator = (const AliNeutralMesonSelection & source)
	97	{
	98	// assignment operator
	99
	100	if(this == &source)return *this;
	101	((TObject *)this)->operator=(source);
	102
	103	fSelect = source.fSelect ;
	104	fM = source.fM ;
	105	fInvMassMaxCut = source.fInvMassMaxCut ;
	106	fInvMassMinCut = source.fInvMassMinCut ;
	107	fAngleMaxParam = source.fAngleMaxParam ;
	108	fMinPt = source.fMinPt ;
109	fDeltaPhiMaxCut = source.fDeltaPhiMaxCut ;
110	fDeltaPhiMinCut = source.fDeltaPhiMinCut ;
111	fRatioMaxCut = source.fRatioMaxCut ;
112	fRatioMinCut = source.fRatioMinCut ;
113	fKeepNeutralMesonHistos = source.fKeepNeutralMesonHistos;
114
115	fhAnglePairNoCut = source. fhAnglePairNoCut ;
116	fhAnglePairCorrelationCut = source. fhAnglePairCorrelationCut ;
117	fhAnglePairOpeningAngleCut = source. fhAnglePairOpeningAngleCut ;
118	fhAnglePairAllCut = source. fhAnglePairAllCut ;
119	fhInvMassPairNoCut = source.fhInvMassPairNoCut ;
120	fhInvMassPairCorrelationCut = source.fhInvMassPairCorrelationCut ;
121	fhInvMassPairOpeningAngleCut = source.fhInvMassPairOpeningAngleCut ;
122	fhInvMassPairAllCut = source.fhInvMassPairAllCut ;
123
124	return *this;
125
126	}
127
128	//____________________________________________________________________________
129	AliNeutralMesonSelection::~AliNeutralMesonSelection()
130	{
131	// Remove all pointers
132
133	delete fhAnglePairNoCut ;
134	delete fhAnglePairCorrelationCut ;
135	delete fhAnglePairOpeningAngleCut ;
136	delete fhAnglePairAllCut ;
137	delete fhInvMassPairNoCut ;
138	delete fhInvMassPairCorrelationCut ;
139	delete fhInvMassPairOpeningAngleCut ;
140	delete fhInvMassPairAllCut ;
141
142	}
143
144
145
146	//________________________________________________________________________
147	TList * AliNeutralMesonSelection::GetCreateOutputObjects()
148	{
149
150	// Create histograms to be saved in output file and
151	// store them in outputContainer
152	TList * outputContainer = new TList() ;
153	outputContainer->SetName("MesonDecayHistos") ;
154
155	fhAnglePairNoCut = new TH2F
156	("AnglePairNoCut",
157	"Angle between all #gamma pair vs E_{#pi^{0}}",200,0,50,200,0,0.2);
158	fhAnglePairNoCut->SetYTitle("Angle (rad)");
159	fhAnglePairNoCut->SetXTitle("E_{ #pi^{0}} (GeV)");
160
161	fhAnglePairOpeningAngleCut = new TH2F
162	("AnglePairOpeningAngleCut",
163	"Angle between all #gamma pair (opening angle + azimuth cut) vs E_{#pi^{0}}"
164	,200,0,50,200,0,0.2);
165	fhAnglePairOpeningAngleCut->SetYTitle("Angle (rad)");
166	fhAnglePairOpeningAngleCut->SetXTitle("E_{ #pi^{0}} (GeV)");
167
168	fhAnglePairAllCut = new TH2F
169	("AnglePairAllCut",
170	"Angle between all #gamma pair (opening angle + inv mass cut+azimuth) vs E_{#pi^{0}}"
171	,200,0,50,200,0,0.2);
172	fhAnglePairAllCut->SetYTitle("Angle (rad)");
173	fhAnglePairAllCut->SetXTitle("E_{ #pi^{0}} (GeV)");
174
175	//
176	fhInvMassPairNoCut = new TH2F
177	("InvMassPairNoCut","Invariant Mass of all #gamma pair vs E_{#pi^{0}}",
178	120,0,120,360,0,0.5);
179	fhInvMassPairNoCut->SetYTitle("Invariant Mass (GeV/c^{2})");
180	fhInvMassPairNoCut->SetXTitle("E_{ #pi^{0}} (GeV)");
181
182	fhInvMassPairOpeningAngleCut = new TH2F
183	("InvMassPairOpeningAngleCut",
184	"Invariant Mass of #gamma pair (angle cut) vs E_{#pi^{0}}",
185	120,0,120,360,0,0.5);
186	fhInvMassPairOpeningAngleCut->SetYTitle("Invariant Mass (GeV/c^{2})");
187	fhInvMassPairOpeningAngleCut->SetXTitle(" E_{#pi^{0}}(GeV)");
188
189	fhInvMassPairAllCut = new TH2F
190	("InvMassPairAllCut",
191	"Invariant Mass of #gamma pair (opening angle+invmass cut) vs E_{#pi^{0}}",
192	120,0,120,360,0,0.5);
193	fhInvMassPairAllCut->SetYTitle("Invariant Mass (GeV/c^{2})");
194	fhInvMassPairAllCut->SetXTitle("E_{#pi^{0}}(GeV)");
195
196	fhAnglePairCorrelationCut = new TH2F
197	("AnglePairCorrelationCut",
198	"Angle between correlated #gamma pair vs E_{#pi^{0}}",200,0,50,200,0,0.2);
199	fhAnglePairCorrelationCut->SetYTitle("Angle (rad)");
200	fhAnglePairCorrelationCut->SetXTitle("E_{ #pi^{0}} (GeV)");
201
202	fhInvMassPairCorrelationCut = new TH2F
203	("InvMassPairCorrelationCut","Invariant Mass of correlated #gamma pair vs E_{#pi^{0}}",
204	120,0,120,360,0,0.5);
205	fhInvMassPairCorrelationCut->SetYTitle("Invariant Mass (GeV/c^{2})");
206	fhInvMassPairCorrelationCut->SetXTitle("E_{ #pi^{0}} (GeV)");
207
208	outputContainer->Add(fhAnglePairNoCut) ;
209	outputContainer->Add(fhAnglePairOpeningAngleCut) ;
210	outputContainer->Add(fhAnglePairAllCut) ;
211
212	outputContainer->Add(fhInvMassPairNoCut) ;
213	outputContainer->Add(fhInvMassPairOpeningAngleCut) ;
214	outputContainer->Add(fhInvMassPairAllCut) ;
215
216	outputContainer->Add(fhAnglePairCorrelationCut) ;
217	outputContainer->Add(fhInvMassPairCorrelationCut) ;
218
219	return outputContainer;
220	}
221
222	//____________________________________________________________________________
223	void AliNeutralMesonSelection::InitParameters()
224	{
225
226	//Initialize the parameters of the analysis.
227	fKeepNeutralMesonHistos = kTRUE ;
228
229	//-------------kHadron, kJetLeadCone-----------------
230	fAngleMaxParam.Set(4) ;
231	fAngleMaxParam.AddAt(0.4,0);//={0.4,-0.25,0.025,-2e-4};
232	fAngleMaxParam.AddAt(-0.25,1) ;
233	fAngleMaxParam.AddAt(0.025,2) ;
234	fAngleMaxParam.AddAt(-2e-4,3) ;
235
236	fInvMassMaxCut = 0.16 ;
237	fInvMassMinCut = 0.11 ;
238
239	fM = 0.1349766;//neutralMeson mass
240
241	fMinPt = 0. ;
242	fDeltaPhiMaxCut = 4.5;
243	fDeltaPhiMinCut = 1.5 ;
244	fRatioMaxCut = 1.0 ;
245	fRatioMinCut = 0.1 ;
246	}
247
248	//__________________________________________________________________________-
3bb2c538	249	Bool_t AliNeutralMesonSelection::IsAngleInWindow(const Float_t angle,const Float_t e) const {
bdcfac30	250	//Check if the opening angle of the candidate pairs is inside
	251	//our selection windowd
	252
	253	Bool_t result = kFALSE;
	254	Double_t max = fAngleMaxParam.At(0)TMath::Exp(fAngleMaxParam.At(1)e)
	255	+fAngleMaxParam.At(2)+fAngleMaxParam.At(3)*e;
	256	Double_t arg = (ee-2fMfM)/(ee);
	257	Double_t min = 100. ;
	258	if(arg>0.)
	259	min = TMath::ACos(arg);
	260
	261	if((angle<max)&&(angle>=min))
	262	result = kTRUE;
	263
	264	return result;
	265	}
	266
	267	//____________________________________________________________________________
3bb2c538	268	Bool_t AliNeutralMesonSelection::CutPtPhi(Double_t ptg, Double_t phig, Double_t pt, Double_t phi) const
bdcfac30	269	{
	270	//Select pair if delta
	271	Bool_t cut = kFALSE ;
	272
	273	if(fSelect == kNoSelectPhiPt) cut = kTRUE ;
	274	else if((phig-phi) > fDeltaPhiMinCut && ((phig-phi) < fDeltaPhiMaxCut)){
	275	//Cut on pt
	276	if((fSelect == kSelectPhiPtRatio && ptg > 0. && pt/ptg > fRatioMinCut && pt/ptg < fRatioMaxCut) \|\|
	277	(fSelect == kSelectPhiMinPt && pt > fMinPt) ) cut = kTRUE ;
	278	}
	279	else cut = kFALSE ;
	280
	281	return cut ;
	282
	283	}
	284
	285	//____________________________________________________________________________
	286	Bool_t AliNeutralMesonSelection::SelectPair(TParticle * pGamma, TLorentzVector gammai, TLorentzVector gammaj)
	287	{
	288
	289	//Search for the neutral pion within selection cuts
	290
	291	Double_t ptg = pGamma->Pt();
	292	Double_t phig = pGamma->Phi() ;
	293	Bool_t goodpair = kFALSE ;
	294
	295	Double_t pt = (gammai+gammaj).Pt();
	296	Double_t phi = (gammai+gammaj).Phi();
	297	if(phi < 0)
	298	phi+=TMath::TwoPi();
	299	Double_t invmass = (gammai+gammaj).M();
	300	Double_t angle = gammaj.Angle(gammai.Vect());
	301	Double_t e = (gammai+gammaj).E();
	302
	303	//Fill histograms with no cuts applied.
	304	fhAnglePairNoCut->Fill(e,angle);
	305	fhInvMassPairNoCut->Fill(e,invmass);
	306
	307	//Cut on phig-phi meson
	308	if(CutPtPhi(ptg, phig, pt, phi)){
	309
	310	fhAnglePairCorrelationCut ->Fill(e,angle);
	311	fhInvMassPairCorrelationCut->Fill(e,invmass);
	312
	313	//Cut on the aperture of the pair
	314	if(IsAngleInWindow(angle,e)){
	315	fhAnglePairOpeningAngleCut ->Fill(e,angle);
	316	fhInvMassPairOpeningAngleCut->Fill(e,invmass);
	317	AliDebug(2,Form("Angle cut: pt %f, phi %f",pt,phi));
	318
	319	//Cut on the invariant mass of the pair
	320	if((invmass>fInvMassMinCut) && (invmass<fInvMassMaxCut)){
	321	fhInvMassPairAllCut ->Fill(e,invmass);
	322	fhAnglePairAllCut ->Fill(e,angle);
	323	goodpair = kTRUE;
	324	AliDebug(2,Form("IM cut: pt %f, phi %f",pt,phi));
	325	}//(invmass>0.125) && (invmass<0.145)
	326	}//Opening angle cut
	327	} // cut on pt and phi
	328
	329
	330	return goodpair;
	331
	332	}
333
334	//__________________________________________________________________
335	void AliNeutralMesonSelection::Print(const Option_t * opt) const
336	{
337
338	//Print some relevant parameters set for the analysis
339	if(! opt)
340	return;
341
342	Info("Print", "%s %s", GetName(), GetTitle() ) ;
343
344	printf("mass : %f \n", fM );
345	printf("Invariant mass limits : %f < m < %f \n", fInvMassMinCut , fInvMassMinCut );
346	printf("Angle selection param: \n");
347	printf("p0 : %f", fAngleMaxParam.At(0));
348	printf("p1 : %f", fAngleMaxParam.At(1));
349	printf("p2 : %f", fAngleMaxParam.At(2));
350	printf("p3 : %f", fAngleMaxParam.At(3));
351
352	printf("pT meson > %f\n", fMinPt) ;
353	printf("Phi gamma-meson < %f\n", fDeltaPhiMaxCut) ;
354	printf("Phi gamma-meson > %f\n", fDeltaPhiMinCut) ;
355	printf("pT meson / pT Gamma < %f\n", fRatioMaxCut) ;
356	printf("pT meson / pT Gamma > %f\n", fRatioMinCut) ;
357	printf("Keep Neutral Meson Histos = %d\n",fKeepNeutralMesonHistos);
358
359	}