[u/mrichter/AliRoot.git] / PWG4 / PartCorrBase / 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: AliNeutralMesonSelection.cxx 27413 2008-07-18 13:28:12Z gconesab $ */

//_________________________________________________________________________
// 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.
//-- Author: Gustavo Conesa (INFN-LNF)

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

ClassImp(AliNeutralMesonSelection)
  
  
//____________________________________________________________________________
  AliNeutralMesonSelection::AliNeutralMesonSelection() : 
    TObject(), fM(0),
    fInvMassMaxCut(0.), fInvMassMinCut(0.),
    fAngleMaxParam(),  fKeepNeutralMesonHistos(0), 
    fhAnglePairNoCut(0), fhAnglePairOpeningAngleCut(0), 
    fhAnglePairAllCut(0), 
    fhInvMassPairNoCut(0), fhInvMassPairOpeningAngleCut(0), 
    fhInvMassPairAllCut(0),
    fHistoNEBins(0),   fHistoEMax(0.),   fHistoEMin(0.),
    fHistoNPtBins(0),  fHistoPtMax(0.),  fHistoPtMin(0.),
    fHistoNAngleBins(0), fHistoAngleMax(0.), fHistoAngleMin(0.),
    fHistoNIMBins(0), fHistoIMMax(0.), fHistoIMMin(0.)
{
  //Default Ctor
  
  //Initialize parameters
  
  // kGammaHadron and kGammaJet 
  fAngleMaxParam.Set(4) ;
  fAngleMaxParam.Reset(0.);
  
  //Initialize parameters
  InitParameters();
}

//____________________________________________________________________________
AliNeutralMesonSelection::AliNeutralMesonSelection(const AliNeutralMesonSelection & g) :   
  TObject(), fM(g.fM),
  fInvMassMaxCut(g.fInvMassMaxCut), fInvMassMinCut(g.fInvMassMinCut),
  fAngleMaxParam(g.fAngleMaxParam),
  fKeepNeutralMesonHistos(g.fKeepNeutralMesonHistos),
  fhAnglePairNoCut(g. fhAnglePairNoCut), 
  fhAnglePairOpeningAngleCut(g. fhAnglePairOpeningAngleCut), 
  fhAnglePairAllCut(g. fhAnglePairAllCut), 
  fhInvMassPairNoCut(g.fhInvMassPairNoCut),  
  fhInvMassPairOpeningAngleCut(g.fhInvMassPairOpeningAngleCut), 
  fhInvMassPairAllCut(g.fhInvMassPairAllCut),
  fHistoNEBins(g.fHistoNEBins),   fHistoEMax(g.fHistoEMax),   fHistoEMin(g.fHistoEMin),
  fHistoNPtBins(g.fHistoNPtBins),   fHistoPtMax(g.fHistoPtMax),   fHistoPtMin(g.fHistoPtMin),
  fHistoNAngleBins(g.fHistoNAngleBins), fHistoAngleMax(g.fHistoAngleMax), fHistoAngleMin(g.fHistoAngleMin),
  fHistoNIMBins(g.fHistoNIMBins), fHistoIMMax(g.fHistoIMMax), fHistoIMMin(g.fHistoIMMin)
{
  // cpy ctor
}

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

  fM = source.fM ;
  fInvMassMaxCut = source.fInvMassMaxCut ; 
  fInvMassMinCut = source.fInvMassMinCut ;
  fAngleMaxParam = source.fAngleMaxParam ;
  fKeepNeutralMesonHistos = source.fKeepNeutralMesonHistos;
 
  fhAnglePairNoCut = source. fhAnglePairNoCut ; 
  fhAnglePairOpeningAngleCut = source. fhAnglePairOpeningAngleCut ; 
  fhAnglePairAllCut = source. fhAnglePairAllCut ; 
  fhInvMassPairNoCut = source.fhInvMassPairNoCut ; 
  fhInvMassPairOpeningAngleCut = source.fhInvMassPairOpeningAngleCut ; 
  fhInvMassPairAllCut = source.fhInvMassPairAllCut ; 
  
  fHistoNEBins = source.fHistoNEBins;   fHistoEMax = source.fHistoEMax;   fHistoEMin = source.fHistoEMin;
  fHistoNPtBins = source.fHistoNPtBins;   fHistoPtMax = source.fHistoPtMax;   fHistoPtMin = source.fHistoPtMin;
  fHistoNAngleBins = source.fHistoNAngleBins; fHistoAngleMax = source.fHistoAngleMax; fHistoAngleMin = source.fHistoAngleMin;
  fHistoNIMBins = source.fHistoNIMBins; fHistoIMMax = source.fHistoIMMax; fHistoIMMin = source.fHistoIMMin;
  
  return *this;
  
}

//____________________________________________________________________________
AliNeutralMesonSelection::~AliNeutralMesonSelection() 
{
  //dtor

  if(!fKeepNeutralMesonHistos){
    //Histograms initialized and filled but not passed to output container
    //delete here, I am not sure this is correct
    
    if(fhAnglePairNoCut) delete fhAnglePairNoCut;
    if(fhAnglePairOpeningAngleCut) delete fhAnglePairOpeningAngleCut; 
    if(fhAnglePairAllCut) delete fhAnglePairAllCut;
    if(fhInvMassPairNoCut) delete fhInvMassPairNoCut;
    if(fhInvMassPairOpeningAngleCut) delete fhInvMassPairOpeningAngleCut;
    if(fhInvMassPairAllCut) delete fhInvMassPairAllCut; 

  }
  
}
//________________________________________________________________________
TList *  AliNeutralMesonSelection::GetCreateOutputObjects()
{  
  // Create histograms to be saved in output file and 
  // store them in outputContainer of the analysis class that calls this class.
  
  TList * outputContainer = new TList() ; 
  outputContainer->SetName("MesonDecayHistos") ; 
  
  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;
}

//____________________________________________________________________________
void AliNeutralMesonSelection::InitParameters()
{
  
  //Initialize the parameters of the analysis.
  fKeepNeutralMesonHistos = kFALSE ;
  
  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, pi0
  
 //Histogrammes settings
  fHistoNEBins = 100 ;
  fHistoEMax   = 50 ;
  fHistoEMin   = 0.  ;  
  
  fHistoNPtBins = 240 ;
  fHistoPtMax   = 120 ;
  fHistoPtMin   = 0.  ;

  fHistoNAngleBins = 200 ;
  fHistoAngleMax   = 0.2;
  fHistoAngleMin   = 0.  ;

  fHistoNIMBins = 300 ;
  fHistoIMMax   = 0.5   ;
  fHistoIMMin   = 0.  ;  
}

//__________________________________________________________________________-
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::SelectPair(TLorentzVector gammai, TLorentzVector gammaj)  
{  
  
  //Search for the neutral pion within selection cuts
  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 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
  
  return goodpair; 
  
}

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

  //Print some relevant parameters set for the analysis
  if(! opt)
    return;
  
  printf("**** Print %s %s ****\n", 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\n", fAngleMaxParam.At(0));
  printf("p1 :     %f\n", fAngleMaxParam.At(1));
  printf("p2 :     %f\n", fAngleMaxParam.At(2));
  printf("p3 :     %f\n", fAngleMaxParam.At(3));

  printf("Keep Neutral Meson Histos = %d\n",fKeepNeutralMesonHistos);
  
  if(fKeepNeutralMesonHistos){
    printf("Histograms: %3.1f < E  < %3.1f,  Nbin = %d\n", fHistoEMin,  fHistoEMax,  fHistoNEBins);
    printf("Histograms: %3.1f < pT < %3.1f,  Nbin = %d\n", fHistoPtMin,  fHistoPtMax,  fHistoNPtBins);
    printf("Histograms: %3.1f < angle < %3.1f, Nbin = %d\n", fHistoAngleMin, fHistoAngleMax, fHistoNAngleBins);
    printf("Histograms: %3.1f < IM < %3.1f, Nbin = %d\n", fHistoIMMin, fHistoIMMax, fHistoNIMBins);
    
  }
  
}
Commit	Line	Data
1c5acb87	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: AliNeutralMesonSelection.cxx 27413 2008-07-18 13:28:12Z gconesab $ */
	16
	17	//_________________________________________________________________________
	18	// Class that contains methods to select candidate pairs to neutral meson
	19	// 2 main selections, invariant mass around pi0 (also any other mass),
	20	// apperture angle to distinguish from combinatorial.
	21	//-- Author: Gustavo Conesa (INFN-LNF)
	22
	23	// --- ROOT system ---
	24	#include <TLorentzVector.h>
	25	#include <TH2.h>
	26	#include <TList.h>
1c5acb87	27
	28	//---- AliRoot system ----
	29	#include "AliNeutralMesonSelection.h"
1c5acb87	30
	31	ClassImp(AliNeutralMesonSelection)
	32
	33
	34	//____________________________________________________________________________
	35	AliNeutralMesonSelection::AliNeutralMesonSelection() :
	36	TObject(), fM(0),
	37	fInvMassMaxCut(0.), fInvMassMinCut(0.),
	38	fAngleMaxParam(), fKeepNeutralMesonHistos(0),
	39	fhAnglePairNoCut(0), fhAnglePairOpeningAngleCut(0),
	40	fhAnglePairAllCut(0),
	41	fhInvMassPairNoCut(0), fhInvMassPairOpeningAngleCut(0),
	42	fhInvMassPairAllCut(0),
477d6cee	43	fHistoNEBins(0), fHistoEMax(0.), fHistoEMin(0.),
	44	fHistoNPtBins(0), fHistoPtMax(0.), fHistoPtMin(0.),
	45	fHistoNAngleBins(0), fHistoAngleMax(0.), fHistoAngleMin(0.),
	46	fHistoNIMBins(0), fHistoIMMax(0.), fHistoIMMin(0.)
	47	{
1c5acb87	48	//Default Ctor
	49
	50	//Initialize parameters
	51
	52	// kGammaHadron and kGammaJet
	53	fAngleMaxParam.Set(4) ;
	54	fAngleMaxParam.Reset(0.);
	55
	56	//Initialize parameters
	57	InitParameters();
	58	}
	59
	60	//____________________________________________________________________________
	61	AliNeutralMesonSelection::AliNeutralMesonSelection(const AliNeutralMesonSelection & g) :
	62	TObject(), fM(g.fM),
	63	fInvMassMaxCut(g.fInvMassMaxCut), fInvMassMinCut(g.fInvMassMinCut),
	64	fAngleMaxParam(g.fAngleMaxParam),
	65	fKeepNeutralMesonHistos(g.fKeepNeutralMesonHistos),
	66	fhAnglePairNoCut(g. fhAnglePairNoCut),
	67	fhAnglePairOpeningAngleCut(g. fhAnglePairOpeningAngleCut),
	68	fhAnglePairAllCut(g. fhAnglePairAllCut),
	69	fhInvMassPairNoCut(g.fhInvMassPairNoCut),
	70	fhInvMassPairOpeningAngleCut(g.fhInvMassPairOpeningAngleCut),
	71	fhInvMassPairAllCut(g.fhInvMassPairAllCut),
	72	fHistoNEBins(g.fHistoNEBins), fHistoEMax(g.fHistoEMax), fHistoEMin(g.fHistoEMin),
	73	fHistoNPtBins(g.fHistoNPtBins), fHistoPtMax(g.fHistoPtMax), fHistoPtMin(g.fHistoPtMin),
	74	fHistoNAngleBins(g.fHistoNAngleBins), fHistoAngleMax(g.fHistoAngleMax), fHistoAngleMin(g.fHistoAngleMin),
	75	fHistoNIMBins(g.fHistoNIMBins), fHistoIMMax(g.fHistoIMMax), fHistoIMMin(g.fHistoIMMin)
	76	{
	77	// cpy ctor
	78	}
	79
	80	//_________________________________________________________________________
	81	AliNeutralMesonSelection & AliNeutralMesonSelection::operator = (const AliNeutralMesonSelection & source)
	82	{
	83	// assignment operator
	84
	85	if(this == &source)return *this;
	86	((TObject *)this)->operator=(source);
	87
	88	fM = source.fM ;
	89	fInvMassMaxCut = source.fInvMassMaxCut ;
	90	fInvMassMinCut = source.fInvMassMinCut ;
	91	fAngleMaxParam = source.fAngleMaxParam ;
	92	fKeepNeutralMesonHistos = source.fKeepNeutralMesonHistos;
	93
	94	fhAnglePairNoCut = source. fhAnglePairNoCut ;
	95	fhAnglePairOpeningAngleCut = source. fhAnglePairOpeningAngleCut ;
	96	fhAnglePairAllCut = source. fhAnglePairAllCut ;
	97	fhInvMassPairNoCut = source.fhInvMassPairNoCut ;
	98	fhInvMassPairOpeningAngleCut = source.fhInvMassPairOpeningAngleCut ;
	99	fhInvMassPairAllCut = source.fhInvMassPairAllCut ;
	100
	101	fHistoNEBins = source.fHistoNEBins; fHistoEMax = source.fHistoEMax; fHistoEMin = source.fHistoEMin;
	102	fHistoNPtBins = source.fHistoNPtBins; fHistoPtMax = source.fHistoPtMax; fHistoPtMin = source.fHistoPtMin;
	103	fHistoNAngleBins = source.fHistoNAngleBins; fHistoAngleMax = source.fHistoAngleMax; fHistoAngleMin = source.fHistoAngleMin;
	104	fHistoNIMBins = source.fHistoNIMBins; fHistoIMMax = source.fHistoIMMax; fHistoIMMin = source.fHistoIMMin;
	105
	106	return *this;
	107
	108	}
	109
	110	//____________________________________________________________________________
	111	AliNeutralMesonSelection::~AliNeutralMesonSelection()
112	{
113	//dtor
114
115	if(!fKeepNeutralMesonHistos){
116	//Histograms initialized and filled but not passed to output container
117	//delete here, I am not sure this is correct
118
119	if(fhAnglePairNoCut) delete fhAnglePairNoCut;
120	if(fhAnglePairOpeningAngleCut) delete fhAnglePairOpeningAngleCut;
121	if(fhAnglePairAllCut) delete fhAnglePairAllCut;
122	if(fhInvMassPairNoCut) delete fhInvMassPairNoCut;
123	if(fhInvMassPairOpeningAngleCut) delete fhInvMassPairOpeningAngleCut;
124	if(fhInvMassPairAllCut) delete fhInvMassPairAllCut;
125
126	}
127
128	}
129	//________________________________________________________________________
130	TList * AliNeutralMesonSelection::GetCreateOutputObjects()
131	{
477d6cee	132	// Create histograms to be saved in output file and
	133	// store them in outputContainer of the analysis class that calls this class.
	134
	135	TList * outputContainer = new TList() ;
	136	outputContainer->SetName("MesonDecayHistos") ;
	137
	138	fhAnglePairNoCut = new TH2F
1c5acb87	139	("AnglePairNoCut",
1c5acb87	140	"Angle between all #gamma pair vs E_{#pi^{0}}",fHistoNEBins,fHistoEMin,fHistoEMax,fHistoNAngleBins,fHistoAngleMin,fHistoAngleMax);
477d6cee	141	fhAnglePairNoCut->SetYTitle("Angle (rad)");
	142	fhAnglePairNoCut->SetXTitle("E_{ #pi^{0}} (GeV)");
	143
	144	fhAnglePairOpeningAngleCut = new TH2F
1c5acb87	145	("AnglePairOpeningAngleCut",
	146	"Angle between all #gamma pair (opening angle + azimuth cut) vs E_{#pi^{0}}"
	147	,fHistoNEBins,fHistoEMin,fHistoEMax,fHistoNAngleBins,fHistoAngleMin,fHistoAngleMax);
477d6cee	148	fhAnglePairOpeningAngleCut->SetYTitle("Angle (rad)");
	149	fhAnglePairOpeningAngleCut->SetXTitle("E_{ #pi^{0}} (GeV)");
	150
	151	fhAnglePairAllCut = new TH2F
1c5acb87	152	("AnglePairAllCut",
	153	"Angle between all #gamma pair (opening angle + inv mass cut+azimuth) vs E_{#pi^{0}}"
	154	,fHistoNEBins,fHistoEMin,fHistoEMax,fHistoNAngleBins,fHistoAngleMin,fHistoAngleMax);
477d6cee	155	fhAnglePairAllCut->SetYTitle("Angle (rad)");
	156	fhAnglePairAllCut->SetXTitle("E_{ #pi^{0}} (GeV)");
	157
	158	//
	159	fhInvMassPairNoCut = new TH2F
1c5acb87	160	("InvMassPairNoCut","Invariant Mass of all #gamma pair vs E_{#pi^{0}}",
1c5acb87	161	fHistoNPtBins,fHistoPtMin,fHistoPtMax,fHistoNIMBins,fHistoIMMin,fHistoIMMax);
477d6cee	162	fhInvMassPairNoCut->SetYTitle("Invariant Mass (GeV/c^{2})");
	163	fhInvMassPairNoCut->SetXTitle("E_{ #pi^{0}} (GeV)");
	164
	165	fhInvMassPairOpeningAngleCut = new TH2F
1c5acb87	166	("InvMassPairOpeningAngleCut",
	167	"Invariant Mass of #gamma pair (angle cut) vs E_{#pi^{0}}",
	168	fHistoNPtBins,fHistoPtMin,fHistoPtMax,fHistoNIMBins,fHistoIMMin,fHistoIMMax);
477d6cee	169	fhInvMassPairOpeningAngleCut->SetYTitle("Invariant Mass (GeV/c^{2})");
	170	fhInvMassPairOpeningAngleCut->SetXTitle(" E_{#pi^{0}}(GeV)");
	171
	172	fhInvMassPairAllCut = new TH2F
1c5acb87	173	("InvMassPairAllCut",
	174	"Invariant Mass of #gamma pair (opening angle+invmass cut) vs E_{#pi^{0}}",
	175	fHistoNPtBins,fHistoPtMin,fHistoPtMax,fHistoNIMBins,fHistoIMMin,fHistoIMMax);
477d6cee	176	fhInvMassPairAllCut->SetYTitle("Invariant Mass (GeV/c^{2})");
	177	fhInvMassPairAllCut->SetXTitle("E_{#pi^{0}}(GeV)");
	178
	179	outputContainer->Add(fhAnglePairNoCut) ;
	180	outputContainer->Add(fhAnglePairOpeningAngleCut) ;
	181	outputContainer->Add(fhAnglePairAllCut) ;
	182
	183	outputContainer->Add(fhInvMassPairNoCut) ;
	184	outputContainer->Add(fhInvMassPairOpeningAngleCut) ;
	185	outputContainer->Add(fhInvMassPairAllCut) ;
	186
	187	return outputContainer;
1c5acb87	188	}
1c5acb87	189
477d6cee	190	//____________________________________________________________________________
1c5acb87	191	void AliNeutralMesonSelection::InitParameters()
1c5acb87	192	{
477d6cee	193
1c5acb87	194	//Initialize the parameters of the analysis.
1c5acb87	195	fKeepNeutralMesonHistos = kFALSE ;
477d6cee	196
1c5acb87	197	fAngleMaxParam.Set(4) ;
	198	fAngleMaxParam.AddAt(0.4,0);//={0.4,-0.25,0.025,-2e-4};
	199	fAngleMaxParam.AddAt(-0.25,1) ;
	200	fAngleMaxParam.AddAt(0.025,2) ;
	201	fAngleMaxParam.AddAt(-2e-4,3) ;
	202
	203	fInvMassMaxCut = 0.16 ;
	204	fInvMassMinCut = 0.11 ;
	205
	206	fM = 0.1349766;//neutralMeson mass, pi0
	207
	208	//Histogrammes settings
	209	fHistoNEBins = 100 ;
	210	fHistoEMax = 50 ;
	211	fHistoEMin = 0. ;
	212
	213	fHistoNPtBins = 240 ;
	214	fHistoPtMax = 120 ;
	215	fHistoPtMin = 0. ;
	216
	217	fHistoNAngleBins = 200 ;
	218	fHistoAngleMax = 0.2;
	219	fHistoAngleMin = 0. ;
	220
	221	fHistoNIMBins = 300 ;
	222	fHistoIMMax = 0.5 ;
	223	fHistoIMMin = 0. ;
	224	}
	225
	226	//__________________________________________________________________________-
	227	Bool_t AliNeutralMesonSelection::IsAngleInWindow(const Float_t angle,const Float_t e) const {
	228	//Check if the opening angle of the candidate pairs is inside
	229	//our selection windowd
	230
	231	Bool_t result = kFALSE;
	232	Double_t max = fAngleMaxParam.At(0)TMath::Exp(fAngleMaxParam.At(1)e)
	233	+fAngleMaxParam.At(2)+fAngleMaxParam.At(3)*e;
	234	Double_t arg = (ee-2fMfM)/(ee);
	235	Double_t min = 100. ;
	236	if(arg>0.)
	237	min = TMath::ACos(arg);
	238
	239	if((angle<max)&&(angle>=min))
	240	result = kTRUE;
	241
	242	return result;
	243	}
	244
	245	//____________________________________________________________________________
	246	Bool_t AliNeutralMesonSelection::SelectPair(TLorentzVector gammai, TLorentzVector gammaj)
	247	{
	248
	249	//Search for the neutral pion within selection cuts
	250	Bool_t goodpair = kFALSE ;
	251
477d6cee	252	// Double_t pt = (gammai+gammaj).Pt();
1c5acb87	253	Double_t phi = (gammai+gammaj).Phi();
	254	if(phi < 0)
	255	phi+=TMath::TwoPi();
	256	Double_t invmass = (gammai+gammaj).M();
	257	Double_t angle = gammaj.Angle(gammai.Vect());
	258	Double_t e = (gammai+gammaj).E();
	259
	260	//Fill histograms with no cuts applied.
	261	fhAnglePairNoCut->Fill(e,angle);
	262	fhInvMassPairNoCut->Fill(e,invmass);
	263
	264	//Cut on the aperture of the pair
	265	if(IsAngleInWindow(angle,e)){
	266	fhAnglePairOpeningAngleCut ->Fill(e,angle);
	267	fhInvMassPairOpeningAngleCut->Fill(e,invmass);
477d6cee	268	//AliDebug(2,Form("Angle cut: pt %f, phi %f",pt,phi));
1c5acb87	269
	270	//Cut on the invariant mass of the pair
	271	if((invmass>fInvMassMinCut) && (invmass<fInvMassMaxCut)){
	272	fhInvMassPairAllCut ->Fill(e,invmass);
	273	fhAnglePairAllCut ->Fill(e,angle);
	274	goodpair = kTRUE;
477d6cee	275	//AliDebug(2,Form("IM cut: pt %f, phi %f",pt,phi));
1c5acb87	276	}//(invmass>0.125) && (invmass<0.145)
	277	}//Opening angle cut
	278
	279	return goodpair;
	280
	281	}
	282
	283	//__________________________________________________________________
	284	void AliNeutralMesonSelection::Print(const Option_t * opt) const
	285	{
	286
	287	//Print some relevant parameters set for the analysis
	288	if(! opt)
	289	return;
	290
a3aebfff	291	printf("** Print %s %s **\n", GetName(), GetTitle() ) ;
1c5acb87	292
	293	printf("mass : %f \n", fM );
	294	printf("Invariant mass limits : %f < m < %f \n", fInvMassMinCut , fInvMassMinCut );
	295	printf("Angle selection param: \n");
fbeaf916	296	printf("p0 : %f\n", fAngleMaxParam.At(0));
	297	printf("p1 : %f\n", fAngleMaxParam.At(1));
	298	printf("p2 : %f\n", fAngleMaxParam.At(2));
	299	printf("p3 : %f\n", fAngleMaxParam.At(3));
1c5acb87	300
	301	printf("Keep Neutral Meson Histos = %d\n",fKeepNeutralMesonHistos);
	302
	303	if(fKeepNeutralMesonHistos){
477d6cee	304	printf("Histograms: %3.1f < E < %3.1f, Nbin = %d\n", fHistoEMin, fHistoEMax, fHistoNEBins);
	305	printf("Histograms: %3.1f < pT < %3.1f, Nbin = %d\n", fHistoPtMin, fHistoPtMax, fHistoNPtBins);
	306	printf("Histograms: %3.1f < angle < %3.1f, Nbin = %d\n", fHistoAngleMin, fHistoAngleMax, fHistoNAngleBins);
	307	printf("Histograms: %3.1f < IM < %3.1f, Nbin = %d\n", fHistoIMMin, fHistoIMMax, fHistoNIMBins);
	308
1c5acb87	309	}
477d6cee	310
1c5acb87	311	}