[u/mrichter/AliRoot.git] / PWG4 / PartCorrDep / AliAnaParticleJetFinderCorrelation.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: AliAnaParticleJetFinderCorrelation.cxx 22232 2007-11-17 16:39:49Z gustavo $ */

//_________________________________________________________________________
// Class for the analysis of particle (direct gamma) -jet (jet found with finder) correlations
//*-- Author: Gustavo Conesa (LNF-INFN) 
//////////////////////////////////////////////////////////////////////////////


// --- ROOT system ---
#include "TH2F.h"
#include "TClonesArray.h"
#include "TClass.h"
//#include "Riostream.h"

//---- AliRoot system ----
#include "AliCaloTrackReader.h"
#include "AliAODJet.h"
#include "AliAnaParticleJetFinderCorrelation.h" 
#include "AliAODPWG4ParticleCorrelation.h"
#include "AliAODTrack.h"
#include "AliAODCaloCluster.h"
#include "AliAODEvent.h"

ClassImp(AliAnaParticleJetFinderCorrelation)
  

//____________________________________________________________________________
  AliAnaParticleJetFinderCorrelation::AliAnaParticleJetFinderCorrelation() : 
    AliAnaPartCorrBaseClass(),  
    fDeltaPhiMaxCut(0.), fDeltaPhiMinCut(0.), fRatioMaxCut(0.),  fRatioMinCut(0.), 
    fConeSize(0), fPtThresholdInCone(0),fUseJetRefTracks(0), fMakeCorrelationInHistoMaker(0), fSelectIsolated(0),
    fhDeltaEta(0), fhDeltaPhi(0), fhDeltaPt(0), fhPtRatio(0), fhPt(0),
    fhFFz(0),fhFFxi(0),fhFFpt(0),fhNTracksInCone(0)
{
  //Default Ctor
  
  //Initialize parameters
  InitParameters();
}
/*
//____________________________________________________________________________
AliAnaParticleJetFinderCorrelation::AliAnaParticleJetFinderCorrelation(const AliAnaParticleJetFinderCorrelation & pjf) :   
  AliAnaPartCorrBaseClass(pjf),  
  fDeltaPhiMaxCut(pjf.fDeltaPhiMaxCut), fDeltaPhiMinCut(pjf.fDeltaPhiMinCut), 
  fRatioMaxCut(pjf.fRatioMaxCut), fRatioMinCut(pjf.fRatioMinCut), 
  fConeSize(pjf.fConeSize), fPtThresholdInCone(pjf.fPtThresholdInCone),   
  fUseJetRefTracks(pjf.fUseJetRefTracks), fMakeCorrelationInHistoMaker(pjf.fMakeCorrelationInHistoMaker),
  fSelectIsolated(pjf.fSelectIsolated),	
  fhDeltaEta(pjf.fhDeltaEta), fhDeltaPhi(pjf.fhDeltaPhi), 
  fhDeltaPt(pjf.fhDeltaPt), fhPtRatio(pjf.fhPtRatio), fhPt(pjf.fhPt), 
  fhFFz(pjf.fhFFz),fhFFxi(pjf.fhFFxi),fhFFpt(pjf.fhFFpt),
  fhNTracksInCone(pjf.fhNTracksInCone)
{
  // cpy ctor
  
}

//_________________________________________________________________________
AliAnaParticleJetFinderCorrelation & AliAnaParticleJetFinderCorrelation::operator = (const AliAnaParticleJetFinderCorrelation & pjf)
{
  // assignment operator
  
  if(this == &pjf)return *this;
  ((AliAnaPartCorrBaseClass *)this)->operator=(pjf);
  
  fDeltaPhiMaxCut    = pjf.fDeltaPhiMaxCut ; 
  fDeltaPhiMinCut    = pjf.fDeltaPhiMinCut ; 
  fRatioMaxCut       = pjf.fRatioMaxCut ;
  fRatioMinCut       = pjf.fRatioMinCut ; 
  fConeSize          = pjf.fConeSize ; 
  fPtThresholdInCone = pjf.fPtThresholdInCone ;
  fUseJetRefTracks   = pjf.fUseJetRefTracks ;
  fMakeCorrelationInHistoMaker = pjf.fMakeCorrelationInHistoMaker ;  
  fSelectIsolated    = pjf.fSelectIsolated ;
  
  //Histograms
  fhDeltaEta = pjf.fhDeltaEta;
  fhDeltaPhi = pjf.fhDeltaPhi;
  fhDeltaPt  = pjf.fhDeltaPt;
  fhPtRatio  = pjf.fhPtRatio;
  fhPt       = pjf.fhPt;
  
  fhFFz   = pjf.fhFFz;
  fhFFxi  = pjf.fhFFxi;
  fhFFpt  = pjf.fhFFpt;
  fhNTracksInCone = pjf.fhNTracksInCone;  
  
  return *this;
  
}
*/
//____________________________________________________________________________
//AliAnaParticleJetFinderCorrelation::~AliAnaParticleJetFinderCorrelation() 
//{
//   // Remove all pointers except analysis output pointers.
// 
//}


//________________________________________________________________________
TList *  AliAnaParticleJetFinderCorrelation::GetCreateOutputObjects()
{  
  // Create histograms to be saved in output file and 
  // store them in fOutputContainer
    
  TList * outputContainer = new TList() ; 
  outputContainer->SetName("ParticleJetFinderHistos") ; 
  
  Int_t nptbins  = GetHistoPtBins();
  //	Int_t nphibins = GetHistoPhiBins();
  //	Int_t netabins = GetHistoEtaBins();
  Float_t ptmax  = GetHistoPtMax();
  //	Float_t phimax = GetHistoPhiMax();
  //	Float_t etamax = GetHistoEtaMax();
  Float_t ptmin  = GetHistoPtMin();
  //	Float_t phimin = GetHistoPhiMin();
//	Float_t etamin = GetHistoEtaMin();	
  
  fhDeltaPhi  = new TH2F("DeltaPhi","#phi_{jet} - #phi_{trigger} vs p_{T trigger}",nptbins,ptmin,ptmax,100,-4,4); 
  fhDeltaPhi->SetYTitle("#Delta #phi");
  fhDeltaPhi->SetXTitle("p_{T trigger} (GeV/c)");
  outputContainer->Add(fhDeltaPhi);
  
  fhDeltaEta  = new TH2F("DeltaEta","#eta_{jet} - #eta_{trigger} vs p_{T trigger}",nptbins,ptmin,ptmax,100,-5,5); 
  fhDeltaEta->SetYTitle("#Delta #eta");
  fhDeltaEta->SetXTitle("p_{T trigger} (GeV/c)");
  outputContainer->Add(fhDeltaEta);
  
  fhDeltaPt  = new TH2F("DeltaPt","p_{T trigger} - #p_{T jet} vs p_{T trigger}",nptbins,ptmin,ptmax,100,-50,50); 
  fhDeltaPt->SetYTitle("#Delta #p_{T}");
  fhDeltaPt->SetXTitle("p_{T trigger} (GeV/c)"); 
  outputContainer->Add(fhDeltaPt);
  
  fhPtRatio  = new TH2F("PtRatio","p_{T jet} / #p_{T trigger} vs p_{T trigger}",nptbins,ptmin,ptmax,200,0,2.); 
  fhPtRatio->SetYTitle("ratio");
  fhPtRatio->SetXTitle("p_{T trigger} (GeV/c)");
  outputContainer->Add(fhPtRatio);
  
  fhPt  = new TH2F("Pt","p_{T jet} vs p_{T trigger}",nptbins,ptmin,ptmax,nptbins,ptmin,ptmax); 
  fhPt->SetYTitle("p_{T jet}(GeV/c)");
  fhPt->SetXTitle("p_{T trigger} (GeV/c)");
  outputContainer->Add(fhPt);
  
  fhFFz  = new TH2F("FFz","z = p_{T i charged}/p_{T trigger} vs p_{T trigger}",nptbins,ptmin,ptmax,200,0.,2);  
  fhFFz->SetYTitle("z");
  fhFFz->SetXTitle("p_{T trigger}");
  outputContainer->Add(fhFFz) ;
	
  fhFFxi  = new TH2F("FFxi","#xi = ln(p_{T trigger}/p_{T i charged}) vs p_{T trigger}", nptbins,ptmin,ptmax,100,0.,10.); 
  fhFFxi->SetYTitle("#xi");
  fhFFxi->SetXTitle("p_{T trigger}");
  outputContainer->Add(fhFFxi) ;
  
  fhFFpt  = new TH2F("FFpt","#xi = p_{T i charged}) vs p_{T trigger}", nptbins,ptmin,ptmax,100,0.,50.); 
  fhFFpt->SetYTitle("p_{T charged hadron}");
  fhFFpt->SetXTitle("p_{T trigger}");
  outputContainer->Add(fhFFpt) ;
  
  fhNTracksInCone  = new TH2F("NTracksInCone","#xi = p_{T i charged}) vs p_{T trigger}", nptbins,ptmin,ptmax,100,0.,50.); 
  fhNTracksInCone->SetYTitle("p_{T charged hadron}");
  fhNTracksInCone->SetXTitle("p_{T trigger}");
  outputContainer->Add(fhNTracksInCone) ;
  
  return outputContainer;

}

//____________________________________________________________________________
void AliAnaParticleJetFinderCorrelation::InitParameters()
{
  
  //Initialize the parameters of the analysis.
  SetInputAODName("PWG4Particle");
  AddToHistogramsName("AnaJetFinderCorr_");

  fDeltaPhiMinCut    = 1.5 ;
  fDeltaPhiMaxCut    = 4.5 ; 
  fRatioMaxCut       = 1.2 ; 
  fRatioMinCut       = 0.3 ; 
  fConeSize          = 0.3 ;
  fPtThresholdInCone = 0.5 ;
  fUseJetRefTracks   = kFALSE ;
  fMakeCorrelationInHistoMaker = kFALSE ;
  fSelectIsolated = kFALSE;
  }

//__________________________________________________________________
Int_t  AliAnaParticleJetFinderCorrelation::SelectJet(AliAODPWG4Particle * particle, AliAODEvent *event) const 
{
  //Returns the index of the jet that is opposite to the particle
  
  Int_t njets = event->GetNJets() ;	
  
  AliAODJet * jet = 0 ;
  Int_t index = -1;
  for(Int_t ijet = 0; ijet < njets ; ijet++){
    jet = event->GetJet(ijet) ;	  
    Float_t dphi  = TMath::Abs(particle->Phi()-jet->Phi());
    Float_t ratio = jet->Pt()/particle->Pt();
    if(GetDebug() > 3)
      printf("AliAnaParticleJetFinderCorrelation::SelectJet() - Jet %d, Ratio pT %2.3f, Delta phi %2.3f\n",ijet,ratio,dphi);	  
    Float_t dphiprev= 10000;
    if((dphi > fDeltaPhiMinCut) && (dphi<fDeltaPhiMaxCut) &&
       (ratio > fRatioMinCut) && (ratio < fRatioMaxCut)  &&
       (TMath::Abs(dphi-3.14159) < TMath::Abs(dphiprev-3.14159))){//In case there is more than one jet select the most opposite.
    dphiprev = dphi;
    index = ijet ;	
  }//Selection cuts
}//AOD Jet loop

return index ;

}

//__________________________________________________________________
void  AliAnaParticleJetFinderCorrelation::MakeAnalysisFillAOD() 
{  
  //Particle-Jet Correlation Analysis, fill AODs

  //Get the event, check if there are AODs available, if not, skip this analysis
  AliAODEvent * event = NULL;
  if(GetReader()->GetOutputEvent()) 
  {
    event = dynamic_cast<AliAODEvent*>(GetReader()->GetOutputEvent()); 
  }
  else if(GetReader()->GetDataType() == AliCaloTrackReader::kAOD) 
  {   
    event = dynamic_cast<AliAODEvent*>(GetReader()->GetInputEvent()); 
  }
  else 
  {
    if(GetDebug() > 3) printf("AliAnaParticleJetFinderCorrelation::MakeAnalysisFillAOD() - There are no jets available for this analysis\n");
    return;
  }
    
  if(!GetInputAODBranch() || !event){
    printf("AliAnaParticleJetFinderCorrelation::MakeAnalysisFillAOD() - No input particles in AOD with name branch < %s > \n",GetInputAODName().Data());
    abort();
  }
  
  if(strcmp(GetInputAODBranch()->GetClass()->GetName(), "AliAODPWG4ParticleCorrelation")){
    printf("AliAnaParticleJetFinderCorrelation::MakeAnalysisFillAOD() - Wrong type of AOD object, change AOD class name in input AOD: It should be <AliAODPWG4ParticleCorrelation> and not <%s> \n",GetInputAODBranch()->GetClass()->GetName());
    abort();
  }
	
  Int_t ntrig =  GetInputAODBranch()->GetEntriesFast() ;  
  if(GetDebug() > 3){
    printf("AliAnaParticleJetFinderCorrelation::MakeAnalysisFillAOD() - Begin jet finder  correlation analysis, fill AODs \n");
    printf("AliAnaParticleJetFinderCorrelation::MakeAnalysisFillAOD() - In particle branch aod entries %d\n", ntrig);
    printf("AliAnaParticleJetFinderCorrelation::MakeAnalysisFillAOD() - In jet      branch aod entries %d\n", event->GetNJets());
  }
  
  //Loop on stored AOD particles, trigger
  for(Int_t iaod = 0; iaod < ntrig ; iaod++){
    AliAODPWG4ParticleCorrelation* particle =  (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
	
    //Correlate with jets
    Int_t ijet = SelectJet(particle,event);
    if(ijet > -1){
      if(GetDebug() > 2) printf ("AliAnaParticleJetFinderCorrelation::MakeAnalysisFillAOD() - Jet with index %d selected \n",ijet);
      AliAODJet *jet = event->GetJet(ijet);   	 
      particle->SetRefJet(jet);	
    }
  } // input aod loop		  
  
  if(GetDebug() > 1 ) printf("AliAnaParticleJetFinderCorrelation::MakeAnalysisFillAOD() - End fill AODs \n");
} 

//__________________________________________________________________
void  AliAnaParticleJetFinderCorrelation::MakeAnalysisFillHistograms() 
{
  //Particle-Jet Correlation Analysis, fill histograms
  
  //Get the event, check if there are AODs available, if not, skip this analysis
  AliAODEvent * event = NULL;
  if(GetReader()->GetOutputEvent()) 
  {
    event = dynamic_cast<AliAODEvent*>(GetReader()->GetOutputEvent()); 
  }
  else if(GetReader()->GetDataType() == AliCaloTrackReader::kAOD) 
  {  
    event = dynamic_cast<AliAODEvent*>(GetReader()->GetInputEvent()); 
  }
  else {
    if(GetDebug() > 3) printf("AliAnaParticleJetFinderCorrelation::MakeAnalysisFillHistograms() - There are no jets available for this analysis\n");
    return;
  }
  
  if(!GetInputAODBranch() || !event){
    printf("AliAnaParticleJetFinderCorrelation::MakeAnalysisFillHistograms() - No input particles in AOD with name branch < %s > \n",GetInputAODName().Data());
    abort();
  }
  
  Int_t ntrig   =  GetInputAODBranch()->GetEntriesFast() ;    
  if(GetDebug() > 1){
    printf("AliAnaParticleJetFinderCorrelation::MakeAnalysisFillHistograms() - Begin jet finder  correlation analysis, fill histograms \n");
    printf("AliAnaParticleJetFinderCorrelation::MakeAnalysisFillHistograms() - In particle branch aod entries %d\n", ntrig);
    printf("AliAnaParticleJetFinderCorrelation::MakeAnalysisFillHistograms() - In jet output branch aod entries %d\n", event->GetNJets());
  }
  
  //Loop on stored AOD particles, trigger
  for(Int_t iaod = 0; iaod < ntrig ; iaod++){
    AliAODPWG4ParticleCorrelation* particlecorr =  (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));	
    
    if(OnlyIsolated() && !particlecorr->IsIsolated()) continue;
    
    //Recover the jet correlated, found previously.
    AliAODJet	* jet = particlecorr->GetJet();
    //If correlation not made before, do it now.
    if(fMakeCorrelationInHistoMaker){
      //Correlate with jets
      Int_t ijet = SelectJet(particlecorr,event);
      if(ijet > -1){
        if(GetDebug() > 2) printf ("AliAnaParticleJetFinderCorrelation::MakeAnalysisFillHistograms() - Jet with index %d selected \n",ijet);
        jet = event->GetJet(ijet);
        particlecorr->SetRefJet(jet);	
      }
    }
    
    if (!jet) continue ;
    
    //Fill Histograms
    
    Double_t ptTrig = particlecorr->Pt();
    Double_t ptJet = jet->Pt();
    Double_t phiTrig = particlecorr->Phi();
    Double_t phiJet = jet->Phi();
    Double_t etaTrig = particlecorr->Eta();
    Double_t etaJet = jet->Eta();
    //printf("pT trigger %2.3f, pT jet %2.3f, Delta phi %2.3f, Delta eta %2.3f, Delta pT %2.3f, ratio %2.3f \n",
    //	ptTrig,ptJet, phiJet-phiTrig, etaJet-etaTrig, ptTrig-ptJet, ptJet/ptTrig);
    fhDeltaPt ->Fill(ptTrig, ptTrig-ptJet);
    fhDeltaPhi->Fill(ptTrig, phiJet-phiTrig);
    fhDeltaEta->Fill(ptTrig, etaJet-etaTrig);
    fhPtRatio ->Fill(ptTrig, ptJet/ptTrig);
    fhPt      ->Fill(ptTrig, ptJet);
    
    //Fragmentation function
    Float_t	 rad = 0, pt = 0, eta = 0, phi = 0;
    Int_t	 npartcone = 0;
    TVector3 p3;
    
    Int_t ntracks =  0;
    if(!fUseJetRefTracks)
      ntracks =GetAODCTS()->GetEntriesFast();
    else //If you want to use jet tracks from JETAN
      ntracks =  (jet->GetRefTracks())->GetEntriesFast();
    
    AliAODTrack* track = 0x0 ;
    for(Int_t ipr = 0;ipr < ntracks ; ipr ++ ){
      if(!fUseJetRefTracks)
        track = (AliAODTrack *) (GetAODCTS()->At(ipr)) ; 
      else //If you want to use jet tracks from JETAN
        track = (AliAODTrack *) ((jet->GetRefTracks())->At(ipr));
      
      p3.SetXYZ(track->Px(),track->Py(),track->Pz());
      pt    = p3.Pt();
      eta  = p3.Eta();
      phi  = p3.Phi() ;
      if(phi < 0) phi+=TMath::TwoPi();
      
      //Check if there is any particle inside cone with pt larger than  fPtThreshold
      rad = TMath::Sqrt((eta-etaJet)*(eta-etaJet)+ (phi-phiJet)*(phi-phiJet));
      if(rad < fConeSize  && pt > fPtThresholdInCone){	
        //printf("charged in jet cone pt %f, phi %f, eta %f, R %f \n",pt,phi,eta,rad);
        npartcone++;
        fhFFz ->Fill(ptTrig, pt/ptTrig);
        fhFFxi->Fill(ptTrig, TMath::Log(ptTrig/pt));
        fhFFpt->Fill(ptTrig, pt);
      }
    }//Tracks
    fhNTracksInCone->Fill(ptTrig, npartcone);
  }//AOD trigger particle loop
  if(GetDebug() > 1) printf("AliAnaParticleJetFinderCorrelation::MakeAnalysisFillHistograms() - End fill histograms \n");
} 


//__________________________________________________________________
void AliAnaParticleJetFinderCorrelation::Print(const Option_t * opt) const
{
  
  //Print some relevant parameters set for the analysis
  if(! opt)
    return;
  
  printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
  AliAnaPartCorrBaseClass::Print(" ");

  printf("Phi trigger-jet        <     %3.2f\n", fDeltaPhiMaxCut) ; 
  printf("Phi trigger-jet        >     %3.2f\n", fDeltaPhiMinCut) ;
  printf("pT Ratio trigger/jet   <     %3.2f\n", fRatioMaxCut) ; 
  printf("pT Ratio trigger/jet   >     %3.2f\n", fRatioMinCut) ;
  printf("fConeSize              =     %3.2f\n", fConeSize) ; 
  printf("fPtThresholdInCone     =     %3.2f\n", fPtThresholdInCone) ;
  printf("fUseJetRefTracks	   =     %d\n",    fUseJetRefTracks) ;
  printf("fMakeCorrelationInHistoMaker	   =     %d\n",    fMakeCorrelationInHistoMaker) ;	
  printf("Isolated Trigger?  %d\n", fSelectIsolated) ;
  
}
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: AliAnaParticleJetFinderCorrelation.cxx 22232 2007-11-17 16:39:49Z gustavo $ */
	16
	17	//_________________________________________________________________________
	18	// Class for the analysis of particle (direct gamma) -jet (jet found with finder) correlations
	19	//*-- Author: Gustavo Conesa (LNF-INFN)
	20	//////////////////////////////////////////////////////////////////////////////
	21
	22
	23	// --- ROOT system ---
	24	#include "TH2F.h"
477d6cee	25	#include "TClonesArray.h"
9415d854	26	#include "TClass.h"
1c5acb87	27	//#include "Riostream.h"
	28
	29	//---- AliRoot system ----
	30	#include "AliCaloTrackReader.h"
	31	#include "AliAODJet.h"
	32	#include "AliAnaParticleJetFinderCorrelation.h"
1c5acb87	33	#include "AliAODPWG4ParticleCorrelation.h"
477d6cee	34	#include "AliAODTrack.h"
	35	#include "AliAODCaloCluster.h"
	36	#include "AliAODEvent.h"
1c5acb87	37
477d6cee	38	ClassImp(AliAnaParticleJetFinderCorrelation)
1c5acb87	39
	40
	41	//____________________________________________________________________________
	42	AliAnaParticleJetFinderCorrelation::AliAnaParticleJetFinderCorrelation() :
477d6cee	43	AliAnaPartCorrBaseClass(),
1c5acb87	44	fDeltaPhiMaxCut(0.), fDeltaPhiMinCut(0.), fRatioMaxCut(0.), fRatioMinCut(0.),
477d6cee	45	fConeSize(0), fPtThresholdInCone(0),fUseJetRefTracks(0), fMakeCorrelationInHistoMaker(0), fSelectIsolated(0),
	46	fhDeltaEta(0), fhDeltaPhi(0), fhDeltaPt(0), fhPtRatio(0), fhPt(0),
	47	fhFFz(0),fhFFxi(0),fhFFpt(0),fhNTracksInCone(0)
1c5acb87	48	{
477d6cee	49	//Default Ctor
	50
	51	//Initialize parameters
	52	InitParameters();
1c5acb87	53	}
78219bac	54	/*
1c5acb87	55	//____________________________________________________________________________
1c5acb87	56	AliAnaParticleJetFinderCorrelation::AliAnaParticleJetFinderCorrelation(const AliAnaParticleJetFinderCorrelation & pjf) :
477d6cee	57	AliAnaPartCorrBaseClass(pjf),
	58	fDeltaPhiMaxCut(pjf.fDeltaPhiMaxCut), fDeltaPhiMinCut(pjf.fDeltaPhiMinCut),
	59	fRatioMaxCut(pjf.fRatioMaxCut), fRatioMinCut(pjf.fRatioMinCut),
	60	fConeSize(pjf.fConeSize), fPtThresholdInCone(pjf.fPtThresholdInCone),
	61	fUseJetRefTracks(pjf.fUseJetRefTracks), fMakeCorrelationInHistoMaker(pjf.fMakeCorrelationInHistoMaker),
	62	fSelectIsolated(pjf.fSelectIsolated),
	63	fhDeltaEta(pjf.fhDeltaEta), fhDeltaPhi(pjf.fhDeltaPhi),
	64	fhDeltaPt(pjf.fhDeltaPt), fhPtRatio(pjf.fhPtRatio), fhPt(pjf.fhPt),
	65	fhFFz(pjf.fhFFz),fhFFxi(pjf.fhFFxi),fhFFpt(pjf.fhFFpt),
	66	fhNTracksInCone(pjf.fhNTracksInCone)
1c5acb87	67	{
477d6cee	68	// cpy ctor
477d6cee	69
1c5acb87	70	}
	71
	72	//_________________________________________________________________________
	73	AliAnaParticleJetFinderCorrelation & AliAnaParticleJetFinderCorrelation::operator = (const AliAnaParticleJetFinderCorrelation & pjf)
	74	{
477d6cee	75	// assignment operator
	76
	77	if(this == &pjf)return *this;
	78	((AliAnaPartCorrBaseClass *)this)->operator=(pjf);
1c5acb87	79
477d6cee	80	fDeltaPhiMaxCut = pjf.fDeltaPhiMaxCut ;
	81	fDeltaPhiMinCut = pjf.fDeltaPhiMinCut ;
	82	fRatioMaxCut = pjf.fRatioMaxCut ;
	83	fRatioMinCut = pjf.fRatioMinCut ;
	84	fConeSize = pjf.fConeSize ;
	85	fPtThresholdInCone = pjf.fPtThresholdInCone ;
	86	fUseJetRefTracks = pjf.fUseJetRefTracks ;
	87	fMakeCorrelationInHistoMaker = pjf.fMakeCorrelationInHistoMaker ;
	88	fSelectIsolated = pjf.fSelectIsolated ;
1c5acb87	89
477d6cee	90	//Histograms
	91	fhDeltaEta = pjf.fhDeltaEta;
	92	fhDeltaPhi = pjf.fhDeltaPhi;
	93	fhDeltaPt = pjf.fhDeltaPt;
	94	fhPtRatio = pjf.fhPtRatio;
	95	fhPt = pjf.fhPt;
	96
	97	fhFFz = pjf.fhFFz;
	98	fhFFxi = pjf.fhFFxi;
	99	fhFFpt = pjf.fhFFpt;
	100	fhNTracksInCone = pjf.fhNTracksInCone;
1c5acb87	101
1c5acb87	102	return *this;
477d6cee	103
1c5acb87	104	}
78219bac	105	*/
1c5acb87	106	//____________________________________________________________________________
	107	//AliAnaParticleJetFinderCorrelation::~AliAnaParticleJetFinderCorrelation()
	108	//{
	109	// // Remove all pointers except analysis output pointers.
	110	//
	111	//}
	112
	113
	114	//________________________________________________________________________
	115	TList * AliAnaParticleJetFinderCorrelation::GetCreateOutputObjects()
	116	{
477d6cee	117	// Create histograms to be saved in output file and
477d6cee	118	// store them in fOutputContainer
a3aebfff	119
477d6cee	120	TList * outputContainer = new TList() ;
	121	outputContainer->SetName("ParticleJetFinderHistos") ;
	122
5a2dbc3c	123	Int_t nptbins = GetHistoPtBins();
	124	// Int_t nphibins = GetHistoPhiBins();
	125	// Int_t netabins = GetHistoEtaBins();
477d6cee	126	Float_t ptmax = GetHistoPtMax();
	127	// Float_t phimax = GetHistoPhiMax();
	128	// Float_t etamax = GetHistoEtaMax();
	129	Float_t ptmin = GetHistoPtMin();
	130	// Float_t phimin = GetHistoPhiMin();
1c5acb87	131	// Float_t etamin = GetHistoEtaMin();
477d6cee	132
4df35693	133	fhDeltaPhi = new TH2F("DeltaPhi","#phi_{jet} - #phi_{trigger} vs p_{T trigger}",nptbins,ptmin,ptmax,100,-4,4);
477d6cee	134	fhDeltaPhi->SetYTitle("#Delta #phi");
	135	fhDeltaPhi->SetXTitle("p_{T trigger} (GeV/c)");
	136	outputContainer->Add(fhDeltaPhi);
	137
4df35693	138	fhDeltaEta = new TH2F("DeltaEta","#eta_{jet} - #eta_{trigger} vs p_{T trigger}",nptbins,ptmin,ptmax,100,-5,5);
477d6cee	139	fhDeltaEta->SetYTitle("#Delta #eta");
	140	fhDeltaEta->SetXTitle("p_{T trigger} (GeV/c)");
	141	outputContainer->Add(fhDeltaEta);
	142
4df35693	143	fhDeltaPt = new TH2F("DeltaPt","p_{T trigger} - #p_{T jet} vs p_{T trigger}",nptbins,ptmin,ptmax,100,-50,50);
477d6cee	144	fhDeltaPt->SetYTitle("#Delta #p_{T}");
	145	fhDeltaPt->SetXTitle("p_{T trigger} (GeV/c)");
	146	outputContainer->Add(fhDeltaPt);
	147
9415d854	148	fhPtRatio = new TH2F("PtRatio","p_{T jet} / #p_{T trigger} vs p_{T trigger}",nptbins,ptmin,ptmax,200,0,2.);
477d6cee	149	fhPtRatio->SetYTitle("ratio");
	150	fhPtRatio->SetXTitle("p_{T trigger} (GeV/c)");
	151	outputContainer->Add(fhPtRatio);
	152
9415d854	153	fhPt = new TH2F("Pt","p_{T jet} vs p_{T trigger}",nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
9415d854	154	fhPt->SetYTitle("p_{T jet}(GeV/c)");
477d6cee	155	fhPt->SetXTitle("p_{T trigger} (GeV/c)");
	156	outputContainer->Add(fhPt);
	157
	158	fhFFz = new TH2F("FFz","z = p_{T i charged}/p_{T trigger} vs p_{T trigger}",nptbins,ptmin,ptmax,200,0.,2);
	159	fhFFz->SetYTitle("z");
	160	fhFFz->SetXTitle("p_{T trigger}");
	161	outputContainer->Add(fhFFz) ;
1c5acb87	162
477d6cee	163	fhFFxi = new TH2F("FFxi","#xi = ln(p_{T trigger}/p_{T i charged}) vs p_{T trigger}", nptbins,ptmin,ptmax,100,0.,10.);
	164	fhFFxi->SetYTitle("#xi");
	165	fhFFxi->SetXTitle("p_{T trigger}");
	166	outputContainer->Add(fhFFxi) ;
	167
4df35693	168	fhFFpt = new TH2F("FFpt","#xi = p_{T i charged}) vs p_{T trigger}", nptbins,ptmin,ptmax,100,0.,50.);
477d6cee	169	fhFFpt->SetYTitle("p_{T charged hadron}");
	170	fhFFpt->SetXTitle("p_{T trigger}");
	171	outputContainer->Add(fhFFpt) ;
	172
4df35693	173	fhNTracksInCone = new TH2F("NTracksInCone","#xi = p_{T i charged}) vs p_{T trigger}", nptbins,ptmin,ptmax,100,0.,50.);
477d6cee	174	fhNTracksInCone->SetYTitle("p_{T charged hadron}");
	175	fhNTracksInCone->SetXTitle("p_{T trigger}");
	176	outputContainer->Add(fhNTracksInCone) ;
	177
	178	return outputContainer;
a3aebfff	179
1c5acb87	180	}
1c5acb87	181
477d6cee	182	//____________________________________________________________________________
1c5acb87	183	void AliAnaParticleJetFinderCorrelation::InitParameters()
1c5acb87	184	{
477d6cee	185
1c5acb87	186	//Initialize the parameters of the analysis.
a3aebfff	187	SetInputAODName("PWG4Particle");
	188	AddToHistogramsName("AnaJetFinderCorr_");
	189
1c5acb87	190	fDeltaPhiMinCut = 1.5 ;
	191	fDeltaPhiMaxCut = 4.5 ;
	192	fRatioMaxCut = 1.2 ;
	193	fRatioMinCut = 0.3 ;
	194	fConeSize = 0.3 ;
	195	fPtThresholdInCone = 0.5 ;
	196	fUseJetRefTracks = kFALSE ;
	197	fMakeCorrelationInHistoMaker = kFALSE ;
477d6cee	198	fSelectIsolated = kFALSE;
1c5acb87	199	}
	200
	201	//__________________________________________________________________
f37fa8d2	202	Int_t AliAnaParticleJetFinderCorrelation::SelectJet(AliAODPWG4Particle * particle, AliAODEvent *event) const
1c5acb87	203	{
477d6cee	204	//Returns the index of the jet that is opposite to the particle
477d6cee	205
f37fa8d2	206	Int_t njets = event->GetNJets() ;
477d6cee	207
76c58de9	208	AliAODJet * jet = 0 ;
477d6cee	209	Int_t index = -1;
477d6cee	210	for(Int_t ijet = 0; ijet < njets ; ijet++){
f37fa8d2	211	jet = event->GetJet(ijet) ;
477d6cee	212	Float_t dphi = TMath::Abs(particle->Phi()-jet->Phi());
	213	Float_t ratio = jet->Pt()/particle->Pt();
	214	if(GetDebug() > 3)
a3aebfff	215	printf("AliAnaParticleJetFinderCorrelation::SelectJet() - Jet %d, Ratio pT %2.3f, Delta phi %2.3f\n",ijet,ratio,dphi);
477d6cee	216	Float_t dphiprev= 10000;
	217	if((dphi > fDeltaPhiMinCut) && (dphi<fDeltaPhiMaxCut) &&
	218	(ratio > fRatioMinCut) && (ratio < fRatioMaxCut) &&
	219	(TMath::Abs(dphi-3.14159) < TMath::Abs(dphiprev-3.14159))){//In case there is more than one jet select the most opposite.
	220	dphiprev = dphi;
	221	index = ijet ;
	222	}//Selection cuts
	223	}//AOD Jet loop
1c5acb87	224
477d6cee	225	return index ;
1c5acb87	226
	227	}
	228
	229	//__________________________________________________________________
	230	void AliAnaParticleJetFinderCorrelation::MakeAnalysisFillAOD()
	231	{
477d6cee	232	//Particle-Jet Correlation Analysis, fill AODs
f37fa8d2	233
	234	//Get the event, check if there are AODs available, if not, skip this analysis
	235	AliAODEvent * event = NULL;
	236	if(GetReader()->GetOutputEvent())
	237	{
	238	event = dynamic_cast<AliAODEvent*>(GetReader()->GetOutputEvent());
	239	}
	240	else if(GetReader()->GetDataType() == AliCaloTrackReader::kAOD)
	241	{
	242	event = dynamic_cast<AliAODEvent*>(GetReader()->GetInputEvent());
	243	}
	244	else
	245	{
	246	if(GetDebug() > 3) printf("AliAnaParticleJetFinderCorrelation::MakeAnalysisFillAOD() - There are no jets available for this analysis\n");
	247	return;
	248	}
	249
f8006433	250	if(!GetInputAODBranch() \|\| !event){
a3aebfff	251	printf("AliAnaParticleJetFinderCorrelation::MakeAnalysisFillAOD() - No input particles in AOD with name branch < %s > \n",GetInputAODName().Data());
477d6cee	252	abort();
477d6cee	253	}
a3aebfff	254
9415d854	255	if(strcmp(GetInputAODBranch()->GetClass()->GetName(), "AliAODPWG4ParticleCorrelation")){
f37fa8d2	256	printf("AliAnaParticleJetFinderCorrelation::MakeAnalysisFillAOD() - Wrong type of AOD object, change AOD class name in input AOD: It should be <AliAODPWG4ParticleCorrelation> and not <%s> \n",GetInputAODBranch()->GetClass()->GetName());
f37fa8d2	257	abort();
9415d854	258	}
9415d854	259
477d6cee	260	Int_t ntrig = GetInputAODBranch()->GetEntriesFast() ;
477d6cee	261	if(GetDebug() > 3){
a3aebfff	262	printf("AliAnaParticleJetFinderCorrelation::MakeAnalysisFillAOD() - Begin jet finder correlation analysis, fill AODs \n");
a3aebfff	263	printf("AliAnaParticleJetFinderCorrelation::MakeAnalysisFillAOD() - In particle branch aod entries %d\n", ntrig);
f37fa8d2	264	printf("AliAnaParticleJetFinderCorrelation::MakeAnalysisFillAOD() - In jet branch aod entries %d\n", event->GetNJets());
477d6cee	265	}
f37fa8d2	266
477d6cee	267	//Loop on stored AOD particles, trigger
	268	for(Int_t iaod = 0; iaod < ntrig ; iaod++){
	269	AliAODPWG4ParticleCorrelation* particle = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
a3aebfff	270
477d6cee	271	//Correlate with jets
f37fa8d2	272	Int_t ijet = SelectJet(particle,event);
477d6cee	273	if(ijet > -1){
a3aebfff	274	if(GetDebug() > 2) printf ("AliAnaParticleJetFinderCorrelation::MakeAnalysisFillAOD() - Jet with index %d selected \n",ijet);
f37fa8d2	275	AliAODJet *jet = event->GetJet(ijet);
477d6cee	276	particle->SetRefJet(jet);
477d6cee	277	}
a3aebfff	278	} // input aod loop
477d6cee	279
a3aebfff	280	if(GetDebug() > 1 ) printf("AliAnaParticleJetFinderCorrelation::MakeAnalysisFillAOD() - End fill AODs \n");
1c5acb87	281	}
	282
	283	//__________________________________________________________________
	284	void AliAnaParticleJetFinderCorrelation::MakeAnalysisFillHistograms()
	285	{
477d6cee	286	//Particle-Jet Correlation Analysis, fill histograms
477d6cee	287
f37fa8d2	288	//Get the event, check if there are AODs available, if not, skip this analysis
	289	AliAODEvent * event = NULL;
	290	if(GetReader()->GetOutputEvent())
	291	{
	292	event = dynamic_cast<AliAODEvent*>(GetReader()->GetOutputEvent());
	293	}
	294	else if(GetReader()->GetDataType() == AliCaloTrackReader::kAOD)
	295	{
	296	event = dynamic_cast<AliAODEvent*>(GetReader()->GetInputEvent());
	297	}
	298	else {
	299	if(GetDebug() > 3) printf("AliAnaParticleJetFinderCorrelation::MakeAnalysisFillHistograms() - There are no jets available for this analysis\n");
	300	return;
	301	}
	302
f8006433	303	if(!GetInputAODBranch() \|\| !event){
a3aebfff	304	printf("AliAnaParticleJetFinderCorrelation::MakeAnalysisFillHistograms() - No input particles in AOD with name branch < %s > \n",GetInputAODName().Data());
477d6cee	305	abort();
477d6cee	306	}
f37fa8d2	307
477d6cee	308	Int_t ntrig = GetInputAODBranch()->GetEntriesFast() ;
477d6cee	309	if(GetDebug() > 1){
a3aebfff	310	printf("AliAnaParticleJetFinderCorrelation::MakeAnalysisFillHistograms() - Begin jet finder correlation analysis, fill histograms \n");
a3aebfff	311	printf("AliAnaParticleJetFinderCorrelation::MakeAnalysisFillHistograms() - In particle branch aod entries %d\n", ntrig);
f37fa8d2	312	printf("AliAnaParticleJetFinderCorrelation::MakeAnalysisFillHistograms() - In jet output branch aod entries %d\n", event->GetNJets());
477d6cee	313	}
	314
	315	//Loop on stored AOD particles, trigger
	316	for(Int_t iaod = 0; iaod < ntrig ; iaod++){
a3aebfff	317	AliAODPWG4ParticleCorrelation* particlecorr = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
477d6cee	318
a3aebfff	319	if(OnlyIsolated() && !particlecorr->IsIsolated()) continue;
477d6cee	320
477d6cee	321	//Recover the jet correlated, found previously.
a3aebfff	322	AliAODJet * jet = particlecorr->GetJet();
477d6cee	323	//If correlation not made before, do it now.
	324	if(fMakeCorrelationInHistoMaker){
	325	//Correlate with jets
f37fa8d2	326	Int_t ijet = SelectJet(particlecorr,event);
477d6cee	327	if(ijet > -1){
f37fa8d2	328	if(GetDebug() > 2) printf ("AliAnaParticleJetFinderCorrelation::MakeAnalysisFillHistograms() - Jet with index %d selected \n",ijet);
	329	jet = event->GetJet(ijet);
	330	particlecorr->SetRefJet(jet);
477d6cee	331	}
	332	}
	333
	334	if (!jet) continue ;
	335
	336	//Fill Histograms
	337
a3aebfff	338	Double_t ptTrig = particlecorr->Pt();
477d6cee	339	Double_t ptJet = jet->Pt();
a3aebfff	340	Double_t phiTrig = particlecorr->Phi();
477d6cee	341	Double_t phiJet = jet->Phi();
a3aebfff	342	Double_t etaTrig = particlecorr->Eta();
477d6cee	343	Double_t etaJet = jet->Eta();
	344	//printf("pT trigger %2.3f, pT jet %2.3f, Delta phi %2.3f, Delta eta %2.3f, Delta pT %2.3f, ratio %2.3f \n",
	345	// ptTrig,ptJet, phiJet-phiTrig, etaJet-etaTrig, ptTrig-ptJet, ptJet/ptTrig);
	346	fhDeltaPt ->Fill(ptTrig, ptTrig-ptJet);
	347	fhDeltaPhi->Fill(ptTrig, phiJet-phiTrig);
	348	fhDeltaEta->Fill(ptTrig, etaJet-etaTrig);
	349	fhPtRatio ->Fill(ptTrig, ptJet/ptTrig);
	350	fhPt ->Fill(ptTrig, ptJet);
	351
	352	//Fragmentation function
	353	Float_t rad = 0, pt = 0, eta = 0, phi = 0;
	354	Int_t npartcone = 0;
	355	TVector3 p3;
477d6cee	356
	357	Int_t ntracks = 0;
	358	if(!fUseJetRefTracks)
	359	ntracks =GetAODCTS()->GetEntriesFast();
	360	else //If you want to use jet tracks from JETAN
	361	ntracks = (jet->GetRefTracks())->GetEntriesFast();
	362
f8006433	363	AliAODTrack* track = 0x0 ;
477d6cee	364	for(Int_t ipr = 0;ipr < ntracks ; ipr ++ ){
477d6cee	365	if(!fUseJetRefTracks)
f37fa8d2	366	track = (AliAODTrack *) (GetAODCTS()->At(ipr)) ;
477d6cee	367	else //If you want to use jet tracks from JETAN
f37fa8d2	368	track = (AliAODTrack *) ((jet->GetRefTracks())->At(ipr));
477d6cee	369
	370	p3.SetXYZ(track->Px(),track->Py(),track->Pz());
	371	pt = p3.Pt();
	372	eta = p3.Eta();
	373	phi = p3.Phi() ;
	374	if(phi < 0) phi+=TMath::TwoPi();
	375
	376	//Check if there is any particle inside cone with pt larger than fPtThreshold
	377	rad = TMath::Sqrt((eta-etaJet)(eta-etaJet)+ (phi-phiJet)(phi-phiJet));
	378	if(rad < fConeSize && pt > fPtThresholdInCone){
f37fa8d2	379	//printf("charged in jet cone pt %f, phi %f, eta %f, R %f \n",pt,phi,eta,rad);
	380	npartcone++;
	381	fhFFz ->Fill(ptTrig, pt/ptTrig);
	382	fhFFxi->Fill(ptTrig, TMath::Log(ptTrig/pt));
	383	fhFFpt->Fill(ptTrig, pt);
477d6cee	384	}
	385	}//Tracks
	386	fhNTracksInCone->Fill(ptTrig, npartcone);
	387	}//AOD trigger particle loop
a3aebfff	388	if(GetDebug() > 1) printf("AliAnaParticleJetFinderCorrelation::MakeAnalysisFillHistograms() - End fill histograms \n");
1c5acb87	389	}
	390
	391
	392	//__________________________________________________________________
	393	void AliAnaParticleJetFinderCorrelation::Print(const Option_t * opt) const
	394	{
1c5acb87	395
477d6cee	396	//Print some relevant parameters set for the analysis
	397	if(! opt)
	398	return;
	399
a3aebfff	400	printf("** Print %s %s **\n", GetName(), GetTitle() ) ;
	401	AliAnaPartCorrBaseClass::Print(" ");
	402
477d6cee	403	printf("Phi trigger-jet < %3.2f\n", fDeltaPhiMaxCut) ;
	404	printf("Phi trigger-jet > %3.2f\n", fDeltaPhiMinCut) ;
	405	printf("pT Ratio trigger/jet < %3.2f\n", fRatioMaxCut) ;
	406	printf("pT Ratio trigger/jet > %3.2f\n", fRatioMinCut) ;
	407	printf("fConeSize = %3.2f\n", fConeSize) ;
	408	printf("fPtThresholdInCone = %3.2f\n", fPtThresholdInCone) ;
	409	printf("fUseJetRefTracks = %d\n", fUseJetRefTracks) ;
	410	printf("fMakeCorrelationInHistoMaker = %d\n", fMakeCorrelationInHistoMaker) ;
	411	printf("Isolated Trigger? %d\n", fSelectIsolated) ;
	412
1c5acb87	413	}
1c5acb87	414