/* History of cvs commits:
*
* $Log$
- * Revision 1.3 2007/03/08 10:24:32 schutz
- * Coding convention
- *
- * Revision 1.2 2007/02/09 18:40:40 schutz
- * B\bNew version from Gustavo
- *
- * Revision 1.1 2007/01/23 17:17:29 schutz
- * New Gamma package
+ * Revision 1.4.4.4 2007/07/26 10:32:09 schutz
+ * new analysis classes in the the new analysis framework
*
*
*/
//_________________________________________________________________________
-// Class for the analysis of gamma correlations (gamma-jet,
-// gamma-hadron and isolation cut.
-// This class contains 3 main methods: one to fill lists of particles (ESDs) comming
-// from the CTS (ITS+TPC) and the calorimeters; the other one tags a candidate
-// cluster as isolated; the last one search in the
-// corresponing calorimeter for the highest energy cluster, identified it as
-// prompt photon;
+// Class for the prompt gamma analysis, isolation cut
//
// Class created from old AliPHOSGammaJet
// (see AliRoot versions previous Release 4-09)
//
//*-- Author: Gustavo Conesa (LNF-INFN)
//////////////////////////////////////////////////////////////////////////////
-
-
-// --- ROOT system ---
-
-#include <TFile.h>
+
+
+// --- ROOT system ---
#include <TParticle.h>
#include <TH2.h>
-#include <TChain.h>
+#include <TList.h>
#include "AliAnaGammaDirect.h"
-#include "AliESD.h"
-#include "AliESDtrack.h"
-#include "AliESDCaloCluster.h"
#include "Riostream.h"
#include "AliLog.h"
-
+
ClassImp(AliAnaGammaDirect)
-
+
//____________________________________________________________________________
- AliAnaGammaDirect::AliAnaGammaDirect(const char *name) :
- AliAnalysisTask(name,""), fChain(0), fESD(0),
- fOutputContainer(new TObjArray(100)),
- fPrintInfo(0), fMinGammaPt(0.),
- fCalorimeter(""), fEMCALPID(0),fPHOSPID(0),
- fEMCALPhotonWeight(0.), fEMCALPi0Weight(0.), fPHOSPhotonWeight(0.),
+ AliAnaGammaDirect::AliAnaGammaDirect() :
+ TObject(),
+ fMinGammaPt(0.),
fConeSize(0.),fPtThreshold(0.),fPtSumThreshold(0),
- fMakeICMethod(0),fhNGamma(0),fhPhiGamma(0),fhEtaGamma(0)
+ fICMethod(0),fhNGamma(0),fhPhiGamma(0),fhEtaGamma(0),
+ //kSeveralIC
+ fNCones(0),fNPtThres(0), fConeSizes(), fPtThresholds(),
+ fhPtThresIsolated(), fhPtSumIsolated()
{
- //Ctor
-
+ //default ctor
+
//Initialize parameters
InitParameters();
- // Input slot #0 works with an Ntuple
- DefineInput(0, TChain::Class());
- // Output slot #0 writes into a TH1 container
- DefineOutput(0, TObjArray::Class()) ;
-
}
-
//____________________________________________________________________________
AliAnaGammaDirect::AliAnaGammaDirect(const AliAnaGammaDirect & g) :
- AliAnalysisTask(g), fChain(g.fChain), fESD(g.fESD),
- fOutputContainer(g. fOutputContainer), fPrintInfo(g.fPrintInfo),
- fMinGammaPt(g.fMinGammaPt), fCalorimeter(g.fCalorimeter),
- fEMCALPID(g.fEMCALPID),fPHOSPID(g.fPHOSPID),
- fEMCALPhotonWeight(g.fEMCALPhotonWeight),
- fEMCALPi0Weight(g.fEMCALPi0Weight),
- fPHOSPhotonWeight(g.fPHOSPhotonWeight),
+ TObject(g),
+ fMinGammaPt(g.fMinGammaPt),
fConeSize(g.fConeSize),
fPtThreshold(g.fPtThreshold),
fPtSumThreshold(g.fPtSumThreshold),
- fMakeICMethod(g.fMakeICMethod),
- fhNGamma(g.fhNGamma),fhPhiGamma(g.fhPhiGamma),fhEtaGamma(g.fhEtaGamma)
+ fICMethod(g.fICMethod),
+ fhNGamma(g.fhNGamma),fhPhiGamma(g.fhPhiGamma),fhEtaGamma(g.fhEtaGamma),
+ //kSeveralIC
+ fNCones(g.fNCones),fNPtThres(g.fNPtThres), fConeSizes(),fPtThresholds(),
+ fhPtThresIsolated(), fhPtSumIsolated()
{
// cpy ctor
- SetName (g.GetName()) ;
- SetTitle(g.GetTitle()) ;
+
+ //kSeveralIC
+ for(Int_t i = 0; i < fNCones ; i++){
+ fConeSizes[i] = g.fConeSizes[i];
+ fhPtSumIsolated[i] = g.fhPtSumIsolated[i];
+ for(Int_t j = 0; j < fNPtThres ; j++)
+ fhPtThresIsolated[i][j] = g.fhPtThresIsolated[i][j];
+ }
+
+ for(Int_t i = 0; i < fNPtThres ; i++)
+ fPtThresholds[i]= g.fPtThresholds[i];
}
//_________________________________________________________________________
AliAnaGammaDirect & AliAnaGammaDirect::operator = (const AliAnaGammaDirect & source)
{
// assignment operator
-
+
if(&source == this) return *this;
- fChain = source.fChain ;
- fESD = source.fESD ;
- fOutputContainer = source.fOutputContainer ;
- fPrintInfo = source.fPrintInfo ;
fMinGammaPt = source.fMinGammaPt ;
- fCalorimeter = source. fCalorimeter ;
- fEMCALPID = source.fEMCALPID ;
- fPHOSPID = source.fPHOSPID ;
- fEMCALPhotonWeight = source. fEMCALPhotonWeight ;
- fEMCALPi0Weight = source.fEMCALPi0Weight ;
- fPHOSPhotonWeight = source.fPHOSPhotonWeight ;
fConeSize = source.fConeSize ;
fPtThreshold = source.fPtThreshold ;
fPtSumThreshold = source.fPtSumThreshold ;
- fMakeICMethod = source.fMakeICMethod ;
+ fICMethod = source.fICMethod ;
fhNGamma = source.fhNGamma ;
fhPhiGamma = source.fhPhiGamma ;
fhEtaGamma = source.fhEtaGamma ;
-
-
+
+ //kSeveralIC
+ fNCones = source.fNCones ;
+ fNPtThres = source.fNPtThres ;
+
+ for(Int_t i = 0; i < fNCones ; i++){
+ fConeSizes[i] = source.fConeSizes[i];
+ fhPtSumIsolated[i] = source.fhPtSumIsolated[i] ;
+ for(Int_t j = 0; j < fNPtThres ; j++)
+ fhPtThresIsolated[i][j] = source.fhPtThresIsolated[i][j] ;
+ }
+
+ for(Int_t i = 0; i < fNPtThres ; i++)
+ fPtThresholds[i]= source.fPtThresholds[i];
+
return *this;
-
+
}
//____________________________________________________________________________
AliAnaGammaDirect::~AliAnaGammaDirect()
{
// Remove all pointers
- fOutputContainer->Clear() ;
- delete fOutputContainer ;
+
delete fhNGamma ;
delete fhPhiGamma ;
delete fhEtaGamma ;
-
-}
-
-//______________________________________________________________________________
-void AliAnaGammaDirect::ConnectInputData(const Option_t*)
-{
- // Initialisation of branch container and histograms
-
- AliInfo(Form("*** Initialization of %s", GetName())) ;
- // Get input data
- fChain = dynamic_cast<TChain *>(GetInputData(0)) ;
- if (!fChain) {
- AliError(Form("Input 0 for %s not found\n", GetName()));
- return ;
- }
+ //kSeveralIC
+ delete [] fhPtThresIsolated ;
+ delete [] fhPtSumIsolated ;
- // One should first check if the branch address was taken by some other task
- char ** address = (char **)GetBranchAddress(0, "ESD");
- if (address) {
- fESD = (AliESD*)(*address);
- } else {
- fESD = new AliESD();
- SetBranchAddress(0, "ESD", &fESD);
- }
-
}
//________________________________________________________________________
-void AliAnaGammaDirect::CreateOutputObjects()
+TList * AliAnaGammaDirect::GetCreateOutputObjects()
{
// Create histograms to be saved in output file and
- // store them in fOutputContainer
-
- fOutputContainer = new TObjArray(3) ;
-
+ // store them in outputContainer
+ TList * outputContainer = new TList() ;
+ outputContainer->SetName("DirectGammaHistos") ;
+
//Histograms of highest gamma identified in Event
- fhNGamma = new TH1F("NGamma","Number of #gamma over PHOS",240,0,120);
+ fhNGamma = new TH1F("NGamma","Number of #gamma over calorimeter",240,0,120);
fhNGamma->SetYTitle("N");
fhNGamma->SetXTitle("p_{T #gamma}(GeV/c)");
- fOutputContainer->Add(fhNGamma) ;
+ outputContainer->Add(fhNGamma) ;
fhPhiGamma = new TH2F
("PhiGamma","#phi_{#gamma}",200,0,120,200,0,7);
fhPhiGamma->SetYTitle("#phi");
fhPhiGamma->SetXTitle("p_{T #gamma} (GeV/c)");
- fOutputContainer->Add(fhPhiGamma) ;
+ outputContainer->Add(fhPhiGamma) ;
fhEtaGamma = new TH2F
("EtaGamma","#phi_{#gamma}",200,0,120,200,-0.8,0.8);
fhEtaGamma->SetYTitle("#eta");
fhEtaGamma->SetXTitle("p_{T #gamma} (GeV/c)");
- fOutputContainer->Add(fhEtaGamma) ;
-
-}
-
-//____________________________________________________________________________
-void AliAnaGammaDirect::CreateParticleList(TClonesArray * pl,
- TClonesArray * plCTS,
- TClonesArray * plEMCAL,
- TClonesArray * plPHOS){
-
- //Create a list of particles from the ESD. These particles have been measured
- //by the Central Tracking system (TPC+ITS), PHOS and EMCAL
-
- Int_t index = pl->GetEntries() ;
- Int_t npar = 0 ;
- Float_t *pid = new Float_t[AliPID::kSPECIESN];
- AliDebug(3,"Fill particle lists");
- if(fPrintInfo)
- AliInfo(Form("fCalorimeter %s",fCalorimeter.Data()));
-
- Double_t v[3] ; //vertex ;
- fESD->GetVertex()->GetXYZ(v) ;
-
- //########### PHOS ##############
- if(fCalorimeter == "PHOS"){
-
- Int_t begphos = fESD->GetFirstPHOSCluster();
- Int_t endphos = fESD->GetFirstPHOSCluster() +
- fESD->GetNumberOfPHOSClusters() ;
- Int_t indexNePHOS = plPHOS->GetEntries() ;
- AliDebug(3,Form("First PHOS particle %d, last particle %d", begphos,endphos));
-
- if(fCalorimeter == "PHOS"){
- for (npar = begphos; npar < endphos; npar++) {//////////////PHOS track loop
- AliESDCaloCluster * clus = fESD->GetCaloCluster(npar) ; // retrieve track from esd
-
- //Create a TParticle to fill the particle list
- TLorentzVector momentum ;
- clus->GetMomentum(momentum);
- TParticle * particle = new TParticle() ;
- //particle->SetMomentum(px,py,pz,en) ;
- particle->SetMomentum(momentum) ;
-
- AliDebug(4,Form("PHOS clusters: pt %f, phi %f, eta %f", particle->Pt(),particle->Phi(),particle->Eta()));
-
- //Select only photons
-
- pid=clus->GetPid();
- //cout<<"pid "<<pid[AliPID::kPhoton]<<endl ;
- if( !fPHOSPID)
- new((*plPHOS)[indexNePHOS++]) TParticle(*particle) ;
- else if( pid[AliPID::kPhoton] > fPHOSPhotonWeight)
- new((*plPHOS)[indexNePHOS++]) TParticle(*particle) ;
- }
- }
+ outputContainer->Add(fhEtaGamma) ;
+
+ if(fICMethod == kSeveralIC){
+ char name[128];
+ char title[128];
+ for(Int_t icone = 0; icone<fNCones; icone++){
+ sprintf(name,"PtSumIsolated_Cone_%d",icone);
+ sprintf(title,"Candidate cone sum p_{T} for cone size %d vs candidate p_{T}",icone);
+ fhPtSumIsolated[icone] = new TH2F(name, title,240,0,120,120,0,10);
+ fhPtSumIsolated[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
+ fhPtSumIsolated[icone]->SetXTitle("p_{T} (GeV/c)");
+ outputContainer->Add(fhPtSumIsolated[icone]) ;
+
+ for(Int_t ipt = 0; ipt<fNPtThres;ipt++){
+ sprintf(name,"PtThresIsol_Cone_%d_Pt%d",icone,ipt);
+ sprintf(title,"Isolated candidate p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
+ fhPtThresIsolated[icone][ipt] = new TH1F(name, title,240,0,120);
+ fhPtThresIsolated[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
+ outputContainer->Add(fhPtThresIsolated[icone][ipt]) ;
+ }//icone loop
+ }//ipt loop
}
- //########## #####################
- //Prepare bool array for EMCAL track matching
-
- // step 1, set the flag in a way that it rejects all not-V1 clusters
- // but at this point all V1 clusters are accepted
-
- Int_t begem = fESD->GetFirstEMCALCluster();
- Int_t endem = fESD->GetFirstEMCALCluster() +
- fESD->GetNumberOfEMCALClusters() ;
-
-// if(endem < begem+12)
-// AliError("Number of pseudoclusters smaller than 12");
- Bool_t *useCluster = new Bool_t[endem+1];
-
-// for (npar = 0; npar < endem; npar++){
-// if(npar < begem+12)
-// useCluster[npar] =kFALSE; //EMCAL Pseudoclusters and PHOS clusters
-// else
-// useCluster[npar] =kTRUE; //EMCAL clusters
-// }
- for (npar = 0; npar < endem; npar++)
- useCluster[npar] =kFALSE; //EMCAL Pseudoclusters and clusters
-
- //########### CTS (TPC+ITS) #####################
- Int_t begtpc = 0 ;
- Int_t endtpc = fESD->GetNumberOfTracks() ;
- Int_t indexCh = plCTS->GetEntries() ;
- AliDebug(3,Form("First CTS particle %d, last particle %d", begtpc,endtpc));
-
- Int_t iemcalMatch = -1 ;
- for (npar = begtpc; npar < endtpc; npar++) {////////////// track loop
- AliESDtrack * track = fESD->GetTrack(npar) ; // retrieve track from esd
+ return outputContainer ;
- // step 2 for EMCAL matching, change the flag for all matched clusters found in tracks
- iemcalMatch = track->GetEMCALcluster();
- if(iemcalMatch > 0) useCluster[iemcalMatch] = kTRUE; // reject matched cluster
-
- //We want tracks fitted in the detectors:
- ULong_t status=AliESDtrack::kTPCrefit;
- status|=AliESDtrack::kITSrefit;
-
- //We want tracks whose PID bit is set:
- // ULong_t status =AliESDtrack::kITSpid;
- // status|=AliESDtrack::kTPCpid;
-
- if ( (track->GetStatus() & status) == status) {//Check if the bits we want are set
- // Do something with the tracks which were successfully
- // re-fitted
- Double_t en = 0; //track ->GetTPCsignal() ;
- Double_t mom[3];
- track->GetPxPyPz(mom) ;
- Double_t px = mom[0];
- Double_t py = mom[1];
- Double_t pz = mom[2]; //Check with TPC people if this is correct.
- Int_t pdg = 11; //Give any charged PDG code, in this case electron.
- //I just want to tag the particle as charged
- TParticle * particle = new TParticle(pdg, 1, -1, -1, -1, -1,
- px, py, pz, en, v[0], v[1], v[2], 0);
-
- //TParticle * particle = new TParticle() ;
- //particle->SetMomentum(px,py,pz,en) ;
-
- new((*plCTS)[indexCh++]) TParticle(*particle) ;
- new((*pl)[index++]) TParticle(*particle) ;
- }
- }
-
- //################ EMCAL ##############
-
- Int_t indexNe = plEMCAL->GetEntries() ;
-
- AliDebug(3,Form("First EMCAL particle %d, last particle %d",begem,endem));
-
- for (npar = begem; npar < endem; npar++) {//////////////EMCAL track loop
- AliESDCaloCluster * clus = fESD->GetCaloCluster(npar) ; // retrieve track from esd
- Int_t clustertype= clus->GetClusterType();
- if(clustertype == AliESDCaloCluster::kClusterv1 && !useCluster[npar] ){
- TLorentzVector momentum ;
- clus->GetMomentum(momentum);
- TParticle * particle = new TParticle() ;
- //particle->SetMomentum(px,py,pz,en) ;
- particle->SetMomentum(momentum) ;
- cout<<"GOOD EMCAL "<<particle->Pt()<<endl;
- pid=clus->GetPid();
- if(fCalorimeter == "EMCAL")
- {
- TParticle * particle = new TParticle() ;
- //particle->SetMomentum(px,py,pz,en) ;
- AliDebug(4,Form("EMCAL clusters: pt %f, phi %f, eta %f", particle->Pt(),particle->Phi(),particle->Eta()));
- if(!fEMCALPID) //Only identified particles
- new((*plEMCAL)[indexNe++]) TParticle(*particle) ;
- else if(pid[AliPID::kPhoton] > fEMCALPhotonWeight)
- new((*plEMCAL)[indexNe++]) TParticle(*particle) ;
- }
- else
- {
- Int_t pdg = 0;
- if(fEMCALPID)
- {
- if( pid[AliPID::kPhoton] > fEMCALPhotonWeight)
- pdg = 22;
- else if( pid[AliPID::kPi0] > fEMCALPi0Weight)
- pdg = 111;
- }
- else
- pdg = 22; //No PID, assume all photons
-
- TParticle * particle = new TParticle(pdg, 1, -1, -1, -1, -1,
- momentum.Px(), momentum.Py(), momentum.Pz(), momentum.E(), v[0], v[1], v[2], 0);
- AliDebug(4,Form("EMCAL clusters: pt %f, phi %f, eta %f", particle->Pt(),particle->Phi(),particle->Eta()));
-
- new((*plEMCAL)[indexNe++]) TParticle(*particle) ;
- new((*pl)[index++]) TParticle(*particle) ;
- }
- }
- }
-
- AliDebug(3,"Particle lists filled");
-
}
-
-
//____________________________________________________________________________
-void AliAnaGammaDirect::Exec(Option_t *)
-{
-
- // Processing of one event
-
- //Get ESDs
- Long64_t entry = fChain->GetReadEntry() ;
-
- if (!fESD) {
- AliError("fESD is not connected to the input!") ;
- return ;
- }
-
- if (fPrintInfo)
- AliInfo(Form("%s ----> Processing event # %lld", (dynamic_cast<TChain *>(fChain))->GetFile()->GetName(), entry)) ;
-
- //CreateTLists with arrays of TParticles. Filled with particles only relevant for the analysis.
-
- TClonesArray * particleList = new TClonesArray("TParticle",1000); // All particles refitted in CTS and detected in EMCAL (jet)
- TClonesArray * plCTS = new TClonesArray("TParticle",1000); // All particles refitted in Central Tracking System (ITS+TPC)
- TClonesArray * plNe = new TClonesArray("TParticle",1000); // All particles measured in Jet Calorimeter (EMCAL)
- TClonesArray * plPHOS = new TClonesArray("TParticle",1000); // All particles measured in PHOS as Gamma calorimeter
- TClonesArray * plEMCAL = new TClonesArray("TParticle",1000); // All particles measured in EMCAL as Gamma calorimeter
-
- TParticle *pGamma = new TParticle(); //It will contain the kinematics of the found prompt gamma
-
- //Fill lists with photons, neutral particles and charged particles
- //look for the highest energy photon in the event inside fCalorimeter
-
- AliDebug(2, "Fill particle lists, get prompt gamma");
-
- //Fill particle lists
- CreateParticleList(particleList, plCTS,plEMCAL,plPHOS);
-
- if(fCalorimeter == "PHOS")
- plNe = plPHOS;
- if(fCalorimeter == "EMCAL")
- plNe = plEMCAL;
-
-
- Bool_t iIsInPHOSorEMCAL = kFALSE ; //To check if Gamma was in any calorimeter
- //Search highest energy prompt gamma in calorimeter
- GetPromptGamma(plNe, plCTS, pGamma, iIsInPHOSorEMCAL) ;
-
- AliDebug(1, Form("Is Gamma in %s? %d",fCalorimeter.Data(),iIsInPHOSorEMCAL));
-
- //If there is any photon candidate in fCalorimeter
- if(iIsInPHOSorEMCAL){
- if (fPrintInfo)
- AliInfo(Form("Prompt Gamma: pt %f, phi %f, eta %f", pGamma->Pt(),pGamma->Phi(),pGamma->Eta())) ;
-
- }//Gamma in Calo
-
- AliDebug(2, "End of analysis, delete pointers");
-
- particleList->Delete() ;
- plCTS->Delete() ;
- plNe->Delete() ;
- plEMCAL->Delete() ;
- plPHOS->Delete() ;
- pGamma->Delete();
-
- PostData(0, fOutputContainer);
-
- delete plNe ;
- delete plCTS ;
- //delete plPHOS ;
- //delete plEMCAL ;
- delete particleList ;
-
- // delete pGamma;
-
-}
-
-
-//____________________________________________________________________________
-void AliAnaGammaDirect::GetPromptGamma(TClonesArray * pl, TClonesArray * plCTS, TParticle *pGamma, Bool_t &Is) const
+void AliAnaGammaDirect::GetPromptGamma(TClonesArray * pl, TClonesArray * plCTS, TParticle *pGamma, Bool_t &found) const
{
//Search for the prompt photon in Calorimeter with pt > fMinGammaPt
index = ipr ;
pt = particle->Pt();
pGamma->SetMomentum(particle->Px(),particle->Py(),particle->Pz(),particle->Energy());
- AliDebug(4,Form("Cluster in calo: pt %f, phi %f, eta %f", pGamma->Pt(),pGamma->Phi(),pGamma->Eta())) ;
- Is = kTRUE;
+ found = kTRUE;
}
}
//Do Isolation?
- if(fMakeICMethod && Is)
+ if( ( fICMethod == kPtIC || fICMethod == kSumPtIC ) && found)
{
Float_t coneptsum = 0 ;
Bool_t icPtThres = kFALSE;
Bool_t icPtSum = kFALSE;
MakeIsolationCut(plCTS,pl, pGamma, index,
icPtThres, icPtSum,coneptsum);
- if(fMakeICMethod == 1) //Pt thres method
- Is = icPtThres ;
- if(fMakeICMethod == 2) //Pt cone sum method
- Is = icPtSum ;
+ if(fICMethod == kPtIC) //Pt thres method
+ found = icPtThres ;
+ if(fICMethod == kSumPtIC) //Pt cone sum method
+ found = icPtSum ;
}
- if(Is){
- AliDebug(3,Form("Cluster with p_{T} larger than %f found in calorimeter ", fMinGammaPt)) ;
- AliDebug(3,Form("Gamma: pt %f, phi %f, eta %f", pGamma->Pt(),pGamma->Phi(),pGamma->Eta())) ;
- //Fill prompt gamma histograms
- fhNGamma ->Fill(pGamma->Pt());
+ if(found){
+ AliDebug(1,Form("Cluster with p_{T} larger than %f found in calorimeter ", fMinGammaPt)) ;
+ AliDebug(1,Form("Gamma: pt %f, phi %f, eta %f", pGamma->Pt(),pGamma->Phi(),pGamma->Eta())) ;
+ //Fill prompt gamma histograms
+ fhNGamma->Fill(pGamma->Pt());
fhPhiGamma->Fill( pGamma->Pt(),pGamma->Phi());
fhEtaGamma->Fill(pGamma->Pt(),pGamma->Eta());
}
{
//Initialize the parameters of the analysis.
- fCalorimeter="PHOS";
- fPrintInfo = kTRUE;
fMinGammaPt = 5. ;
//Fill particle lists when PID is ok
- fEMCALPID = kFALSE;
- fPHOSPID = kFALSE;
- fEMCALPhotonWeight = 0.5 ;
- fEMCALPi0Weight = 0.5 ;
- fPHOSPhotonWeight = 0.8 ;
fConeSize = 0.2 ;
fPtThreshold = 2.0;
fPtSumThreshold = 1.;
- fMakeICMethod = 1; // 0 don't isolate, 1 pt thresh method, 2 cone pt sum method
+ fICMethod = kNoIC; // 0 don't isolate, 1 pt thresh method, 2 cone pt sum method
+
+ //-----------kSeveralIC-----------------
+ fNCones = 4 ;
+ fNPtThres = 4 ;
+ fConeSizes[0] = 0.1; fConeSizes[1] = 0.2; fConeSizes[2] = 0.3; fConeSizes[3] = 0.4;
+ fPtThresholds[0]=1.; fPtThresholds[1]=2.; fPtThresholds[2]=3.; fPtThresholds[3]=4.;
+
}
//__________________________________________________________________
Float_t phi = -100. ;
Float_t rad = -100 ;
Int_t n = 0 ;
- TParticle * particle = new TParticle();
+ TParticle * particle = new TParticle;
coneptsum = 0;
icmpt = kFALSE;
}
-void AliAnaGammaDirect::Print(const Option_t * opt) const
+//__________________________________________________________________
+void AliAnaGammaDirect::MakeSeveralICAnalysis(TClonesArray * plCalo, TClonesArray * plCTS)
{
+ //Isolation Cut Analysis for both methods and different pt cuts and cones
+ if (fICMethod != kSeveralIC)
+ AliFatal("Remember to set in config file: directGamma->SetICMethod(kSeveralIC)");
+
+ for(Int_t ipr = 0; ipr < plCalo->GetEntries() ; ipr ++ ){
+ TParticle * pCandidate = dynamic_cast<TParticle *>(plCalo->At(ipr)) ;
+
+ if(pCandidate->Pt() > fMinGammaPt){
+
+ Bool_t icPtThres = kFALSE;
+ Bool_t icPtSum = kFALSE;
+
+ Float_t ptC = pCandidate->Pt() ;
+
+ fhNGamma->Fill(ptC);
+ fhPhiGamma->Fill( ptC,pCandidate->Phi());
+ fhEtaGamma->Fill(ptC,pCandidate->Eta());
+
+ for(Int_t icone = 0; icone<fNCones; icone++){
+ fConeSize=fConeSizes[icone] ;
+ Float_t coneptsum = 0 ;
+ for(Int_t ipt = 0; ipt<fNPtThres;ipt++){
+ fPtThreshold=fPtThresholds[ipt] ;
+ MakeIsolationCut(plCTS,plCalo, pCandidate, ipr, icPtThres, icPtSum,coneptsum);
+ AliDebug(1,Form("Candidate pt %f, pt in cone %f, Isolated? ICPt %d, ICSum %d",
+ pCandidate->Pt(), coneptsum, icPtThres, icPtSum));
+
+ fhPtThresIsolated[icone][ipt]->Fill(ptC);
+ }//pt thresh loop
+ fhPtSumIsolated[icone]->Fill(ptC,coneptsum) ;
+ }//cone size loop
+ }//min pt candidate
+ }//candidate loop
+}
+
+void AliAnaGammaDirect::Print(const Option_t * opt) const
+{
+
//Print some relevant parameters set for the analysis
if(! opt)
return;
-
+
Info("Print", "%s %s", GetName(), GetTitle() ) ;
- printf("IC method = %d\n", fMakeICMethod) ;
+
+ printf("Min Gamma pT = %f\n", fMinGammaPt) ;
+ printf("IC method = %d\n", fICMethod) ;
printf("Cone Size = %f\n", fConeSize) ;
- printf("pT threshold = %f\n", fPtThreshold) ;
- printf("pT sum threshold = %f\n", fPtSumThreshold) ;
- printf("Min Gamma pT = %f\n", fMinGammaPt) ;
- printf("Calorimeter = %s\n", fCalorimeter.Data()) ;
-}
-
-void AliAnaGammaDirect::Terminate(Option_t *)
-{
- // The Terminate() function is the last function to be called during
- // a query. It always runs on the client, it can be used to present
- // the results graphically or save the results to file.
-
+ if(fICMethod == kPtIC) printf("pT threshold = %f\n", fPtThreshold) ;
+ if(fICMethod == kSumPtIC) printf("pT sum threshold = %f\n", fPtSumThreshold) ;
+
+ if(fICMethod == kSeveralIC){
+ printf("N Cone Sizes = %d\n", fNCones) ;
+ printf("N pT thresholds = %d\n", fNPtThres) ;
+ printf("Cone Sizes = \n") ;
+ for(Int_t i = 0; i < fNCones; i++)
+ printf(" %f;", fConeSizes[i]) ;
+ printf(" \n") ;
+ for(Int_t i = 0; i < fNPtThres; i++)
+ printf(" %f;", fPtThresholds[i]) ;
+ }
-}
+ printf(" \n") ;
+
+}