// fSecondInputFileName(""),fSecondInputFirstEvent(0),
// fAODCTSNormalInputEntries(0), fAODEMCALNormalInputEntries(0),
// fAODPHOSNormalInputEntries(0),
- fTrackStatus(0), fESDtrackCuts(0), fTrackMult(0), fTrackMultEtaCut(0.9),
+ fTrackStatus(0), fESDtrackCuts(0), fTrackMult(0), fTrackMultEtaCut(0.8),
fReadStack(kFALSE), fReadAODMCParticles(kFALSE),
fDeltaAODFileName("deltaAODPartCorr.root"),fFiredTriggerClassName(""),
fAnaLED(kFALSE),fTaskName(""),fCaloUtils(0x0),
fMixedEvent(NULL), fNMixedEvent(1), fVertex(NULL),
- fWriteOutputDeltaAOD(kFALSE),fOldAOD(kFALSE){
+ fWriteOutputDeltaAOD(kFALSE),fOldAOD(kFALSE),fCaloFilterPatch(kFALSE),
+ fEMCALClustersListName(""),fZvtxCut(0.)
+{
//Ctor
//Initialize parameters
}
if(fAODEMCAL){
- if(fDataType!=kMC)fAODEMCAL->Clear() ;
+ if(fDataType!=kMC)fAODEMCAL->Clear("C") ;
else fAODEMCAL->Delete() ;
delete fAODEMCAL ;
}
if(fAODPHOS){
- if(fDataType!=kMC)fAODPHOS->Clear() ;
+ if(fDataType!=kMC)fAODPHOS->Clear("C") ;
else fAODPHOS->Delete() ;
delete fAODPHOS ;
}
- if(fEMCALCells){
- fEMCALCells->Clear() ;
- delete fEMCALCells ;
- }
-
- if(fPHOSCells){
- fPHOSCells->Clear() ;
- delete fPHOSCells ;
- }
+// if(fEMCALCells){
+// delete fEMCALCells ;
+// }
+//
+// if(fPHOSCells){
+// delete fPHOSCells ;
+// }
if(fVertex){
for (Int_t i = 0; i < fNMixedEvent; i++) {
//_________________________________________________________________________
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 = 0;
- 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);
+ // 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 = 0;
+ 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;
- }
- }
+ 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;
+ }
+ }
if(jet) delete jet;
- }
-
- return kTRUE ;
-
+ }
+
+ return kTRUE ;
+
}
//____________________________________________________________________________
//____________________________________________________________________________
AliAODMCHeader* AliCaloTrackReader::GetAODMCHeader(Int_t input) const {
- //Return MC header in AOD. Do it for the corresponding input event.
+ //Return MC header in AOD. Do it for the corresponding input event.
AliAODMCHeader *mch = NULL;
- if(fDataType == kAOD){
- //Normal input AOD
- if(input == 0) {
+ if(fDataType == kAOD){
+ //Normal input AOD
+ if(input == 0) {
mch = (AliAODMCHeader*)((AliAODEvent*)fInputEvent)->FindListObject("mcheader");
}
-// //Second input AOD
-// else if(input == 1){
-// mch = (AliAODMCHeader*) fSecondInputAODEvent->FindListObject("mcheader");
-// }
- else {
- printf("AliCaloTrackReader::GetAODMCHeader() - wrong AOD input index, %d\n",input);
- }
- }
- else {
- printf("AliCaloTrackReader::GetAODMCHeader() - Input are not AODs\n");
- }
+ // //Second input AOD
+ // else if(input == 1){
+ // mch = (AliAODMCHeader*) fSecondInputAODEvent->FindListObject("mcheader");
+ // }
+ else {
+ printf("AliCaloTrackReader::GetAODMCHeader() - wrong AOD input index, %d\n",input);
+ }
+ }
+ else {
+ printf("AliCaloTrackReader::GetAODMCHeader() - Input are not AODs\n");
+ }
return mch;
}
//_______________________________________________________________
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;
- }
-
+ //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");
fESDtrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010();
+ fV0ADC[0] = 0; fV0ADC[1] = 0;
+ fV0Mul[0] = 0; fV0Mul[1] = 0;
+
+ fZvtxCut = 10.;
}
//Fill Vertex array
FillVertexArray();
+ //Reject events with Z vertex too large, only for SE analysis, if not, cut on the analysis code
+ if(!GetMixedEvent() && TMath::Abs(fVertex[0][2]) > fZvtxCut) return kFALSE;
//Fill the arrays with cluster/tracks/cells data
if(fFillEMCALCells)
if(fFillPHOSCells)
FillInputPHOSCells();
- if(fFillCTS)
+ if(fFillCTS){
FillInputCTS();
+ //Accept events with at least one track
+ if(fTrackMult == 0) return kFALSE;
+ }
+
+ //In case of data produced with calo filter, some information stored in non usual places
+ if(IsCaloFilterPatchOn()){
+ fTrackMult = (Int_t) ((AliAODHeader*)fInputEvent->GetHeader())->GetCentrality();
+ //printf("Track multiplicity %d \n",fTrackMult);
+ if(fTrackMult == 0) return kFALSE;
+ }
+
if(fFillEMCAL)
FillInputEMCAL();
if(fFillPHOS)
FillInputPHOS();
+ FillInputVZERO();
return kTRUE ;
}
// 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();
+ if(fAODEMCAL) fAODEMCAL -> Clear("C");
+ if(fAODPHOS) fAODPHOS -> Clear("C");
+// if(fEMCALCells) fEMCALCells -> Clear("");
+// if(fPHOSCells) fPHOSCells -> Clear("");
+
+ fV0ADC[0] = 0; fV0ADC[1] = 0;
+ fV0Mul[0] = 0; fV0Mul[1] = 0;
+
}
//____________________________________________________________________________
}
if (!fMixedEvent) { //Single event analysis
-
- if(fDataType!=kMC)fInputEvent->GetPrimaryVertex()->GetXYZ(fVertex[0]);
- else {
+ if(fDataType!=kMC){
+
+ if(fInputEvent->GetPrimaryVertex()){
+ fInputEvent->GetPrimaryVertex()->GetXYZ(fVertex[0]);
+ }
+ else {
+ printf("AliCaloTrackReader::FillVertexArray() - NULL primary vertex\n");
+ fVertex[0][0]=0.; fVertex[0][1]=0.; fVertex[0][2]=0.;
+ }//Primary vertex pointer do not exist
+
+ } else {//MC read event
fVertex[0][0]=0.; fVertex[0][1]=0.; fVertex[0][2]=0.;
}
}
-
//____________________________________________________________________________
void AliCaloTrackReader::FillInputCTS() {
//Return array with Central Tracking System (CTS) tracks
if(fDataType==kESD && !fESDtrackCuts->AcceptTrack((AliESDtrack*)track)) continue;
+ // Track filter selection
+ //if (fTrackFilter) {
+ // selectInfo = fTrackFilter->IsSelected(esdTrack);
+ // if (!selectInfo && !(esd->GetPrimaryVertex())->UsesTrack(esdTrack->GetID())) continue;
+ // }
+
//Count the tracks in eta < 0.9
//printf("Eta %f cut %f\n",TMath::Abs(track->Eta()),fTrackMultEtaCut);
if(TMath::Abs(track->Eta())< fTrackMultEtaCut) fTrackMult++;
//
}
- //____________________________________________________________________________
+//____________________________________________________________________________
+void AliCaloTrackReader::FillInputEMCALAlgorithm(AliVCluster * clus, const Int_t iclus) {
+ //Fill the EMCAL data in the array, do it
+
+ Int_t vindex = 0 ;
+ if (fMixedEvent)
+ vindex = fMixedEvent->EventIndexForCaloCluster(iclus);
+
+ //Check if the cluster contains any bad channel and if close to calorimeter borders
+ if(GetCaloUtils()->ClusterContainsBadChannel("EMCAL",clus->GetCellsAbsId(), clus->GetNCells()))
+ return;
+ if(!GetCaloUtils()->CheckCellFiducialRegion(clus, (AliVCaloCells*)fInputEvent->GetEMCALCells(), fInputEvent, vindex))
+ return;
+
+ TLorentzVector momentum ;
+
+ clus->GetMomentum(momentum, fVertex[vindex]);
+
+ if(fEMCALPtMin < momentum.Pt()){
+
+ if(fCheckFidCut && !fFiducialCut->IsInFiducialCut(momentum,"EMCAL"))
+ return;
+
+ if(fDebug > 2 && momentum.E() > 0.1)
+ printf("AliCaloTrackReader::FillInputEMCAL() - Selected clusters E %3.2f, pt %3.2f, phi %3.2f, eta %3.2f\n",
+ momentum.E(),momentum.Pt(),momentum.Phi()*TMath::RadToDeg(),momentum.Eta());
+
+ //Float_t pos[3];
+ //clus->GetPosition(pos);
+ //printf("Before Corrections: e %f, x %f, y %f, z %f\n",clus->E(),pos[0],pos[1],pos[2]);
+
+ //Recalibrate the cluster energy
+ if(GetCaloUtils()->IsRecalibrationOn()) {
+ Float_t energy = GetCaloUtils()->RecalibrateClusterEnergy(clus, GetEMCALCells());
+ clus->SetE(energy);
+ //printf("Recalibrated Energy %f\n",clus->E());
+ GetCaloUtils()->RecalculateClusterShowerShapeParameters(GetEMCALCells(),clus);
+ GetCaloUtils()->RecalculateClusterPID(clus);
+
+ }
+
+ //Recalculate distance to bad channels, if new list of bad channels provided
+ GetCaloUtils()->RecalculateClusterDistanceToBadChannel(GetEMCALCells(),clus);
+
+ //Recalculate cluster position
+ if(GetCaloUtils()->IsRecalculationOfClusterPositionOn()){
+ GetCaloUtils()->RecalculateClusterPosition(GetEMCALCells(),clus);
+ //clus->GetPosition(pos);
+ //printf("After Corrections: e %f, x %f, y %f, z %f\n",clus->E(),pos[0],pos[1],pos[2]);
+ }
+
+ //Correct non linearity
+ if(GetCaloUtils()->IsCorrectionOfClusterEnergyOn()){
+ GetCaloUtils()->CorrectClusterEnergy(clus) ;
+ //printf("Linearity Corrected Energy %f\n",clus->E());
+ }
+
+ if (fMixedEvent)
+ clus->SetID(iclus) ;
+
+ fAODEMCAL->Add(clus);
+ }
+}
+
+//____________________________________________________________________________
void AliCaloTrackReader::FillInputEMCAL() {
//Return array with EMCAL clusters in aod format
if(fDebug > 2 ) printf("AliCaloTrackReader::FillInputEMCAL()\n");
//Loop to select clusters in fiducial cut and fill container with aodClusters
- Int_t nclusters = fInputEvent->GetNumberOfCaloClusters();
- for (Int_t iclus = 0; iclus < nclusters; iclus++) {
- AliVCluster * clus = 0;
- if ( (clus = fInputEvent->GetCaloCluster(iclus)) ) {
- if (IsEMCALCluster(clus)){
-
- //Check if the cluster contains any bad channel and if close to calorimeter borders
- Int_t vindex = 0 ;
- if (fMixedEvent)
- vindex = fMixedEvent->EventIndexForCaloCluster(iclus);
-
- if(GetCaloUtils()->ClusterContainsBadChannel("EMCAL",clus->GetCellsAbsId(), clus->GetNCells()))
- continue;
- if(!GetCaloUtils()->CheckCellFiducialRegion(clus, (AliVCaloCells*)fInputEvent->GetEMCALCells(), fInputEvent, vindex))
- continue;
-
- TLorentzVector momentum ;
-
- clus->GetMomentum(momentum, fVertex[vindex]);
-
- if(fEMCALPtMin < momentum.Pt()){
-
- if(fCheckFidCut && !fFiducialCut->IsInFiducialCut(momentum,"EMCAL"))
- continue;
-
- if(fDebug > 2 && momentum.E() > 0.1)
- printf("AliCaloTrackReader::FillInputEMCAL() - Selected clusters E %3.2f, pt %3.2f, phi %3.2f, eta %3.2f\n",
- momentum.E(),momentum.Pt(),momentum.Phi()*TMath::RadToDeg(),momentum.Eta());
-
- //Float_t pos[3];
- //clus->GetPosition(pos);
- //printf("Before Corrections: e %f, x %f, y %f, z %f\n",clus->E(),pos[0],pos[1],pos[2]);
-
- //Recalibrate the cluster energy
- if(GetCaloUtils()->IsRecalibrationOn()) {
- Float_t energy = GetCaloUtils()->RecalibrateClusterEnergy(clus, GetEMCALCells());
- clus->SetE(energy);
- //printf("Recalibrated Energy %f\n",clus->E());
- GetCaloUtils()->RecalculateClusterShowerShapeParameters(GetEMCALCells(),clus);
- GetCaloUtils()->RecalculateClusterPID(clus);
-
- }
-
- //Correct non linearity
- if(GetCaloUtils()->IsCorrectionOfClusterEnergyOn()){
- GetCaloUtils()->CorrectClusterEnergy(clus) ;
- //printf("Linearity Corrected Energy %f\n",clus->E());
- }
-
- //Recalculate cluster position
- if(GetCaloUtils()->IsRecalculationOfClusterPositionOn()){
- GetCaloUtils()->RecalculateClusterPosition(GetEMCALCells(),clus);
- //clus->GetPosition(pos);
- //printf("After Corrections: e %f, x %f, y %f, z %f\n",clus->E(),pos[0],pos[1],pos[2]);
- }
-
- if (fMixedEvent) {
- clus->SetID(iclus) ;
- }
-
- fAODEMCAL->Add(clus);
-
- }//Pt and Fiducial cut passed.
- }//EMCAL cluster
- }// cluster exists
- }// cluster loop
+ if(fEMCALClustersListName==""){
+ Int_t nclusters = fInputEvent->GetNumberOfCaloClusters();
+ for (Int_t iclus = 0; iclus < nclusters; iclus++) {
+ AliVCluster * clus = 0;
+ if ( (clus = fInputEvent->GetCaloCluster(iclus)) ) {
+ if (IsEMCALCluster(clus)){
+ FillInputEMCALAlgorithm(clus, iclus);
+ }//EMCAL cluster
+ }// cluster exists
+ }// cluster loop
+ }//Get the clusters from the input event
+ else {
+ TClonesArray * clusterList = dynamic_cast<TClonesArray*> (fOutputEvent->FindListObject(fEMCALClustersListName));
+ if(!clusterList){
+ printf("AliCaloTrackReader::FillInputEMCAL() - Wrong name of list with clusters? <%s>\n",fEMCALClustersListName.Data());
+ return;
+ }
+ Int_t nclusters = clusterList->GetEntriesFast();
+ for (Int_t iclus = 0; iclus < nclusters; iclus++) {
+ AliVCluster * clus = dynamic_cast<AliVCluster*> (clusterList->At(iclus));
+ //printf("E %f\n",clus->E());
+ if (clus) FillInputEMCALAlgorithm(clus, iclus);
+ else printf("AliCaloTrackReader::FillInputEMCAL() - Null cluster in list!\n");
+ }// cluster loop
+ }
//Recalculate track matching
GetCaloUtils()->RecalculateClusterTrackMatching(fInputEvent);
}
+
//____________________________________________________________________________
Bool_t AliCaloTrackReader::IsEMCALCluster(AliVCluster* cluster) const {
// Check if it is a cluster from EMCAL. For old AODs cluster type has
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff