From 3bfc473233cf8cefc5139f3f53757b8292cf06a9 Mon Sep 17 00:00:00 2001 From: gconesab Date: Tue, 27 Sep 2011 17:04:08 +0000 Subject: [PATCH] First implementation of the time calibration for analysis in AliEMCALRecoUtils All the other classes addapt to this change AliAnaCalorimeterQA: Some minor fixes, put back the track matching related histograms AliCalTrackReader: Add a swich to recalculate the clusters or not, even if the AliEMCALRecoUtils switch are on, avoid multiple aplication of calibration parameters in case of using clusterization task. AliAnalysisTaskEMCALClusterizer: Possibility to add matched tracks in AODCaloCluster added, this needs a patch in AliAODCaloCluster not committed yet. --- EMCAL/AliEMCALRecoUtils.cxx | 277 +++- EMCAL/AliEMCALRecoUtils.h | 79 +- PWG4/CaloCalib/AliAnalysisTaskCaloFilter.cxx | 2 +- .../AliAnalysisTaskEMCALClusterize.cxx | 27 +- PWG4/PartCorrBase/AliCaloTrackReader.cxx | 121 +- PWG4/PartCorrBase/AliCaloTrackReader.h | 9 +- PWG4/PartCorrDep/AliAnaCalorimeterQA.cxx | 1269 +++++++++-------- PWG4/PartCorrDep/AliAnaCalorimeterQA.h | 42 +- 8 files changed, 1078 insertions(+), 748 deletions(-) diff --git a/EMCAL/AliEMCALRecoUtils.cxx b/EMCAL/AliEMCALRecoUtils.cxx index 89037772710..7a22a372f9f 100644 --- a/EMCAL/AliEMCALRecoUtils.cxx +++ b/EMCAL/AliEMCALRecoUtils.cxx @@ -57,7 +57,7 @@ #include "AliEMCALRecoUtils.h" #include "AliEMCALGeometry.h" #include "AliEMCALTrack.h" -#include "AliEMCALCalibTimeDepCorrection.h" +#include "AliEMCALCalibTimeDepCorrection.h" // Run dependent #include "AliEMCALPIDUtils.h" ClassImp(AliEMCALRecoUtils) @@ -67,8 +67,9 @@ AliEMCALRecoUtils::AliEMCALRecoUtils(): fParticleType(kPhoton), fPosAlgo(kUnchanged), fW0(4.), fNonLinearityFunction(kNoCorrection), fNonLinearThreshold(30), fSmearClusterEnergy(kFALSE), fRandom(), - fRecalibration(kFALSE), fEMCALRecalibrationFactors(), - fUseTimeCorrectionFactors(kFALSE), fTimeCorrectionFactorsSet(kFALSE), + fCellsRecalibrated(kFALSE), fRecalibration(kFALSE), fEMCALRecalibrationFactors(), + fTimeRecalibration(kFALSE), fEMCALTimeRecalibrationFactors(), + fUseRunCorrectionFactors(kFALSE), fRunCorrectionFactorsSet(kFALSE), fRemoveBadChannels(kFALSE), fRecalDistToBadChannels(kFALSE), fEMCALBadChannelMap(), fNCellsFromEMCALBorder(0), fNoEMCALBorderAtEta0(kTRUE), fRejectExoticCluster(kFALSE), fPIDUtils(), fAODFilterMask(32), @@ -132,8 +133,10 @@ AliEMCALRecoUtils::AliEMCALRecoUtils(const AliEMCALRecoUtils & reco) fParticleType(reco.fParticleType), fPosAlgo(reco.fPosAlgo), fW0(reco.fW0), fNonLinearityFunction(reco.fNonLinearityFunction), fNonLinearThreshold(reco.fNonLinearThreshold), fSmearClusterEnergy(reco.fSmearClusterEnergy), fRandom(), + fCellsRecalibrated(reco.fCellsRecalibrated), fRecalibration(reco.fRecalibration), fEMCALRecalibrationFactors(reco.fEMCALRecalibrationFactors), - fUseTimeCorrectionFactors(reco.fUseTimeCorrectionFactors), fTimeCorrectionFactorsSet(reco.fTimeCorrectionFactorsSet), + fTimeRecalibration(reco.fTimeRecalibration), fEMCALTimeRecalibrationFactors(reco.fEMCALTimeRecalibrationFactors), + fUseRunCorrectionFactors(reco.fUseRunCorrectionFactors), fRunCorrectionFactorsSet(reco.fRunCorrectionFactorsSet), fRemoveBadChannels(reco.fRemoveBadChannels), fRecalDistToBadChannels(reco.fRecalDistToBadChannels), fEMCALBadChannelMap(reco.fEMCALBadChannelMap), fNCellsFromEMCALBorder(reco.fNCellsFromEMCALBorder), fNoEMCALBorderAtEta0(reco.fNoEMCALBorderAtEta0), @@ -184,10 +187,15 @@ AliEMCALRecoUtils & AliEMCALRecoUtils::operator = (const AliEMCALRecoUtils & rec fNonLinearThreshold = reco.fNonLinearThreshold; fSmearClusterEnergy = reco.fSmearClusterEnergy; + fCellsRecalibrated = reco.fCellsRecalibrated; fRecalibration = reco.fRecalibration; fEMCALRecalibrationFactors = reco.fEMCALRecalibrationFactors; - fUseTimeCorrectionFactors = reco.fUseTimeCorrectionFactors; - fTimeCorrectionFactorsSet = reco.fTimeCorrectionFactorsSet; + + fTimeRecalibration = reco.fTimeRecalibration; + fEMCALTimeRecalibrationFactors = reco.fEMCALTimeRecalibrationFactors; + + fUseRunCorrectionFactors = reco.fUseRunCorrectionFactors; + fRunCorrectionFactorsSet = reco.fRunCorrectionFactorsSet; fRemoveBadChannels = reco.fRemoveBadChannels; fRecalDistToBadChannels = reco.fRecalDistToBadChannels; @@ -285,6 +293,11 @@ AliEMCALRecoUtils::~AliEMCALRecoUtils() delete fEMCALRecalibrationFactors; } + if(fEMCALTimeRecalibrationFactors) { + fEMCALTimeRecalibrationFactors->Clear(); + delete fEMCALTimeRecalibrationFactors; + } + if(fEMCALBadChannelMap) { fEMCALBadChannelMap->Clear(); delete fEMCALBadChannelMap; @@ -716,6 +729,31 @@ void AliEMCALRecoUtils::InitEMCALRecalibrationFactors(){ TH1::AddDirectory(oldStatus); } +//________________________________________________________________ +void AliEMCALRecoUtils::InitEMCALTimeRecalibrationFactors(){ + //Init EMCAL recalibration factors + AliDebug(2,"AliCalorimeterUtils::InitEMCALRecalibrationFactors()"); + //In order to avoid rewriting the same histograms + Bool_t oldStatus = TH1::AddDirectoryStatus(); + TH1::AddDirectory(kFALSE); + + fEMCALTimeRecalibrationFactors = new TObjArray(4); + for (int i = 0; i < 4; i++) + fEMCALTimeRecalibrationFactors->Add(new TH1F(Form("hAllTimeAvBC%d",i), + Form("hAllTimeAvBC%d",i), + 11521,0.,11521) ); + //Init the histograms with 1 + for (Int_t bc = 0; bc < 4; bc++) { + + SetEMCALChannelTimeRecalibrationFactor(bc,0.); + } + + fEMCALTimeRecalibrationFactors->SetOwner(kTRUE); + fEMCALTimeRecalibrationFactors->Compress(); + + //In order to avoid rewriting the same histograms + TH1::AddDirectory(oldStatus); +} //________________________________________________________________ void AliEMCALRecoUtils::InitEMCALBadChannelStatusMap(){ @@ -730,9 +768,7 @@ void AliEMCALRecoUtils::InitEMCALBadChannelStatusMap(){ for (int i = 0; i < 10; i++) { fEMCALBadChannelMap->Add(new TH2I(Form("EMCALBadChannelMap_Mod%d",i),Form("EMCALBadChannelMap_Mod%d",i), 48, 0, 48, 24, 0, 24)); } - - //delete hTemp; - + fEMCALBadChannelMap->SetOwner(kTRUE); fEMCALBadChannelMap->Compress(); @@ -741,8 +777,10 @@ void AliEMCALRecoUtils::InitEMCALBadChannelStatusMap(){ } //________________________________________________________________ -void AliEMCALRecoUtils::RecalibrateClusterEnergy(AliEMCALGeometry* geom, AliVCluster * cluster, AliVCaloCells * cells){ - // Recalibrate the cluster energy, considering the recalibration map and the energy of the cells that compose the cluster. +void AliEMCALRecoUtils::RecalibrateClusterEnergy(AliEMCALGeometry* geom, AliVCluster * cluster, AliVCaloCells * cells, const Int_t bc){ + // Recalibrate the cluster energy and Time, considering the recalibration map + // and the energy of the cells and time that compose the cluster. + // bc= bunch crossing number returned by esdevent->GetBunchCrossNumber(); if(!cluster){ AliInfo("Cluster pointer null!"); @@ -756,34 +794,160 @@ void AliEMCALRecoUtils::RecalibrateClusterEnergy(AliEMCALGeometry* geom, AliVClu //Initialize some used variables Float_t energy = 0; - Int_t absId = -1; - Int_t icol = -1, irow = -1, imod=1; + Int_t absId =-1; + Int_t icol =-1, irow =-1, imod=1; Float_t factor = 1, frac = 0; - + Int_t absIdMax = -1; + Float_t emax = 0; + //Loop on the cells, get the cell amplitude and recalibration factor, multiply and and to the new energy for(Int_t icell = 0; icell < ncells; icell++){ absId = index[icell]; frac = fraction[icell]; if(frac < 1e-5) frac = 1; //in case of EMCAL, this is set as 0 since unfolding is off - Int_t iTower = -1, iIphi = -1, iIeta = -1; - geom->GetCellIndex(absId,imod,iTower,iIphi,iIeta); - if(fEMCALRecalibrationFactors->GetEntries() <= imod) continue; - geom->GetCellPhiEtaIndexInSModule(imod,iTower,iIphi, iIeta,irow,icol); - factor = GetEMCALChannelRecalibrationFactor(imod,icol,irow); - AliDebug(2,Form("AliEMCALRecoUtils::RecalibrateClusterEnergy - recalibrate cell: module %d, col %d, row %d, cell fraction %f,recalibration factor %f, cell energy %f\n", - imod,icol,irow,frac,factor,cells->GetCellAmplitude(absId))); - + + if(!fCellsRecalibrated && IsRecalibrationOn()){ + + // Energy + Int_t iTower = -1, iIphi = -1, iIeta = -1; + geom->GetCellIndex(absId,imod,iTower,iIphi,iIeta); + if(fEMCALRecalibrationFactors->GetEntries() <= imod) continue; + geom->GetCellPhiEtaIndexInSModule(imod,iTower,iIphi, iIeta,irow,icol); + factor = GetEMCALChannelRecalibrationFactor(imod,icol,irow); + + AliDebug(2,Form("AliEMCALRecoUtils::RecalibrateClusterEnergy - recalibrate cell: module %d, col %d, row %d, cell fraction %f,recalibration factor %f, cell energy %f\n", + imod,icol,irow,frac,factor,cells->GetCellAmplitude(absId))); + + } + energy += cells->GetCellAmplitude(absId)*factor*frac; + + if(emax < cells->GetCellAmplitude(absId)*factor*frac){ + emax = cells->GetCellAmplitude(absId)*factor*frac; + absIdMax = absId; + } + } - - AliDebug(2,Form("AliEMCALRecoUtils::RecalibrateClusterEnergy - Energy before %f, after %f\n",cluster->E(),energy)); - - cluster->SetE(energy); - + cluster->SetE(energy); + + AliDebug(2,Form("AliEMCALRecoUtils::RecalibrateClusterEnergy - Energy before %f, after %f\n",cluster->E(),energy)); + + // Recalculate time of cluster only for ESDs + if(!strcmp("AliESDCaloCluster",Form("%s",cluster->ClassName()))){ + + // Time + Double_t weightedTime = 0; + Double_t weight = 0; + Double_t weightTot = 0; + Double_t maxcellTime = 0; + for(Int_t icell = 0; icell < ncells; icell++){ + absId = index[icell]; + frac = fraction[icell]; + if(frac < 1e-5) frac = 1; //in case of EMCAL, this is set as 0 since unfolding is off + + Double_t celltime = cells->GetCellTime(absId); + RecalibrateCellTime(absId, bc, celltime); + if(absId == absIdMax) maxcellTime = celltime; + + if(!fCellsRecalibrated){ + + Int_t iTower = -1, iIphi = -1, iIeta = -1; + geom->GetCellIndex(absId,imod,iTower,iIphi,iIeta); + if(fEMCALRecalibrationFactors->GetEntries() <= imod) continue; + geom->GetCellPhiEtaIndexInSModule(imod,iTower,iIphi, iIeta,irow,icol); + factor = GetEMCALChannelRecalibrationFactor(imod,icol,irow); + + AliDebug(2,Form("AliEMCALRecoUtils::RecalibrateClusterEnergy - recalibrate cell: module %d, col %d, row %d, cell fraction %f,recalibration factor %f, cell energy %f\n", + imod,icol,irow,frac,factor,cells->GetCellTime(absId))); + + } + + weight = GetCellWeight(cells->GetCellAmplitude(absId)*factor*frac , energy ); + weightTot += weight; + weightedTime += celltime * weight; + + } + + if(weightTot > 0) + cluster->SetTOF(weightedTime/weightTot); + else + cluster->SetTOF(maxcellTime); + + } +} + +//________________________________________________________________ +void AliEMCALRecoUtils::RecalibrateCells(AliEMCALGeometry* geom, AliVCaloCells * cells, Int_t bc){ + // Recalibrate the cells time and energy, considering the recalibration map and the energy + // of the cells that compose the cluster. + // bc= bunch crossing number returned by esdevent->GetBunchCrossNumber(); + + if(!IsRecalibrationOn()) return; + + if(!cells){ + AliInfo("Cells pointer null!"); + return; + } + + fCellsRecalibrated = kTRUE; + + Int_t absId =-1; + Int_t icol =-1, irow =-1, imod = 1; + Int_t iTower =-1, iIeta =-1, iIphi =-1; + + Int_t nEMcell = cells->GetNumberOfCells() ; + + for (Int_t iCell = 0; iCell < nEMcell; iCell++) { + + absId = cells->GetCellNumber(iCell); + + // Energy + Float_t factor = 1; + if(IsRecalibrationOn()){ + geom->GetCellIndex(absId,imod,iTower,iIphi,iIeta); + if(fEMCALRecalibrationFactors->GetEntries() <= imod) continue; + geom->GetCellPhiEtaIndexInSModule(imod,iTower,iIphi, iIeta,irow,icol); + factor = GetEMCALChannelRecalibrationFactor(imod,icol,irow); + } + + Float_t cellE = cells->GetAmplitude(iCell) * factor ; + + //Time + Double_t celltime = cells->GetCellTime(absId); + RecalibrateCellTime(absId, bc, celltime); + + //Set new values + cells->SetCell(iCell,cells->GetCellNumber(iCell),cellE, celltime); + + } + } +//________________________________________________________________ +void AliEMCALRecoUtils::RecalibrateCellTime(const Int_t absId, const Int_t bc, Double_t & celltime){ + // Recalibrate time of cell with absID considering the recalibration map + // bc= bunch crossing number returned by esdevent->GetBunchCrossNumber(); + + if(!fCellsRecalibrated && IsTimeRecalibrationOn()){ +// printf("cell time org %g, ",celltime); + Double_t timeBCoffset = 0.; + if( bc%4 ==0 || bc%4==1) timeBCoffset = 100.*1.e-9; //in ns + + Double_t celloffset = GetEMCALChannelTimeRecalibrationFactor(bc%4,absId)*1.e-9; + +// printf("absId %d, time %f bc %d-%d: bc0 %f, bc1 %f, bc2 %f, bc3 %f \n", absId, celltime*1.e9,bc, bc%4, +// GetEMCALChannelTimeRecalibrationFactor(0,absId),GetEMCALChannelTimeRecalibrationFactor(1,absId), +// GetEMCALChannelTimeRecalibrationFactor(2,absId),GetEMCALChannelTimeRecalibrationFactor(3,absId)); + + celltime -= timeBCoffset ; + celltime -= celloffset ; +// printf("new %g\n",celltime); + } + +} + //__________________________________________________ void AliEMCALRecoUtils::RecalculateClusterPosition(AliEMCALGeometry *geom, AliVCaloCells* cells, AliVCluster* clu) { @@ -825,20 +989,26 @@ void AliEMCALRecoUtils::RecalculateClusterPositionFromTowerGlobal(AliEMCALGeomet //printf("** Cluster energy %f, ncells %d, depth %f\n",clEnergy,clu->GetNCells(),depth); for (Int_t iDig=0; iDig< clu->GetNCells(); iDig++) { - absId = clu->GetCellAbsId(iDig); - fraction = clu->GetCellAmplitudeFraction(iDig); - if(fraction < 1e-4) fraction = 1.; // in case unfolding is off - geom->GetCellIndex(absId,iSM,iTower,iIphi,iIeta); - geom->GetCellPhiEtaIndexInSModule(iSM,iTower,iIphi, iIeta,iphi,ieta); - if(IsRecalibrationOn()) { - recalFactor = GetEMCALChannelRecalibrationFactor(iSM,ieta,iphi); + absId = clu->GetCellAbsId(iDig); + fraction = clu->GetCellAmplitudeFraction(iDig); + if(fraction < 1e-4) fraction = 1.; // in case unfolding is off + + if(!fCellsRecalibrated){ + + geom->GetCellIndex(absId,iSM,iTower,iIphi,iIeta); + geom->GetCellPhiEtaIndexInSModule(iSM,iTower,iIphi, iIeta,iphi,ieta); + + if(IsRecalibrationOn()) { + recalFactor = GetEMCALChannelRecalibrationFactor(iSM,ieta,iphi); + } } + eCell = cells->GetCellAmplitude(absId)*fraction*recalFactor; weight = GetCellWeight(eCell,clEnergy); - //printf("cell energy %f, weight %f\n",eCell,weight); totalWeight += weight; + geom->RelPosCellInSModule(absId,depth,pLocal[0],pLocal[1],pLocal[2]); //printf("pLocal (%f,%f,%f), SM %d, absId %d\n",pLocal[0],pLocal[1],pLocal[2],iSupModMax,absId); geom->GetGlobal(pLocal,pGlobal,iSupModMax); @@ -906,17 +1076,24 @@ void AliEMCALRecoUtils::RecalculateClusterPositionFromTowerIndex(AliEMCALGeometr absId = clu->GetCellAbsId(iDig); fraction = clu->GetCellAmplitudeFraction(iDig); if(fraction < 1e-4) fraction = 1.; // in case unfolding is off - geom->GetCellIndex(absId,iSupMod,iTower,iIphi,iIeta); - geom->GetCellPhiEtaIndexInSModule(iSupMod,iTower,iIphi, iIeta,iphi,ieta); - + if (iDig==0) startingSM = iSupMod; else if(iSupMod != startingSM) areInSameSM = kFALSE; eCell = cells->GetCellAmplitude(absId); - if(IsRecalibrationOn()) { - recalFactor = GetEMCALChannelRecalibrationFactor(iSupMod,ieta,iphi); + geom->GetCellIndex(absId,iSupMod,iTower,iIphi,iIeta); + geom->GetCellPhiEtaIndexInSModule(iSupMod,iTower,iIphi, iIeta,iphi,ieta); + + if(!fCellsRecalibrated){ + + if(IsRecalibrationOn()) { + + recalFactor = GetEMCALChannelRecalibrationFactor(iSupMod,ieta,iphi); + + } } + eCell = cells->GetCellAmplitude(absId)*fraction*recalFactor; weight = GetCellWeight(eCell,clEnergy); @@ -1091,9 +1268,15 @@ void AliEMCALRecoUtils::RecalculateClusterShowerShapeParameters(AliEMCALGeometry //Get the cell energy, if recalibration is on, apply factors fraction = cluster->GetCellAmplitudeFraction(iDigit); if(fraction < 1e-4) fraction = 1.; // in case unfolding is off - if(IsRecalibrationOn()) { - recalFactor = GetEMCALChannelRecalibrationFactor(iSupMod,ieta,iphi); + + if(!fCellsRecalibrated){ + + if(IsRecalibrationOn()) { + recalFactor = GetEMCALChannelRecalibrationFactor(iSupMod,ieta,iphi); + } + } + eCell = cells->GetCellAmplitude(cluster->GetCellAbsId(iDigit))*fraction*recalFactor; if(cluster->E() > 0 && eCell > 0){ @@ -1751,14 +1934,14 @@ void AliEMCALRecoUtils::Print(const Option_t *) const } //_____________________________________________________________________ -void AliEMCALRecoUtils::SetTimeDependentCorrections(Int_t runnumber){ +void AliEMCALRecoUtils::SetRunDependentCorrections(Int_t runnumber){ //Get EMCAL time dependent corrections from file and put them in the recalibration histograms //Do it only once and only if it is requested - if(!fUseTimeCorrectionFactors) return; - if(fTimeCorrectionFactorsSet) return; + if(!fUseRunCorrectionFactors) return; + if(fRunCorrectionFactorsSet) return; - printf("AliEMCALRecoUtils::GetTimeDependentCorrections() - Get Correction Factors for Run number %d\n",runnumber); + AliInfo(Form("AliEMCALRecoUtils::GetRunDependentCorrections() - Get Correction Factors for Run number %d\n",runnumber)); AliEMCALCalibTimeDepCorrection *corr = new AliEMCALCalibTimeDepCorrection(); corr->ReadRootInfo(Form("CorrectionFiles/Run%d_Correction.root",runnumber)); @@ -1776,6 +1959,6 @@ void AliEMCALRecoUtils::SetTimeDependentCorrections(Int_t runnumber){ } } } - fTimeCorrectionFactorsSet = kTRUE; + fRunCorrectionFactorsSet = kTRUE; } diff --git a/EMCAL/AliEMCALRecoUtils.h b/EMCAL/AliEMCALRecoUtils.h index f2292ef5f83..21b2cb207e1 100644 --- a/EMCAL/AliEMCALRecoUtils.h +++ b/EMCAL/AliEMCALRecoUtils.h @@ -129,10 +129,11 @@ public: else { AliInfo(Form("Index %d larger than 2, do nothing\n",i)) ; } } //----------------------------------------------------- - //Recalibration + // Energy Recalibration //----------------------------------------------------- - - void RecalibrateClusterEnergy(AliEMCALGeometry* geom, AliVCluster* cluster, AliVCaloCells * cells) ; + + void RecalibrateCells(AliEMCALGeometry* geom, AliVCaloCells * cells, Int_t bc) ; // Energy and Time + void RecalibrateClusterEnergy(AliEMCALGeometry* geom, AliVCluster* cluster, AliVCaloCells * cells, const Int_t bc=0) ; // Energy and time Bool_t IsRecalibrationOn() const { return fRecalibration ; } void SwitchOffRecalibration() { fRecalibration = kFALSE ; } @@ -140,27 +141,52 @@ public: if(!fEMCALRecalibrationFactors)InitEMCALRecalibrationFactors() ; } void InitEMCALRecalibrationFactors() ; - //Recalibrate channels with time dependent corrections - void SwitchOffTimeDepCorrection() { fUseTimeCorrectionFactors = kFALSE ; } - void SwitchOnTimeDepCorrection() { fUseTimeCorrectionFactors = kTRUE ; - SwitchOnRecalibration() ; } - void SetTimeDependentCorrections(Int_t runnumber); - + TH2F * GetEMCALChannelRecalibrationFactors(Int_t iSM) const { return (TH2F*)fEMCALRecalibrationFactors->At(iSM) ; } + void SetEMCALChannelRecalibrationFactors(TObjArray *map) { fEMCALRecalibrationFactors = map ; } + void SetEMCALChannelRecalibrationFactors(Int_t iSM , TH2F* h) { fEMCALRecalibrationFactors->AddAt(h,iSM) ; } + Float_t GetEMCALChannelRecalibrationFactor(Int_t iSM , Int_t iCol, Int_t iRow) const { - if(fEMCALRecalibrationFactors) - return (Float_t) ((TH2F*)fEMCALRecalibrationFactors->At(iSM))->GetBinContent(iCol,iRow); - else return 1 ; } + if(fEMCALRecalibrationFactors) + return (Float_t) ((TH2F*)fEMCALRecalibrationFactors->At(iSM))->GetBinContent(iCol,iRow); + else return 1 ; } void SetEMCALChannelRecalibrationFactor(Int_t iSM , Int_t iCol, Int_t iRow, Double_t c = 1) { - if(!fEMCALRecalibrationFactors) InitEMCALRecalibrationFactors() ; - ((TH2F*)fEMCALRecalibrationFactors->At(iSM))->SetBinContent(iCol,iRow,c) ; } + if(!fEMCALRecalibrationFactors) InitEMCALRecalibrationFactors() ; + ((TH2F*)fEMCALRecalibrationFactors->At(iSM))->SetBinContent(iCol,iRow,c) ; } + + //Recalibrate channels energy with run dependent corrections + void SwitchOffRunDepCorrection() { fUseRunCorrectionFactors = kFALSE ; } + void SwitchOnRunDepCorrection() { fUseRunCorrectionFactors = kTRUE ; + SwitchOnRecalibration() ; } + void SetRunDependentCorrections(Int_t runnumber); + + //----------------------------------------------------- + // Time Recalibration + //----------------------------------------------------- + + void RecalibrateCellTime(const Int_t absId, const Int_t bc, Double_t & time); + + Bool_t IsTimeRecalibrationOn() const { return fTimeRecalibration ; } + void SwitchOffTimeRecalibration() { fTimeRecalibration = kFALSE ; } + void SwitchOnTimeRecalibration() { fTimeRecalibration = kTRUE ; + if(!fEMCALTimeRecalibrationFactors)InitEMCALTimeRecalibrationFactors() ; } + void InitEMCALTimeRecalibrationFactors() ; + + Float_t GetEMCALChannelTimeRecalibrationFactor(Int_t bc, Int_t absID) const { + if(fEMCALTimeRecalibrationFactors) + return (Float_t) ((TH1F*)fEMCALTimeRecalibrationFactors->At(bc))->GetBinContent(absID); + else return 1 ; } + + void SetEMCALChannelTimeRecalibrationFactor(Int_t bc,Int_t absID, Double_t c = 1) { + if(!fEMCALTimeRecalibrationFactors) InitEMCALTimeRecalibrationFactors() ; + ((TH1F*)fEMCALTimeRecalibrationFactors->At(bc))->SetBinContent(absID,c) ; } + + TH1F * GetEMCALChannelTimeRecalibrationFactors(Int_t bc) const { return (TH1F*)fEMCALTimeRecalibrationFactors->At(bc) ; } + void SetEMCALChannelTimeRecalibrationFactors(TObjArray *map) { fEMCALTimeRecalibrationFactors = map ; } + void SetEMCALChannelTimeRecalibrationFactors(Int_t bc , TH1F* h) { fEMCALTimeRecalibrationFactors->AddAt(h,bc) ; } - TH2F * GetEMCALChannelRecalibrationFactors(Int_t iSM) const { return (TH2F*)fEMCALRecalibrationFactors->At(iSM) ; } - void SetEMCALChannelRecalibrationFactors(TObjArray *map) { fEMCALRecalibrationFactors = map ; } - void SetEMCALChannelRecalibrationFactors(Int_t iSM , TH2F* h) { fEMCALRecalibrationFactors->AddAt(h,iSM) ; } - //----------------------------------------------------- - //Modules fiducial region, remove clusters in borders + // Modules fiducial region, remove clusters in borders //----------------------------------------------------- Bool_t CheckCellFiducialRegion(AliEMCALGeometry* geom, AliVCluster* cluster, AliVCaloCells* cells) ; @@ -312,13 +338,18 @@ private: Float_t fSmearClusterParam[3]; // Smearing parameters TRandom3 fRandom; // Random generator - // Recalibration + // Energy Recalibration + Bool_t fCellsRecalibrated; // Internal bool to check if cells (time/energy) where recalibrated and not recalibrate them when recalculating different things Bool_t fRecalibration; // Switch on or off the recalibration TObjArray* fEMCALRecalibrationFactors; // Array of histograms with map of recalibration factors, EMCAL - // Recalibrate with run dependent corrections - Bool_t fUseTimeCorrectionFactors; // Use Time Dependent Correction - Bool_t fTimeCorrectionFactorsSet; // Time Correction set at leat once + // Time Recalibration + Bool_t fTimeRecalibration; // Switch on or off the time recalibration + TObjArray* fEMCALTimeRecalibrationFactors; // Array of histograms with map of time recalibration factors, EMCAL + + // Recalibrate with run dependent corrections, energy + Bool_t fUseRunCorrectionFactors; // Use Run Dependent Correction + Bool_t fRunCorrectionFactorsSet; // Run Correction set at leat once // Bad Channels Bool_t fRemoveBadChannels; // Check the channel status provided and remove clusters with bad channels @@ -363,7 +394,7 @@ private: Float_t fCutMaxDCAToVertexZ; // Track-to-vertex cut in max absolute distance in z-plane Bool_t fCutDCAToVertex2D; // If true a 2D DCA cut is made. Tracks are accepted if sqrt((DCAXY / fCutMaxDCAToVertexXY)^2 + (DCAZ / fCutMaxDCAToVertexZ)^2) < 1 AND sqrt((DCAXY / fCutMinDCAToVertexXY)^2 + (DCAZ / fCutMinDCAToVertexZ)^2) > 1 - ClassDef(AliEMCALRecoUtils, 13) + ClassDef(AliEMCALRecoUtils, 14) }; diff --git a/PWG4/CaloCalib/AliAnalysisTaskCaloFilter.cxx b/PWG4/CaloCalib/AliAnalysisTaskCaloFilter.cxx index fbccca26623..1b2f2b61137 100644 --- a/PWG4/CaloCalib/AliAnalysisTaskCaloFilter.cxx +++ b/PWG4/CaloCalib/AliAnalysisTaskCaloFilter.cxx @@ -344,7 +344,7 @@ void AliAnalysisTaskCaloFilter::UserExec(Option_t */*option*/) }//Load matrices from Data //Recover time dependent corrections, put then in recalibration histograms. Do it once - fEMCALRecoUtils->SetTimeDependentCorrections(InputEvent()->GetRunNumber()); + fEMCALRecoUtils->SetRunDependentCorrections(InputEvent()->GetRunNumber()); }//first event diff --git a/PWG4/CaloCalib/AliAnalysisTaskEMCALClusterize.cxx b/PWG4/CaloCalib/AliAnalysisTaskEMCALClusterize.cxx index ed234ae258f..f1d5c1ed0bd 100644 --- a/PWG4/CaloCalib/AliAnalysisTaskEMCALClusterize.cxx +++ b/PWG4/CaloCalib/AliAnalysisTaskEMCALClusterize.cxx @@ -459,7 +459,7 @@ void AliAnalysisTaskEMCALClusterize::UserExec(Option_t *) }//Load matrices from Data //Recover time dependent corrections, put then in recalibration histograms. Do it once - fRecoUtils->SetTimeDependentCorrections(InputEvent()->GetRunNumber()); + fRecoUtils->SetRunDependentCorrections(InputEvent()->GetRunNumber()); }//first event @@ -583,6 +583,8 @@ void AliAnalysisTaskEMCALClusterize::UserExec(Option_t *) for(Int_t icell=0; icell < clus->GetNCells(); icell++ ){ fCellLabels[index[icell]] = label; fCellSecondLabels[index[icell]] = label2; + //printf("Clusterizer in : TOF %g\n",clus->GetTOF()*1.e9); + fCellTime[icell] = clus->GetTOF(); //printf("1) absID %d, label[0] %d label[1] %d\n",index[icell], fCellLabels[index[icell]],fCellSecondLabels[index[icell]]); @@ -593,14 +595,16 @@ void AliAnalysisTaskEMCALClusterize::UserExec(Option_t *) // Create digits //................. - Int_t idigit = 0; - Int_t id = -1; - Float_t amp = -1; - Float_t time = -1; + Int_t idigit = 0; + Int_t id = -1; + Float_t amp = -1; + Double_t time = -1; TTree *digitsTree = new TTree("digitstree","digitstree"); digitsTree->Branch("EMCAL","TClonesArray", &fDigitsArr, 32000); + Int_t bc = InputEvent()->GetBunchCrossNumber(); + for (Int_t icell = 0; icell < cells->GetNumberOfCells(); icell++) { if (cells->GetCell(icell, cellNumber, cellAmplitude, cellTime) != kTRUE) @@ -629,10 +633,17 @@ void AliAnalysisTaskEMCALClusterize::UserExec(Option_t *) //printf("CalibFactor %f times %f for id %d\n",fRecoUtils->GetEMCALChannelRecalibrationFactor(imod,ieta,iphi),amp,id); amp *=fRecoUtils->GetEMCALChannelRecalibrationFactor(imod,ieta,iphi); } - + // In case of AOD analysis cell time is 0, approximate replacing by time of the cluster the digit belongs. if (time*1e9 < 1.) time = fCellTime[id]; + // Recalibrate time + fRecoUtils->RecalibrateCellTime(id,bc,time); + +// printf("Clusterizer: Id %d, Time org %e, Time new %e; Amp org %f, Amp new %f\n", +// id, cells->GetTime(icell),time, cells->GetAmplitude(icell),amp); + + //Create the digit, put a fake primary deposited energy to trick the clusterizer when checking the most likely primary new((*fDigitsArr)[idigit]) AliEMCALDigit( fCellLabels[id], fCellLabels[id],id, amp, time,AliEMCALDigit::kHG,idigit, 0, 0, 1); //if(fCellLabels[id]>=0)printf("2) Digit cell %d, label %d\n",id,fCellLabels[id]) ; @@ -705,8 +716,8 @@ void AliAnalysisTaskEMCALClusterize::UserExec(Option_t *) AliAODCaloCluster *newCluster = (AliAODCaloCluster *) fCaloClusterArr->At(i); //if(Entry()==0) Info("UserExec","newCluster E %f\n", newCluster->E()); - //Add matched track, if any, only with ESDs - if(esdevent && fDoTrackMatching){ + //Add matched track + if(fDoTrackMatching){ Int_t trackIndex = fRecoUtils->GetMatchedTrackIndex(i); if(trackIndex >= 0){ newCluster->AddTrackMatched(event->GetTrack(trackIndex)); diff --git a/PWG4/PartCorrBase/AliCaloTrackReader.cxx b/PWG4/PartCorrBase/AliCaloTrackReader.cxx index 4a3beeb24ac..045a21ad2c8 100755 --- a/PWG4/PartCorrBase/AliCaloTrackReader.cxx +++ b/PWG4/PartCorrBase/AliCaloTrackReader.cxx @@ -51,25 +51,31 @@ ClassImp(AliCaloTrackReader) //____________________________________________________________________________ AliCaloTrackReader::AliCaloTrackReader() : - TObject(), fEventNumber(-1), //fCurrentFileName(""), - fDataType(0), fDebug(0), - fFiducialCut(0x0), fCheckFidCut(kFALSE), fComparePtHardAndJetPt(kFALSE), fPtHardAndJetPtFactor(7), - fCTSPtMin(0), fEMCALPtMin(0),fPHOSPtMin(0), - fCTSPtMax(1000), fEMCALPtMax(1000),fPHOSPtMax(1000), fAODBranchList(new TList ), + TObject(), fEventNumber(-1), //fCurrentFileName(""), + fDataType(0), fDebug(0), + fFiducialCut(0x0), fCheckFidCut(kFALSE), + fComparePtHardAndJetPt(0), fPtHardAndJetPtFactor(7), + fCTSPtMin(0), fEMCALPtMin(0), fPHOSPtMin(0), + fCTSPtMax(1000), fEMCALPtMax(1000), fPHOSPtMax(1000), + fAODBranchList(new TList ), fCTSTracks(new TObjArray()), fEMCALClusters(new TObjArray()), fPHOSClusters(new TObjArray()), - fEMCALCells(0x0), fPHOSCells(0x0), - fInputEvent(0x0), fOutputEvent(0x0),fMC(0x0), - fFillCTS(0),fFillEMCAL(0),fFillPHOS(0), - fFillEMCALCells(0),fFillPHOSCells(0), fSelectEmbeddedClusters(kFALSE), - fTrackStatus(0), fTrackFilterMask(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),fCaloFilterPatch(kFALSE), - fEMCALClustersListName(""),fZvtxCut(0.), - fDoEventSelection(kFALSE), fDoV0ANDEventSelection(kFALSE), fUseEventsWithPrimaryVertex(kFALSE), - fTriggerAnalysis (new AliTriggerAnalysis), fCentralityClass("V0M"),fCentralityOpt(10), + fEMCALCells(0x0), fPHOSCells(0x0), + fInputEvent(0x0), fOutputEvent(0x0),fMC(0x0), + fFillCTS(0), fFillEMCAL(0), fFillPHOS(0), + fFillEMCALCells(0), fFillPHOSCells(0), + fRecalculateClusters(kFALSE),fSelectEmbeddedClusters(kFALSE), + fTrackStatus(0), fTrackFilterMask(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), fCaloFilterPatch(kFALSE), + fEMCALClustersListName(""), fZvtxCut(0.), + fDoEventSelection(kFALSE), fDoV0ANDEventSelection(kFALSE), fUseEventsWithPrimaryVertex(kFALSE), + fTriggerAnalysis (new AliTriggerAnalysis), + fCentralityClass("V0M"), fCentralityOpt(10), fEventPlaneMethod("Q") { @@ -322,7 +328,8 @@ void AliCaloTrackReader::Print(const Option_t * opt) const printf("Track filter mask (AODs) = %d\n", (Int_t) fTrackFilterMask) ; printf("Track Mult Eta Cut = %d\n", (Int_t) fTrackMultEtaCut) ; printf("Write delta AOD = %d\n", fWriteOutputDeltaAOD) ; - + printf("Recalculate Clusters = %d\n", fRecalculateClusters) ; + if(fComparePtHardAndJetPt) printf("Compare jet pt and pt hard to accept event, factor = %2.2f",fPtHardAndJetPtFactor); @@ -681,12 +688,12 @@ void AliCaloTrackReader::FillInputEMCALAlgorithm(AliVCluster * clus, const Int_t //Reject clusters with bad channels, close to borders and exotic; if(!GetCaloUtils()->GetEMCALRecoUtils()->IsGoodCluster(clus,GetCaloUtils()->GetEMCALGeometry(),GetEMCALCells())) return; -// //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; -// + // //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; + // //Mask all cells in collumns facing ALICE thick material if requested if(GetCaloUtils()->GetNMaskCellColumns()){ @@ -721,35 +728,37 @@ void AliCaloTrackReader::FillInputEMCALAlgorithm(AliVCluster * clus, const Int_t //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]); + if(fRecalculateClusters){ + //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()); - } - - //In case of MC analysis, to match resolution/calibration in real data - Float_t rdmEnergy = GetCaloUtils()->GetEMCALRecoUtils()->SmearClusterEnergy(clus); - // printf("\t Energy %f, smeared %f\n", clus->E(),rdmEnergy); - clus->SetE(rdmEnergy); + + //In case of MC analysis, to match resolution/calibration in real data + Float_t rdmEnergy = GetCaloUtils()->GetEMCALRecoUtils()->SmearClusterEnergy(clus); + // printf("\t Energy %f, smeared %f\n", clus->E(),rdmEnergy); + clus->SetE(rdmEnergy); + } if (fMixedEvent) clus->SetID(iclus) ; @@ -764,6 +773,12 @@ void AliCaloTrackReader::FillInputEMCAL() { if(fDebug > 2 ) printf("AliCaloTrackReader::FillInputEMCAL()\n"); + // First recalibrate cells, time or energy + // if(GetCaloUtils()->IsRecalibrationOn()) + // GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCells(GetCaloUtils()->GetEMCALGeometry(), + // GetEMCALCells(), + // fInputEvent->GetBunchCrossNumber()); + //Loop to select clusters in fiducial cut and fill container with aodClusters if(fEMCALClustersListName==""){ Int_t nclusters = fInputEvent->GetNumberOfCaloClusters(); @@ -845,10 +860,14 @@ void AliCaloTrackReader::FillInputPHOS() { printf("AliCaloTrackReader::FillInputPHOS() - 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()); + if(fRecalculateClusters){ + //Recalibrate the cluster energy - if(GetCaloUtils()->IsRecalibrationOn()) { - Float_t energy = GetCaloUtils()->RecalibrateClusterEnergy(clus, (AliAODCaloCells*)GetPHOSCells()); - clus->SetE(energy); + if(GetCaloUtils()->IsRecalibrationOn()) { + Float_t energy = GetCaloUtils()->RecalibrateClusterEnergy(clus, (AliAODCaloCells*)GetPHOSCells()); + clus->SetE(energy); + } + } if (fMixedEvent) { diff --git a/PWG4/PartCorrBase/AliCaloTrackReader.h b/PWG4/PartCorrBase/AliCaloTrackReader.h index 2a2d8436339..e62e9cf67b9 100755 --- a/PWG4/PartCorrBase/AliCaloTrackReader.h +++ b/PWG4/PartCorrBase/AliCaloTrackReader.h @@ -156,6 +156,10 @@ public: void SwitchOnPHOSCells() { fFillPHOSCells = kTRUE ; } void SwitchOffPHOSCells() { fFillPHOSCells = kFALSE ; } + Bool_t AreClustersRecalculated() const { return fRecalculateClusters ; } + void SwitchOnClusterRecalculation() { fRecalculateClusters = kTRUE ; } + void SwitchOffClusterRecalculation() { fRecalculateClusters = kFALSE ; } + Bool_t IsEmbeddedClusterSelectionOn() const { return fSelectEmbeddedClusters ; } void SwitchOnEmbeddedClustersSelection() { fSelectEmbeddedClusters = kTRUE ; } void SwitchOffEmbeddedClustersSelection() { fSelectEmbeddedClusters = kFALSE ; } @@ -353,7 +357,8 @@ public: Bool_t fFillPHOS; // use data from PHOS Bool_t fFillEMCALCells; // use data from EMCAL Bool_t fFillPHOSCells; // use data from PHOS - Bool_t fSelectEmbeddedClusters; // Use only simulated clusters that come from embedding. + Bool_t fRecalculateClusters; // Correct clusters, recalculate them if recalibration parameters is given + Bool_t fSelectEmbeddedClusters; // Use only simulated clusters that come from embedding. ULong_t fTrackStatus ; // Track selection bit, select tracks refitted in TPC, ITS ... ULong_t fTrackFilterMask ; // Track selection bit, for AODs (any difference with track status?) @@ -396,7 +401,7 @@ public: Int_t fCentralityBin[2]; // Minimum and maximum value of the centrality for the analysis TString fEventPlaneMethod; // Name of event plane method, by default "Q" - ClassDef(AliCaloTrackReader,32) + ClassDef(AliCaloTrackReader,33) } ; diff --git a/PWG4/PartCorrDep/AliAnaCalorimeterQA.cxx b/PWG4/PartCorrDep/AliAnaCalorimeterQA.cxx index b653d963c65..9d1be340be1 100755 --- a/PWG4/PartCorrDep/AliAnaCalorimeterQA.cxx +++ b/PWG4/PartCorrDep/AliAnaCalorimeterQA.cxx @@ -144,12 +144,12 @@ fhGenMCE(), fhGenMCEtaPhi(), fhGenMCAccE(), fhGenMCAccEtaPhi(), //matched MC -fhEMVxyz(0), fhEMR(0), fhHaVxyz(0), fhHaR(0) -//fh1pOverE(0), fh1dR(0), fh2EledEdx(0), fh2MatchdEdx(0), -//fhMCEle1pOverE(0), fhMCEle1dR(0), fhMCEle2MatchdEdx(0), -//fhMCChHad1pOverE(0), fhMCChHad1dR(0), fhMCChHad2MatchdEdx(0), -//fhMCNeutral1pOverE(0), fhMCNeutral1dR(0), fhMCNeutral2MatchdEdx(0),fh1pOverER02(0), -//fhMCEle1pOverER02(0), fhMCChHad1pOverER02(0), fhMCNeutral1pOverER02(0) +fhEMVxyz(0), fhEMR(0), fhHaVxyz(0), fhHaR(0), +fh1pOverE(0), fh1dR(0), fh2EledEdx(0), fh2MatchdEdx(0), +fhMCEle1pOverE(0), fhMCEle1dR(0), fhMCEle2MatchdEdx(0), +fhMCChHad1pOverE(0), fhMCChHad1dR(0), fhMCChHad2MatchdEdx(0), +fhMCNeutral1pOverE(0), fhMCNeutral1dR(0), fhMCNeutral2MatchdEdx(0),fh1pOverER02(0), +fhMCEle1pOverER02(0), fhMCChHad1pOverER02(0), fhMCNeutral1pOverER02(0) { //Default Ctor @@ -199,9 +199,9 @@ TList * AliAnaCalorimeterQA::GetCreateOutputObjects() Int_t netabins = GetHistoEtaBins(); Float_t etamax = GetHistoEtaMax(); Float_t etamin = GetHistoEtaMin(); Int_t nmassbins = GetHistoMassBins(); Float_t massmax = GetHistoMassMax(); Float_t massmin = GetHistoMassMin(); Int_t nasymbins = GetHistoAsymmetryBins(); Float_t asymmax = GetHistoAsymmetryMax(); Float_t asymmin = GetHistoAsymmetryMin(); - //Int_t nPoverEbins = GetHistoPOverEBins(); Float_t pOverEmax = GetHistoPOverEMax(); Float_t pOverEmin = GetHistoPOverEMin(); - //Int_t ndedxbins = GetHistodEdxBins(); Float_t dedxmax = GetHistodEdxMax(); Float_t dedxmin = GetHistodEdxMin(); - //Int_t ndRbins = GetHistodRBins(); Float_t dRmax = GetHistodRMax(); Float_t dRmin = GetHistodRMin(); + Int_t nPoverEbins = GetHistoPOverEBins(); Float_t pOverEmax = GetHistoPOverEMax(); Float_t pOverEmin = GetHistoPOverEMin(); + Int_t ndedxbins = GetHistodEdxBins(); Float_t dedxmax = GetHistodEdxMax(); Float_t dedxmin = GetHistodEdxMin(); + Int_t ndRbins = GetHistodRBins(); Float_t dRmax = GetHistodRMax(); Float_t dRmin = GetHistodRMin(); Int_t ntimebins = GetHistoTimeBins(); Float_t timemax = GetHistoTimeMax(); Float_t timemin = GetHistoTimeMin(); Int_t nclbins = GetHistoNClustersBins(); Int_t nclmax = GetHistoNClustersMax(); Int_t nclmin = GetHistoNClustersMin(); Int_t ncebins = GetHistoNCellsBins(); Int_t ncemax = GetHistoNCellsMax(); Int_t ncemin = GetHistoNCellsMin(); @@ -548,29 +548,29 @@ TList * AliAnaCalorimeterQA::GetCreateOutputObjects() outputContainer->Add(fhEtaPhiECharged); } -// fh1pOverE = new TH2F("h1pOverE","TRACK matches p/E",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax); -// fh1pOverE->SetYTitle("p/E"); -// fh1pOverE->SetXTitle("p_{T} (GeV/c)"); -// outputContainer->Add(fh1pOverE); -// -// fh1dR = new TH1F("h1dR","TRACK matches dR",ndRbins,dRmin,dRmax); -// fh1dR->SetXTitle("#Delta R (rad)"); -// outputContainer->Add(fh1dR) ; -// -// fh2MatchdEdx = new TH2F("h2MatchdEdx","dE/dx vs. p for all matches",nptbins,ptmin,ptmax,ndedxbins,dedxmin,dedxmax); -// fh2MatchdEdx->SetXTitle("p (GeV/c)"); -// fh2MatchdEdx->SetYTitle(""); -// outputContainer->Add(fh2MatchdEdx); -// -// fh2EledEdx = new TH2F("h2EledEdx","dE/dx vs. p for electrons",nptbins,ptmin,ptmax,ndedxbins,dedxmin,dedxmax); -// fh2EledEdx->SetXTitle("p (GeV/c)"); -// fh2EledEdx->SetYTitle(""); -// outputContainer->Add(fh2EledEdx) ; -// -// fh1pOverER02 = new TH2F("h1pOverER02","TRACK matches p/E, all",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax); -// fh1pOverER02->SetYTitle("p/E"); -// fh1pOverER02->SetXTitle("p_{T} (GeV/c)"); -// outputContainer->Add(fh1pOverER02); + fh1pOverE = new TH2F("h1pOverE","TRACK matches p/E",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax); + fh1pOverE->SetYTitle("p/E"); + fh1pOverE->SetXTitle("p_{T} (GeV/c)"); + outputContainer->Add(fh1pOverE); + + fh1dR = new TH1F("h1dR","TRACK matches dR",ndRbins,dRmin,dRmax); + fh1dR->SetXTitle("#Delta R (rad)"); + outputContainer->Add(fh1dR) ; + + fh2MatchdEdx = new TH2F("h2MatchdEdx","dE/dx vs. p for all matches",nptbins,ptmin,ptmax,ndedxbins,dedxmin,dedxmax); + fh2MatchdEdx->SetXTitle("p (GeV/c)"); + fh2MatchdEdx->SetYTitle(""); + outputContainer->Add(fh2MatchdEdx); + + fh2EledEdx = new TH2F("h2EledEdx","dE/dx vs. p for electrons",nptbins,ptmin,ptmax,ndedxbins,dedxmin,dedxmax); + fh2EledEdx->SetXTitle("p (GeV/c)"); + fh2EledEdx->SetYTitle(""); + outputContainer->Add(fh2EledEdx) ; + + fh1pOverER02 = new TH2F("h1pOverER02","TRACK matches p/E, all",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax); + fh1pOverER02->SetYTitle("p/E"); + fh1pOverER02->SetXTitle("p_{T} (GeV/c)"); + outputContainer->Add(fh1pOverER02); } if(fFillAllPi0Histo){ @@ -1113,62 +1113,62 @@ TList * AliAnaCalorimeterQA::GetCreateOutputObjects() //Track Matching -// fhMCEle1pOverE = new TH2F("hMCEle1pOverE","TRACK matches p/E, MC electrons",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax); -// fhMCEle1pOverE->SetYTitle("p/E"); -// fhMCEle1pOverE->SetXTitle("p_{T} (GeV/c)"); -// outputContainer->Add(fhMCEle1pOverE); -// -// fhMCEle1dR = new TH1F("hMCEle1dR","TRACK matches dR, MC electrons",ndRbins,dRmin,dRmax); -// fhMCEle1dR->SetXTitle("#Delta R (rad)"); -// outputContainer->Add(fhMCEle1dR) ; -// -// fhMCEle2MatchdEdx = new TH2F("hMCEle2MatchdEdx","dE/dx vs. p for all matches, MC electrons",nptbins,ptmin,ptmax,ndedxbins,dedxmin,dedxmax); -// fhMCEle2MatchdEdx->SetXTitle("p (GeV/c)"); -// fhMCEle2MatchdEdx->SetYTitle(""); -// outputContainer->Add(fhMCEle2MatchdEdx); -// -// fhMCChHad1pOverE = new TH2F("hMCChHad1pOverE","TRACK matches p/E, MC charged hadrons",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax); -// fhMCChHad1pOverE->SetYTitle("p/E"); -// fhMCChHad1pOverE->SetXTitle("p_{T} (GeV/c)"); -// outputContainer->Add(fhMCChHad1pOverE); -// -// fhMCChHad1dR = new TH1F("hMCChHad1dR","TRACK matches dR, MC charged hadrons",ndRbins,dRmin,dRmax); -// fhMCChHad1dR->SetXTitle("#Delta R (rad)"); -// outputContainer->Add(fhMCChHad1dR) ; -// -// fhMCChHad2MatchdEdx = new TH2F("hMCChHad2MatchdEdx","dE/dx vs. p for all matches, MC charged hadrons",nptbins,ptmin,ptmax,ndedxbins,dedxmin,dedxmax); -// fhMCChHad2MatchdEdx->SetXTitle("p (GeV/c)"); -// fhMCChHad2MatchdEdx->SetYTitle(""); -// outputContainer->Add(fhMCChHad2MatchdEdx); -// -// fhMCNeutral1pOverE = new TH2F("hMCNeutral1pOverE","TRACK matches p/E, MC neutrals",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax); -// fhMCNeutral1pOverE->SetYTitle("p/E"); -// fhMCNeutral1pOverE->SetXTitle("p_{T} (GeV/c)"); -// outputContainer->Add(fhMCNeutral1pOverE); -// -// fhMCNeutral1dR = new TH1F("hMCNeutral1dR","TRACK matches dR, MC neutrals",ndRbins,dRmin,dRmax); -// fhMCNeutral1dR->SetXTitle("#Delta R (rad)"); -// outputContainer->Add(fhMCNeutral1dR) ; -// -// fhMCNeutral2MatchdEdx = new TH2F("hMCNeutral2MatchdEdx","dE/dx vs. p for all matches, MC neutrals",nptbins,ptmin,ptmax,ndedxbins,dedxmin,dedxmax); -// fhMCNeutral2MatchdEdx->SetXTitle("p (GeV/c)"); -// fhMCNeutral2MatchdEdx->SetYTitle(""); -// outputContainer->Add(fhMCNeutral2MatchdEdx); -// -// fhMCEle1pOverER02 = new TH2F("hMCEle1pOverER02","TRACK matches p/E, MC electrons",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax); -// fhMCEle1pOverER02->SetYTitle("p/E"); -// fhMCEle1pOverER02->SetXTitle("p_{T} (GeV/c)"); -// outputContainer->Add(fhMCEle1pOverER02); -// -// fhMCChHad1pOverER02 = new TH2F("hMCChHad1pOverER02","TRACK matches p/E, MC charged hadrons",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax); -// fhMCChHad1pOverER02->SetYTitle("p/E"); -// fhMCChHad1pOverER02->SetXTitle("p_{T} (GeV/c)"); -// outputContainer->Add(fhMCChHad1pOverER02); -// -// fhMCNeutral1pOverER02 = new TH2F("hMCNeutral1pOverER02","TRACK matches p/E, MC neutrals",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax); -// fhMCNeutral1pOverER02->SetYTitle("p/E"); -// fhMCNeutral1pOverER02->SetXTitle("p_{T} (GeV/c)"); -// outputContainer->Add(fhMCNeutral1pOverER02); + fhMCEle1pOverE = new TH2F("hMCEle1pOverE","TRACK matches p/E, MC electrons",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax); + fhMCEle1pOverE->SetYTitle("p/E"); + fhMCEle1pOverE->SetXTitle("p_{T} (GeV/c)"); + outputContainer->Add(fhMCEle1pOverE); + + fhMCEle1dR = new TH1F("hMCEle1dR","TRACK matches dR, MC electrons",ndRbins,dRmin,dRmax); + fhMCEle1dR->SetXTitle("#Delta R (rad)"); + outputContainer->Add(fhMCEle1dR) ; + + fhMCEle2MatchdEdx = new TH2F("hMCEle2MatchdEdx","dE/dx vs. p for all matches, MC electrons",nptbins,ptmin,ptmax,ndedxbins,dedxmin,dedxmax); + fhMCEle2MatchdEdx->SetXTitle("p (GeV/c)"); + fhMCEle2MatchdEdx->SetYTitle(""); + outputContainer->Add(fhMCEle2MatchdEdx); + + fhMCChHad1pOverE = new TH2F("hMCChHad1pOverE","TRACK matches p/E, MC charged hadrons",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax); + fhMCChHad1pOverE->SetYTitle("p/E"); + fhMCChHad1pOverE->SetXTitle("p_{T} (GeV/c)"); + outputContainer->Add(fhMCChHad1pOverE); + + fhMCChHad1dR = new TH1F("hMCChHad1dR","TRACK matches dR, MC charged hadrons",ndRbins,dRmin,dRmax); + fhMCChHad1dR->SetXTitle("#Delta R (rad)"); + outputContainer->Add(fhMCChHad1dR) ; + + fhMCChHad2MatchdEdx = new TH2F("hMCChHad2MatchdEdx","dE/dx vs. p for all matches, MC charged hadrons",nptbins,ptmin,ptmax,ndedxbins,dedxmin,dedxmax); + fhMCChHad2MatchdEdx->SetXTitle("p (GeV/c)"); + fhMCChHad2MatchdEdx->SetYTitle(""); + outputContainer->Add(fhMCChHad2MatchdEdx); + + fhMCNeutral1pOverE = new TH2F("hMCNeutral1pOverE","TRACK matches p/E, MC neutrals",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax); + fhMCNeutral1pOverE->SetYTitle("p/E"); + fhMCNeutral1pOverE->SetXTitle("p_{T} (GeV/c)"); + outputContainer->Add(fhMCNeutral1pOverE); + + fhMCNeutral1dR = new TH1F("hMCNeutral1dR","TRACK matches dR, MC neutrals",ndRbins,dRmin,dRmax); + fhMCNeutral1dR->SetXTitle("#Delta R (rad)"); + outputContainer->Add(fhMCNeutral1dR) ; + + fhMCNeutral2MatchdEdx = new TH2F("hMCNeutral2MatchdEdx","dE/dx vs. p for all matches, MC neutrals",nptbins,ptmin,ptmax,ndedxbins,dedxmin,dedxmax); + fhMCNeutral2MatchdEdx->SetXTitle("p (GeV/c)"); + fhMCNeutral2MatchdEdx->SetYTitle(""); + outputContainer->Add(fhMCNeutral2MatchdEdx); + + fhMCEle1pOverER02 = new TH2F("hMCEle1pOverER02","TRACK matches p/E, MC electrons",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax); + fhMCEle1pOverER02->SetYTitle("p/E"); + fhMCEle1pOverER02->SetXTitle("p_{T} (GeV/c)"); + outputContainer->Add(fhMCEle1pOverER02); + + fhMCChHad1pOverER02 = new TH2F("hMCChHad1pOverER02","TRACK matches p/E, MC charged hadrons",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax); + fhMCChHad1pOverER02->SetYTitle("p/E"); + fhMCChHad1pOverER02->SetXTitle("p_{T} (GeV/c)"); + outputContainer->Add(fhMCChHad1pOverER02); + + fhMCNeutral1pOverER02 = new TH2F("hMCNeutral1pOverER02","TRACK matches p/E, MC neutrals",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax); + fhMCNeutral1pOverER02->SetYTitle("p/E"); + fhMCNeutral1pOverER02->SetXTitle("p_{T} (GeV/c)"); + outputContainer->Add(fhMCNeutral1pOverER02); } // for(Int_t i = 0; i < outputContainer->GetEntries() ; i++) @@ -1228,22 +1228,25 @@ void AliAnaCalorimeterQA::Print(const Option_t * opt) const void AliAnaCalorimeterQA::MakeAnalysisFillHistograms() { //Fill Calorimeter QA histograms + + TLorentzVector mom ; TLorentzVector mom2 ; TObjArray * caloClusters = NULL; - Int_t nLabel = 0 ; - Int_t *labels = 0x0; - Int_t nCaloClusters = 0; - Int_t nCaloClustersAccepted = 0; - Int_t nCaloCellsPerCluster = 0; - Bool_t matched = kFALSE; - Int_t nModule = -1; - + Int_t nLabel = 0 ; + Int_t *labels = 0x0; + Int_t nCaloClusters = 0 ; + Int_t nCaloClustersAccepted = 0 ; + Int_t nCaloCellsPerCluster = 0 ; + Bool_t matched = kFALSE; + Int_t nModule =-1 ; + Int_t bc = GetReader()->GetInputEvent()->GetBunchCrossNumber(); + //Get vertex for photon momentum calculation and event selection Double_t v[3] = {0,0,0}; //vertex ; GetReader()->GetVertex(v); if (TMath::Abs(v[2]) > GetZvertexCut()) return ; - + //Play with the MC stack if available //Get the MC arrays and do some checks if(IsDataMC()){ @@ -1255,7 +1258,7 @@ void AliAnaCalorimeterQA::MakeAnalysisFillHistograms() //Fill some pure MC histograms, only primaries. for(Int_t i=0 ; iGetNprimary(); i++){//Only primary particles, for all MC transport put GetNtrack() TParticle *primary = GetMCStack()->Particle(i) ; - //printf("i %d, %s: status = %d, primary? %d\n",i, primary->GetName(), primary->GetStatusCode(), primary->IsPrimary()); + if (primary->GetStatusCode() > 11) continue; //Working for PYTHIA and simple generators, check for HERWIG primary->Momentum(mom); MCHistograms(mom,TMath::Abs(primary->GetPdgCode())); @@ -1269,11 +1272,9 @@ void AliAnaCalorimeterQA::MakeAnalysisFillHistograms() //Fill some pure MC histograms, only primaries. for(Int_t i=0 ; i < (GetReader()->GetAODMCParticles(0))->GetEntriesFast(); i++){ AliAODMCParticle *aodprimary = (AliAODMCParticle*) (GetReader()->GetAODMCParticles(0))->At(i) ; - //printf("i %d, %s: primary? %d physical primary? %d, flag %d\n", - // i,(TDatabasePDG::Instance()->GetParticle(aodprimary->GetPdgCode()))->GetName(), - // aodprimary->IsPrimary(), aodprimary->IsPhysicalPrimary(), aodprimary->GetFlag()); + if (!aodprimary->IsPrimary()) continue; //accept all which is not MC transport generated. Don't know how to avoid partons - //aodprimary->Momentum(mom); + mom.SetPxPyPzE(aodprimary->Px(),aodprimary->Py(),aodprimary->Pz(),aodprimary->E()); MCHistograms(mom,TMath::Abs(aodprimary->GetPdgCode())); } //primary loop @@ -1281,436 +1282,503 @@ void AliAnaCalorimeterQA::MakeAnalysisFillHistograms() } }// is data and MC - //Get List with CaloClusters if (fCalorimeter == "PHOS") caloClusters = GetPHOSClusters(); else if (fCalorimeter == "EMCAL") caloClusters = GetEMCALClusters(); else AliFatal(Form("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - Wrong calorimeter name <%s>, END\n", fCalorimeter.Data())); - - if(!caloClusters) { - AliFatal(Form("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - No CaloClusters available\n")); + //Get List with CaloCells + AliVCaloCells * cell = 0x0; + if(fCalorimeter == "PHOS") cell = GetPHOSCells(); + else cell = GetEMCALCells(); + + if(!caloClusters || !cell) { + AliFatal(Form("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - No CaloClusters or CaloCells available\n")); + return; // trick coverity } - else{ + + //---------------------------------------------------------- + //Correlate Calorimeters and V0 and track Multiplicity + //---------------------------------------------------------- + + if(fCorrelate) Correlate(); + + //---------------------------------------------------------- + // CALOCLUSTERS + //---------------------------------------------------------- + + nCaloClusters = caloClusters->GetEntriesFast() ; + Int_t *nClustersInModule = new Int_t[fNModules]; + for(Int_t imod = 0; imod < fNModules; imod++ ) nClustersInModule[imod] = 0; + + if(GetDebug() > 0) + printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - In %s there are %d clusters \n", fCalorimeter.Data(), nCaloClusters); + + AliVTrack * track = 0x0; + Float_t pos[3] ; + Double_t tof = 0; + + //Loop over CaloClusters + //if(nCaloClusters > 0)printf("QA : Vertex Cut passed %f, cut %f, entries %d, %s\n",v[2], 40., nCaloClusters, fCalorimeter.Data()); + for(Int_t iclus = 0; iclus < nCaloClusters; iclus++){ - //---------------------------------------------------------- - //Correlate Calorimeters and V0 and track Multiplicity - //---------------------------------------------------------- - - if(fCorrelate) Correlate(); + if(GetDebug() > 0) printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - cluster: %d/%d, data %d \n", + iclus+1,nCaloClusters,GetReader()->GetDataType()); - //---------------------------------------------------------- - // CALOCLUSTERS - //---------------------------------------------------------- + AliVCluster* clus = (AliVCluster*)caloClusters->At(iclus); + + //Get cluster kinematics + clus->GetPosition(pos); + clus->GetMomentum(mom,v); - nCaloClusters = caloClusters->GetEntriesFast() ; - Int_t *nClustersInModule = new Int_t[fNModules]; - for(Int_t imod = 0; imod < fNModules; imod++ ) nClustersInModule[imod] = 0; + //Check only certain regions + Bool_t in = kTRUE; + if(IsFiducialCutOn()) in = GetFiducialCut()->IsInFiducialCut(mom,fCalorimeter) ; + if(!in) continue; + //MC labels + nLabel = clus->GetNLabels(); + labels = clus->GetLabels(); - if(GetDebug() > 0) - printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - In %s there are %d clusters \n", fCalorimeter.Data(), nCaloClusters); + //Cells per cluster + nCaloCellsPerCluster = clus->GetNCells(); + //if(mom.E() > 10 && nCaloCellsPerCluster == 1 ) printf("%s:************** E = %f ********** ncells = %d\n",fCalorimeter.Data(), mom.E(),nCaloCellsPerCluster); - //AliVTrack * track = 0x0; - Float_t pos[3] ; - Double_t tof = 0; - //Loop over CaloClusters - //if(nCaloClusters > 0)printf("QA : Vertex Cut passed %f, cut %f, entries %d, %s\n",v[2], 40., nCaloClusters, fCalorimeter.Data()); - for(Int_t iclus = 0; iclus < nCaloClusters; iclus++){ - - if(GetDebug() > 0) printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - cluster: %d/%d, data %d \n", - iclus+1,nCaloClusters,GetReader()->GetDataType()); - - AliVCluster* clus = (AliVCluster*)caloClusters->At(iclus); - AliVCaloCells * cell = 0x0; - if(fCalorimeter == "PHOS") cell = GetPHOSCells(); - else cell = GetEMCALCells(); - - //Get cluster kinematics - clus->GetPosition(pos); - clus->GetMomentum(mom,v); - tof = clus->GetTOF()*1e9; - if(tof < fTimeCutMin || tof > fTimeCutMax) continue; - - //Check only certain regions - Bool_t in = kTRUE; - if(IsFiducialCutOn()) in = GetFiducialCut()->IsInFiducialCut(mom,fCalorimeter) ; - if(!in) continue; - - //MC labels - nLabel = clus->GetNLabels(); - labels = clus->GetLabels(); - - //Cells per cluster - nCaloCellsPerCluster = clus->GetNCells(); - //if(mom.E() > 10 && nCaloCellsPerCluster == 1 ) printf("%s:************** E = %f ********** ncells = %d\n",fCalorimeter.Data(), mom.E(),nCaloCellsPerCluster); - - // Cluster mathed with track? - matched = GetCaloPID()->IsTrackMatched(clus,GetCaloUtils()); - - - //====================== - //Cells in cluster - //====================== - - //Get list of contributors - UShort_t * indexList = clus->GetCellsAbsId() ; - - if(fFillAllPosHisto) FillCellPositionHistograms(nCaloCellsPerCluster,indexList,pos,mom.E()); - - // Get the fraction of the cluster energy that carries the cell with highest energy - Float_t maxCellFraction = 0.; - Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(cell, clus,maxCellFraction); - Int_t smMax =0; Int_t ietaaMax=-1; Int_t iphiiMax = 0; Int_t rcuMax = 0; - smMax = GetModuleNumberCellIndexes(absIdMax,fCalorimeter, ietaaMax, iphiiMax, rcuMax); - Int_t dIeta = 0; - Int_t dIphi = 0; - Double_t tmax = cell->GetCellTime(absIdMax)*1e9; - //Float_t emax = cell->GetCellAmplitude(absIdMax); - - if (clus->E() < 2.){ - fhLambda0vsClusterMaxCellDiffE0->Fill(clus->GetM02(), maxCellFraction); - fhNCellsvsClusterMaxCellDiffE0 ->Fill(nCaloCellsPerCluster,maxCellFraction); - } - else if(clus->E() < 6.){ - fhLambda0vsClusterMaxCellDiffE2->Fill(clus->GetM02(), maxCellFraction); - fhNCellsvsClusterMaxCellDiffE2 ->Fill(nCaloCellsPerCluster,maxCellFraction); - } - else{ - fhLambda0vsClusterMaxCellDiffE6->Fill(clus->GetM02(), maxCellFraction); - fhNCellsvsClusterMaxCellDiffE6 ->Fill(nCaloCellsPerCluster,maxCellFraction); + // Cluster mathed with track? + matched = GetCaloPID()->IsTrackMatched(clus,GetCaloUtils()); + + + //====================== + //Cells in cluster + //====================== + + //Get list of contributors + UShort_t * indexList = clus->GetCellsAbsId() ; + + if(fFillAllPosHisto) FillCellPositionHistograms(nCaloCellsPerCluster,indexList,pos,mom.E()); + + // Get the fraction of the cluster energy that carries the cell with highest energy + Float_t maxCellFraction = 0.; + Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(cell, clus,maxCellFraction); + Int_t smMax =0; Int_t ietaaMax=-1; Int_t iphiiMax = 0; Int_t rcuMax = 0; + smMax = GetModuleNumberCellIndexes(absIdMax,fCalorimeter, ietaaMax, iphiiMax, rcuMax); + Int_t dIeta = 0; + Int_t dIphi = 0; + + if (clus->E() < 2.){ + fhLambda0vsClusterMaxCellDiffE0->Fill(clus->GetM02(), maxCellFraction); + fhNCellsvsClusterMaxCellDiffE0 ->Fill(nCaloCellsPerCluster,maxCellFraction); + } + else if(clus->E() < 6.){ + fhLambda0vsClusterMaxCellDiffE2->Fill(clus->GetM02(), maxCellFraction); + fhNCellsvsClusterMaxCellDiffE2 ->Fill(nCaloCellsPerCluster,maxCellFraction); + } + else{ + fhLambda0vsClusterMaxCellDiffE6->Fill(clus->GetM02(), maxCellFraction); + fhNCellsvsClusterMaxCellDiffE6 ->Fill(nCaloCellsPerCluster,maxCellFraction); + } + + fhNCellsPerClusterNoCut ->Fill(clus->E(), nCaloCellsPerCluster); + nModule = GetModuleNumber(clus); + if(nModule >=0 && nModule < fNModules) fhNCellsPerClusterModNoCut[nModule]->Fill(clus->E(), nCaloCellsPerCluster); + + fhClusterMaxCellDiffNoCut->Fill(clus->E(),maxCellFraction); + + // Look inside the cluster + // Calculate average time of cells in cluster and weighted average + Float_t rawEnergy = 0; + Float_t factor = 1; + Double_t time = 0; + Double_t tmax = cell->GetCellTime(absIdMax); + Double_t averTime = 0; + Double_t weightedTime = 0; + Double_t weight = 0; + Double_t averTimeNoMax = 0; + Double_t weightedTimeNoMax = 0; + Double_t weightNoMax = 0; + + for (Int_t ipos = 0; ipos < nCaloCellsPerCluster; ipos++) { + Int_t icol = -1; + Int_t irow = -1; + Int_t iRCU = -1; + Int_t id = indexList[ipos]; + Int_t nModuleCell = GetModuleNumberCellIndexes(id,fCalorimeter, icol, irow, iRCU); + + //Recalibrate energy + if(GetCaloUtils()->IsRecalibrationOn()){ + if(fCalorimeter == "PHOS") { + factor = GetCaloUtils()->GetPHOSChannelRecalibrationFactor(nModuleCell,icol,irow); + } + else { + factor = GetCaloUtils()->GetEMCALChannelRecalibrationFactor(nModuleCell,icol,irow); + + } } - fhNCellsPerClusterNoCut ->Fill(clus->E(), nCaloCellsPerCluster); - nModule = GetModuleNumber(clus); - if(nModule >=0 && nModule < fNModules) fhNCellsPerClusterModNoCut[nModule]->Fill(clus->E(), nCaloCellsPerCluster); + //Recalibrate time + time = cell->GetCellTime(id); + if(fCalorimeter=="EMCAL" && GetCaloUtils()->GetEMCALRecoUtils()->IsTimeRecalibrationOn()){ + GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTime(id,bc,time); + if(id==absIdMax) tmax = time; + } - fhClusterMaxCellDiffNoCut->Fill(clus->E(),maxCellFraction); +// printf("QA cluster : Id %d, Time org %e, Time new %e; Amp org %f, Amp new %f\n", +// id, cell->GetCellTime(id),time, cell->GetCellAmplitude(id),cell->GetCellAmplitude(id)*factor); - // Calculate average time of cells in cluster and weighted average - Double_t averTime = 0; - Double_t weightedTime = 0; - Double_t weight = 0; - Double_t averTimeNoMax = 0; - Double_t weightedTimeNoMax = 0; - Double_t weightNoMax = 0; + // Get the energy of the cluster without the non linearity correction + rawEnergy += cell->GetCellAmplitude(id)*factor; - //if(GetReader()->GetDataType()==AliCaloTrackReader::kESD) { - Float_t rawEnergy = 0; - for (Int_t ipos = 0; ipos < nCaloCellsPerCluster; ipos++) - rawEnergy += cell->GetCellAmplitude(indexList[ipos]); + Double_t w = GetCaloUtils()->GetEMCALRecoUtils()->GetCellWeight(cell->GetCellAmplitude(id),rawEnergy); + averTime += time*1e9; + weightedTime += time*1e9 * w; + weight += w; + if(id != absIdMax){ + averTimeNoMax += time*1e9; + weightedTimeNoMax += time*1e9 * w; + weightNoMax += w; + } + } + + tmax*=1.e9; + averTime /= nCaloCellsPerCluster; + if(nCaloCellsPerCluster > 1 ) averTimeNoMax /= (nCaloCellsPerCluster-1); + + if(weight > 0 ) tof = weightedTime / weight; + else { + + tof = clus->GetTOF()*1e9; + + printf("AliAnaCalorimeterQA:: Null weight! Investigate: E %f GeV, ncells %d, time max cell %f ns, average time %f ns, absIdMax %d, SM %d\n", + rawEnergy,nCaloCellsPerCluster, tmax, averTime,absIdMax,GetModuleNumber(clus)); + } + + //printf("QA: E %f, tof org %g, tof new %g, ncells %d\n",clus->E(),clus->GetTOF()*1.e9,tof, clus->GetNCells()); + + if(weightNoMax > 0 ) weightedTimeNoMax /= weightNoMax; + + //Cut on time of clusters + if(tof < fTimeCutMin || tof > fTimeCutMax){ + printf("Remove cluster with TOF %f\n",tof); + continue; + } + + //====================== + //Bad clusters selection + //====================== + //Check bad clusters if rejection was not on + + Bool_t badCluster = kFALSE; + if(fCalorimeter=="EMCAL" && !GetCaloUtils()->GetEMCALRecoUtils()->IsRejectExoticCluster()){ + //Bad clusters histograms + Float_t minNCells = 1; + if(clus->E() > 7) minNCells = TMath::Max(1,TMath::FloorNint(1 + TMath::Log(clus->E() - 7 )*1.7 )); + if(nCaloCellsPerCluster <= minNCells) { + //if(clus->GetM02() < 0.05) { - for (Int_t ipos = 0; ipos < nCaloCellsPerCluster; ipos++) { - Double_t w = TMath::Max( 0., 4.5 + TMath::Log(cell->GetCellAmplitude(indexList[ipos]) / rawEnergy )); - averTime += cell->GetCellTime(indexList[ipos])*1e9; - weightedTime += cell->GetCellTime(indexList[ipos])*1e9 * w; - weight += w; - if(indexList[ipos]!=absIdMax){ - averTimeNoMax += cell->GetCellTime(indexList[ipos])*1e9; - weightedTimeNoMax += cell->GetCellTime(indexList[ipos])*1e9 * w; - weightNoMax += w; - } - } + badCluster = kTRUE; - averTime /= nCaloCellsPerCluster; - if(nCaloCellsPerCluster > 1 ) averTimeNoMax /= (nCaloCellsPerCluster-1); - - if(weight > 0 ) weightedTime /= weight; - else { - printf("AliAnaCalorimeterQA:: Null weight! Investigate: E %f GeV, ncells %d, time max cell %f ns, average time %f ns, absIdMax %d, SM %d\n", - rawEnergy,nCaloCellsPerCluster, tmax, averTime,absIdMax,GetModuleNumber(clus)); - } + fhBadClusterEnergy ->Fill(clus->E()); + fhBadClusterTimeEnergy ->Fill(clus->E(),tof); + fhBadClusterMaxCellDiff->Fill(clus->E(),maxCellFraction); - if(weightNoMax > 0 ) weightedTimeNoMax /= weightNoMax; - - //} // only possible in ESDs - - //====================== - //Bad clusters selection - //====================== - //Check bad clusters if rejection was not on - - Bool_t badCluster = kFALSE; - if(fCalorimeter=="EMCAL" && !GetCaloUtils()->GetEMCALRecoUtils()->IsRejectExoticCluster()){ - //Bad clusters histograms - Float_t minNCells = 1; - if(clus->E() > 7) minNCells = TMath::Max(1,TMath::FloorNint(1 + TMath::Log(clus->E() - 7 )*1.7 )); - if(nCaloCellsPerCluster <= minNCells) { - //if(clus->GetM02() < 0.05) { + if(clus->GetM02() > 0 || TMath::Abs(clus->GetM20()) > 0 || clus->GetDispersion() > 0){ - badCluster = kTRUE; + fhBadClusterL0 ->Fill(clus->E(),clus->GetM02()); + fhBadClusterL1 ->Fill(clus->E(),clus->GetM20()); + fhBadClusterD ->Fill(clus->E(),clus->GetDispersion()); - fhBadClusterEnergy ->Fill(clus->E()); - fhBadClusterTimeEnergy ->Fill(clus->E(),tof); - fhBadClusterMaxCellDiff->Fill(clus->E(),maxCellFraction); - - if(clus->GetM02() > 0 || TMath::Abs(clus->GetM20()) > 0 || clus->GetDispersion() > 0){ - - fhBadClusterL0 ->Fill(clus->E(),clus->GetM02()); - fhBadClusterL1 ->Fill(clus->E(),clus->GetM20()); - fhBadClusterD ->Fill(clus->E(),clus->GetDispersion()); - - } + } + + //Clusters in event time difference + + for(Int_t iclus2 = 0; iclus2 < nCaloClusters; iclus2++ ){ - //Clusters in event time difference + AliVCluster* clus2 = (AliVCluster*)caloClusters->At(iclus2); - for(Int_t iclus2 = 0; iclus2 < nCaloClusters; iclus2++ ){ - - AliVCluster* clus2 = (AliVCluster*)caloClusters->At(iclus2); - - if(clus->GetID()==clus2->GetID()) continue; - - if(clus->GetM02() > 0.01) - fhBadClusterPairDiffTimeE ->Fill(clus->E(), tof-clus2->GetTOF()*1.e9); - - } // loop + if(clus->GetID()==clus2->GetID()) continue; - // Max cell compared to other cells in cluster - if(GetReader()->GetDataType()==AliCaloTrackReader::kESD) { - fhBadClusterMaxCellDiffAverageTime->Fill(clus->E(),tmax-averTime); - fhBadClusterMaxCellDiffWeightTime ->Fill(clus->E(),tmax-weightedTime); - fhBadClusterDiffWeightAverTime ->Fill(clus->E(),weightedTime-averTime); - fhBadClusterMaxCellDiffAverageNoMaxTime->Fill(clus->E(),tmax-averTimeNoMax); - fhBadClusterMaxCellDiffWeightNoMaxTime ->Fill(clus->E(),tmax-weightedTimeNoMax); - //printf("E %f, tmax %f, aver %f, weigh %f, averNoMax %f, weightNoMax %f\n", - // clus->E(),tmax,averTime,weightedTime,averTimeNoMax,weightedTimeNoMax); - } - - if( weight > 0 ) fhBadClusterNoMaxCellWeight ->Fill(clus->E(),weightNoMax/weight); - //if( weight > 0 && weightNoMax > 0.0)printf("weight %f, weightNoMax %f\n",weight,weightNoMax); - - for (Int_t ipos = 0; ipos < nCaloCellsPerCluster; ipos++) { - Int_t absId = indexList[ipos]; - if(absId!=absIdMax){ - Float_t frac = cell->GetCellAmplitude(absId)/cell->GetCellAmplitude(absIdMax); - fhBadClusterMaxCellCloseCellRatio->Fill(mom.E(),frac); - fhBadClusterMaxCellCloseCellDiff ->Fill(mom.E(),cell->GetCellAmplitude(absIdMax)-cell->GetCellAmplitude(absId)); - if(GetReader()->GetDataType()==AliCaloTrackReader::kESD) { - Float_t diff = (tmax-cell->GetCellTime(absId)*1e9); - fhBadCellTimeSpreadRespectToCellMax->Fill(mom.E(), diff); + if(clus->GetM02() > 0.01 && clus2->GetM02() > 0.01) { + Double_t tof2 = clus2->GetTOF(); + // approximation in case of time recalibration + if(fCalorimeter=="EMCAL" && GetCaloUtils()->GetEMCALRecoUtils()->IsTimeRecalibrationOn()){ + Float_t maxCellFraction2 = -1; + Int_t absIdMax2 = GetCaloUtils()->GetMaxEnergyCell(cell, clus2,maxCellFraction2); + tof2 = cell->GetCellTime(absIdMax2); + GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTime(absIdMax2,bc,tof2); + } + fhBadClusterPairDiffTimeE ->Fill(clus->E(), tof-tof2*1.e9); + } + } // loop + + // Max cell compared to other cells in cluster + if(GetReader()->GetDataType()==AliCaloTrackReader::kESD) { + fhBadClusterMaxCellDiffAverageTime ->Fill(clus->E(),tmax-averTime); + fhBadClusterMaxCellDiffWeightTime ->Fill(clus->E(),tmax-weightedTime); + fhBadClusterDiffWeightAverTime ->Fill(clus->E(),weightedTime-averTime); + fhBadClusterMaxCellDiffAverageNoMaxTime->Fill(clus->E(),tmax-averTimeNoMax); + fhBadClusterMaxCellDiffWeightNoMaxTime ->Fill(clus->E(),tmax-weightedTimeNoMax); + //printf("E %f, tmax %f, aver %f, weigh %f, averNoMax %f, weightNoMax %f\n", + // clus->E(),tmax,averTime,weightedTime,averTimeNoMax,weightedTimeNoMax); + } + + if( weight > 0 ) fhBadClusterNoMaxCellWeight ->Fill(clus->E(),weightNoMax/weight); + //if( weight > 0 && weightNoMax > 0.0)printf("weight %f, weightNoMax %f\n",weight,weightNoMax); + + for (Int_t ipos = 0; ipos < nCaloCellsPerCluster; ipos++) { + Int_t absId = indexList[ipos]; + if(absId!=absIdMax){ + Float_t frac = cell->GetCellAmplitude(absId)/cell->GetCellAmplitude(absIdMax); + fhBadClusterMaxCellCloseCellRatio->Fill(mom.E(),frac); + fhBadClusterMaxCellCloseCellDiff ->Fill(mom.E(),cell->GetCellAmplitude(absIdMax)-cell->GetCellAmplitude(absId)); + time = cell->GetCellTime(absId); + if(GetReader()->GetDataType()==AliCaloTrackReader::kESD) { + + if(fCalorimeter=="EMCAL" && GetCaloUtils()->GetEMCALRecoUtils()->IsTimeRecalibrationOn()){ + GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTime(absId,bc,time); } - }// Not max - }//loop - }//Bad cluster - } + Float_t diff = (tmax-time*1e9); + fhBadCellTimeSpreadRespectToCellMax->Fill(mom.E(), diff); + + } // ESD + }// Not max + }//loop + }//Bad cluster + } + + if(!badCluster){ + + fhClusterMaxCellDiff->Fill(clus->E(),maxCellFraction); + fhClusterTimeEnergy ->Fill(mom.E(),tof); - if(!badCluster){ - - fhClusterMaxCellDiff->Fill(clus->E(),maxCellFraction); - fhClusterTimeEnergy ->Fill(mom.E(),tof); + //Clusters in event time difference + for(Int_t iclus2 = 0; iclus2 < nCaloClusters; iclus2++ ){ - //Clusters in event time difference - for(Int_t iclus2 = 0; iclus2 < nCaloClusters; iclus2++ ){ + AliVCluster* clus2 = (AliVCluster*) caloClusters->At(iclus2); + + if(clus->GetID()==clus2->GetID()) continue; + + if(clus->GetM02() > 0.01 && clus2->GetM02() > 0.01) { - AliVCluster* clus2 = (AliVCluster*) caloClusters->At(iclus2); + Double_t tof2 = clus2->GetTOF(); + // approximation in case of time recalibration + if(fCalorimeter=="EMCAL" && GetCaloUtils()->GetEMCALRecoUtils()->IsTimeRecalibrationOn()){ + Float_t maxCellFraction2 = -1; + Int_t absIdMax2 = GetCaloUtils()->GetMaxEnergyCell(cell, clus2,maxCellFraction2); + tof2 = cell->GetCellTime(absIdMax2); + GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTime(absIdMax2,bc,tof2); + } - if(clus->GetID()==clus2->GetID()) continue; + fhClusterPairDiffTimeE ->Fill(clus->E(), tof-tof2*1.e9); + } + } + + if(nCaloCellsPerCluster > 1){ + + // check time of cells respect to max energy cell + + if(GetReader()->GetDataType()==AliCaloTrackReader::kESD) { + fhClusterMaxCellDiffAverageTime ->Fill(clus->E(),tmax-averTime); + fhClusterMaxCellDiffWeightTime ->Fill(clus->E(),tmax-weightedTime); + fhClusterDiffWeightAverTime ->Fill(clus->E(), weightedTime-averTime); + fhClusterMaxCellDiffAverageNoMaxTime->Fill(clus->E(),tmax-averTimeNoMax); + fhClusterMaxCellDiffWeightNoMaxTime ->Fill(clus->E(),tmax-weightedTimeNoMax); - if(clus->GetM02() > 0.01) { - fhClusterPairDiffTimeE ->Fill(clus->E(), tof-clus2->GetTOF()*1.e9); - } - } + } - if(nCaloCellsPerCluster > 1){ + if( weight > 0 )fhClusterNoMaxCellWeight ->Fill(clus->E(),weightNoMax/weight); + + for (Int_t ipos = 0; ipos < nCaloCellsPerCluster; ipos++) { - // check time of cells respect to max energy cell + Int_t absId = indexList[ipos]; + if(absId == absIdMax) continue; + + Float_t frac = cell->GetCellAmplitude(absId)/cell->GetCellAmplitude(absIdMax); + fhClusterMaxCellCloseCellRatio->Fill(mom.E(),frac); + fhClusterMaxCellCloseCellDiff ->Fill(mom.E(),cell->GetCellAmplitude(absIdMax)-cell->GetCellAmplitude(absId)); if(GetReader()->GetDataType()==AliCaloTrackReader::kESD) { - fhClusterMaxCellDiffAverageTime->Fill(clus->E(),tmax-averTime); - fhClusterMaxCellDiffWeightTime ->Fill(clus->E(),tmax-weightedTime); - fhClusterDiffWeightAverTime ->Fill(clus->E(), weightedTime-averTime); - fhClusterMaxCellDiffAverageNoMaxTime->Fill(clus->E(),tmax-averTimeNoMax); - fhClusterMaxCellDiffWeightNoMaxTime ->Fill(clus->E(),tmax-weightedTimeNoMax); - + + time = cell->GetCellTime(absId); + if(fCalorimeter=="EMCAL" && GetCaloUtils()->GetEMCALRecoUtils()->IsTimeRecalibrationOn()){ + GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTime(absId,bc,time); + } + + Float_t diff = (tmax-time*1.0e9); + fhCellTimeSpreadRespectToCellMax->Fill(mom.E(), diff); + if(TMath::Abs(TMath::Abs(diff) > 100) && mom.E() > 1 ) fhCellIdCellLargeTimeSpread->Fill(absId); } - if( weight > 0 )fhClusterNoMaxCellWeight ->Fill(clus->E(),weightNoMax/weight); - - for (Int_t ipos = 0; ipos < nCaloCellsPerCluster; ipos++) { + Int_t sm =0; Int_t ietaa=-1; Int_t iphii = 0; Int_t rcu = 0; + sm = GetModuleNumberCellIndexes(absId,fCalorimeter, ietaa, iphii, rcu); + if(dIphi < TMath::Abs(iphii-iphiiMax)) dIphi = TMath::Abs(iphii-iphiiMax); + if(smMax==sm){ + if(dIeta < TMath::Abs(ietaa-ietaaMax)) dIeta = TMath::Abs(ietaa-ietaaMax); + } + else { + Int_t ietaaShift = ietaa; + Int_t ietaaMaxShift = ietaaMax; + if (ietaa > ietaaMax) ietaaMaxShift+=48; + else ietaaShift +=48; + if(dIeta < TMath::Abs(ietaaShift-ietaaMaxShift)) dIeta = TMath::Abs(ietaaShift-ietaaMaxShift); + } + + //if(TMath::Abs(clus->GetM20()) < 0.0001 && clus->GetNCells() > 3){ + // printf("Good : E %f, mcells %d, l0 %f, l1 %f, d %f, cell max t %f, cluster TOF %f, sm %d, icol %d, irow %d; Max icol %d, irow %d \n", + // clus->E(), clus->GetNCells(),clus->GetM02(), clus->GetM20(), clus->GetDispersion(),tmax, tof,sm,ietaa,iphii, ietaaMax, iphiiMax); + //} + + + }// fill cell-cluster histogram loop + + if(nCaloCellsPerCluster > 3){ + + if (mom.E() < 2 ) fhDeltaIEtaDeltaIPhiE0[matched]->Fill(dIeta,dIphi); + else if(mom.E() < 6 ) fhDeltaIEtaDeltaIPhiE2[matched]->Fill(dIeta,dIphi); + else fhDeltaIEtaDeltaIPhiE6[matched]->Fill(dIeta,dIphi); + + Float_t dIA = 1.*(dIphi-dIeta)/(dIeta+dIphi); + fhDeltaIA[matched]->Fill(mom.E(),dIA); + + if(mom.E() > 0.5){ - Int_t absId = indexList[ipos]; - if(absId == absIdMax) continue; + fhDeltaIAL0[matched]->Fill(clus->GetM02(),dIA); + fhDeltaIAL1[matched]->Fill(clus->GetM20(),dIA); + fhDeltaIANCells[matched]->Fill(nCaloCellsPerCluster,dIA); - Float_t frac = cell->GetCellAmplitude(absId)/cell->GetCellAmplitude(absIdMax); - fhClusterMaxCellCloseCellRatio->Fill(mom.E(),frac); - fhClusterMaxCellCloseCellDiff ->Fill(mom.E(),cell->GetCellAmplitude(absIdMax)-cell->GetCellAmplitude(absId)); - - if(GetReader()->GetDataType()==AliCaloTrackReader::kESD) { - - Float_t diff = (tmax-cell->GetCellTime(absId)*1e9); - fhCellTimeSpreadRespectToCellMax->Fill(mom.E(), diff); - if(TMath::Abs(TMath::Abs(diff) > 100) && mom.E() > 1 ) fhCellIdCellLargeTimeSpread->Fill(absId); + } + + if(IsDataMC()){ + Int_t tag = GetMCAnalysisUtils()->CheckOrigin(labels,nLabel, GetReader(),0); + if( (GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPhoton) || + GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPi0) || + GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCEta) ) && + !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion) ){ + fhDeltaIAMC[0]->Fill(mom.E(),dIA);//Pure Photon } - - Int_t sm =0; Int_t ietaa=-1; Int_t iphii = 0; Int_t rcu = 0; - sm = GetModuleNumberCellIndexes(absId,fCalorimeter, ietaa, iphii, rcu); - if(dIphi < TMath::Abs(iphii-iphiiMax)) dIphi = TMath::Abs(iphii-iphiiMax); - if(smMax==sm){ - if(dIeta < TMath::Abs(ietaa-ietaaMax)) dIeta = TMath::Abs(ietaa-ietaaMax); + else if ( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCElectron) && + !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion) ){ + fhDeltaIAMC[1]->Fill(mom.E(),dIA);//Pure electron } - else { - Int_t ietaaShift = ietaa; - Int_t ietaaMaxShift = ietaaMax; - if (ietaa > ietaaMax) ietaaMaxShift+=48; - else ietaaShift +=48; - if(dIeta < TMath::Abs(ietaaShift-ietaaMaxShift)) dIeta = TMath::Abs(ietaaShift-ietaaMaxShift); + else if ( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion) ){ + fhDeltaIAMC[2]->Fill(mom.E(),dIA);//Converted cluster } - - //if(TMath::Abs(clus->GetM20()) < 0.0001 && clus->GetNCells() > 3){ - // printf("Good : E %f, mcells %d, l0 %f, l1 %f, d %f, cell max t %f, cluster TOF %f, sm %d, icol %d, irow %d; Max icol %d, irow %d \n", - // clus->E(), clus->GetNCells(),clus->GetM02(), clus->GetM20(), clus->GetDispersion(),tmax, tof,sm,ietaa,iphii, ietaaMax, iphiiMax); - //} - - - }// fill cell-cluster histogram loop - - if(nCaloCellsPerCluster > 3){ - - if (mom.E() < 2 ) fhDeltaIEtaDeltaIPhiE0[matched]->Fill(dIeta,dIphi); - else if(mom.E() < 6 ) fhDeltaIEtaDeltaIPhiE2[matched]->Fill(dIeta,dIphi); - else fhDeltaIEtaDeltaIPhiE6[matched]->Fill(dIeta,dIphi); - - Float_t dIA = 1.*(dIphi-dIeta)/(dIeta+dIphi); - fhDeltaIA[matched]->Fill(mom.E(),dIA); - - if(mom.E() > 0.5){ - - fhDeltaIAL0[matched]->Fill(clus->GetM02(),dIA); - fhDeltaIAL1[matched]->Fill(clus->GetM20(),dIA); - fhDeltaIANCells[matched]->Fill(nCaloCellsPerCluster,dIA); - + else{ + fhDeltaIAMC[3]->Fill(mom.E(),dIA);//Hadrons } - if(IsDataMC()){ - Int_t tag = GetMCAnalysisUtils()->CheckOrigin(labels,nLabel, GetReader(),0); - if( (GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPhoton) || - GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPi0) || - GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCEta) ) && - !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion) ){ - fhDeltaIAMC[0]->Fill(mom.E(),dIA);//Pure Photon - } - else if ( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCElectron) && - !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion) ){ - fhDeltaIAMC[1]->Fill(mom.E(),dIA);//Pure electron - } - else if ( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion) ){ - fhDeltaIAMC[2]->Fill(mom.E(),dIA);//Converted cluster - } - else{ - fhDeltaIAMC[3]->Fill(mom.E(),dIA);//Hadrons - } - - } - }// 4 cells at least in cluster for size study - - }//check time and energy of cells respect to max energy cell if cluster of more than 1 cell - + } + }// 4 cells at least in cluster for size study - //Get module of cluster - - nCaloClustersAccepted++; - if(nModule >=0 && nModule < fNModules) { - if (fCalorimeter=="EMCAL" && mom.E() > 2*fEMCALCellAmpMin) nClustersInModule[nModule]++; - else if(fCalorimeter=="PHOS" && mom.E() > 2*fPHOSCellAmpMin ) nClustersInModule[nModule]++; + }//check time and energy of cells respect to max energy cell if cluster of more than 1 cell + + + //Get module of cluster + + nCaloClustersAccepted++; + if(nModule >=0 && nModule < fNModules) { + if (fCalorimeter=="EMCAL" && mom.E() > 2*fEMCALCellAmpMin) nClustersInModule[nModule]++; + else if(fCalorimeter=="PHOS" && mom.E() > 2*fPHOSCellAmpMin ) nClustersInModule[nModule]++; + } + + //----------------------------------------------------------- + //Fill histograms related to single cluster or track matching + //----------------------------------------------------------- + + // Get Track matched + if(matched){ + + if(!strcmp("AliESDCaloCluster",Form("%s",clus->ClassName()))) + { + track = dynamic_cast (GetReader()->GetInputEvent()->GetTrack(clus->GetTrackMatchedIndex())); } + else //AOD + { + track = dynamic_cast (clus->GetTrackMatched(0)); + } + } + + ClusterHistograms(mom, pos, nCaloCellsPerCluster, nModule, matched,track, labels, nLabel); + + + //----------------------------------------------------------- + //Invariant mass + //----------------------------------------------------------- + if(fFillAllPi0Histo){ + if(GetDebug()>1) printf("Invariant mass \n"); - //----------------------------------------------------------- - //Fill histograms related to single cluster or track matching - //----------------------------------------------------------- - - ClusterHistograms(mom, pos, nCaloCellsPerCluster, nModule, matched,0x0/*track*/, labels, nLabel); - - - //----------------------------------------------------------- - //Invariant mass - //----------------------------------------------------------- - if(fFillAllPi0Histo){ - if(GetDebug()>1) printf("Invariant mass \n"); - - //do not do for bad vertex - // Float_t fZvtxCut = 40. ; - if(v[2]<-GetZvertexCut() || v[2]> GetZvertexCut()) continue ; //Event can not be used (vertex, centrality,... cuts not fulfilled) - - Int_t nModule2 = -1; - if (nCaloClusters > 1 && nCaloCellsPerCluster > 1) { - for(Int_t jclus = iclus + 1 ; jclus < nCaloClusters ; jclus++) { - AliVCluster* clus2 = (AliVCluster*)caloClusters->At(jclus); - if( clus2->GetNCells() > 1 ){ - - //Get cluster kinematics - clus2->GetMomentum(mom2,v); - - //Check only certain regions - Bool_t in2 = kTRUE; - if(IsFiducialCutOn()) in2 = GetFiducialCut()->IsInFiducialCut(mom2,fCalorimeter) ; - if(!in2) continue; - - //Get module of cluster - nModule2 = GetModuleNumber(clus2); - - //Fill invariant mass histograms - - //All modules - fhIM ->Fill((mom+mom2).Pt(),(mom+mom2).M()); - - //Single module - if(nModule == nModule2 && nModule >=0 && nModule < fNModules) - fhIMMod[nModule]->Fill((mom+mom2).Pt(),(mom+mom2).M()); - - //Asymetry histograms - fhAsym->Fill((mom+mom2).Pt(),TMath::Abs((mom.E()-mom2.E())/(mom.E()+mom2.E()))); - } - }// 2nd cluster loop - } // At least one cell in cluster and one cluster in the array - }//Fill Pi0 + //do not do for bad vertex + // Float_t fZvtxCut = 40. ; + if(v[2]<-GetZvertexCut() || v[2]> GetZvertexCut()) continue ; //Event can not be used (vertex, centrality,... cuts not fulfilled) - }//good cluster + Int_t nModule2 = -1; + if (nCaloClusters > 1 && nCaloCellsPerCluster > 1) { + for(Int_t jclus = iclus + 1 ; jclus < nCaloClusters ; jclus++) { + AliVCluster* clus2 = (AliVCluster*)caloClusters->At(jclus); + if( clus2->GetNCells() > 1 ){ + + //Get cluster kinematics + clus2->GetMomentum(mom2,v); + + //Check only certain regions + Bool_t in2 = kTRUE; + if(IsFiducialCutOn()) in2 = GetFiducialCut()->IsInFiducialCut(mom2,fCalorimeter) ; + if(!in2) continue; + + //Get module of cluster + nModule2 = GetModuleNumber(clus2); + + //Fill invariant mass histograms + + //All modules + fhIM ->Fill((mom+mom2).Pt(),(mom+mom2).M()); + + //Single module + if(nModule == nModule2 && nModule >=0 && nModule < fNModules) + fhIMMod[nModule]->Fill((mom+mom2).Pt(),(mom+mom2).M()); + + //Asymetry histograms + fhAsym->Fill((mom+mom2).Pt(),TMath::Abs((mom.E()-mom2.E())/(mom.E()+mom2.E()))); + } + }// 2nd cluster loop + } // At least one cell in cluster and one cluster in the array + }//Fill Pi0 - }//cluster loop + }//good cluster - //Number of clusters histograms - if(nCaloClustersAccepted > 0) fhNClusters->Fill(nCaloClustersAccepted); - // Number of clusters per module - for(Int_t imod = 0; imod < fNModules; imod++ ){ - if(GetDebug() > 1) - printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - module %d calo %s clusters %d\n", imod, fCalorimeter.Data(), nClustersInModule[imod]); - fhNClustersMod->Fill(nClustersInModule[imod],imod); - } - delete [] nClustersInModule; - //delete caloClusters; - }// calo clusters array exists + }//cluster loop + + //Number of clusters histograms + if(nCaloClustersAccepted > 0) fhNClusters->Fill(nCaloClustersAccepted); + // Number of clusters per module + for(Int_t imod = 0; imod < fNModules; imod++ ){ + if(GetDebug() > 1) + printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - module %d calo %s clusters %d\n", imod, fCalorimeter.Data(), nClustersInModule[imod]); + fhNClustersMod->Fill(nClustersInModule[imod],imod); + } + delete [] nClustersInModule; + //delete caloClusters; //---------------------------------------------------------- // CALOCELLS //---------------------------------------------------------- - - AliVCaloCells * cell = 0x0; - Int_t ncells = 0; - if(fCalorimeter == "PHOS") - cell = GetPHOSCells(); - else - cell = GetEMCALCells(); - - if(!cell){ - AliFatal(Form("No %s CELLS available for analysis",fCalorimeter.Data())); - return; // just to trick coverity - } + + Int_t ncells = cell->GetNumberOfCells(); if(GetDebug() > 0) - printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - %s cell entries %d\n", fCalorimeter.Data(), cell->GetNumberOfCells()); + printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - %s cell entries %d\n", fCalorimeter.Data(), ncells ); //Init arrays and used variables Int_t *nCellsInModule = new Int_t[fNModules]; for(Int_t imod = 0; imod < fNModules; imod++ ) nCellsInModule[imod] = 0; - Int_t icol = -1; - Int_t irow = -1; - Int_t iRCU = -1; - Float_t amp = 0.; - Float_t time = 0.; - Int_t id = -1; - Float_t recalF = 1.; - + Int_t icol = -1; + Int_t irow = -1; + Int_t iRCU = -1; + Float_t amp = 0.; + Double_t time = 0.; + Int_t id = -1; + Float_t recalF = 1.; + for (Int_t iCell = 0; iCell < cell->GetNumberOfCells(); iCell++) { if(GetDebug() > 2) printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - Cell : amp %f, absId %d \n", cell->GetAmplitude(iCell), cell->GetCellNumber(iCell)); @@ -1733,27 +1801,40 @@ void AliAnaCalorimeterQA::MakeAnalysisFillHistograms() } } // use bad channel map + amp = cell->GetAmplitude(iCell)*recalF; + time = cell->GetTime(iCell); + id = cell->GetCellNumber(iCell); + //Get Recalibration factor if set if (GetCaloUtils()->IsRecalibrationOn()) { - if(fCalorimeter == "PHOS") recalF = GetCaloUtils()->GetPHOSChannelRecalibrationFactor(nModule,icol,irow); - else recalF = GetCaloUtils()->GetEMCALChannelRecalibrationFactor(nModule,icol,irow); + if(fCalorimeter == "PHOS") { + amp *= GetCaloUtils()->GetPHOSChannelRecalibrationFactor(nModule,icol,irow); + } + else { + amp *= GetCaloUtils()->GetEMCALChannelRecalibrationFactor(nModule,icol,irow); + GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTime(id,bc,time); +// printf("QA : Id %d, Time org %e, Time new %e; Amp org %f, Amp new %f\n", +// id, cell->GetTime(iCell),time, cell->GetAmplitude(iCell),amp); + } //if(fCalorimeter == "PHOS")printf("Recalibration factor (sm,row,col)=(%d,%d,%d) - %f\n",nModule,icol,irow,recalF); } - amp = cell->GetAmplitude(iCell)*recalF; - time = cell->GetTime(iCell)*1e9;//transform time to ns + //Transform time to ns + time *= 1.0e9; //Remove noisy channels, only possible in ESDs if(GetReader()->GetDataType() == AliCaloTrackReader::kESD){ - if(time < fTimeCutMin || time > fTimeCutMax) continue; + if(time < fTimeCutMin || time > fTimeCutMax){ + printf("Remove cluster with TOF %f\n",time); + continue; + } } //if(amp > 3 && fCalorimeter=="EMCAL") printf("Amp = %f, time = %f, (mod, col, row)= (%d,%d,%d)\n", // amp,time,nModule,icol,irow); - id = cell->GetCellNumber(iCell); fhAmplitude->Fill(amp); fhAmpId ->Fill(amp,id); - + if ((fCalorimeter=="EMCAL" && amp > fEMCALCellAmpMin) || (fCalorimeter=="PHOS" && amp > fPHOSCellAmpMin)) { @@ -1768,7 +1849,7 @@ void AliAnaCalorimeterQA::MakeAnalysisFillHistograms() else { irows = irow + fNMaxRows * fNModules; } - + fhGridCellsMod ->Fill(icols,irows); fhGridCellsEMod->Fill(icols,irows,amp); @@ -1778,9 +1859,9 @@ void AliAnaCalorimeterQA::MakeAnalysisFillHistograms() fhTimeId ->Fill(time,id); fhTimeAmp ->Fill(amp,time); fhGridCellsTimeMod->Fill(icols,irows,time); - + fhTimeAmpPerRCU [nModule*fNRCU+iRCU]->Fill(amp, time); - + } } @@ -1845,7 +1926,7 @@ void AliAnaCalorimeterQA::MakeAnalysisFillHistograms() //_____________________________________________________________________________________________ void AliAnaCalorimeterQA::ClusterHistograms(const TLorentzVector mom, Float_t *pos, const Int_t nCaloCellsPerCluster,const Int_t nModule, - const Bool_t matched, const AliVTrack * /*track*/, + const Bool_t matched, const AliVTrack * track, const Int_t * labels, const Int_t nLabels){ //Fill CaloCluster related histograms @@ -2158,118 +2239,118 @@ void AliAnaCalorimeterQA::ClusterHistograms(const TLorentzVector mom, } //printf("track index %d ntracks %d\n", esd->GetNumberOfTracks()); -// //Study the track and matched cluster if track exists. -// if(!track) return; -// Double_t emcpos[3] = {0.,0.,0.}; -// Double_t emcmom[3] = {0.,0.,0.}; -// Double_t radius = 441.0; //[cm] EMCAL radius +13cm -// Double_t bfield = 0.; -// Double_t tphi = 0; -// Double_t teta = 0; -// Double_t tmom = 0; -// Double_t tpt = 0; -// Double_t tmom2 = 0; -// Double_t tpcSignal = 0; -// Bool_t okpos = kFALSE; -// Bool_t okmom = kFALSE; -// Bool_t okout = kFALSE; -// Int_t nITS = 0; -// Int_t nTPC = 0; -// -// //In case of ESDs get the parameters in this way -// if(GetReader()->GetDataType()==AliCaloTrackReader::kESD) { -// if (track->GetOuterParam() ) { -// okout = kTRUE; -// -// bfield = GetReader()->GetInputEvent()->GetMagneticField(); -// okpos = track->GetOuterParam()->GetXYZAt(radius,bfield,emcpos); -// okmom = track->GetOuterParam()->GetPxPyPzAt(radius,bfield,emcmom); -// if(!(okpos && okmom)) return; -// -// TVector3 position(emcpos[0],emcpos[1],emcpos[2]); -// TVector3 momentum(emcmom[0],emcmom[1],emcmom[2]); -// tphi = position.Phi(); -// teta = position.Eta(); -// tmom = momentum.Mag(); -// -// tpt = track->Pt(); -// tmom2 = track->P(); -// tpcSignal = track->GetTPCsignal(); -// -// nITS = track->GetNcls(0); -// nTPC = track->GetNcls(1); -// }//Outer param available -// }// ESDs -// else if(GetReader()->GetDataType()==AliCaloTrackReader::kAOD) { -// AliAODPid* pid = (AliAODPid*) ((AliAODTrack *) track)->GetDetPid(); -// if (pid) { -// okout = kTRUE; -// pid->GetEMCALPosition(emcpos); -// pid->GetEMCALMomentum(emcmom); -// -// TVector3 position(emcpos[0],emcpos[1],emcpos[2]); -// TVector3 momentum(emcmom[0],emcmom[1],emcmom[2]); -// tphi = position.Phi(); -// teta = position.Eta(); -// tmom = momentum.Mag(); -// -// tpt = track->Pt(); -// tmom2 = track->P(); -// tpcSignal = pid->GetTPCsignal(); -// -// }//pid -// }//AODs -// -// if(okout){ -// printf("okout\n"); -// Double_t deta = teta - eta; -// Double_t dphi = tphi - phi; -// if(dphi > TMath::Pi()) dphi -= 2*TMath::Pi(); -// if(dphi < -TMath::Pi()) dphi += 2*TMath::Pi(); -// Double_t dR = sqrt(dphi*dphi + deta*deta); -// -// Double_t pOverE = tmom/e; -// -// fh1pOverE->Fill(tpt, pOverE); -// if(dR < 0.02) fh1pOverER02->Fill(tpt,pOverE); -// -// fh1dR->Fill(dR); -// fh2MatchdEdx->Fill(tmom2,tpcSignal); -// -// if(IsDataMC() && primary){ -// Int_t pdg = primary->GetPdgCode(); -// Double_t charge = TDatabasePDG::Instance()->GetParticle(pdg)->Charge(); -// -// if(TMath::Abs(pdg) == 11){ -// fhMCEle1pOverE->Fill(tpt,pOverE); -// fhMCEle1dR->Fill(dR); -// fhMCEle2MatchdEdx->Fill(tmom2,tpcSignal); -// if(dR < 0.02) fhMCEle1pOverER02->Fill(tpt,pOverE); -// } -// else if(charge!=0){ -// fhMCChHad1pOverE->Fill(tpt,pOverE); -// fhMCChHad1dR->Fill(dR); -// fhMCChHad2MatchdEdx->Fill(tmom2,tpcSignal); -// if(dR < 0.02) fhMCChHad1pOverER02->Fill(tpt,pOverE); -// } -// else if(charge == 0){ -// fhMCNeutral1pOverE->Fill(tpt,pOverE); -// fhMCNeutral1dR->Fill(dR); -// fhMCNeutral2MatchdEdx->Fill(tmom2,tpcSignal); -// if(dR < 0.02) fhMCNeutral1pOverER02->Fill(tpt,pOverE); -// } -// }//DataMC -// -// if(dR < 0.02 && pOverE > 0.5 && pOverE < 1.5 -// && nCaloCellsPerCluster > 1 && nITS > 3 && nTPC > 20) { -// fh2EledEdx->Fill(tmom2,tpcSignal); -// } -// } -// else{//no ESD external param or AODPid -// -// if(GetDebug() >= 0) printf("No ESD external param or AliAODPid \n"); -// -// }//No out params + //Study the track and matched cluster if track exists. + if(!track) return; + Double_t emcpos[3] = {0.,0.,0.}; + Double_t emcmom[3] = {0.,0.,0.}; + Double_t radius = 441.0; //[cm] EMCAL radius +13cm + Double_t bfield = 0.; + Double_t tphi = 0; + Double_t teta = 0; + Double_t tmom = 0; + Double_t tpt = 0; + Double_t tmom2 = 0; + Double_t tpcSignal = 0; + Bool_t okpos = kFALSE; + Bool_t okmom = kFALSE; + Bool_t okout = kFALSE; + Int_t nITS = 0; + Int_t nTPC = 0; + + //In case of ESDs get the parameters in this way + if(GetReader()->GetDataType()==AliCaloTrackReader::kESD) { + if (track->GetOuterParam() ) { + okout = kTRUE; + + bfield = GetReader()->GetInputEvent()->GetMagneticField(); + okpos = track->GetOuterParam()->GetXYZAt(radius,bfield,emcpos); + okmom = track->GetOuterParam()->GetPxPyPzAt(radius,bfield,emcmom); + if(!(okpos && okmom)) return; + + TVector3 position(emcpos[0],emcpos[1],emcpos[2]); + TVector3 momentum(emcmom[0],emcmom[1],emcmom[2]); + tphi = position.Phi(); + teta = position.Eta(); + tmom = momentum.Mag(); + + tpt = track->Pt(); + tmom2 = track->P(); + tpcSignal = track->GetTPCsignal(); + + nITS = track->GetNcls(0); + nTPC = track->GetNcls(1); + }//Outer param available + }// ESDs + else if(GetReader()->GetDataType()==AliCaloTrackReader::kAOD) { + AliAODPid* pid = (AliAODPid*) ((AliAODTrack *) track)->GetDetPid(); + if (pid) { + okout = kTRUE; + pid->GetEMCALPosition(emcpos); + pid->GetEMCALMomentum(emcmom); + + TVector3 position(emcpos[0],emcpos[1],emcpos[2]); + TVector3 momentum(emcmom[0],emcmom[1],emcmom[2]); + tphi = position.Phi(); + teta = position.Eta(); + tmom = momentum.Mag(); + + tpt = track->Pt(); + tmom2 = track->P(); + tpcSignal = pid->GetTPCsignal(); + + }//pid + }//AODs + + if(okout){ + //printf("okout\n"); + Double_t deta = teta - eta; + Double_t dphi = tphi - phi; + if(dphi > TMath::Pi()) dphi -= 2*TMath::Pi(); + if(dphi < -TMath::Pi()) dphi += 2*TMath::Pi(); + Double_t dR = sqrt(dphi*dphi + deta*deta); + + Double_t pOverE = tmom/e; + + fh1pOverE->Fill(tpt, pOverE); + if(dR < 0.02) fh1pOverER02->Fill(tpt,pOverE); + + fh1dR->Fill(dR); + fh2MatchdEdx->Fill(tmom2,tpcSignal); + + if(IsDataMC() && primary){ + Int_t pdg = primary->GetPdgCode(); + Double_t charge = TDatabasePDG::Instance()->GetParticle(pdg)->Charge(); + + if(TMath::Abs(pdg) == 11){ + fhMCEle1pOverE->Fill(tpt,pOverE); + fhMCEle1dR->Fill(dR); + fhMCEle2MatchdEdx->Fill(tmom2,tpcSignal); + if(dR < 0.02) fhMCEle1pOverER02->Fill(tpt,pOverE); + } + else if(charge!=0){ + fhMCChHad1pOverE->Fill(tpt,pOverE); + fhMCChHad1dR->Fill(dR); + fhMCChHad2MatchdEdx->Fill(tmom2,tpcSignal); + if(dR < 0.02) fhMCChHad1pOverER02->Fill(tpt,pOverE); + } + else if(charge == 0){ + fhMCNeutral1pOverE->Fill(tpt,pOverE); + fhMCNeutral1dR->Fill(dR); + fhMCNeutral2MatchdEdx->Fill(tmom2,tpcSignal); + if(dR < 0.02) fhMCNeutral1pOverER02->Fill(tpt,pOverE); + } + }//DataMC + + if(dR < 0.02 && pOverE > 0.5 && pOverE < 1.5 + && nCaloCellsPerCluster > 1 && nITS > 3 && nTPC > 20) { + fh2EledEdx->Fill(tmom2,tpcSignal); + } + } + else{//no ESD external param or AODPid + + if(GetDebug() >= 0) printf("No ESD external param or AliAODPid \n"); + + }//No out params }//matched clusters with tracks }// Clusters diff --git a/PWG4/PartCorrDep/AliAnaCalorimeterQA.h b/PWG4/PartCorrDep/AliAnaCalorimeterQA.h index 612a790d077..74f0132adcc 100755 --- a/PWG4/PartCorrDep/AliAnaCalorimeterQA.h +++ b/PWG4/PartCorrDep/AliAnaCalorimeterQA.h @@ -291,27 +291,27 @@ public: TH2F * fhHaR ; //! Hadron distance to vertex vs rec energy //Histograms for MC track-matching -// TH2F * fh1pOverE; //! p/E for track-cluster matches -// TH1F * fh1dR; //! distance between projected track and cluster -// TH2F * fh2EledEdx; //! dE/dx vs. momentum for electron candidates -// TH2F * fh2MatchdEdx; //! dE/dx vs. momentum for all matches -// -// TH2F * fhMCEle1pOverE; //! p/E for track-cluster matches, MC electrons -// TH1F * fhMCEle1dR; //! distance between projected track and cluster, MC electrons -// TH2F * fhMCEle2MatchdEdx; //! dE/dx vs. momentum for all matches, MC electrons -// -// TH2F * fhMCChHad1pOverE; //! p/E for track-cluster matches, MC charged hadrons -// TH1F * fhMCChHad1dR; //! distance between projected track and cluster, MC charged hadrons -// TH2F * fhMCChHad2MatchdEdx; //! dE/dx vs. momentum for all matches, MC charged -// -// TH2F * fhMCNeutral1pOverE; //! p/E for track-cluster matches, MC neutral -// TH1F * fhMCNeutral1dR; //! distance between projected track and cluster, MC neutral -// TH2F * fhMCNeutral2MatchdEdx; //! dE/dx vs. momentum for all matches, MC neutral -// -// TH2F * fh1pOverER02; //! p/E for track-cluster matches, dR > 0.2 -// TH2F * fhMCEle1pOverER02; //! p/E for track-cluster matches, dR > 0.2, MC electrons -// TH2F * fhMCChHad1pOverER02; //! p/E for track-cluster matches, dR > 0.2, MC charged hadrons -// TH2F * fhMCNeutral1pOverER02; //! p/E for track-cluster matches, dR > 0.2, MC neutral + TH2F * fh1pOverE; //! p/E for track-cluster matches + TH1F * fh1dR; //! distance between projected track and cluster + TH2F * fh2EledEdx; //! dE/dx vs. momentum for electron candidates + TH2F * fh2MatchdEdx; //! dE/dx vs. momentum for all matches + + TH2F * fhMCEle1pOverE; //! p/E for track-cluster matches, MC electrons + TH1F * fhMCEle1dR; //! distance between projected track and cluster, MC electrons + TH2F * fhMCEle2MatchdEdx; //! dE/dx vs. momentum for all matches, MC electrons + + TH2F * fhMCChHad1pOverE; //! p/E for track-cluster matches, MC charged hadrons + TH1F * fhMCChHad1dR; //! distance between projected track and cluster, MC charged hadrons + TH2F * fhMCChHad2MatchdEdx; //! dE/dx vs. momentum for all matches, MC charged + + TH2F * fhMCNeutral1pOverE; //! p/E for track-cluster matches, MC neutral + TH1F * fhMCNeutral1dR; //! distance between projected track and cluster, MC neutral + TH2F * fhMCNeutral2MatchdEdx; //! dE/dx vs. momentum for all matches, MC neutral + + TH2F * fh1pOverER02; //! p/E for track-cluster matches, dR > 0.2 + TH2F * fhMCEle1pOverER02; //! p/E for track-cluster matches, dR > 0.2, MC electrons + TH2F * fhMCChHad1pOverER02; //! p/E for track-cluster matches, dR > 0.2, MC charged hadrons + TH2F * fhMCNeutral1pOverER02; //! p/E for track-cluster matches, dR > 0.2, MC neutral ClassDef(AliAnaCalorimeterQA,19) } ; -- 2.43.0