[u/mrichter/AliRoot.git] / PWG4 / PartCorrBase / AliCaloTrackReader.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:  $ */

//_________________________________________________________________________
// Base class for reading data: MonteCarlo, ESD or AOD, of PHOS EMCAL and 
// Central Barrel Tracking detectors (CTS).
// Not all MC particles/tracks/clusters are kept, some kinematical/fidutial restrictions are done.
// Mother class of : AliCaloTrackESDReader: Fills ESD data in 3 TObjArrays (PHOS, EMCAL, CTS)
//                 : AliCaloTrackMCReader: Fills Kinematics data in 3 TObjArrays (PHOS, EMCAL, CTS)
//                 : AliCaloTrackAODReader: Fills AOD data in 3 TObjArrays (PHOS, EMCAL, CTS) 
//                
//-- Author: Gustavo Conesa (LNF-INFN) 
//////////////////////////////////////////////////////////////////////////////


// --- ROOT system ---
#include <TObjArray.h>
#include "TFile.h"

//---- ANALYSIS system ----
#include "AliCaloTrackReader.h"
#include "AliFidutialCut.h"
#include "AliMCEvent.h"
#include "AliAODMCHeader.h"
#include "AliGenPythiaEventHeader.h"
#include "AliAODMCParticle.h"

ClassImp(AliCaloTrackReader)
  
  
//____________________________________________________________________________
  AliCaloTrackReader::AliCaloTrackReader() : 
    TObject(), fEventNumber(-1), fCurrentFileName(""),fDataType(0), fDebug(0), 
    fFidutialCut(0x0), fComparePtHardAndJetPt(kFALSE), fPtHardAndJetPtFactor(7),
    fCTSPtMin(0), fEMCALPtMin(0),fPHOSPtMin(0),
    fAODCTS(new TObjArray()), fAODEMCAL(new TObjArray()), fAODPHOS(new TObjArray()),
    fEMCALCells(0x0), fPHOSCells(0x0),
    fInputEvent(0x0), fOutputEvent(0x0),fMC(0x0),
    fFillCTS(0),fFillEMCAL(0),fFillPHOS(0),
    fFillEMCALCells(0),fFillPHOSCells(0), 
	fSecondInputAODTree(0x0), fSecondInputAODEvent(0x0),
	fSecondInputFileName(""),fSecondInputFirstEvent(0), 
	fAODCTSNormalInputEntries(0), fAODEMCALNormalInputEntries(0), 
    fAODPHOSNormalInputEntries(0), fTrackStatus(0), 
	fReadStack(kFALSE), fReadAODMCParticles(kFALSE)
{
  //Ctor
  
  //Initialize parameters
  InitParameters();
}

//____________________________________________________________________________
AliCaloTrackReader::AliCaloTrackReader(const AliCaloTrackReader & g) :   
  TObject(g), fEventNumber(g.fEventNumber), fCurrentFileName(g.fCurrentFileName), 
  fDataType(g.fDataType), fDebug(g.fDebug),
  fFidutialCut(g.fFidutialCut),
  fComparePtHardAndJetPt(g.fComparePtHardAndJetPt),
  fPtHardAndJetPtFactor(g.fPtHardAndJetPtFactor),
  fCTSPtMin(g.fCTSPtMin), fEMCALPtMin(g.fEMCALPtMin),fPHOSPtMin(g.fPHOSPtMin), 
  fAODCTS(new TObjArray(*g.fAODCTS)),  
  fAODEMCAL(new TObjArray(*g.fAODEMCAL)),
  fAODPHOS(new TObjArray(*g.fAODPHOS)),
  fEMCALCells(new TNamed(*g.fEMCALCells)),
  fPHOSCells(new TNamed(*g.fPHOSCells)),
  fInputEvent(g.fInputEvent), fOutputEvent(g.fOutputEvent), fMC(g.fMC),
  fFillCTS(g.fFillCTS),fFillEMCAL(g.fFillEMCAL),fFillPHOS(g.fFillPHOS),
  fFillEMCALCells(g.fFillEMCALCells),fFillPHOSCells(g.fFillPHOSCells),
  fSecondInputAODTree(g.fSecondInputAODTree), 
  fSecondInputAODEvent(g.fSecondInputAODEvent),
  fSecondInputFileName(g.fSecondInputFileName), 
  fSecondInputFirstEvent(g.fSecondInputFirstEvent),
  fAODCTSNormalInputEntries(g.fAODCTSNormalInputEntries), 
  fAODEMCALNormalInputEntries(g.fAODEMCALNormalInputEntries), 
  fAODPHOSNormalInputEntries(g.fAODPHOSNormalInputEntries),
  fTrackStatus(g.fTrackStatus),
  fReadStack(g.fReadStack), fReadAODMCParticles(g.fReadAODMCParticles)
{
  // cpy ctor
  
}

//_________________________________________________________________________
AliCaloTrackReader & AliCaloTrackReader::operator = (const AliCaloTrackReader & source)
{
  // assignment operator
  
  if(&source == this) return *this;
  
  fDataType    = source.fDataType ;
  fDebug       = source.fDebug ;
  fEventNumber = source.fEventNumber ;
  fCurrentFileName = source.fCurrentFileName ;
  fFidutialCut = source.fFidutialCut;
	
  fComparePtHardAndJetPt = source.fComparePtHardAndJetPt;
  fPtHardAndJetPtFactor  = source.fPtHardAndJetPtFactor;
	
  fCTSPtMin    = source.fCTSPtMin ;
  fEMCALPtMin  = source.fEMCALPtMin ;
  fPHOSPtMin   = source.fPHOSPtMin ; 
  
  fAODCTS     = new TObjArray(*source.fAODCTS) ;
  fAODEMCAL   = new TObjArray(*source.fAODEMCAL) ;
  fAODPHOS    = new TObjArray(*source.fAODPHOS) ;
  fEMCALCells = new TNamed(*source.fEMCALCells) ;
  fPHOSCells  = new TNamed(*source.fPHOSCells) ;

  fInputEvent  = source.fInputEvent;
  fOutputEvent = source.fOutputEvent;
  fMC          = source.fMC;
  
  fFillCTS        = source.fFillCTS;
  fFillEMCAL      = source.fFillEMCAL;
  fFillPHOS       = source.fFillPHOS;
  fFillEMCALCells = source.fFillEMCALCells;
  fFillPHOSCells  = source.fFillPHOSCells;

  fSecondInputAODTree    = source.fSecondInputAODTree;
  fSecondInputAODEvent   = source.fSecondInputAODEvent;
  fSecondInputFileName   = source.fSecondInputFileName;
  fSecondInputFirstEvent = source.fSecondInputFirstEvent;

  fAODCTSNormalInputEntries   = source.fAODCTSNormalInputEntries; 
  fAODEMCALNormalInputEntries = source.fAODEMCALNormalInputEntries; 
  fAODPHOSNormalInputEntries  = source.fAODPHOSNormalInputEntries;
	
  fTrackStatus        = source.fTrackStatus;
  fReadStack          = source.fReadStack;
  fReadAODMCParticles = source.fReadAODMCParticles;	
	
  return *this;
  
}

//_________________________________
AliCaloTrackReader::~AliCaloTrackReader() {
  //Dtor
  
  if(fFidutialCut) delete fFidutialCut ;
	
  if(fAODCTS){
    fAODCTS->Clear() ; 
    delete fAODCTS ;
  }
  
  if(fAODEMCAL){
    fAODEMCAL->Clear() ; 
    delete fAODEMCAL ;
  }
  
  if(fAODPHOS){
    fAODPHOS->Clear() ; 
    delete fAODPHOS ;
  }
  
  if(fEMCALCells){
    fEMCALCells->Clear() ; 
    delete fEMCALCells ;
  }
  
  if(fPHOSCells){
    fPHOSCells->Clear() ; 
    delete fPHOSCells ;
  }

  if(fInputEvent)  delete fInputEvent ;
  if(fOutputEvent) delete fOutputEvent ;
  if(fMC)          delete fMC ;  
	
  if(fSecondInputAODTree){
	fSecondInputAODTree->Clear();
	delete fSecondInputAODTree;
  }
	
  if(fSecondInputAODEvent) delete fSecondInputAODEvent ;

}


//_________________________________________________________________________
Bool_t AliCaloTrackReader::ComparePtHardAndJetPt(){
	// Check the event, if the requested ptHard is much larger than the jet pT, then there is a problem.
	// Only for PYTHIA.
	if(!fReadStack) return kTRUE; //Information not filtered to AOD
	
	if(!strcmp(GetGenEventHeader()->ClassName(), "AliGenPythiaEventHeader")){
		TParticle * jet =  new TParticle;
		AliGenPythiaEventHeader* pygeh= (AliGenPythiaEventHeader*) GetGenEventHeader();
		Int_t nTriggerJets =  pygeh->NTriggerJets();
		Float_t ptHard = pygeh->GetPtHard();
		
		//if(fDebug > 1) printf("AliMCAnalysisUtils::PythiaEventHeader: Njets: %d, pT Hard %f\n",nTriggerJets, ptHard);
	    Float_t tmpjet[]={0,0,0,0};
		for(Int_t ijet = 0; ijet< nTriggerJets; ijet++){
			pygeh->TriggerJet(ijet, tmpjet);
			jet = new TParticle(94, 21, -1, -1, -1, -1, tmpjet[0],tmpjet[1],tmpjet[2],tmpjet[3], 0,0,0,0);
			//Compare jet pT and pt Hard
			//if(fDebug > 1) printf("AliMCAnalysisUtils:: %d pycell jet pT %f\n",ijet, jet->Pt());
			if(jet->Pt() > fPtHardAndJetPtFactor * ptHard) {
				printf("AliMCAnalysisUtils::PythiaEventHeader: Njets: %d, pT Hard %2.2f, pycell jet pT %2.2f, rejection factor %1.1f\n",
					   nTriggerJets, ptHard, jet->Pt(), fPtHardAndJetPtFactor);
				return kFALSE;
			}
		}
	}
	
	return kTRUE ;
	
}

//____________________________________________________________________________
AliStack* AliCaloTrackReader::GetStack() const {
  //Return pointer to stack
  if(fMC)
    return fMC->Stack();
  else{
    if(fDebug > 1) printf("AliCaloTrackReader::GetStack() - Stack is not available\n"); 
    return 0x0 ;
  }
}

//____________________________________________________________________________
AliHeader* AliCaloTrackReader::GetHeader() const {
  //Return pointer to header
  if(fMC)
    return fMC->Header();
  else{
    printf("AliCaloTrackReader::Header is not available\n"); 
    return 0x0 ;
  }
}
//____________________________________________________________________________
AliGenEventHeader* AliCaloTrackReader::GetGenEventHeader() const {
  //Return pointer to Generated event header
  if(fMC)
    return fMC->GenEventHeader();
  else{
    printf("AliCaloTrackReader::GenEventHeader is not available\n"); 
    return 0x0 ;
  }
}

//____________________________________________________________________________
TClonesArray* AliCaloTrackReader::GetAODMCParticles(Int_t input) const {
	//Return list of particles in AOD. Do it for the corresponding input event.
	if(fDataType == kAOD){
	 //Normal input AOD
	 if(input == 0) return (TClonesArray*)((AliAODEvent*)fInputEvent)->FindListObject("mcparticles");
	  //Second input AOD
	 else if(input == 1 && fSecondInputAODEvent) return (TClonesArray*) fSecondInputAODEvent->FindListObject("mcparticles");	
	 else {
	     printf("AliCaloTrackReader::GetAODMCParticles() - wrong AOD input index? %d, or non existing tree? \n",input); 
		 return 0x0;
	 }
	}
	else {
		printf("AliCaloTrackReader::GetAODMCParticles() - Input are not AODs\n"); 
		return 0x0;
	}
}

//____________________________________________________________________________
AliAODMCHeader* AliCaloTrackReader::GetAODMCHeader(Int_t input) const {
	//Return MC header in AOD. Do it for the corresponding input event.
	if(fDataType == kAOD){
		//Normal input AOD
		if(input == 0) return (AliAODMCHeader*)((AliAODEvent*)fInputEvent)->FindListObject("mcheader");
		//Second input AOD
		else if(input == 1) return  (AliAODMCHeader*) fSecondInputAODEvent->FindListObject("mcheader");	
		else {
			printf("AliCaloTrackReader::GetAODMCHeader() - wrong AOD input index, %d\n",input);
			return 0x0;
		}
	}
	else {
		printf("AliCaloTrackReader::GetAODMCHeader() - Input are not AODs\n");
		return 0x0;
	}
}

//_______________________________________________________________
void AliCaloTrackReader::Init()
{
	//Init reader. Method to be called in AliAnaPartCorrMaker
	
	//Get the file with second input events if the filename is given
	//Get the tree and connect the AODEvent. Only with AODs

	if(fReadStack && fReadAODMCParticles){
		printf("AliCaloTrackReader::Init() - Cannot access stack and mcparticles at the same time, change them \n");
		fReadStack = kFALSE;
		fReadAODMCParticles = kFALSE;
	}
	
	if(fSecondInputFileName!=""){
		if(fDataType == kAOD){
			TFile * input2   = new TFile(fSecondInputFileName,"read");
			printf("AliCaloTrackReader::Init() - Second input file opened: %s, size %d \n", input2->GetName(), (Int_t) input2->GetSize());
			fSecondInputAODTree = (TTree*) input2->Get("aodTree");
			if(fSecondInputAODTree) printf("AliCaloTrackReader::Init() - Second input tree opened: %s, entries %d \n", 
										   fSecondInputAODTree->GetName(), (Int_t) fSecondInputAODTree->GetEntries());
			else{
			 printf("AliCaloTrackReader::Init() - Second input tree not available, STOP \n");
			 abort();
			}
			fSecondInputAODEvent = new AliAODEvent;
			fSecondInputAODEvent->ReadFromTree(fSecondInputAODTree);
			if(fSecondInputFirstEvent >= fSecondInputAODTree->GetEntriesFast()){
				printf("AliCaloTrackReader::Init() - Requested first event of second input %d, is larger than number of events %d, STOP\n", 
					   fSecondInputFirstEvent, (Int_t) fSecondInputAODTree->GetEntriesFast());
				abort();
			}
		}
		else printf("AliCaloTrackReader::Init() - Second input not added, reader is not AOD\n");
	}
	
	
}
//_______________________________________________________________
void AliCaloTrackReader::InitParameters()
{
 
  //Initialize the parameters of the analysis.
  fDataType   = kESD ;
  fCTSPtMin   = 0.2 ;
  fEMCALPtMin = 0.2 ;
  fPHOSPtMin  = 0.2 ;

  fFillEMCAL      = kTRUE;
  fFillPHOS       = kTRUE;
  fFillCTS        = kTRUE;
  fFillEMCALCells = kFALSE;
  fFillPHOSCells  = kFALSE;

  fFidutialCut           = new AliFidutialCut();
  fSecondInputFileName   = "" ;
  fSecondInputFirstEvent = 0 ;
  fReadStack             = kFALSE; // Check in the constructor of the other readers if it was set or in the configuration file
  fReadAODMCParticles    = kFALSE; // Check in the constructor of the other readers if it was set or in the configuration file
	
}


//________________________________________________________________
void AliCaloTrackReader::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("Data type      : %d\n", fDataType) ;
  printf("CTS Min pT     : %2.1f GeV/c\n", fCTSPtMin) ;
  printf("EMCAL Min pT   : %2.1f GeV/c\n", fEMCALPtMin) ;
  printf("PHOS Min pT    : %2.1f GeV/c\n", fPHOSPtMin) ;
  printf("Use CTS         =     %d\n", fFillCTS) ;
  printf("Use EMCAL       =     %d\n", fFillEMCAL) ;
  printf("Use PHOS        =     %d\n", fFillPHOS) ;
  printf("Use EMCAL Cells =     %d\n", fFillEMCALCells) ;
  printf("Use PHOS  Cells =     %d\n", fFillPHOSCells) ;
  printf("Track status    =     %d\n", (Int_t) fTrackStatus) ;
  if(fComparePtHardAndJetPt)
	  printf("Compare jet pt and pt hard to accept event, factor = %2.2f",fPtHardAndJetPtFactor);
	
  if(fSecondInputFileName!="") {
	  printf("Second Input File Name     =     %s\n", fSecondInputFileName.Data()) ;
	  printf("Second Input First Event   =     %d\n", fSecondInputFirstEvent) ;
  }
	
	 printf("Read Kine from, stack? %d, AOD ? %d \n", fReadStack, fReadAODMCParticles) ;
	
  printf("    \n") ;
} 

//___________________________________________________
Bool_t AliCaloTrackReader::FillInputEvent(const Int_t iEntry, const char * currentFileName) {
  //Fill the event counter and input lists that are needed, called by the analysis maker.
  
  fEventNumber = iEntry;
  fCurrentFileName = TString(currentFileName);
	
  if((fDataType != kAOD) && ((fOutputEvent->GetCaloClusters())->GetEntriesFast()!=0 ||(fOutputEvent->GetTracks())->GetEntriesFast()!=0)){
    printf("AliCaloTrackReader::AODCaloClusters or AODTracks already filled by the filter, do not use the ESD reader, use the AOD reader, STOP\n");
    abort();
  }
	
  //In case of analysis of events with jets, skip those with jet pt > 5 pt hard	
  if(fComparePtHardAndJetPt && GetStack()) {
		if(!ComparePtHardAndJetPt()) return kFALSE ;
  }

  //In case of mixing events with other AOD file	
  if(fDataType == kAOD && fSecondInputAODTree){
	 
	 if(fDebug > 1) 
		 printf("AliCaloTrackReader::FillInputEvent() - Get event %d from second input AOD file \n", iEntry+fSecondInputFirstEvent);
	 if(fSecondInputAODTree->GetEntriesFast() <= iEntry+fSecondInputFirstEvent) {
		 if(fSecondInputAODTree->GetEntriesFast() == iEntry+fSecondInputFirstEvent) 
			 printf("AliCaloTrackReader::FillInputEvent() - Skip events from event %d, no more events in second AOD file \n", iEntry);
		 return kFALSE;
	 }
	  
	 //Get the Event
	 Int_t nbytes = fSecondInputAODTree->GetEvent(iEntry+fSecondInputFirstEvent);
	 if ( nbytes == 0 ) {//If nothing in AOD
		 printf("AliCaloTrackReader::FillInputEvent() - Nothing in Second AOD input, STOP\n");
		 abort() ; 
	 }
	  
  }
  //if(iEntry > 10) return kFALSE;
  if(fFillCTS)   FillInputCTS();
  if(fFillEMCAL) FillInputEMCAL();
  if(fFillPHOS)  FillInputPHOS();
  if(fFillEMCALCells) FillInputEMCALCells();
  if(fFillPHOSCells)  FillInputPHOSCells();

  return kTRUE ;
}

//__________________________________________________
void AliCaloTrackReader::ResetLists() {
  //  Reset lists, called by the analysis maker 

  if(fAODCTS)   fAODCTS -> Clear();
  if(fAODEMCAL) fAODEMCAL -> Clear();
  if(fAODPHOS)  fAODPHOS -> Clear();
  if(fEMCALCells) fEMCALCells -> Clear();
  if(fPHOSCells)  fPHOSCells -> Clear();

}
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: $ */
	16
	17	//_________________________________________________________________________
	18	// Base class for reading data: MonteCarlo, ESD or AOD, of PHOS EMCAL and
	19	// Central Barrel Tracking detectors (CTS).
	20	// Not all MC particles/tracks/clusters are kept, some kinematical/fidutial restrictions are done.
591cc579	21	// Mother class of : AliCaloTrackESDReader: Fills ESD data in 3 TObjArrays (PHOS, EMCAL, CTS)
	22	// : AliCaloTrackMCReader: Fills Kinematics data in 3 TObjArrays (PHOS, EMCAL, CTS)
	23	// : AliCaloTrackAODReader: Fills AOD data in 3 TObjArrays (PHOS, EMCAL, CTS)
1c5acb87	24	//
	25	//-- Author: Gustavo Conesa (LNF-INFN)
	26	//////////////////////////////////////////////////////////////////////////////
	27
	28
	29	// --- ROOT system ---
591cc579	30	#include <TObjArray.h>
591cc579	31	#include "TFile.h"
1c5acb87	32
	33	//---- ANALYSIS system ----
	34	#include "AliCaloTrackReader.h"
477d6cee	35	#include "AliFidutialCut.h"
477d6cee	36	#include "AliMCEvent.h"
591cc579	37	#include "AliAODMCHeader.h"
	38	#include "AliGenPythiaEventHeader.h"
	39	#include "AliAODMCParticle.h"
1c5acb87	40
	41	ClassImp(AliCaloTrackReader)
	42
	43
	44	//____________________________________________________________________________
	45	AliCaloTrackReader::AliCaloTrackReader() :
a79a2424	46	TObject(), fEventNumber(-1), fCurrentFileName(""),fDataType(0), fDebug(0),
591cc579	47	fFidutialCut(0x0), fComparePtHardAndJetPt(kFALSE), fPtHardAndJetPtFactor(7),
1c5acb87	48	fCTSPtMin(0), fEMCALPtMin(0),fPHOSPtMin(0),
591cc579	49	fAODCTS(new TObjArray()), fAODEMCAL(new TObjArray()), fAODPHOS(new TObjArray()),
1c5acb87	50	fEMCALCells(0x0), fPHOSCells(0x0),
477d6cee	51	fInputEvent(0x0), fOutputEvent(0x0),fMC(0x0),
1c5acb87	52	fFillCTS(0),fFillEMCAL(0),fFillPHOS(0),
591cc579	53	fFillEMCALCells(0),fFillPHOSCells(0),
	54	fSecondInputAODTree(0x0), fSecondInputAODEvent(0x0),
	55	fSecondInputFileName(""),fSecondInputFirstEvent(0),
	56	fAODCTSNormalInputEntries(0), fAODEMCALNormalInputEntries(0),
	57	fAODPHOSNormalInputEntries(0), fTrackStatus(0),
	58	fReadStack(kFALSE), fReadAODMCParticles(kFALSE)
1c5acb87	59	{
	60	//Ctor
	61
	62	//Initialize parameters
	63	InitParameters();
	64	}
	65
	66	//____________________________________________________________________________
	67	AliCaloTrackReader::AliCaloTrackReader(const AliCaloTrackReader & g) :
a79a2424	68	TObject(g), fEventNumber(g.fEventNumber), fCurrentFileName(g.fCurrentFileName),
a79a2424	69	fDataType(g.fDataType), fDebug(g.fDebug),
6639984f	70	fFidutialCut(g.fFidutialCut),
591cc579	71	fComparePtHardAndJetPt(g.fComparePtHardAndJetPt),
591cc579	72	fPtHardAndJetPtFactor(g.fPtHardAndJetPtFactor),
1c5acb87	73	fCTSPtMin(g.fCTSPtMin), fEMCALPtMin(g.fEMCALPtMin),fPHOSPtMin(g.fPHOSPtMin),
591cc579	74	fAODCTS(new TObjArray(*g.fAODCTS)),
	75	fAODEMCAL(new TObjArray(*g.fAODEMCAL)),
	76	fAODPHOS(new TObjArray(*g.fAODPHOS)),
1c5acb87	77	fEMCALCells(new TNamed(*g.fEMCALCells)),
1c5acb87	78	fPHOSCells(new TNamed(*g.fPHOSCells)),
477d6cee	79	fInputEvent(g.fInputEvent), fOutputEvent(g.fOutputEvent), fMC(g.fMC),
1c5acb87	80	fFillCTS(g.fFillCTS),fFillEMCAL(g.fFillEMCAL),fFillPHOS(g.fFillPHOS),
591cc579	81	fFillEMCALCells(g.fFillEMCALCells),fFillPHOSCells(g.fFillPHOSCells),
	82	fSecondInputAODTree(g.fSecondInputAODTree),
	83	fSecondInputAODEvent(g.fSecondInputAODEvent),
	84	fSecondInputFileName(g.fSecondInputFileName),
	85	fSecondInputFirstEvent(g.fSecondInputFirstEvent),
	86	fAODCTSNormalInputEntries(g.fAODCTSNormalInputEntries),
	87	fAODEMCALNormalInputEntries(g.fAODEMCALNormalInputEntries),
	88	fAODPHOSNormalInputEntries(g.fAODPHOSNormalInputEntries),
	89	fTrackStatus(g.fTrackStatus),
	90	fReadStack(g.fReadStack), fReadAODMCParticles(g.fReadAODMCParticles)
1c5acb87	91	{
	92	// cpy ctor
	93
	94	}
	95
	96	//_________________________________________________________________________
	97	AliCaloTrackReader & AliCaloTrackReader::operator = (const AliCaloTrackReader & source)
	98	{
	99	// assignment operator
	100
	101	if(&source == this) return *this;
	102
6639984f	103	fDataType = source.fDataType ;
	104	fDebug = source.fDebug ;
	105	fEventNumber = source.fEventNumber ;
a79a2424	106	fCurrentFileName = source.fCurrentFileName ;
1c5acb87	107	fFidutialCut = source.fFidutialCut;
591cc579	108
29b2ceec	109	fComparePtHardAndJetPt = source.fComparePtHardAndJetPt;
591cc579	110	fPtHardAndJetPtFactor = source.fPtHardAndJetPtFactor;
591cc579	111
29b2ceec	112	fCTSPtMin = source.fCTSPtMin ;
	113	fEMCALPtMin = source.fEMCALPtMin ;
	114	fPHOSPtMin = source.fPHOSPtMin ;
1c5acb87	115
591cc579	116	fAODCTS = new TObjArray(*source.fAODCTS) ;
	117	fAODEMCAL = new TObjArray(*source.fAODEMCAL) ;
	118	fAODPHOS = new TObjArray(*source.fAODPHOS) ;
1c5acb87	119	fEMCALCells = new TNamed(*source.fEMCALCells) ;
6639984f	120	fPHOSCells = new TNamed(*source.fPHOSCells) ;
1c5acb87	121
477d6cee	122	fInputEvent = source.fInputEvent;
	123	fOutputEvent = source.fOutputEvent;
	124	fMC = source.fMC;
1c5acb87	125
6639984f	126	fFillCTS = source.fFillCTS;
	127	fFillEMCAL = source.fFillEMCAL;
	128	fFillPHOS = source.fFillPHOS;
1c5acb87	129	fFillEMCALCells = source.fFillEMCALCells;
6639984f	130	fFillPHOSCells = source.fFillPHOSCells;
1c5acb87	131
591cc579	132	fSecondInputAODTree = source.fSecondInputAODTree;
	133	fSecondInputAODEvent = source.fSecondInputAODEvent;
	134	fSecondInputFileName = source.fSecondInputFileName;
	135	fSecondInputFirstEvent = source.fSecondInputFirstEvent;
	136
	137	fAODCTSNormalInputEntries = source.fAODCTSNormalInputEntries;
	138	fAODEMCALNormalInputEntries = source.fAODEMCALNormalInputEntries;
	139	fAODPHOSNormalInputEntries = source.fAODPHOSNormalInputEntries;
	140
	141	fTrackStatus = source.fTrackStatus;
	142	fReadStack = source.fReadStack;
	143	fReadAODMCParticles = source.fReadAODMCParticles;
	144
1c5acb87	145	return *this;
	146
	147	}
	148
	149	//_________________________________
	150	AliCaloTrackReader::~AliCaloTrackReader() {
	151	//Dtor
	152
	153	if(fFidutialCut) delete fFidutialCut ;
29b2ceec	154
1c5acb87	155	if(fAODCTS){
	156	fAODCTS->Clear() ;
	157	delete fAODCTS ;
	158	}
	159
	160	if(fAODEMCAL){
	161	fAODEMCAL->Clear() ;
	162	delete fAODEMCAL ;
	163	}
	164
	165	if(fAODPHOS){
	166	fAODPHOS->Clear() ;
	167	delete fAODPHOS ;
	168	}
	169
	170	if(fEMCALCells){
	171	fEMCALCells->Clear() ;
	172	delete fEMCALCells ;
	173	}
	174
	175	if(fPHOSCells){
	176	fPHOSCells->Clear() ;
	177	delete fPHOSCells ;
	178	}
	179
477d6cee	180	if(fInputEvent) delete fInputEvent ;
477d6cee	181	if(fOutputEvent) delete fOutputEvent ;
591cc579	182	if(fMC) delete fMC ;
	183
	184	if(fSecondInputAODTree){
	185	fSecondInputAODTree->Clear();
	186	delete fSecondInputAODTree;
	187	}
	188
	189	if(fSecondInputAODEvent) delete fSecondInputAODEvent ;
	190
1c5acb87	191	}
1c5acb87	192
477d6cee	193
591cc579	194	//_________________________________________________________________________
	195	Bool_t AliCaloTrackReader::ComparePtHardAndJetPt(){
	196	// Check the event, if the requested ptHard is much larger than the jet pT, then there is a problem.
	197	// Only for PYTHIA.
	198	if(!fReadStack) return kTRUE; //Information not filtered to AOD
	199
	200	if(!strcmp(GetGenEventHeader()->ClassName(), "AliGenPythiaEventHeader")){
	201	TParticle * jet = new TParticle;
	202	AliGenPythiaEventHeader* pygeh= (AliGenPythiaEventHeader*) GetGenEventHeader();
	203	Int_t nTriggerJets = pygeh->NTriggerJets();
	204	Float_t ptHard = pygeh->GetPtHard();
	205
	206	//if(fDebug > 1) printf("AliMCAnalysisUtils::PythiaEventHeader: Njets: %d, pT Hard %f\n",nTriggerJets, ptHard);
	207	Float_t tmpjet[]={0,0,0,0};
	208	for(Int_t ijet = 0; ijet< nTriggerJets; ijet++){
	209	pygeh->TriggerJet(ijet, tmpjet);
	210	jet = new TParticle(94, 21, -1, -1, -1, -1, tmpjet[0],tmpjet[1],tmpjet[2],tmpjet[3], 0,0,0,0);
	211	//Compare jet pT and pt Hard
	212	//if(fDebug > 1) printf("AliMCAnalysisUtils:: %d pycell jet pT %f\n",ijet, jet->Pt());
	213	if(jet->Pt() > fPtHardAndJetPtFactor * ptHard) {
	214	printf("AliMCAnalysisUtils::PythiaEventHeader: Njets: %d, pT Hard %2.2f, pycell jet pT %2.2f, rejection factor %1.1f\n",
	215	nTriggerJets, ptHard, jet->Pt(), fPtHardAndJetPtFactor);
	216	return kFALSE;
	217	}
	218	}
	219	}
	220
	221	return kTRUE ;
	222
	223	}
	224
1c5acb87	225	//____________________________________________________________________________
	226	AliStack* AliCaloTrackReader::GetStack() const {
	227	//Return pointer to stack
	228	if(fMC)
	229	return fMC->Stack();
	230	else{
477d6cee	231	if(fDebug > 1) printf("AliCaloTrackReader::GetStack() - Stack is not available\n");
1c5acb87	232	return 0x0 ;
	233	}
	234	}
	235
	236	//____________________________________________________________________________
	237	AliHeader* AliCaloTrackReader::GetHeader() const {
	238	//Return pointer to header
	239	if(fMC)
	240	return fMC->Header();
	241	else{
477d6cee	242	printf("AliCaloTrackReader::Header is not available\n");
1c5acb87	243	return 0x0 ;
	244	}
	245	}
	246	//____________________________________________________________________________
	247	AliGenEventHeader* AliCaloTrackReader::GetGenEventHeader() const {
	248	//Return pointer to Generated event header
	249	if(fMC)
	250	return fMC->GenEventHeader();
	251	else{
477d6cee	252	printf("AliCaloTrackReader::GenEventHeader is not available\n");
1c5acb87	253	return 0x0 ;
	254	}
	255	}
	256
591cc579	257	//____________________________________________________________________________
	258	TClonesArray* AliCaloTrackReader::GetAODMCParticles(Int_t input) const {
	259	//Return list of particles in AOD. Do it for the corresponding input event.
	260	if(fDataType == kAOD){
	261	//Normal input AOD
	262	if(input == 0) return (TClonesArray)((AliAODEvent)fInputEvent)->FindListObject("mcparticles");
	263	//Second input AOD
	264	else if(input == 1 && fSecondInputAODEvent) return (TClonesArray*) fSecondInputAODEvent->FindListObject("mcparticles");
	265	else {
	266	printf("AliCaloTrackReader::GetAODMCParticles() - wrong AOD input index? %d, or non existing tree? \n",input);
	267	return 0x0;
	268	}
	269	}
	270	else {
	271	printf("AliCaloTrackReader::GetAODMCParticles() - Input are not AODs\n");
	272	return 0x0;
	273	}
	274	}
	275
	276	//____________________________________________________________________________
	277	AliAODMCHeader* AliCaloTrackReader::GetAODMCHeader(Int_t input) const {
	278	//Return MC header in AOD. Do it for the corresponding input event.
	279	if(fDataType == kAOD){
	280	//Normal input AOD
	281	if(input == 0) return (AliAODMCHeader)((AliAODEvent)fInputEvent)->FindListObject("mcheader");
	282	//Second input AOD
	283	else if(input == 1) return (AliAODMCHeader*) fSecondInputAODEvent->FindListObject("mcheader");
	284	else {
	285	printf("AliCaloTrackReader::GetAODMCHeader() - wrong AOD input index, %d\n",input);
	286	return 0x0;
	287	}
	288	}
	289	else {
	290	printf("AliCaloTrackReader::GetAODMCHeader() - Input are not AODs\n");
	291	return 0x0;
	292	}
	293	}
	294
	295	//_______________________________________________________________
	296	void AliCaloTrackReader::Init()
	297	{
	298	//Init reader. Method to be called in AliAnaPartCorrMaker
	299
	300	//Get the file with second input events if the filename is given
	301	//Get the tree and connect the AODEvent. Only with AODs
f1f0bd84	302
591cc579	303	if(fReadStack && fReadAODMCParticles){
	304	printf("AliCaloTrackReader::Init() - Cannot access stack and mcparticles at the same time, change them \n");
	305	fReadStack = kFALSE;
	306	fReadAODMCParticles = kFALSE;
	307	}
	308
591cc579	309	if(fSecondInputFileName!=""){
	310	if(fDataType == kAOD){
	311	TFile * input2 = new TFile(fSecondInputFileName,"read");
ddcb0d54	312	printf("AliCaloTrackReader::Init() - Second input file opened: %s, size %d \n", input2->GetName(), (Int_t) input2->GetSize());
591cc579	313	fSecondInputAODTree = (TTree*) input2->Get("aodTree");
ddcb0d54	314	if(fSecondInputAODTree) printf("AliCaloTrackReader::Init() - Second input tree opened: %s, entries %d \n",
	315	fSecondInputAODTree->GetName(), (Int_t) fSecondInputAODTree->GetEntries());
	316	else{
	317	printf("AliCaloTrackReader::Init() - Second input tree not available, STOP \n");
	318	abort();
	319	}
591cc579	320	fSecondInputAODEvent = new AliAODEvent;
	321	fSecondInputAODEvent->ReadFromTree(fSecondInputAODTree);
	322	if(fSecondInputFirstEvent >= fSecondInputAODTree->GetEntriesFast()){
	323	printf("AliCaloTrackReader::Init() - Requested first event of second input %d, is larger than number of events %d, STOP\n",
	324	fSecondInputFirstEvent, (Int_t) fSecondInputAODTree->GetEntriesFast());
	325	abort();
	326	}
	327	}
	328	else printf("AliCaloTrackReader::Init() - Second input not added, reader is not AOD\n");
	329	}
	330
	331
	332	}
1c5acb87	333	//_______________________________________________________________
	334	void AliCaloTrackReader::InitParameters()
	335	{
	336
	337	//Initialize the parameters of the analysis.
591cc579	338	fDataType = kESD ;
1c5acb87	339	fCTSPtMin = 0.2 ;
29b2ceec	340	fEMCALPtMin = 0.2 ;
29b2ceec	341	fPHOSPtMin = 0.2 ;
1c5acb87	342
591cc579	343	fFillEMCAL = kTRUE;
	344	fFillPHOS = kTRUE;
	345	fFillCTS = kTRUE;
1c5acb87	346	fFillEMCALCells = kFALSE;
591cc579	347	fFillPHOSCells = kFALSE;
1c5acb87	348
591cc579	349	fFidutialCut = new AliFidutialCut();
	350	fSecondInputFileName = "" ;
	351	fSecondInputFirstEvent = 0 ;
	352	fReadStack = kFALSE; // Check in the constructor of the other readers if it was set or in the configuration file
	353	fReadAODMCParticles = kFALSE; // Check in the constructor of the other readers if it was set or in the configuration file
	354
1c5acb87	355	}
	356
	357
	358	//________________________________________________________________
	359	void AliCaloTrackReader::Print(const Option_t * opt) const
	360	{
	361
	362	//Print some relevant parameters set for the analysis
	363	if(! opt)
	364	return;
	365
	366	printf("*** Print: %s %s ****\n", GetName(), GetTitle() ) ;
	367	printf("Data type : %d\n", fDataType) ;
	368	printf("CTS Min pT : %2.1f GeV/c\n", fCTSPtMin) ;
	369	printf("EMCAL Min pT : %2.1f GeV/c\n", fEMCALPtMin) ;
	370	printf("PHOS Min pT : %2.1f GeV/c\n", fPHOSPtMin) ;
	371	printf("Use CTS = %d\n", fFillCTS) ;
	372	printf("Use EMCAL = %d\n", fFillEMCAL) ;
	373	printf("Use PHOS = %d\n", fFillPHOS) ;
	374	printf("Use EMCAL Cells = %d\n", fFillEMCALCells) ;
	375	printf("Use PHOS Cells = %d\n", fFillPHOSCells) ;
591cc579	376	printf("Track status = %d\n", (Int_t) fTrackStatus) ;
	377	if(fComparePtHardAndJetPt)
	378	printf("Compare jet pt and pt hard to accept event, factor = %2.2f",fPtHardAndJetPtFactor);
	379
	380	if(fSecondInputFileName!="") {
	381	printf("Second Input File Name = %s\n", fSecondInputFileName.Data()) ;
	382	printf("Second Input First Event = %d\n", fSecondInputFirstEvent) ;
	383	}
	384
	385	printf("Read Kine from, stack? %d, AOD ? %d \n", fReadStack, fReadAODMCParticles) ;
29b2ceec	386
1c5acb87	387	printf(" \n") ;
	388	}
	389
	390	//___________________________________________________
29b2ceec	391	Bool_t AliCaloTrackReader::FillInputEvent(const Int_t iEntry, const char * currentFileName) {
6639984f	392	//Fill the event counter and input lists that are needed, called by the analysis maker.
	393
	394	fEventNumber = iEntry;
a79a2424	395	fCurrentFileName = TString(currentFileName);
a79a2424	396
477d6cee	397	if((fDataType != kAOD) && ((fOutputEvent->GetCaloClusters())->GetEntriesFast()!=0 \|\|(fOutputEvent->GetTracks())->GetEntriesFast()!=0)){
29b2ceec	398	printf("AliCaloTrackReader::AODCaloClusters or AODTracks already filled by the filter, do not use the ESD reader, use the AOD reader, STOP\n");
477d6cee	399	abort();
477d6cee	400	}
29b2ceec	401
	402	//In case of analysis of events with jets, skip those with jet pt > 5 pt hard
	403	if(fComparePtHardAndJetPt && GetStack()) {
591cc579	404	if(!ComparePtHardAndJetPt()) return kFALSE ;
29b2ceec	405	}
591cc579	406
	407	//In case of mixing events with other AOD file
	408	if(fDataType == kAOD && fSecondInputAODTree){
	409
	410	if(fDebug > 1)
	411	printf("AliCaloTrackReader::FillInputEvent() - Get event %d from second input AOD file \n", iEntry+fSecondInputFirstEvent);
	412	if(fSecondInputAODTree->GetEntriesFast() <= iEntry+fSecondInputFirstEvent) {
	413	if(fSecondInputAODTree->GetEntriesFast() == iEntry+fSecondInputFirstEvent)
	414	printf("AliCaloTrackReader::FillInputEvent() - Skip events from event %d, no more events in second AOD file \n", iEntry);
	415	return kFALSE;
	416	}
	417
	418	//Get the Event
	419	Int_t nbytes = fSecondInputAODTree->GetEvent(iEntry+fSecondInputFirstEvent);
	420	if ( nbytes == 0 ) {//If nothing in AOD
	421	printf("AliCaloTrackReader::FillInputEvent() - Nothing in Second AOD input, STOP\n");
	422	abort() ;
	423	}
	424
	425	}
	426	//if(iEntry > 10) return kFALSE;
477d6cee	427	if(fFillCTS) FillInputCTS();
1c5acb87	428	if(fFillEMCAL) FillInputEMCAL();
477d6cee	429	if(fFillPHOS) FillInputPHOS();
1c5acb87	430	if(fFillEMCALCells) FillInputEMCALCells();
477d6cee	431	if(fFillPHOSCells) FillInputPHOSCells();
1c5acb87	432
29b2ceec	433	return kTRUE ;
1c5acb87	434	}
	435
	436	//__________________________________________________
	437	void AliCaloTrackReader::ResetLists() {
	438	// Reset lists, called by the analysis maker
	439
477d6cee	440	if(fAODCTS) fAODCTS -> Clear();
1c5acb87	441	if(fAODEMCAL) fAODEMCAL -> Clear();
477d6cee	442	if(fAODPHOS) fAODPHOS -> Clear();
1c5acb87	443	if(fEMCALCells) fEMCALCells -> Clear();
477d6cee	444	if(fPHOSCells) fPHOSCells -> Clear();
1c5acb87	445
1c5acb87	446	}