* provided "as is" without express or implied warranty. *
**************************************************************************/
-
///////////////////////////////////////////////////////////////////////////////
// //
-// EMCAL tender, apply corrections to EMCAl clusters //
-// and do track matching. //
+// EMCAL tender, apply corrections to EMCAL clusters and do track matching. //
+// //
// Author: Deepa Thomas (Utrecht University) //
+// Later mods/rewrite: Jiri Kral (University of Jyvaskyla) //
+// S. Aiola, C. Loizides : Make it work for AODs //
// //
///////////////////////////////////////////////////////////////////////////////
#include <TObjArray.h>
#include <TClonesArray.h>
#include <TGeoGlobalMagField.h>
-#include <AliLog.h>
-#include <AliESDEvent.h>
-#include <AliAnalysisManager.h>
-#include <AliTender.h>
-#include "AliOADBContainer.h"
-#include "AliCDBManager.h"
-#include "AliCDBStorage.h"
-#include "AliCDBEntry.h"
-#include "AliMagF.h"
-#include "AliESDCaloCluster.h"
-#include "AliEMCALTenderSupply.h"
-#include "AliEMCALGeometry.h"
-#include "AliEMCALRecoUtils.h"
-#include "AliEMCALClusterizer.h"
-#include "AliEMCALRecParam.h"
-#include "AliEMCALRecPoint.h"
#include "AliEMCALAfterBurnerUF.h"
+#include "AliEMCALClusterizer.h"
#include "AliEMCALClusterizerNxN.h"
#include "AliEMCALClusterizerv1.h"
#include "AliEMCALClusterizerv2.h"
#include "AliEMCALDigit.h"
+#include "AliEMCALGeometry.h"
+#include "AliEMCALRecParam.h"
+#include "AliEMCALRecParam.h"
+#include "AliEMCALRecPoint.h"
+#include "AliEMCALRecoUtils.h"
+#include "AliEMCALTenderSupply.h"
+#include "AliESDCaloCluster.h"
+#include "AliMagF.h"
+#include "AliOADBContainer.h"
+#include "AliAODEvent.h"
+#include "AliAnalysisManager.h"
+#include "AliESDEvent.h"
+#include "AliLog.h"
+#include "AliTender.h"
+#include "AliAODMCParticle.h"
ClassImp(AliEMCALTenderSupply)
AliEMCALTenderSupply::AliEMCALTenderSupply() :
- AliTenderSupply()
- ,fEMCALGeo(0x0)
- ,fEMCALGeoName("EMCAL_COMPLETEV1")
- ,fEMCALRecoUtils(0)
- ,fConfigName("")
- ,fDebugLevel(0)
- ,fNonLinearFunc(AliEMCALRecoUtils::kNoCorrection)
- ,fNonLinearThreshold(30)
- ,fReCalibCluster(kFALSE)
- ,fReCalibCell(kFALSE)
- ,fRecalClusPos(kFALSE)
- ,fFiducial(kFALSE)
- ,fNCellsFromEMCALBorder(1)
- ,fRecalDistToBadChannels(kFALSE)
- ,fInputTree(0)
- ,fInputFile(0)
- ,fFilepass(0)
- ,fMass(0.139)
- ,fStep(50)
- ,fCutEtaPhiSum(kTRUE)
- ,fCutEtaPhiSeparate(kFALSE)
- ,fRcut(0.05)
- ,fEtacut(0.025)
- ,fPhicut(0.05)
- ,fBasePath(".")
- ,fReClusterize(kFALSE)
- ,fClusterizer(0)
- ,fGeomMatrixSet(kFALSE)
- ,fLoadGeomMatrices(kFALSE)
- ,fRecParam(0)
- ,fOCDBpath(" ")
- ,fUnfolder(0)
- ,fDigitsArr(0)
- ,fClusterArr(0)
+AliTenderSupply()
+,fTask(0)
+,fRun(0)
+,fEMCALGeo(0x0)
+,fEMCALGeoName("")
+,fEMCALRecoUtils(0)
+,fConfigName("")
+,fDebugLevel(0)
+,fNonLinearFunc(-1)
+,fNonLinearThreshold(-1)
+,fReCalibCluster(kFALSE)
+,fUpdateCell(kFALSE)
+,fCalibrateEnergy(kFALSE)
+,fCalibrateTime(kFALSE)
+,fDoNonLinearity(kFALSE)
+,fBadCellRemove(kFALSE)
+,fRejectExoticCells(kFALSE)
+,fRejectExoticClusters(kFALSE)
+,fClusterBadChannelCheck(kFALSE)
+,fRecalClusPos(kFALSE)
+,fFiducial(kFALSE)
+,fNCellsFromEMCALBorder(-1)
+,fRecalDistToBadChannels(kFALSE)
+,fRecalShowerShape(kFALSE)
+,fInputTree(0)
+,fInputFile(0)
+,fGetPassFromFileName(kTRUE)
+,fFilepass(0)
+,fMass(-1)
+,fStep(-1)
+,fCutEtaPhiSum(kTRUE)
+,fCutEtaPhiSeparate(kFALSE)
+,fRcut(-1)
+,fEtacut(-1)
+,fPhicut(-1)
+,fBasePath("")
+,fReClusterize(kFALSE)
+,fClusterizer(0)
+,fGeomMatrixSet(kFALSE)
+,fLoadGeomMatrices(kFALSE)
+,fRecParam(0x0)
+,fDoTrackMatch(kFALSE)
+,fDoUpdateOnly(kFALSE)
+,fUnfolder(0)
+,fDigitsArr(0)
+,fClusterArr(0)
+,fMisalignSurvey(kdefault)
+,fExoticCellFraction(-1)
+,fExoticCellDiffTime(-1)
+,fExoticCellMinAmplitude(-1)
+,fSetCellMCLabelFromCluster(0)
+,fTempClusterArr(0)
+,fRemapMCLabelForAODs(0)
+
{
// Default constructor.
+
+ for(Int_t i = 0; i < 12; i++) fEMCALMatrix[i] = 0 ;
+ for(Int_t j = 0; j < 12672; j++) fOrgClusterCellId[j] = -1;
}
//_____________________________________________________
AliEMCALTenderSupply::AliEMCALTenderSupply(const char *name, const AliTender *tender) :
- AliTenderSupply(name,tender)
- ,fEMCALGeo(0x0)
- ,fEMCALGeoName("EMCAL_COMPLETEV1")
- ,fEMCALRecoUtils(0)
- ,fConfigName("")
- ,fDebugLevel(0)
- ,fNonLinearFunc(AliEMCALRecoUtils::kNoCorrection)
- ,fNonLinearThreshold(30)
- ,fReCalibCluster(kFALSE)
- ,fReCalibCell(kFALSE)
- ,fRecalClusPos(kFALSE)
- ,fFiducial(kFALSE)
- ,fNCellsFromEMCALBorder(1)
- ,fRecalDistToBadChannels(kFALSE)
- ,fInputTree(0)
- ,fInputFile(0)
- ,fFilepass(0)
- ,fMass(0.139)
- ,fStep(50)
- ,fCutEtaPhiSum(kTRUE)
- ,fCutEtaPhiSeparate(kFALSE)
- ,fRcut(0.05)
- ,fEtacut(0.025)
- ,fPhicut(0.05)
- ,fBasePath(".")
- ,fReClusterize(kFALSE)
- ,fClusterizer(0)
- ,fGeomMatrixSet(kFALSE)
- ,fLoadGeomMatrices(kFALSE)
- ,fRecParam(0)
- ,fOCDBpath(" ")
- ,fUnfolder(0)
- ,fDigitsArr(0)
- ,fClusterArr(0)
+AliTenderSupply(name,tender)
+,fTask(0)
+,fRun(0)
+,fEMCALGeo(0x0)
+,fEMCALGeoName("")
+,fEMCALRecoUtils(0)
+,fConfigName("")
+,fDebugLevel(0)
+,fNonLinearFunc(-1)
+,fNonLinearThreshold(-1)
+,fReCalibCluster(kFALSE)
+,fUpdateCell(kFALSE)
+,fCalibrateEnergy(kFALSE)
+,fCalibrateTime(kFALSE)
+,fDoNonLinearity(kFALSE)
+,fBadCellRemove(kFALSE)
+,fRejectExoticCells(kFALSE)
+,fRejectExoticClusters(kFALSE)
+,fClusterBadChannelCheck(kFALSE)
+,fRecalClusPos(kFALSE)
+,fFiducial(kFALSE)
+,fNCellsFromEMCALBorder(-1)
+,fRecalDistToBadChannels(kFALSE)
+,fRecalShowerShape(kFALSE)
+,fInputTree(0)
+,fInputFile(0)
+,fGetPassFromFileName(kTRUE)
+,fFilepass("")
+,fMass(-1)
+,fStep(-1)
+,fCutEtaPhiSum(kTRUE)
+,fCutEtaPhiSeparate(kFALSE)
+,fRcut(-1)
+,fEtacut(-1)
+,fPhicut(-1)
+,fBasePath("")
+,fReClusterize(kFALSE)
+,fClusterizer(0)
+,fGeomMatrixSet(kFALSE)
+,fLoadGeomMatrices(kFALSE)
+,fRecParam(0x0)
+,fDoTrackMatch(kFALSE)
+,fDoUpdateOnly(kFALSE)
+,fUnfolder(0)
+,fDigitsArr(0)
+,fClusterArr(0)
+,fMisalignSurvey(kdefault)
+,fExoticCellFraction(-1)
+,fExoticCellDiffTime(-1)
+,fExoticCellMinAmplitude(-1)
+,fSetCellMCLabelFromCluster(0)
+,fTempClusterArr(0)
+,fRemapMCLabelForAODs(0)
+
{
// Named constructor
+
+ for(Int_t i = 0; i < 12; i++) fEMCALMatrix[i] = 0 ;
+ for(Int_t j = 0; j < 12672; j++) fOrgClusterCellId[j] = -1;
+}
- for(Int_t i = 0; i < 10; i++) fEMCALMatrix[i] = 0 ;
-
- fRecParam = new AliEMCALRecParam;
- fEMCALRecoUtils = new AliEMCALRecoUtils();
- fDigitsArr = new TClonesArray("AliEMCALDigit",1000);
+//_____________________________________________________
+AliEMCALTenderSupply::AliEMCALTenderSupply(const char *name, AliAnalysisTaskSE *task) :
+AliTenderSupply(name)
+,fTask(task)
+,fRun(0)
+,fEMCALGeo(0x0)
+,fEMCALGeoName("")
+,fEMCALRecoUtils(0)
+,fConfigName("")
+,fDebugLevel(0)
+,fNonLinearFunc(-1)
+,fNonLinearThreshold(-1)
+,fReCalibCluster(kFALSE)
+,fUpdateCell(kFALSE)
+,fCalibrateEnergy(kFALSE)
+,fCalibrateTime(kFALSE)
+,fDoNonLinearity(kFALSE)
+,fBadCellRemove(kFALSE)
+,fRejectExoticCells(kFALSE)
+,fRejectExoticClusters(kFALSE)
+,fClusterBadChannelCheck(kFALSE)
+,fRecalClusPos(kFALSE)
+,fFiducial(kFALSE)
+,fNCellsFromEMCALBorder(-1)
+,fRecalDistToBadChannels(kFALSE)
+,fRecalShowerShape(kFALSE)
+,fInputTree(0)
+,fInputFile(0)
+,fGetPassFromFileName(kTRUE)
+,fFilepass("")
+,fMass(-1)
+,fStep(-1)
+,fCutEtaPhiSum(kTRUE)
+,fCutEtaPhiSeparate(kFALSE)
+,fRcut(-1)
+,fEtacut(-1)
+,fPhicut(-1)
+,fBasePath("")
+,fReClusterize(kFALSE)
+,fClusterizer(0)
+,fGeomMatrixSet(kFALSE)
+,fLoadGeomMatrices(kFALSE)
+,fRecParam(0x0)
+,fDoTrackMatch(kFALSE)
+,fDoUpdateOnly(kFALSE)
+,fUnfolder(0)
+,fDigitsArr(0)
+,fClusterArr(0)
+,fMisalignSurvey(kdefault)
+,fExoticCellFraction(-1)
+,fExoticCellDiffTime(-1)
+,fExoticCellMinAmplitude(-1)
+,fSetCellMCLabelFromCluster(0)
+,fTempClusterArr(0)
+,fRemapMCLabelForAODs(0)
+{
+ // Named constructor.
+
+ for(Int_t i = 0; i < 12; i++) fEMCALMatrix[i] = 0 ;
+ for(Int_t j = 0; j < 12672; j++) fOrgClusterCellId[j] = -1;
}
//_____________________________________________________
AliEMCALTenderSupply::~AliEMCALTenderSupply()
{
//Destructor
-
- delete fEMCALRecoUtils;
- delete fClusterizer;
- delete fUnfolder;
- if (fDigitsArr){
- fDigitsArr->Clear("C");
- delete fDigitsArr;
+
+ if (!AliAnalysisManager::GetAnalysisManager()->IsProofMode())
+ {
+ delete fEMCALRecoUtils;
+ delete fRecParam;
+ delete fUnfolder;
+ if (!fClusterizer)
+ {
+ fDigitsArr->Clear("C");
+ delete fDigitsArr;
+ } else
+ {
+ delete fClusterizer;
+ fDigitsArr = 0;
+ }
}
}
+//_____________________________________________________
+void AliEMCALTenderSupply::SetDefaults()
+{
+ // Set default settings.
+
+ SwitchOnBadCellRemove();
+ SwitchOnExoticCellRemove();
+ SwitchOnCalibrateEnergy();
+ SwitchOnCalibrateTime();
+ SwitchOnUpdateCell();
+ SwitchOnReclustering();
+ SwitchOnClusterBadChannelCheck();
+ SwitchOnClusterExoticChannelCheck();
+ SwitchOnTrackMatch();
+}
+
+//_____________________________________________________
+Bool_t AliEMCALTenderSupply::RunChanged() const
+{
+ // Get run number.
+
+ return (fTender && fTender->RunChanged()) || (fTask && fRun != fTask->InputEvent()->GetRunNumber());
+}
+
//_____________________________________________________
void AliEMCALTenderSupply::Init()
{
// Initialise EMCAL tender.
if (fDebugLevel>0)
- AliInfo("Init EMCAL Tender supply");
-
- if(fConfigName.Length()>0 && gROOT->LoadMacro(fConfigName) >=0) {
+ AliWarning("Init EMCAL Tender supply");
+
+ if (fConfigName.Length()>0 && gROOT->LoadMacro(fConfigName) >=0) {
AliDebug(1, Form("Loading settings from macro %s", fConfigName.Data()));
AliEMCALTenderSupply *tender = (AliEMCALTenderSupply*)gInterpreter->ProcessLine("ConfigEMCALTenderSupply()");
fDebugLevel = tender->fDebugLevel;
fEMCALGeoName = tender->fEMCALGeoName;
- delete fEMCALRecoUtils;
fEMCALRecoUtils = tender->fEMCALRecoUtils;
fConfigName = tender->fConfigName;
fNonLinearFunc = tender->fNonLinearFunc;
fNonLinearThreshold = tender->fNonLinearThreshold;
fReCalibCluster = tender->fReCalibCluster;
- fReCalibCell = tender->fReCalibCell;
+ fUpdateCell = tender->fUpdateCell;
fRecalClusPos = tender->fRecalClusPos;
- fFiducial = tender->fFiducial;
+ fCalibrateEnergy = tender->fCalibrateEnergy;
+ fCalibrateTime = tender->fCalibrateTime;
+ fFiducial = tender->fFiducial;
fNCellsFromEMCALBorder = tender->fNCellsFromEMCALBorder;
fRecalDistToBadChannels = tender->fRecalDistToBadChannels;
+ fRecalShowerShape = tender->fRecalShowerShape;
+ fClusterBadChannelCheck = tender->fClusterBadChannelCheck;
+ fBadCellRemove = tender->fBadCellRemove;
+ fRejectExoticCells = tender->fRejectExoticCells;
+ fRejectExoticClusters = tender->fRejectExoticClusters;
fMass = tender->fMass;
fStep = tender->fStep;
+ fCutEtaPhiSum = tender->fCutEtaPhiSum;
+ fCutEtaPhiSeparate = tender->fCutEtaPhiSeparate;
fRcut = tender->fRcut;
fEtacut = tender->fEtacut;
fPhicut = tender->fPhicut;
fReClusterize = tender->fReClusterize;
fLoadGeomMatrices = tender->fLoadGeomMatrices;
fRecParam = tender->fRecParam;
- fOCDBpath = tender->fOCDBpath;
- for(Int_t i = 0; i < 10; i++)
+ fDoNonLinearity = tender->fDoNonLinearity;
+ fDoTrackMatch = tender->fDoTrackMatch;
+ fDoUpdateOnly = tender->fDoUpdateOnly;
+ fMisalignSurvey = tender->fMisalignSurvey;
+ fExoticCellFraction = tender->fExoticCellFraction;
+ fExoticCellDiffTime = tender->fExoticCellDiffTime;
+ fExoticCellMinAmplitude = tender->fExoticCellMinAmplitude;
+
+ for(Int_t i = 0; i < 12; i++)
fEMCALMatrix[i] = tender->fEMCALMatrix[i] ;
}
+
+ if (fDebugLevel>0){
+ AliInfo( "Emcal Tender settings: ======================================" );
+ AliInfo( "------------ Switches --------------------------" );
+ AliInfo( Form( "BadCellRemove : %d", fBadCellRemove ));
+ AliInfo( Form( "ExoticCellRemove : %d", fRejectExoticCells ));
+ AliInfo( Form( "CalibrateEnergy : %d", fCalibrateEnergy ));
+ AliInfo( Form( "CalibrateTime : %d", fCalibrateTime ));
+ AliInfo( Form( "UpdateCell : %d", fUpdateCell ));
+ AliInfo( Form( "DoUpdateOnly : %d", fDoUpdateOnly ));
+ AliInfo( Form( "Reclustering : %d", fReClusterize ));
+ AliInfo( Form( "ClusterBadChannelCheck : %d", fClusterBadChannelCheck ));
+ AliInfo( Form( "ClusterExoticChannelCheck : %d", fRejectExoticClusters ));
+ AliInfo( Form( "CellFiducialRegion : %d", fFiducial ));
+ AliInfo( Form( "ReCalibrateCluster : %d", fReCalibCluster ));
+ AliInfo( Form( "RecalculateClusPos : %d", fRecalClusPos ));
+ AliInfo( Form( "RecalShowerShape : %d", fRecalShowerShape ));
+ AliInfo( Form( "NonLinearityCorrection : %d", fDoNonLinearity ));
+ AliInfo( Form( "RecalDistBadChannel : %d", fRecalDistToBadChannels ));
+ AliInfo( Form( "TrackMatch : %d", fDoTrackMatch ));
+ AliInfo( "------------ Variables -------------------------" );
+ AliInfo( Form( "DebugLevel : %d", fDebugLevel ));
+ AliInfo( Form( "BasePath : %s", fBasePath.Data() ));
+ AliInfo( Form( "ConfigFileName : %s", fConfigName.Data() ));
+ AliInfo( Form( "EMCALGeometryName : %s", fEMCALGeoName.Data() ));
+ AliInfo( Form( "NonLinearityFunction : %d", fNonLinearFunc ));
+ AliInfo( Form( "NonLinearityThreshold : %d", fNonLinearThreshold ));
+ AliInfo( Form( "MisalignmentMatrixSurvey : %d", fMisalignSurvey ));
+ AliInfo( Form( "NumberOfCellsFromEMCALBorder : %d", fNCellsFromEMCALBorder ));
+ AliInfo( Form( "RCut : %f", fRcut ));
+ AliInfo( Form( "Mass : %f", fMass ));
+ AliInfo( Form( "Step : %f", fStep ));
+ AliInfo( Form( "EtaCut : %f", fEtacut ));
+ AliInfo( Form( "PhiCut : %f", fPhicut ));
+ AliInfo( Form( "ExoticCellFraction : %f", fExoticCellFraction ));
+ AliInfo( Form( "ExoticCellDiffTime : %f", fExoticCellDiffTime ));
+ AliInfo( Form( "ExoticCellMinAmplitude : %f", fExoticCellMinAmplitude ));
+ AliInfo( "=============================================================" );
+ }
+
+ // init reco utils
+
+ if(!fEMCALRecoUtils)
+ fEMCALRecoUtils = new AliEMCALRecoUtils;
+
+ // init geometry if requested
+ if (fEMCALGeoName.Length()>0)
+ fEMCALGeo = AliEMCALGeometry::GetInstance(fEMCALGeoName) ;
- // Init goemetry
- fEMCALGeo = AliEMCALGeometry::GetInstance(fEMCALGeoName) ;
+ // digits array
+ fDigitsArr = new TClonesArray("AliEMCALDigit",1000);
- // Initialising Non linearity parameters
- fEMCALRecoUtils->SetNonLinearityThreshold(fNonLinearThreshold);
- fEMCALRecoUtils->SetNonLinearityFunction(fNonLinearFunc);
+ // initialising non-linearity parameters
+ if( fNonLinearThreshold != -1 )
+ fEMCALRecoUtils->SetNonLinearityThreshold(fNonLinearThreshold);
+ if( fNonLinearFunc != -1 )
+ fEMCALRecoUtils->SetNonLinearityFunction(fNonLinearFunc);
- // Setting mass, step size and residual cut
- fEMCALRecoUtils->SetMass(fMass);
- fEMCALRecoUtils->SetStep(fStep);
+ // missalignment function
+ fEMCALRecoUtils->SetPositionAlgorithm(AliEMCALRecoUtils::kPosTowerGlobal);
+
+ // fiducial cut
+ // do not do the eta0 fiducial cut
+ if( fNCellsFromEMCALBorder != -1 )
+ fEMCALRecoUtils->SetNumberOfCellsFromEMCALBorder(fNCellsFromEMCALBorder);
+ fEMCALRecoUtils->SwitchOnNoFiducialBorderInEMCALEta0();
+
+ // exotic cell rejection
+ if( fExoticCellFraction != -1 )
+ fEMCALRecoUtils->SetExoticCellFractionCut( fExoticCellFraction );
+ if( fExoticCellDiffTime != -1 )
+ fEMCALRecoUtils->SetExoticCellDiffTimeCut( fExoticCellDiffTime );
+ if( fExoticCellMinAmplitude != -1 )
+ fEMCALRecoUtils->SetExoticCellMinAmplitudeCut( fExoticCellMinAmplitude );
+
+ // setting track matching parameters ... mass, step size and residual cut
+ if( fMass != -1 )
+ fEMCALRecoUtils->SetMass(fMass);
+ if( fStep != -1 )
+ fEMCALRecoUtils->SetStep(fStep);
+
+ // spatial cut based on separate eta/phi or common processing
if(fCutEtaPhiSum){
fEMCALRecoUtils->SwitchOnCutEtaPhiSum();
- fEMCALRecoUtils->SetCutR(fRcut);
+ if( fRcut != -1 )
+ fEMCALRecoUtils->SetCutR(fRcut);
} else if (fCutEtaPhiSeparate) {
fEMCALRecoUtils->SwitchOnCutEtaPhiSeparate();
- fEMCALRecoUtils->SetCutEta(fEtacut);
- fEMCALRecoUtils->SetCutPhi(fPhicut);
+ if( fEtacut != -1 )
+ fEMCALRecoUtils->SetCutEta(fEtacut);
+ if( fPhicut != -1 )
+ fEMCALRecoUtils->SetCutPhi(fPhicut);
}
}
+//_____________________________________________________
+AliVEvent* AliEMCALTenderSupply::GetEvent()
+{
+ // Return the event pointer.
+
+ if (fTender) {
+ return fTender->GetEvent();
+ }
+ else if (fTask) {
+ return fTask->InputEvent();
+ }
+
+ return 0;
+}
+
//_____________________________________________________
void AliEMCALTenderSupply::ProcessEvent()
{
// Event loop.
+
+ AliVEvent *event = GetEvent();
- AliESDEvent *event = fTender->GetEvent();
if (!event) {
- AliError("ESD event ptr = 0, returning");
+ AliError("Event ptr = 0, returning");
return;
}
-
- // Initialising parameters once per run number
- if(fTender->RunChanged()){
- GetPass();
- if (!InitBadChannels()) {
- AliError("InitBadChannels returned false, returning");
- return;
- }
- if (fRecalClusPos || fReClusterize) {
- if (!InitMisalignMatrix()) {
- AliError("InitMisalignmentMatrix returned false, returning");
+
+ // Initialising parameters once per run number
+ if (RunChanged()) {
+
+ AliWarning( "Run changed, initializing parameters" );
+ fRun = event->GetRunNumber();
+
+ // init geometry if not already done
+ if (fEMCALGeoName.Length()==0) {
+ fEMCALGeoName = "EMCAL_FIRSTYEARV1";
+ if (fRun>139517) {
+ fEMCALGeoName = "EMCAL_COMPLETEV1";
+ }
+ if (fRun>170593) {
+ fEMCALGeoName = "EMCAL_COMPLETE12SMV1";
+ }
+ fEMCALGeo = AliEMCALGeometry::GetInstance(fEMCALGeoName);
+ if (!fEMCALGeo) {
+ AliFatal(Form("Can not create geometry: %s", fEMCALGeoName.Data()));
return;
}
+ }
+
+ // get pass
+ if (fGetPassFromFileName)
+ GetPass();
+
+ // define what recalib parameters are needed for various switches
+ // this is based on implementation in AliEMCALRecoUtils
+ Bool_t needRecoParam = fReClusterize;
+ Bool_t needBadChannels = fBadCellRemove | fClusterBadChannelCheck | fRecalDistToBadChannels | fReClusterize;
+ Bool_t needRecalib = fCalibrateEnergy | fReClusterize;
+ Bool_t needTimecalib = fCalibrateTime | fReClusterize;
+ Bool_t needMisalign = fRecalClusPos | fReClusterize;
+ Bool_t needClusterizer = fReClusterize;
+
+ // initiate reco params with some defaults
+ // will not overwrite, if those have been provided by user
+ if( needRecoParam ){
+ Int_t initRC = InitRecParam();
+
+ if( initRC == 0 )
+ AliInfo("Defaults reco params loaded.");
+ if( initRC > 1 )
+ AliWarning("User defined reco params.");
+ }
+
+ // init bad channels
+ if( needBadChannels ){
+ Int_t fInitBC = InitBadChannels();
+ if (fInitBC==0)
+ AliError("InitBadChannels returned false, returning");
+ if (fInitBC==1)
+ AliWarning("InitBadChannels OK");
+ if (fInitBC>1)
+ AliWarning(Form("No external hot channel set: %d - %s", event->GetRunNumber(), fFilepass.Data()));
}
- if (fReCalibCluster || fReCalibCell) {
- if (!InitRecalib()) {
+
+ // init recalibration factors
+ if( needRecalib )
+ {
+ Int_t fInitRecalib = InitRecalib();
+ if (fInitRecalib==0)
AliError("InitRecalib returned false, returning");
- return;
- }
+ if (fInitRecalib==1)
+ AliWarning("InitRecalib OK");
+ if (fInitRecalib>1)
+ AliWarning(Form("No recalibration available: %d - %s", event->GetRunNumber(), fFilepass.Data()));
+
+ Int_t fInitRunDepRecalib = InitRunDepRecalib();
+ if (fInitRunDepRecalib==0)
+ AliError("InitrunDepRecalib returned false, returning");
+ if (fInitRunDepRecalib==1)
+ AliWarning("InitRecalib OK");
+ if (fInitRunDepRecalib>1)
+ AliWarning(Form("No Temperature recalibration available: %d - %s", event->GetRunNumber(), fFilepass.Data()));
+
}
- if (fReClusterize) {
- if (!InitClusterization()) {
- AliError("InitClusterization returned false, returning");
- return;
- }
+
+
+
+ // init time calibration
+ if( needTimecalib ){
+ Int_t initTC = InitTimeCalibration();
+ if ( !initTC )
+ AliError("InitTimeCalibration returned false, returning");
+ if (initTC==1)
+ AliWarning("InitTimeCalib OK");
+ if( initTC > 1 )
+ AliWarning(Form("No external time calibration set: %d - %s", event->GetRunNumber(), fFilepass.Data()));
}
+ // init misalignment matrix
+ if( needMisalign ) {
+ if (!InitMisalignMatrix())
+ AliError("InitMisalignmentMatrix returned false, returning");
+ else
+ AliWarning("InitMisalignMatrix OK");
+ }
+
+ // init clusterizer
+ if( needClusterizer ) {
+ if (!InitClusterization())
+ AliError("InitClusterization returned false, returning");
+ else
+ AliWarning("InitClusterization OK");
+ }
+
if(fDebugLevel>1)
fEMCALRecoUtils->Print("");
}
+
+ // disable implied switches -------------------------------------------------
+ // AliEMCALRecoUtils or clusterizer functions alredy include some
+ // recalibration so based on those implied callibration te switches
+ // here are disabled to avoid duplication
+
+ // clusterizer does cluster energy recalibration, position recomputation
+ // and shower shape
+ if( fReClusterize ){
+ fReCalibCluster = kFALSE;
+ fRecalClusPos = kFALSE;
+ fRecalShowerShape = kFALSE;
+ }
+
+ // CONFIGURE THE RECO UTILS -------------------------------------------------
+ // configure the reco utils
+ // this option does energy recalibration
+ if( fCalibrateEnergy )
+ fEMCALRecoUtils->SwitchOnRecalibration();
+ else
+ fEMCALRecoUtils->SwitchOffRecalibration();
+
+ // allows time calibration
+ if( fCalibrateTime )
+ fEMCALRecoUtils->SwitchOnTimeRecalibration();
+ else
+ fEMCALRecoUtils->SwitchOffTimeRecalibration();
+
+ // allows to zero bad cells
+ if( fBadCellRemove )
+ fEMCALRecoUtils->SwitchOnBadChannelsRemoval();
+ else
+ fEMCALRecoUtils->SwitchOffBadChannelsRemoval();
+
+ // distance to bad channel recalibration
+ if( fRecalDistToBadChannels )
+ fEMCALRecoUtils->SwitchOnDistToBadChannelRecalculation();
+ else
+ fEMCALRecoUtils->SwitchOffDistToBadChannelRecalculation();
+
+ // exclude exotic cells
+ if( fRejectExoticCells )
+ fEMCALRecoUtils->SwitchOnRejectExoticCell();
+ else
+ fEMCALRecoUtils->SwitchOffRejectExoticCell();
+
+ // exclude clusters with exotic cells
+ if( fRejectExoticClusters )
+ fEMCALRecoUtils->SwitchOnRejectExoticCluster();
+ else
+ fEMCALRecoUtils->SwitchOffRejectExoticCluster();
+
+ // TODO: not implemented switches
+ // SwitchOnClusterEnergySmearing
+ // SwitchOnRunDepCorrection
+ // START PROCESSING ---------------------------------------------------------
// Test if cells present
- AliESDCaloCells *cells= event->GetEMCALCells();
- if (cells->GetNumberOfCells()<=0) {
- AliWarning(Form("Number of EMCAL cells = %d, returning", cells->GetNumberOfCells()));
+ AliVCaloCells *cells= event->GetEMCALCells();
+ if (cells->GetNumberOfCells()<=0)
+ {
+ if(fDebugLevel>1)
+ AliWarning(Form("Number of EMCAL cells = %d, returning", cells->GetNumberOfCells()));
return;
}
+
+ if (fDebugLevel>2)
+ AliInfo(Form("Re-calibrate cluster %d\n",fReCalibCluster));
- // Recalibrate cells
- if (fReCalibCell)
- RecalibrateCells();
+ // mark the cells not recalibrated in case of selected
+ // time, energy recalibration or bad channel removal
+ if( fCalibrateEnergy || fCalibrateTime || fBadCellRemove )
+ fEMCALRecoUtils->ResetCellsCalibrated();
+
+ // CELL RECALIBRATION -------------------------------------------------------
+ // cell objects will be updated
+ // the cell calibrations are also processed locally any time those are needed
+ // in case that the cell objects are not to be updated here for later use
+ if( fUpdateCell )
+ {
+ // do the update
+ // includes exotic cell check in the UpdateCells function - is not provided
+ // by the reco utils
+ UpdateCells();
+
+ // switch off recalibrations so those are not done multiple times
+ // this is just for safety, the recalibrated flag of cell object
+ // should not allow for farther processing anyways
+ fEMCALRecoUtils->SwitchOffRecalibration();
+ fEMCALRecoUtils->SwitchOffTimeRecalibration();
+
+ if (fDoUpdateOnly)
+ return;
+ }
- // Reclusterize
- if(fReClusterize){
+ // RECLUSTERIZATION ---------------------------------------------------------
+ if (fReClusterize)
+ {
FillDigitsArray();
Clusterize();
UpdateClusters();
AliWarning(Form("No cluster array, number of cells in event = %d, returning", cells->GetNumberOfCells()));
return;
}
+
+ // PROCESS SINGLE CLUSTER RECALIBRATION -------------------------------------
+ // now go through clusters one by one and process separate correction
+ // as those were defined or not
Int_t nclusters = clusArr->GetEntriesFast();
- for (Int_t icluster=0; icluster < nclusters; ++icluster) {
+ for (Int_t icluster=0; icluster < nclusters; ++icluster)
+ {
AliVCluster *clust = static_cast<AliVCluster*>(clusArr->At(icluster));
if (!clust)
continue;
if (!clust->IsEMCAL())
continue;
- if (fEMCALRecoUtils->ClusterContainsBadChannel(fEMCALGeo, clust->GetCellsAbsId(), clust->GetNCells())) {
- delete clusArr->RemoveAt(icluster);
- continue; //todo is it really needed to remove it? Or should we flag it?
+ // REMOVE CLUSTERS WITH BAD CELLS -----------------------------
+ if( fClusterBadChannelCheck )
+ {
+ // careful, the the ClusterContainsBadChannel is dependent on
+ // SwitchOnBadChannelsRemoval, switching it ON automatically
+ // and returning to original value after processing
+ Bool_t badRemoval = fEMCALRecoUtils->IsBadChannelsRemovalSwitchedOn();
+ fEMCALRecoUtils->SwitchOnBadChannelsRemoval();
+
+ Bool_t badResult = fEMCALRecoUtils->ClusterContainsBadChannel(fEMCALGeo, clust->GetCellsAbsId(), clust->GetNCells());
+
+ // switch the bad channels removal back
+ if( ! badRemoval )
+ fEMCALRecoUtils->SwitchOffBadChannelsRemoval();
+
+ if( badResult )
+ {
+ delete clusArr->RemoveAt(icluster);
+ continue; //TODO is it really needed to remove it? Or should we flag it?
+ }
}
-
- if (fFiducial){
+
+ // REMOVE EXOTIC CLUSTERS -------------------------------------
+ // does process local cell recalibration energy and time without replacing
+ // the global cell values, in case of no cell recalib done yet
+ if( fRejectExoticClusters )
+ {
+ // careful, the IsExoticCluster is dependent on
+ // SwitchOnRejectExoticCell, switching it ON automatically
+ // and returning to original value after processing
+ Bool_t exRemoval = fEMCALRecoUtils->IsRejectExoticCell();
+ fEMCALRecoUtils->SwitchOnRejectExoticCell();
+
+ // get bunch crossing
+ Int_t bunchCrossNo = event->GetBunchCrossNumber();
+
+ Bool_t exResult = fEMCALRecoUtils->IsExoticCluster(clust, cells, bunchCrossNo );
+
+ // switch the exotic channels removal back
+ if( ! exRemoval )
+ fEMCALRecoUtils->SwitchOffRejectExoticCell();
+
+ if( exResult )
+ {
+ delete clusArr->RemoveAt(icluster);
+ continue; //TODO is it really needed to remove it? Or should we flag it?
+ }
+ }
+
+ // FIDUCIAL CUT -----------------------------------------------
+ if (fFiducial)
+ {
+ // depends on SetNumberOfCellsFromEMCALBorder
+ // SwitchOnNoFiducialBorderInEMCALEta0
if (!fEMCALRecoUtils->CheckCellFiducialRegion(fEMCALGeo, clust, cells)){
delete clusArr->RemoveAt(icluster);
- continue; // todo it would be nice to store the distance
+ continue; //TODO it would be nice to store the distance
}
}
-
- fEMCALRecoUtils->CorrectClusterEnergyLinearity(clust);
- if(fRecalDistToBadChannels)
- fEMCALRecoUtils->RecalculateClusterDistanceToBadChannel(fEMCALGeo, cells, clust);
- if(fReCalibCluster)
+
+ // CLUSTER ENERGY ---------------------------------------------
+ // does process local cell recalibration energy and time without replacing
+ // the global cell values, in case of no cell recalib done yet
+ if( fReCalibCluster ) {
fEMCALRecoUtils->RecalibrateClusterEnergy(fEMCALGeo, clust, cells);
- if(fRecalClusPos)
+ if (clust->E() < 1e-9) {
+ delete clusArr->RemoveAt(icluster);
+ continue;
+ }
+ }
+
+ // CLUSTER POSITION -------------------------------------------
+ // does process local cell energy recalibration, if enabled and cells
+ // not calibrated yet
+ if( fRecalClusPos )
fEMCALRecoUtils->RecalculateClusterPosition(fEMCALGeo, cells, clust);
- }
- clusArr->Compress();
-
- // Track matching
- if (!TGeoGlobalMagField::Instance()->GetField()) {
- const AliESDRun *erun = event->GetESDRun();
- AliMagF *field = AliMagF::CreateFieldMap(erun->GetCurrentL3(),
- erun->GetCurrentDip(),
- AliMagF::kConvLHC,
- kFALSE,
- erun->GetBeamEnergy(),
- erun->GetBeamType());
- TGeoGlobalMagField::Instance()->SetField(field);
- }
- fEMCALRecoUtils->FindMatches(event,0x0,fEMCALGeo);
- Int_t nTracks = event->GetNumberOfTracks();
- if (nTracks>0) {
- SetClusterMatchedToTrack(event);
- SetTracksMatchedToCluster(event);
- }
-}
-
-//_____________________________________________________
-void AliEMCALTenderSupply::SetClusterMatchedToTrack(AliESDEvent *event)
-{
- // Checks if tracks are matched to EMC clusters and set the matched EMCAL cluster index to ESD track.
-
- Int_t nTracks = event->GetNumberOfTracks();
- for (Int_t iTrack = 0; iTrack < nTracks; ++iTrack) {
- AliESDtrack* track = event->GetTrack(iTrack);
- if (!track) {
- AliWarning(Form("Could not receive track %d", iTrack));
- continue;
+
+ // SHOWER SHAPE -----------------------------------------------
+ if( fRecalShowerShape )
+ fEMCALRecoUtils->RecalculateClusterShowerShapeParameters(fEMCALGeo, cells, clust);
+
+ // NONLINEARITY -----------------------------------------------
+ if( fDoNonLinearity )
+ {
+ Float_t correctedEnergy = fEMCALRecoUtils->CorrectClusterEnergyLinearity(clust);
+ clust->SetE(correctedEnergy);
}
- Int_t matchClusIndex = fEMCALRecoUtils->GetMatchedClusterIndex(iTrack);
- track->SetEMCALcluster(matchClusIndex); //sets -1 if track not matched within residual
- if(matchClusIndex != -1)
- track->SetStatus(AliESDtrack::kEMCALmatch);
- }
- if (fDebugLevel>2)
- AliInfo("Track matched to closest cluster");
-}
-
-//_____________________________________________________
-void AliEMCALTenderSupply::SetTracksMatchedToCluster(AliESDEvent *event)
-{
- // Checks if EMC clusters are matched to ESD track.
- // Adds track indexes of all the tracks matched to a cluster withing residuals in ESDCalocluster.
- for (Int_t iClus=0; iClus < event->GetNumberOfCaloClusters(); ++iClus) {
- AliESDCaloCluster *cluster = event->GetCaloCluster(iClus);
- if (!cluster->IsEMCAL())
- continue;
+ // DISTANCE TO BAD CHANNELS -----------------------------------
+ if( fRecalDistToBadChannels )
+ fEMCALRecoUtils->RecalculateClusterDistanceToBadChannel(fEMCALGeo, cells, clust);
+ }
- Int_t nTracks = event->GetNumberOfTracks();
- TArrayI arrayTrackMatched(nTracks);
+ clusArr->Compress();
- // Get the closest track matched to the cluster
- Int_t nMatched = 0;
- Int_t matchTrackIndex = fEMCALRecoUtils->GetMatchedTrackIndex(iClus);
- if (matchTrackIndex != -1) {
- arrayTrackMatched[nMatched] = matchTrackIndex;
- nMatched++;
- }
+ if (!fDoTrackMatch)
+ return;
- // Get all other tracks matched to the cluster
- for(Int_t iTrk=0; iTrk<nTracks; ++iTrk) {
- AliESDtrack* track = event->GetTrack(iTrk);
- if(iTrk == matchTrackIndex) continue;
- if(track->GetEMCALcluster() == iClus){
- arrayTrackMatched[nMatched] = iTrk;
- ++nMatched;
- }
+ // TRACK MATCHING -----------------------------------------------------------
+ if (!TGeoGlobalMagField::Instance()->GetField()) {
+ AliESDEvent *esd = dynamic_cast<AliESDEvent*>(event);
+ if (esd)
+ esd->InitMagneticField();
+ else {
+ AliAODEvent *aod = dynamic_cast<AliAODEvent*>(event);
+ Double_t curSol = 30000*aod->GetMagneticField()/5.00668;
+ Double_t curDip = 6000 *aod->GetMuonMagFieldScale();
+ AliMagF *field = AliMagF::CreateFieldMap(curSol,curDip);
+ TGeoGlobalMagField::Instance()->SetField(field);
}
-
- arrayTrackMatched.Set(nMatched);
- cluster->AddTracksMatched(arrayTrackMatched);
-
- Float_t eta= -999, phi = -999;
- if (matchTrackIndex != -1)
- fEMCALRecoUtils->GetMatchedResiduals(iClus, eta, phi);
- cluster->SetTrackDistance(phi, eta);
}
- if (fDebugLevel>2)
- AliInfo("Cluster matched to tracks");
+ fEMCALRecoUtils->FindMatches(event,0x0,fEMCALGeo);
+ fEMCALRecoUtils->SetClusterMatchedToTrack(event);
+ fEMCALRecoUtils->SetTracksMatchedToCluster(event);
}
//_____________________________________________________
{
// Initialising misalignment matrices
- AliESDEvent *event = fTender->GetEvent();
+ AliVEvent *event = GetEvent();
+
if (!event)
return kFALSE;
-
- if (fGeomMatrixSet) {
- AliInfo("Misalignment matrix already set");
+
+ if (fGeomMatrixSet)
+ {
+ AliInfo("Misalignment matrix already set");
return kTRUE;
}
-
+
if (fDebugLevel>0)
- AliInfo("Initialising misalignment matrix");
-
- fEMCALRecoUtils->SetPositionAlgorithm(AliEMCALRecoUtils::kPosTowerGlobal);
-
+ AliInfo("Initialising misalignment matrix");
+
if (fLoadGeomMatrices) {
- for(Int_t mod=0; mod < fEMCALGeo->GetNumberOfSuperModules(); ++mod) {
+ for(Int_t mod=0; mod < fEMCALGeo->GetNumberOfSuperModules(); ++mod)
+ {
if (fEMCALMatrix[mod]){
if(fDebugLevel > 2)
fEMCALMatrix[mod]->Print();
fGeomMatrixSet = kTRUE;
return kTRUE;
}
-
+
Int_t runGM = event->GetRunNumber();
- if (runGM <= 140000) { //2010 data
- Double_t rotationMatrix[4][9] = {{-0.014587, -0.999892, -0.002031, 0.999892, -0.014591, 0.001979, -0.002009, -0.002002, 0.999996},
- {-0.014587, 0.999892, 0.002031, 0.999892, 0.014591, -0.001979, -0.002009, 0.002002, -0.999996},
- {-0.345864, -0.938278, -0.003412, 0.938276, -0.345874, 0.003010, -0.004004, -0.002161, 0.999990},
- {-0.345861, 0.938280, 0.003412, 0.938276, 0.345874, -0.003010, -0.004004, 0.002161, -0.999990}};
-
- Double_t translationMatrix[4][3] = {{0.351659, 447.576446, 176.269742},
- {1.062577, 446.893974, -173.728870},
- {-154.213287, 419.306156, 176.753692},
- {-153.018950, 418.623681, -173.243605}};
-
- for(Int_t mod=0; mod < (fEMCALGeo->GetEMCGeometry())->GetNumberOfSuperModules(); ++mod) {
- //if(DebugLevel() > 1) fEMCALMatrix[mod]->Print();
- fEMCALMatrix[mod] = new TGeoHMatrix();
- fEMCALMatrix[mod]->SetRotation(rotationMatrix[mod]);
- fEMCALMatrix[mod]->SetTranslation(translationMatrix[mod]);
- fEMCALGeo->SetMisalMatrix(fEMCALMatrix[mod],mod);
- }
- } else if (runGM>140000 && runGM <148531 && (fFilepass == "pass1")) { // 2011 LHC11a pass1 data
- AliOADBContainer emcalgeoCont(Form("emcal2011"));
- emcalgeoCont.InitFromFile(Form("%s/BetaGood.root",fBasePath.Data()),Form("AliEMCALgeo"));
-
- TObjArray *mobj=(TObjArray*)emcalgeoCont.GetObject(100,"survey11byS");
- for(Int_t mod=0; mod < (fEMCALGeo->GetEMCGeometry())->GetNumberOfSuperModules(); mod++){
- fEMCALMatrix[mod] = (TGeoHMatrix*) mobj->At(mod);
- fEMCALGeo->SetMisalMatrix(fEMCALMatrix[mod],mod);
- fEMCALMatrix[mod]->Print();
- }
- } else {
- AliInfo(Form("No external misalignment matrix set: %d - %s, loading from ESD event", runGM, fFilepass.Data()));
- for(Int_t mod=0; mod < (fEMCALGeo->GetEMCGeometry())->GetNumberOfSuperModules(); mod++){
- if(fDebugLevel > 2)
- event->GetEMCALMatrix(mod)->Print();
- if(event->GetEMCALMatrix(mod))
- fEMCALGeo->SetMisalMatrix(event->GetEMCALMatrix(mod),mod) ;
+ TObjArray *mobj = 0;
+
+ if(fMisalignSurvey == kdefault)
+ { //take default alignment corresponding to run no
+ AliOADBContainer emcalgeoCont(Form("emcal"));
+ emcalgeoCont.InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALlocal2master.root",Form("AliEMCALgeo"));
+ mobj=(TObjArray*)emcalgeoCont.GetObject(runGM,"EmcalMatrices");
+ }
+
+ if(fMisalignSurvey == kSurveybyS)
+ { //take alignment at sector level
+ if (runGM <= 140000) { //2010 data
+ AliOADBContainer emcalgeoCont(Form("emcal2010"));
+ emcalgeoCont.InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALlocal2master.root",Form("AliEMCALgeo"));
+ mobj=(TObjArray*)emcalgeoCont.GetObject(100,"survey10");
+ }
+ else if (runGM>140000)
+ { // 2011 LHC11a pass1 data
+ AliOADBContainer emcalgeoCont(Form("emcal2011"));
+ emcalgeoCont.InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALlocal2master.root",Form("AliEMCALgeo"));
+ mobj=(TObjArray*)emcalgeoCont.GetObject(100,"survey11byS");
+ }
+ }
+
+ if(fMisalignSurvey == kSurveybyM)
+ { //take alignment at module level
+ if (runGM <= 140000) { //2010 data
+ AliOADBContainer emcalgeoCont(Form("emcal2010"));
+ emcalgeoCont.InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALlocal2master.root",Form("AliEMCALgeo"));
+ mobj=(TObjArray*)emcalgeoCont.GetObject(100,"survey10");
}
- fGeomMatrixSet=kTRUE;
+ else if (runGM>140000)
+ { // 2011 LHC11a pass1 data
+ AliOADBContainer emcalgeoCont(Form("emcal2011"));
+ emcalgeoCont.InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALlocal2master.root",Form("AliEMCALgeo"));
+ mobj=(TObjArray*)emcalgeoCont.GetObject(100,"survey11byM");
+ }
+ }
+
+ if(!mobj)
+ {
+ AliFatal("Geometry matrix array not found");
+ return kFALSE;
+ }
+
+ for(Int_t mod=0; mod < (fEMCALGeo->GetEMCGeometry())->GetNumberOfSuperModules(); mod++)
+ {
+ fEMCALMatrix[mod] = (TGeoHMatrix*) mobj->At(mod);
+ fEMCALGeo->SetMisalMatrix(fEMCALMatrix[mod],mod);
+ fEMCALMatrix[mod]->Print();
}
+
return kTRUE;
}
//_____________________________________________________
-Bool_t AliEMCALTenderSupply::InitBadChannels()
+Int_t AliEMCALTenderSupply::InitBadChannels()
{
// Initialising bad channel maps
- AliESDEvent *event = fTender->GetEvent();
- if (!event)
- return kFALSE;
+ AliVEvent *event = GetEvent();
+ if (!event)
+ return 0;
+
if (fDebugLevel>0)
AliInfo("Initialising Bad channel map");
+
+ // init default maps first
+ if( !fEMCALRecoUtils->GetEMCALBadChannelStatusMapArray() )
+ fEMCALRecoUtils->InitEMCALBadChannelStatusMap() ;
+
+ Int_t runBC = event->GetRunNumber();
+
+ AliOADBContainer *contBC = new AliOADBContainer("");
+ if (fBasePath!="")
+ { //if fBasePath specified in the ->SetBasePath()
+ if (fDebugLevel>0) AliInfo(Form("Loading Bad Channels OADB from given path %s",fBasePath.Data()));
+
+ TFile *fbad=new TFile(Form("%s/EMCALBadChannels.root",fBasePath.Data()),"read");
+ if (!fbad || fbad->IsZombie())
+ {
+ AliFatal(Form("EMCALBadChannels.root was not found in the path provided: %s",fBasePath.Data()));
+ return 0;
+ }
+
+ if (fbad) delete fbad;
+
+ contBC->InitFromFile(Form("%s/EMCALBadChannels.root",fBasePath.Data()),"AliEMCALBadChannels");
+ }
+ else
+ { // Else choose the one in the $ALICE_ROOT directory
+ if (fDebugLevel>0) AliInfo("Loading Bad Channels OADB from $ALICE_ROOT/OADB/EMCAL");
+
+ TFile *fbad=new TFile("$ALICE_ROOT/OADB/EMCAL/EMCALBadChannels.root","read");
+ if (!fbad || fbad->IsZombie())
+ {
+ AliFatal("$ALICE_ROOT/OADB/EMCAL/EMCALBadChannels.root was not found");
+ return 0;
+ }
+
+ if (fbad) delete fbad;
+
+ contBC->InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALBadChannels.root","AliEMCALBadChannels");
+ }
+
+ TObjArray *arrayBC=(TObjArray*)contBC->GetObject(runBC);
+ if (!arrayBC)
+ {
+ AliError(Form("No external hot channel set for run number: %d", runBC));
+ return 2;
+ }
- if (fFiducial){
- fEMCALRecoUtils->SetNumberOfCellsFromEMCALBorder(fNCellsFromEMCALBorder);
- fEMCALRecoUtils->SwitchOnNoFiducialBorderInEMCALEta0();
+ Int_t sms = fEMCALGeo->GetEMCGeometry()->GetNumberOfSuperModules();
+ for (Int_t i=0; i<sms; ++i)
+ {
+ TH2I *h = fEMCALRecoUtils->GetEMCALChannelStatusMap(i);
+ if (h)
+ delete h;
+ h=(TH2I*)arrayBC->FindObject(Form("EMCALBadChannelMap_Mod%d",i));
+
+ if (!h)
+ {
+ AliError(Form("Can not get EMCALBadChannelMap_Mod%d",i));
+ continue;
+ }
+ h->SetDirectory(0);
+ fEMCALRecoUtils->SetEMCALChannelStatusMap(i,h);
}
+ return 1;
+}
- fEMCALRecoUtils->SwitchOnBadChannelsRemoval();
- if (fRecalDistToBadChannels)
- fEMCALRecoUtils->SwitchOnDistToBadChannelRecalculation();
+//_____________________________________________________
+Int_t AliEMCALTenderSupply::InitRecalib()
+{
+ // Initialising recalibration factors.
+
+ AliVEvent *event = GetEvent();
- Int_t runBC = event->GetRunNumber();
- if (runBC <=140000) { //2010
- TFile *fbad = new TFile(Form("%s/BadChannels.root",fBasePath.Data()),"read");
- if (fbad->IsZombie()) {
- TString fPath = TString(fBasePath.Data());
- if (fPath.Contains("alien")) {
- if (fDebugLevel>1)
- AliInfo("Connecting to alien to get BadChannels.root");
- delete fbad;
- fbad = TFile::Open(Form("%s/BadChannelsDB/BadChannels.root",fBasePath.Data()));
- }
+ if (!event)
+ return 0;
+
+ if (fDebugLevel>0)
+ AliInfo("Initialising recalibration factors");
+
+ // init default maps first
+ if( !fEMCALRecoUtils->GetEMCALRecalibrationFactorsArray() )
+ fEMCALRecoUtils->InitEMCALRecalibrationFactors() ;
+
+ Int_t runRC = event->GetRunNumber();
+
+ AliOADBContainer *contRF=new AliOADBContainer("");
+ if (fBasePath!="")
+ { //if fBasePath specified in the ->SetBasePath()
+ if (fDebugLevel>0) AliInfo(Form("Loading Recalib OADB from given path %s",fBasePath.Data()));
+
+ TFile *fRecalib= new TFile(Form("%s/EMCALRecalib.root",fBasePath.Data()),"read");
+ if (!fRecalib || fRecalib->IsZombie())
+ {
+ AliFatal(Form("EMCALRecalib.root not found in %s",fBasePath.Data()));
+ return 0;
}
- for (Int_t i=0; i<4; ++i) {
- TH2I *h = fEMCALRecoUtils->GetEMCALChannelStatusMap(i);
- if (h)
- delete h;
- h = (TH2I*)fbad->Get(Form("EMCALBadChannelMap_Mod%d",i));
- if (!h) {
- AliError(Form("Can not get EMCALBadChannelMap_Mod%d",i));
- continue;
- }
- h->SetDirectory(0);
- fEMCALRecoUtils->SetEMCALChannelStatusMap(i,h);
+
+ if (fRecalib) delete fRecalib;
+
+ contRF->InitFromFile(Form("%s/EMCALRecalib.root",fBasePath.Data()),"AliEMCALRecalib");
+ }
+ else
+ { // Else choose the one in the $ALICE_ROOT directory
+ if (fDebugLevel>0) AliInfo("Loading Recalib OADB from $ALICE_ROOT/OADB/EMCAL");
+
+ TFile *fRecalib= new TFile("$ALICE_ROOT/OADB/EMCAL/EMCALRecalib.root","read");
+ if (!fRecalib || fRecalib->IsZombie())
+ {
+ AliFatal("$ALICE_ROOT/OADB/EMCAL/EMCALRecalib.root was not found");
+ return 0;
}
- delete fbad;
- return kTRUE;
+
+ if (fRecalib) delete fRecalib;
+
+ contRF->InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALRecalib.root","AliEMCALRecalib");
}
- //2011
- Int_t nSupMod=-1, nModule=-1, nIphi=-1, nIeta=-1, iphi=-1, ieta=-1;
- if (runBC>=144871 && runBC<=146860){ // LHC11a 2.76 TeV pp
-
- if (fFilepass == "pass1") { // pass1
- const Int_t nTowers=89;
- Int_t hotChannels[nTowers]={74, 103, 152, 320, 321, 322, 323, 324, 325, 326, 327, 328, 329, 330, 331, 332, 333, 334, 335, 368, 369, 370, 371, 372, 373, 374,375, 376, 377, 378, 379, 380, 381, 382, 383, 917, 1275, 1288, 1519, 1595, 1860, 1967, 2022, 2026, 2047, 2117, 2298, 2540, 2776, 3135, 3764, 6095, 6111, 6481, 6592, 6800, 6801, 6802, 6803, 6804, 6805, 6806, 6807, 6808, 6809, 6810, 6811, 6812, 6813, 6814, 6815, 7371, 7425, 7430, 7457, 7491, 7709, 8352, 8353, 8356, 8357, 8808, 8810, 8812, 8814, 9056, 9769, 9815, 9837};
- for (Int_t i=0; i<nTowers; ++i) {
- fEMCALGeo->GetCellIndex(hotChannels[i],nSupMod,nModule,nIphi,nIeta);
- fEMCALGeo->GetCellPhiEtaIndexInSModule(nSupMod,nModule,nIphi,nIeta,iphi,ieta);
- fEMCALRecoUtils->SetEMCALChannelStatus(nSupMod, ieta, iphi);
- }
- } else if (fFilepass == "pass2") { // pass2
- const Int_t nTowers=24;
- Int_t hotChannels[nTowers]= {74, 103, 152, 917, 1059, 1175, 1276, 1288, 1376, 1382, 1595, 2022, 2026, 2210, 2540, 2778, 2793, 3135, 3764, 5767, 6481, 7371, 7878, 9769};
- for(Int_t i=0; i<nTowers; ++i) {
- fEMCALGeo->GetCellIndex(hotChannels[i],nSupMod,nModule,nIphi,nIeta);
- fEMCALGeo->GetCellPhiEtaIndexInSModule(nSupMod,nModule,nIphi,nIeta,iphi,ieta);
- fEMCALRecoUtils->SetEMCALChannelStatus(nSupMod, ieta, iphi);
- }
- }
- return kTRUE;
+ TObjArray *recal=(TObjArray*)contRF->GetObject(runRC);
+ if (!recal)
+ {
+ AliError(Form("No Objects for run: %d",runRC));
+ return 2;
+ }
+
+ TObjArray *recalpass=(TObjArray*)recal->FindObject(fFilepass);
+ if (!recalpass)
+ {
+ AliError(Form("No Objects for run: %d - %s",runRC,fFilepass.Data()));
+ return 2;
}
- if (runBC>=151636 && runBC<=155384) { //LHC11c : 7TeV by Rongrong
- const Int_t nTowers=8;
- Int_t hotChannels[nTowers]={917, 2115, 2123, 2540, 6481, 9815, 10113, 10115};
- for(Int_t i=0; i<nTowers; i++) {
- fEMCALGeo->GetCellIndex(hotChannels[i],nSupMod,nModule,nIphi,nIeta);
- fEMCALGeo->GetCellPhiEtaIndexInSModule(nSupMod,nModule,nIphi,nIeta,iphi,ieta);
- fEMCALRecoUtils->SetEMCALChannelStatus(nSupMod, ieta, iphi);
- }
- return kTRUE;
+ TObjArray *recalib=(TObjArray*)recalpass->FindObject("Recalib");
+ if (!recalib)
+ {
+ AliError(Form("No Recalib histos found for %d - %s",runRC,fFilepass.Data()));
+ return 2;
}
- AliInfo(Form("No external hot channel set: %d - %s", runBC, fFilepass.Data()));
- return kTRUE;
+ if (fDebugLevel>0) recalib->Print();
+
+ Int_t sms = fEMCALGeo->GetEMCGeometry()->GetNumberOfSuperModules();
+ for (Int_t i=0; i<sms; ++i)
+ {
+ TH2F *h = fEMCALRecoUtils->GetEMCALChannelRecalibrationFactors(i);
+ if (h)
+ delete h;
+ h = (TH2F*)recalib->FindObject(Form("EMCALRecalFactors_SM%d",i));
+ if (!h)
+ {
+ AliError(Form("Could not load EMCALRecalFactors_SM%d",i));
+ continue;
+ }
+ h->SetDirectory(0);
+ fEMCALRecoUtils->SetEMCALChannelRecalibrationFactors(i,h);
+ }
+ return 1;
}
//_____________________________________________________
-Bool_t AliEMCALTenderSupply::InitRecalib()
+Int_t AliEMCALTenderSupply::InitRunDepRecalib()
{
// Initialising recalibration factors.
-
- AliESDEvent *event = fTender->GetEvent();
+
+ AliVEvent *event = GetEvent();
+
if (!event)
- return kFALSE;
-
+ return 0;
+
if (fDebugLevel>0)
AliInfo("Initialising recalibration factors");
-
- fEMCALRecoUtils->SwitchOnRecalibration();
+
+ // init default maps first
+ if( !fEMCALRecoUtils->GetEMCALRecalibrationFactorsArray() )
+ fEMCALRecoUtils->InitEMCALRecalibrationFactors() ;
+
Int_t runRC = event->GetRunNumber();
+
+ AliOADBContainer *contRF=new AliOADBContainer("");
+ if (fBasePath!="")
+ { //if fBasePath specified in the ->SetBasePath()
+ if (fDebugLevel>0) AliInfo(Form("Loading Recalib OADB from given path %s",fBasePath.Data()));
+
+ TFile *fRunDepRecalib= new TFile(Form("%s/EMCALTemperatureCorrCalib.root",fBasePath.Data()),"read");
+ if (!fRunDepRecalib || fRunDepRecalib->IsZombie())
+ {
+ AliFatal(Form("EMCALTemperatureCorrCalib.root not found in %s",fBasePath.Data()));
+ return 0;
+ }
+
+ if (fRunDepRecalib) delete fRunDepRecalib;
+
+ contRF->InitFromFile(Form("%s/EMCALTemperatureCorrCalib.root",fBasePath.Data()),"AliEMCALRunDepTempCalibCorrections");
+ }
+ else
+ { // Else choose the one in the $ALICE_ROOT directory
+ if (fDebugLevel>0) AliInfo("Loading Recalib OADB from $ALICE_ROOT/OADB/EMCAL");
+
+ TFile *fRunDepRecalib= new TFile("$ALICE_ROOT/OADB/EMCAL/EMCALTemperatureCorrCalib.root","read");
+ if (!fRunDepRecalib || fRunDepRecalib->IsZombie())
+ {
+ AliFatal("$ALICE_ROOT/OADB/EMCAL/EMCALTemperatureCorrCalib.root was not found");
+ return 0;
+ }
+
+ if (fRunDepRecalib) delete fRunDepRecalib;
+
+ contRF->InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALTemperatureCorrCalib.root","AliEMCALRunDepTempCalibCorrections");
+ }
+
+ TH1S *rundeprecal=(TH1S*)contRF->GetObject(runRC);
+ if (!rundeprecal)
+ {
+ AliWarning(Form("No TemperatureCorrCalib Objects for run: %d",runRC));
+ // let's get the closest runnumber instead then..
+ Int_t lower = 0;
+ Int_t ic = 0;
+ Int_t maxEntry = contRF->GetNumberOfEntries();
+
+ while ( (ic < maxEntry) && (contRF->UpperLimit(ic) < runRC) ) {
+ lower = ic;
+ ic++;
+ }
- TFile *fRecalib = 0;
- if (runRC <= 140000) { // 2010
- fRecalib = new TFile(Form("%s/RecalibrationFactors.root",fBasePath.Data()),"read");
- if (fRecalib->IsZombie()) {
- TString fPath = TString(fBasePath.Data());
- if(fPath.Contains("alien")) {
- if (fDebugLevel>1)
- AliInfo("Connecting to alien to get RecalibrationFactors.root");
-
- if( (runRC >= 114737 && runRC <= 117223 && (fFilepass = "pass2")) || //LHC10b pass2
- (runRC >= 118503 && runRC <= 121040 && (fFilepass = "pass2")) || //LHC10c pass2
- (runRC >= 118503 && runRC <= 121040 && (fFilepass = "pass3")) || //LHC10c pass3
- (runRC >= 122195 && runRC <= 126437 && (fFilepass = "pass1")) || //LHC10d pass1
- (runRC >= 127712 && runRC <= 130850 && (fFilepass = "pass1")) || //LHC10e pass1
- (runRC >= 133004 && runRC < 134657 && (fFilepass = "pass1"))) //LHC10f pass1<134657
- {
- fRecalib = TFile::Open(Form("%s/RecalDB/summer_december_2010/RecalibrationFactors.root",fBasePath.Data()));
- }
- else if((runRC >= 122195 && runRC <= 126437 && (fFilepass = "pass2")) || //LHC10d pass2
- (runRC >= 134657 && runRC <= 135029 && (fFilepass = "pass1")) || //LHC10f pass1 >= 134657
- (runRC >= 135654 && runRC <= 136377 && (fFilepass = "pass1")) || //LHC10g pass1
- (runRC >= 136851 && runRC < 137231 && (fFilepass = "pass1"))) //LHC10h pass1 untill christmas
- {
- fRecalib = TFile::Open(Form("%s/RecalDB/december2010/RecalibrationFactors.root",fBasePath.Data()));
- }
- else if(runRC >= 137231 && runRC <= 139517 && (fFilepass = "pass1")) //LHC10h pass1 from christmas
- {
- fRecalib = TFile::Open(Form("%s/RecalDB/summer2010/RecalibrationFactors.root",fBasePath.Data()));
- }
- else {
- AliError(Form("Run number or pass number not found: %d - %s; RECALIBRATION NOT APPLIED", runRC, fFilepass.Data()));
- return kTRUE;
- }
- }
+ Int_t closest = lower;
+ if ( (ic<maxEntry) &&
+ (contRF->LowerLimit(ic)-runRC) < (runRC - contRF->UpperLimit(lower)) ) {
+ closest = ic;
}
- } else if (runRC > 140000) { // 2011
- fRecalib = new TFile(Form("%s/RecalibrationFactors.root",fBasePath.Data()),"read");
- if (fRecalib->IsZombie()){
- TString fPath = TString(fBasePath.Data());
- if(fPath.Contains("alien")) {
- if (fDebugLevel>1)
- AliInfo("Connecting to alien to get RecalibrationFactors.root");
- fRecalib = TFile::Open(Form("%s/RecalDB/2011_v0/RecalibrationFactors.root",fBasePath.Data()));
- if (!fRecalib)
- AliError("Recalibration file not found");
- return kFALSE;
- }
+
+ AliWarning(Form("TemperatureCorrCalib Objects found closest id %d from run: %d", closest, contRF->LowerLimit(closest)));
+ rundeprecal = (TH1S*) contRF->GetObjectByIndex(closest);
+ }
+
+ if (fDebugLevel>0) rundeprecal->Print();
+
+ Int_t nSM = fEMCALGeo->GetEMCGeometry()->GetNumberOfSuperModules();
+
+ for (Int_t ism=0; ism<nSM; ++ism)
+ {
+ for (Int_t icol=0; icol<48; ++icol)
+ {
+ for (Int_t irow=0; irow<24; ++irow)
+ {
+ Float_t factor = fEMCALRecoUtils->GetEMCALChannelRecalibrationFactor(ism,icol,irow);
+
+ Int_t absID = fEMCALGeo->GetAbsCellIdFromCellIndexes(ism, irow, icol); // original calibration factor
+ factor *= rundeprecal->GetBinContent(absID) / 10000. ; // correction dependent on T
+
+ fEMCALRecoUtils->SetEMCALChannelRecalibrationFactor(ism,icol,irow,factor);
+ } // columns
+ } // rows
+ } // SM loop
+
+ return 1;
+}
+
+
+//_____________________________________________________
+Int_t AliEMCALTenderSupply::InitTimeCalibration()
+{
+ // Initialising bad channel maps
+ AliVEvent *event = GetEvent();
+
+ if (!event)
+ return 0;
+
+ if (fDebugLevel>0)
+ AliInfo("Initialising time calibration map");
+
+ // init default maps first
+ if ( !fEMCALRecoUtils->GetEMCALTimeRecalibrationFactorsArray() )
+ fEMCALRecoUtils->InitEMCALTimeRecalibrationFactors() ;
+
+ Int_t runBC = event->GetRunNumber();
+
+ AliOADBContainer *contBC = new AliOADBContainer("");
+ if (fBasePath!="")
+ { //if fBasePath specified in the ->SetBasePath()
+ if (fDebugLevel>0) AliInfo(Form("Loading time calibration OADB from given path %s",fBasePath.Data()));
+
+ TFile *fbad=new TFile(Form("%s/EMCALTimeCalib.root",fBasePath.Data()),"read");
+ if (!fbad || fbad->IsZombie())
+ {
+ AliFatal(Form("EMCALTimeCalib.root was not found in the path provided: %s",fBasePath.Data()));
+ return 0;
+ }
+
+ if (fbad) delete fbad;
+
+ contBC->InitFromFile(Form("%s/EMCALTimeCalib.root",fBasePath.Data()),"AliEMCALTimeCalib");
+ }
+ else
+ { // Else choose the one in the $ALICE_ROOT directory
+ if (fDebugLevel>0) AliInfo("Loading time calibration OADB from $ALICE_ROOT/OADB/EMCAL");
+
+ TFile *fbad=new TFile("$ALICE_ROOT/OADB/EMCAL/EMCALTimeCalib.root","read");
+ if (!fbad || fbad->IsZombie())
+ {
+ AliFatal("$ALICE_ROOT/OADB/EMCAL/EMCALTimeCalib.root was not found");
+ return 0;
+ }
+
+ if (fbad) delete fbad;
+
+ contBC->InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALTimeCalib.root","AliEMCALTimeCalib");
+ }
+
+ TObjArray *arrayBC=(TObjArray*)contBC->GetObject(runBC);
+ if (!arrayBC)
+ {
+ AliError(Form("No external time calibration set for run number: %d", runBC));
+ return 2;
+ }
+
+ // Here, it looks for a specific pass
+ TObjArray *arrayBCpass=(TObjArray*)arrayBC->FindObject(fFilepass);
+ if (!arrayBCpass)
+ {
+ AliError(Form("No external time calibration set for: %d -%s", runBC,fFilepass.Data()));
+ return 2;
+ }
+
+ if (fDebugLevel>0) arrayBCpass->Print();
+
+ for( Int_t i = 0; i < 4; i++ )
+ {
+ TH1F *h = fEMCALRecoUtils->GetEMCALChannelTimeRecalibrationFactors( i );
+ if( h )
+ delete h;
+
+ h = (TH1F*)arrayBCpass->FindObject(Form("hAllTimeAvBC%d",i));
+
+ if (!h)
+ {
+ AliError(Form("Can not get hAllTimeAvBC%d",i));
+ continue;
}
+ h->SetDirectory(0);
+ fEMCALRecoUtils->SetEMCALChannelTimeRecalibrationFactors(i,h);
+ }
+ return 1;
+}
- Int_t sms = fEMCALGeo->GetEMCGeometry()->GetNumberOfSuperModules();
- for (Int_t i=0; i<sms; ++i) {
- TH2F *h = fEMCALRecoUtils->GetEMCALChannelRecalibrationFactors(i);
- if (h)
- delete h;
- h = (TH2F*)fRecalib->Get(Form("EMCALRecalFactors_SM%d",i));
- if (!h) {
- AliError(Form("Could not load EMCALRecalFactors_SM%d",i));
- continue;
- }
- h->SetDirectory(0);
- fEMCALRecoUtils->SetEMCALChannelRecalibrationFactors(i,h);
+//_____________________________________________________
+void AliEMCALTenderSupply::UpdateCells()
+{
+ //Remove bad cells from the cell list
+ //Recalibrate energy and time cells
+ //This is required for later reclusterization
+
+ AliVEvent *event = GetEvent();
+
+ if(!event) return ;
+
+ AliVCaloCells *cells = event->GetEMCALCells();
+ Int_t bunchCrossNo = event->GetBunchCrossNumber();
+
+ fEMCALRecoUtils->RecalibrateCells(cells, bunchCrossNo);
+
+ // remove exotic cells - loop through cells and zero the exotic ones
+ // just like with bad cell rejection in reco utils (inside RecalibrateCells)
+ if( fRejectExoticCells )
+ {
+ Short_t absId =-1;
+ Double_t ecell = 0;
+ Double_t tcell = 0;
+ Double_t efrac = 0;
+ Int_t mclabel = -1;
+ Bool_t isExot = kFALSE;
+
+ // loop through cells
+ Int_t nEMcell = cells->GetNumberOfCells() ;
+ for (Int_t iCell = 0; iCell < nEMcell; iCell++)
+ {
+ cells->GetCell( iCell, absId, ecell, tcell, mclabel, efrac );
+
+ isExot = fEMCALRecoUtils->IsExoticCell( absId, cells, bunchCrossNo );
+ // zero if exotic
+ if( isExot )
+ cells->SetCell( iCell, absId, 0.0, -1.0, mclabel, efrac );
+ } // cell loop
+ } // reject exotic cells
+
+ cells->Sort();
+}
+
+//_____________________________________________________
+TString AliEMCALTenderSupply::GetBeamType()
+{
+ // Get beam type : pp-AA-pA
+ // ESDs have it directly, AODs get it from hardcoded run number ranges
+
+ AliVEvent *event = GetEvent();
+
+ if (!event) {
+ AliError("Couldn't retrieve event!");
+ return "";
+ }
+
+ TString beamType;
+
+ AliESDEvent *esd = dynamic_cast<AliESDEvent*>(event);
+ if (esd) {
+ const AliESDRun *run = esd->GetESDRun();
+ beamType = run->GetBeamType();
+ }
+ else
+ {
+ Int_t runNumber = event->GetRunNumber();
+ if ((runNumber >= 136851 && runNumber <= 139517) // LHC10h
+ || (runNumber >= 166529 && runNumber <= 170593)) // LHC11h
+ {
+ beamType = "A-A";
+ }
+ else
+ {
+ beamType = "p-p";
}
}
- delete fRecalib;
- return kTRUE;
+
+ return beamType;
}
//_____________________________________________________
-void AliEMCALTenderSupply::RecalibrateCells()
+Int_t AliEMCALTenderSupply::InitRecParam()
{
- // Recalibrate cells.
-
- AliESDCaloCells *cells = fTender->GetEvent()->GetEMCALCells();
-
- Int_t nEMCcell = cells->GetNumberOfCells();
- for(Int_t icell=0; icell<nEMCcell; ++icell) {
- Int_t imod = -1, iphi =-1, ieta=-1,iTower = -1, iIphi = -1, iIeta = -1;
- fEMCALGeo->GetCellIndex(cells->GetCellNumber(icell),imod,iTower,iIphi,iIeta);
- fEMCALGeo->GetCellPhiEtaIndexInSModule(imod,iTower,iIphi, iIeta,iphi,ieta);
- Double_t calibFactor = fEMCALRecoUtils->GetEMCALChannelRecalibrationFactor(imod,ieta,iphi);
- cells->SetCell(icell,cells->GetCellNumber(icell),cells->GetAmplitude(icell)*calibFactor,cells->GetTime(icell));
- }
+ // Init reco params if not yet exist (probably shipped by the user already)
+
+ if( fRecParam != 0 )
+ return 2;
+
+ TString beamType = GetBeamType();
+
+ // set some default reco params
+ fRecParam = new AliEMCALRecParam();
+ fRecParam->SetClusteringThreshold(0.100);
+ fRecParam->SetMinECut(0.050);
+ fRecParam->SetTimeCut(250);
+ fRecParam->SetTimeMin(425);
+ fRecParam->SetTimeMax(825);
+ if ( beamType == "A-A") {
+ fRecParam->SetClusterizerFlag(AliEMCALRecParam::kClusterizerv2);
+ }
+ else {
+ fRecParam->SetClusterizerFlag(AliEMCALRecParam::kClusterizerv1);
+ }
+
+ return 0;
}
//_____________________________________________________
Bool_t AliEMCALTenderSupply::InitClusterization()
{
// Initialing clusterization/unfolding algorithm and set all the needed parameters.
+
+ AliVEvent *event = GetEvent();
- AliESDEvent *event=fTender->GetEvent();
if (!event)
return kFALSE;
-
+
if (fDebugLevel>0)
- AliInfo(Form("Initialising reclustering parameters: Clusterizer-%d",fRecParam->GetClusterizerFlag()));
-
+ AliInfo(Form("Initialising reclustering parameters: Clusterizer type: %d",fRecParam->GetClusterizerFlag()));
+
//---setup clusterizer
- delete fClusterizer;
+ if (fClusterizer) {
+ // avoid deleting fDigitsArr
+ fClusterizer->SetDigitsArr(0);
+ delete fClusterizer;
+ fClusterizer = 0;
+ }
if (fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerv1)
fClusterizer = new AliEMCALClusterizerv1 (fEMCALGeo);
- else if(fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerv2)
+ else if (fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerv2)
fClusterizer = new AliEMCALClusterizerv2(fEMCALGeo);
- else if(fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerNxN) {
+ else if (fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerNxN)
+ {
AliEMCALClusterizerNxN *clusterizer = new AliEMCALClusterizerNxN(fEMCALGeo);
clusterizer->SetNRowDiff(fRecParam->GetNRowDiff());
clusterizer->SetNColDiff(fRecParam->GetNColDiff());
fClusterizer = clusterizer;
- } else {
+ }
+ else
+ {
AliFatal(Form("Clusterizer < %d > not available", fRecParam->GetClusterizerFlag()));
return kFALSE;
}
-
+
// Set the clustering parameters
fClusterizer->SetECAClusteringThreshold( fRecParam->GetClusteringThreshold() );
fClusterizer->SetECALogWeight ( fRecParam->GetW0() );
fClusterizer->SetTimeMax ( fRecParam->GetTimeMax() );
fClusterizer->SetInputCalibrated ( kTRUE );
fClusterizer->SetJustClusters ( kTRUE );
-
+
// In case of unfolding after clusterization is requested, set the corresponding parameters
- if (fRecParam->GetUnfold()) {
- for (Int_t i = 0; i < 8; ++i) {
+ if (fRecParam->GetUnfold())
+ {
+ for (Int_t i = 0; i < 8; ++i)
+ {
fClusterizer->SetSSPars(i, fRecParam->GetSSPars(i));
}
- for (Int_t i = 0; i < 3; ++i) {
+ for (Int_t i = 0; i < 3; ++i)
+ {
fClusterizer->SetPar5 (i, fRecParam->GetPar5(i));
fClusterizer->SetPar6 (i, fRecParam->GetPar6(i));
}
fClusterizer->InitClusterUnfolding();
}
-
+
fClusterizer->SetDigitsArr(fDigitsArr);
fClusterizer->SetOutput(0);
fClusterArr = const_cast<TObjArray *>(fClusterizer->GetRecPoints());
void AliEMCALTenderSupply::FillDigitsArray()
{
// Fill digits from cells to a TClonesArray.
+
+ AliVEvent *event = GetEvent();
- AliESDEvent *event = fTender->GetEvent();
- if (!event)
+ if (!event)
return;
+
+ // In case of MC productions done before aliroot tag v5-02-Rev09
+ // assing the cluster label to all the cells belonging to this cluster
+ // very rough
+ Int_t cellLabels[12672];
+ if(fSetCellMCLabelFromCluster)
+ {
+ for (Int_t i = 0; i < 12672; i++)
+ {
+ cellLabels [i] = 0 ;
+ fOrgClusterCellId[i] =-1 ;
+ }
+
+ Int_t nClusters = event->GetNumberOfCaloClusters();
+ for (Int_t i = 0; i < nClusters; i++)
+ {
+ AliVCluster *clus = event->GetCaloCluster(i);
+
+ if(!clus) continue;
+
+ if(!clus->IsEMCAL()) continue ;
+
+ Int_t label = clus->GetLabel();
+ UShort_t * index = clus->GetCellsAbsId() ;
+
+ for(Int_t icell=0; icell < clus->GetNCells(); icell++ )
+ {
+ cellLabels [index[icell]] = label;
+ fOrgClusterCellId[index[icell]] = i ; // index of the original cluster
+ }
+ }// cluster loop
+ }
fDigitsArr->Clear("C");
- AliESDCaloCells *cells = event->GetEMCALCells();
+ AliVCaloCells *cells = event->GetEMCALCells();
Int_t ncells = cells->GetNumberOfCells();
- for (Int_t icell = 0, idigit = 0; icell < ncells; ++icell) {
- Double_t cellAmplitude=0, cellTime=0;
- Short_t cellNumber=0;
- if (cells->GetCell(icell, cellNumber, cellAmplitude, cellTime) != kTRUE)
+ for (Int_t icell = 0, idigit = 0; icell < ncells; ++icell)
+ {
+ Double_t cellAmplitude=0, cellTime=0, efrac = 0;
+ Short_t cellNumber=0;
+ Int_t mcLabel=-1;
+
+ if (cells->GetCell(icell, cellNumber, cellAmplitude, cellTime, mcLabel, efrac) != kTRUE)
break;
- AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(fDigitsArr->New(idigit));
- digit->SetId(cellNumber);
- digit->SetTime(cellTime);
- digit->SetTimeR(cellTime);
- digit->SetIndexInList(idigit);
- digit->SetType(AliEMCALDigit::kHG);
- digit->SetAmplitude(cellAmplitude);
+
+ // Do not add if energy already too low (some cells set to 0 if bad channels)
+ if (cellAmplitude < fRecParam->GetMinECut())
+ continue;
+
+ // If requested, do not include exotic cells
+ if (fEMCALRecoUtils->IsExoticCell(cellNumber,cells,event->GetBunchCrossNumber()))
+ continue;
+
+ if ( fSetCellMCLabelFromCluster ) mcLabel = cellLabels[cellNumber];
+ else if( fRemapMCLabelForAODs ) RemapMCLabelForAODs(mcLabel);
+
+ if (mcLabel > 0 && efrac < 1e-6) efrac = 1;
+
+ new((*fDigitsArr)[idigit]) AliEMCALDigit(mcLabel, mcLabel, cellNumber,
+ (Float_t)cellAmplitude, (Float_t)cellTime,
+ AliEMCALDigit::kHG,idigit, 0, 0, efrac*cellAmplitude);
idigit++;
}
}
void AliEMCALTenderSupply::Clusterize()
{
// Clusterize.
-
+
fClusterizer->Digits2Clusters("");
}
void AliEMCALTenderSupply::UpdateClusters()
{
// Update ESD cluster list.
+
+ AliVEvent *event = GetEvent();
- AliESDEvent *event = fTender->GetEvent();
if (!event)
return;
-
+
TClonesArray *clus = dynamic_cast<TClonesArray*>(event->FindListObject("caloClusters"));
if (!clus)
clus = dynamic_cast<TClonesArray*>(event->FindListObject("CaloClusters"));
- if (!clus) {
- AliError(" Null calo clusters array, returning");
+ if (!clus)
+ {
+ AliError(" Null pointer to calo clusters array, returning");
return;
}
+
+ // Before destroying the orignal list, assign to the rec points the MC labels
+ // of the original clusters, if requested
+ if( fSetCellMCLabelFromCluster == 2 ) SetClustersMCLabelFromOriginalClusters() ;
Int_t nents = clus->GetEntriesFast();
- for (Int_t i=0; i < nents; ++i) {
- AliESDCaloCluster *c = static_cast<AliESDCaloCluster*>(clus->At(i));
+ for (Int_t i=0; i < nents; ++i)
+ {
+ AliVCluster *c = dynamic_cast<AliVCluster*>(clus->At(i));
if (!c)
continue;
if (c->IsEMCAL())
+ {
delete clus->RemoveAt(i);
+ }
}
+
clus->Compress();
+
RecPoints2Clusters(clus);
}
//_____________________________________________________
void AliEMCALTenderSupply::RecPoints2Clusters(TClonesArray *clus)
{
- // Convert AliEMCALRecoPoints to AliESDCaloClusters.
+ // Convert AliEMCALRecoPoints to AliESDCaloClusters/AliAODCaloClusters.
// Cluster energy, global position, cells and their amplitude fractions are restored.
- AliESDEvent *event = fTender->GetEvent();
+ AliVEvent *event = GetEvent();
+
if (!event)
return;
-
+
Int_t ncls = fClusterArr->GetEntriesFast();
- for(Int_t i=0, nout=clus->GetEntriesFast(); i < ncls; ++i) {
+ for(Int_t i=0, nout=clus->GetEntriesFast(); i < ncls; ++i)
+ {
AliEMCALRecPoint *recpoint = static_cast<AliEMCALRecPoint*>(fClusterArr->At(i));
- Int_t ncells_true = 0;
+ Int_t ncellsTrue = 0;
const Int_t ncells = recpoint->GetMultiplicity();
UShort_t absIds[ncells];
Double32_t ratios[ncells];
Int_t *dlist = recpoint->GetDigitsList();
Float_t *elist = recpoint->GetEnergiesList();
- for (Int_t c = 0; c < ncells; ++c) {
+ for (Int_t c = 0; c < ncells; ++c)
+ {
AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(fDigitsArr->At(dlist[c]));
- absIds[ncells_true] = digit->GetId();
- ratios[ncells_true] = elist[c]/digit->GetAmplitude();
- if (ratios[ncells_true] < 0.001)
+ absIds[ncellsTrue] = digit->GetId();
+ ratios[ncellsTrue] = elist[c]/digit->GetAmplitude();
+ if (ratios[ncellsTrue] < 0.001)
continue;
- ++ncells_true;
+ ++ncellsTrue;
}
- if (ncells_true < 1) {
+ if (ncellsTrue < 1)
+ {
AliWarning("Skipping cluster with no cells");
continue;
}
Float_t g[3];
gpos.GetXYZ(g);
- AliESDCaloCluster *c = static_cast<AliESDCaloCluster*>(clus->New(nout++));
+ AliVCluster *c = static_cast<AliVCluster*>(clus->New(nout++));
+ c->SetID(nout-1);
c->SetType(AliVCluster::kEMCALClusterv1);
c->SetE(recpoint->GetEnergy());
c->SetPosition(g);
- c->SetNCells(ncells_true);
+ c->SetNCells(ncellsTrue);
c->SetDispersion(recpoint->GetDispersion());
c->SetEmcCpvDistance(-1); //not yet implemented
c->SetChi2(-1); //not yet implemented
recpoint->GetElipsAxis(elipAxis);
c->SetM02(elipAxis[0]*elipAxis[0]) ;
c->SetM20(elipAxis[1]*elipAxis[1]) ;
- AliESDCaloCluster *cesd = static_cast<AliESDCaloCluster*>(c);
- cesd->SetCellsAbsId(absIds);
- cesd->SetCellsAmplitudeFraction(ratios);
+ c->SetCellsAbsId(absIds);
+ c->SetCellsAmplitudeFraction(ratios);
+
+ //MC labels
+ Int_t parentMult = 0;
+ Int_t *parentList = recpoint->GetParents(parentMult);
+ if (parentMult > 0) c->SetLabel(parentList, parentMult);
+
+ }
+}
+
+//___________________________________________________________
+void AliEMCALTenderSupply::RemapMCLabelForAODs(Int_t & label)
+{
+ // MC label for Cells not remapped after ESD filtering, do it here.
+
+ if(label < 0) return ;
+
+ AliAODEvent * evt = dynamic_cast<AliAODEvent*> (GetEvent()) ;
+ if(!evt) return ;
+
+ TClonesArray * arr = dynamic_cast<TClonesArray*>(evt->FindListObject("mcparticles")) ;
+ if(!arr) return ;
+
+ if(label < arr->GetEntriesFast())
+ {
+ AliAODMCParticle * particle = dynamic_cast<AliAODMCParticle *>(arr->At(label));
+ if(!particle) return ;
+
+ if(label == particle->Label()) return ; // label already OK
+ }
+
+ // loop on the particles list and check if there is one with the same label
+ for(Int_t ind = 0; ind < arr->GetEntriesFast(); ind++ )
+ {
+ AliAODMCParticle * particle = dynamic_cast<AliAODMCParticle *>(arr->At(ind));
+ if(!particle) continue ;
+
+ if(label == particle->Label())
+ {
+ label = ind;
+ return;
+ }
}
+
+ label = -1;
+
+}
+
+//_____________________________________________________________________________________________
+void AliEMCALTenderSupply::SetClustersMCLabelFromOriginalClusters()
+{
+ // Get the original clusters that contribute to the new rec point cluster,
+ // assign the labels of such clusters to the new rec point cluster.
+ // Only approximatedly valid when output are V1 clusters, or input/output clusterizer
+ // are the same handle with care
+ // Copy from same method in AliAnalysisTaskEMCALClusterize, but here modify the recpoint and
+ // not the output calocluster
+
+ Int_t ncls = fClusterArr->GetEntriesFast();
+ for(Int_t irp=0; irp < ncls; ++irp)
+ {
+ TArrayI clArray(300) ; //Weird if more than a few clusters are in the origin ...
+ clArray.Reset();
+ Int_t nClu = 0;
+ Int_t nLabTotOrg = 0;
+ Float_t emax = -1;
+ Int_t idMax = -1;
+
+ AliEMCALRecPoint *clus = static_cast<AliEMCALRecPoint*>(fClusterArr->At(irp));
+
+ //Find the clusters that originally had the cells
+ const Int_t ncells = clus->GetMultiplicity();
+ Int_t *digList = clus->GetDigitsList();
+
+ for ( Int_t iLoopCell = 0 ; iLoopCell < ncells ; iLoopCell++ )
+ {
+ AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(fDigitsArr->At(digList[iLoopCell]));
+ Int_t idCell = digit->GetId();
+
+ if(idCell>=0)
+ {
+ Int_t idCluster = fOrgClusterCellId[idCell];
+ Bool_t set = kTRUE;
+ for(Int_t icl =0; icl < nClu; icl++)
+ {
+ if(((Int_t)clArray.GetAt(icl))==-1 ) continue;
+ if( idCluster == ((Int_t)clArray.GetAt(icl)) ) set = kFALSE;
+ }
+ if( set && idCluster >= 0)
+ {
+ clArray.SetAt(idCluster,nClu++);
+ nLabTotOrg+=(GetEvent()->GetCaloCluster(idCluster))->GetNLabels();
+
+ //Search highest E cluster
+ AliVCluster * clOrg = GetEvent()->GetCaloCluster(idCluster);
+ if(emax < clOrg->E())
+ {
+ emax = clOrg->E();
+ idMax = idCluster;
+ }
+ }
+ }
+ }// cell loop
+
+ // Put the first in the list the cluster with highest energy
+ if(idMax != ((Int_t)clArray.GetAt(0))) // Max not at first position
+ {
+ Int_t maxIndex = -1;
+ Int_t firstCluster = ((Int_t)clArray.GetAt(0));
+ for ( Int_t iLoopCluster = 0 ; iLoopCluster < nClu ; iLoopCluster++ )
+ {
+ if(idMax == ((Int_t)clArray.GetAt(iLoopCluster))) maxIndex = iLoopCluster;
+ }
+
+ if(firstCluster >=0 && idMax >=0)
+ {
+ clArray.SetAt(idMax,0);
+ clArray.SetAt(firstCluster,maxIndex);
+ }
+ }
+
+ // Get the labels list in the original clusters, assign all to the new cluster
+ TArrayI clMCArray(nLabTotOrg) ;
+ clMCArray.Reset();
+
+ Int_t nLabTot = 0;
+ for ( Int_t iLoopCluster = 0 ; iLoopCluster < nClu ; iLoopCluster++ )
+ {
+ Int_t idCluster = (Int_t) clArray.GetAt(iLoopCluster);
+ AliVCluster * clOrg = GetEvent()->GetCaloCluster(idCluster);
+ Int_t nLab = clOrg->GetNLabels();
+
+ for ( Int_t iLab = 0 ; iLab < nLab ; iLab++ )
+ {
+ Int_t lab = clOrg->GetLabelAt(iLab) ;
+ if(lab>=0)
+ {
+ Bool_t set = kTRUE;
+ for(Int_t iLabTot =0; iLabTot < nLabTot; iLabTot++)
+ {
+ if( lab == ((Int_t)clMCArray.GetAt(iLabTot)) ) set = kFALSE;
+ }
+ if( set ) clMCArray.SetAt(lab,nLabTot++);
+ }
+ }
+ }// cluster loop
+
+ // Set the final list of labels to rec point
+
+ Int_t *labels = new Int_t[nLabTot];
+ for(Int_t il = 0; il < nLabTot; il++) labels[il] = clMCArray.GetArray()[il];
+ clus->SetParents(nLabTot,labels);
+
+ } // rec point array
}
+
//_____________________________________________________
void AliEMCALTenderSupply::GetPass()
{
- // Get passx from filename.
-
+ // Get passx from filename
+
AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
fInputTree = mgr->GetTree();
-
- if (!fInputTree) {
- AliError("Ptr to tree = 0, returning");
+
+ if (!fInputTree)
+ {
+ AliError("Pointer to tree = 0, returning");
return;
}
-
+
fInputFile = fInputTree->GetCurrentFile();
if (!fInputFile) {
- AliError("Null input file, returning");
+ AliError("Null pointer input file, returning");
return;
}
-
+
TString fname(fInputFile->GetName());
- if (fname.Contains("pass1")) fFilepass = TString("pass1");
- else if(fname.Contains("pass2")) fFilepass = TString("pass2");
- else if(fname.Contains("pass3")) fFilepass = TString("pass3");
- else {
+ if (fname.Contains("pass1")) fFilepass = TString("pass1");
+ else if (fname.Contains("pass2")) fFilepass = TString("pass2");
+ else if (fname.Contains("pass3")) fFilepass = TString("pass3");
+ else if (fname.Contains("pass4")) fFilepass = TString("pass4");
+ else
+ {
AliError(Form("Pass number string not found: %s", fname.Data()));
return;
}
}
-