// Base class for the clusterization algorithm (pure abstract)
//*--
//*-- Author: Yves Schutz SUBATECH
-// Modif:
-// August 2002 Yves Schutz: clone PHOS as closely as possible and intoduction
-// of new IO (à la PHOS)
+//
+// Clusterization mother class. Contains common methods/data members of different
+// clusterizers. GCB 2010
//////////////////////////////////////////////////////////////////////////////
// --- ROOT system ---
#include "TClonesArray.h"
#include "TTree.h"
+#include <TFile.h>
+class TFolder;
+#include <TMath.h>
+#include <TMinuit.h>
+#include <TTree.h>
+class TSystem;
+#include <TBenchmark.h>
+#include <TBrowser.h>
+#include <TROOT.h>
// --- Standard library ---
-
+#include <cassert>
// --- AliRoot header files ---
#include "AliEMCALClusterizer.h"
+#include "AliEMCALReconstructor.h"
+#include "AliRunLoader.h"
+#include "AliRun.h"
#include "AliLog.h"
+#include "AliEMCAL.h"
+#include "AliEMCALRecPoint.h"
+#include "AliEMCALRecParam.h"
+#include "AliEMCALGeometry.h"
+#include "AliEMCALRecParam.h"
+#include "AliCDBManager.h"
+#include "AliCaloCalibPedestal.h"
+#include "AliEMCALCalibData.h"
+class AliCDBStorage;
+#include "AliCDBEntry.h"
ClassImp(AliEMCALClusterizer)
AliEMCALClusterizer::AliEMCALClusterizer():
fDigitsArr(NULL),
fTreeR(NULL),
- fRecPoints(NULL)
+ fRecPoints(NULL),
+ fGeom(NULL),
+ fCalibData(NULL),
+ fCaloPed(NULL),
+ fADCchannelECA(0.),fADCpedestalECA(0.),
+ fTimeMin(-1.),fTimeMax(1.),fTimeCut(1.),
+ fDefaultInit(kFALSE),fToUnfold(kFALSE),
+ fNumberOfECAClusters(0), fECAClusteringThreshold(0.),
+ fECALocMaxCut(0.),fECAW0(0.),fMinECut(0.)
{
// ctor
+
+ Init();
+
+}
+
+//____________________________________________________________________________
+AliEMCALClusterizer::AliEMCALClusterizer(AliEMCALGeometry* geometry):
+ fDigitsArr(NULL),
+ fTreeR(NULL),
+ fRecPoints(NULL),
+ fGeom(geometry),
+ fCalibData(NULL),
+ fCaloPed(NULL),
+ fADCchannelECA(0.),fADCpedestalECA(0.),
+ fTimeMin(-1.),fTimeMax(1.),fTimeCut(1.),
+ fDefaultInit(kFALSE),fToUnfold(kFALSE),
+ fNumberOfECAClusters(0), fECAClusteringThreshold(0.),
+ fECALocMaxCut(0.),fECAW0(0.),fMinECut(0.)
+{
+ // ctor with the indication of the file where header Tree and digits Tree are stored
+ // use this contructor to avoid usage of Init() which uses runloader
+ // change needed by HLT - MP
+
+ // Note for the future: the use on runloader should be avoided or optional at least
+ // another way is to make Init virtual and protected at least such that the deriving classes can overload
+ // Init() ;
+ //
+
+ if (!fGeom)
+ {
+ AliFatal("Geometry not initialized.");
+ }
+
}
+//____________________________________________________________________________
+AliEMCALClusterizer::AliEMCALClusterizer(AliEMCALGeometry* geometry, AliEMCALCalibData * calib, AliCaloCalibPedestal * caloped):
+ fDigitsArr(NULL),
+ fTreeR(NULL),
+ fRecPoints(NULL),
+ fGeom(geometry),
+ fCalibData(calib),
+ fCaloPed(caloped),
+ fADCchannelECA(0.),fADCpedestalECA(0.),
+ fTimeMin(-1.),fTimeMax(1.),fTimeCut(1.),
+ fDefaultInit(kFALSE),fToUnfold(kFALSE),
+ fNumberOfECAClusters(0), fECAClusteringThreshold(0.),
+ fECALocMaxCut(0.),fECAW0(0.),fMinECut(0.)
+{
+ // ctor, geometry and calibration are initialized elsewhere.
+
+ if (!fGeom)
+ AliFatal("Geometry not initialized.");
+
+}
+
+
//____________________________________________________________________________
AliEMCALClusterizer::~AliEMCALClusterizer()
{
}
}
+//____________________________________________________________________________
+Float_t AliEMCALClusterizer::Calibrate(const Float_t amp, const Float_t time, const Int_t absId)
+{
+
+ // Convert digitized amplitude into energy.
+ // Calibration parameters are taken from calibration data base for raw data,
+ // or from digitizer parameters for simulated data.
+ if(fCalibData){
+
+ if (fGeom==0)
+ AliFatal("Did not get geometry from EMCALLoader") ;
+
+ Int_t iSupMod = -1;
+ Int_t nModule = -1;
+ Int_t nIphi = -1;
+ Int_t nIeta = -1;
+ Int_t iphi = -1;
+ Int_t ieta = -1;
+
+ Bool_t bCell = fGeom->GetCellIndex(absId, iSupMod, nModule, nIphi, nIeta) ;
+ if(!bCell) {
+ fGeom->PrintGeometry();
+ Error("Calibrate()"," Wrong cell id number : %i", absId);
+ assert(0);
+ }
+
+ fGeom->GetCellPhiEtaIndexInSModule(iSupMod,nModule,nIphi, nIeta,iphi,ieta);
+
+ // Check if channel is bad (dead or hot), in this case return 0.
+ // Gustavo: 15-12-09 In case of RAW data this selection is already done, but not in simulation.
+ // for the moment keep it here but remember to do the selection at the sdigitizer level
+ // and remove it from here
+ Int_t channelStatus = (Int_t)(fCaloPed->GetDeadMap(iSupMod))->GetBinContent(ieta,iphi);
+ if(channelStatus == AliCaloCalibPedestal::kHot || channelStatus == AliCaloCalibPedestal::kDead) {
+ AliDebug(2,Form("Tower from SM %d, ieta %d, iphi %d is BAD : status %d !!!",iSupMod,ieta,iphi, channelStatus));
+ return 0;
+ }
+ //Check if time is too large or too small, indication of a noisy channel, remove in this case
+ if(time > fTimeMax || time < fTimeMin) return 0;
+
+ fADCchannelECA = fCalibData->GetADCchannel (iSupMod,ieta,iphi);
+ fADCpedestalECA = fCalibData->GetADCpedestal(iSupMod,ieta,iphi);
+
+ return -fADCpedestalECA + amp * fADCchannelECA ;
+
+ }
+ else //Return energy with default parameters if calibration is not available
+ return -fADCpedestalECA + amp * fADCchannelECA ;
+
+}
+
+//____________________________________________________________________________
+void AliEMCALClusterizer::GetCalibrationParameters()
+{
+ // Set calibration parameters:
+ // if calibration database exists, they are read from database,
+ // otherwise, they are taken from digitizer.
+ //
+ // It is a user responsilibity to open CDB before reconstruction,
+ // for example:
+ // AliCDBStorage* storage = AliCDBManager::Instance()->GetStorage("local://CalibDB");
+
+ //Check if calibration is stored in data base
+
+ if(!fCalibData)
+ {
+ AliCDBEntry *entry = (AliCDBEntry*)
+ AliCDBManager::Instance()->Get("EMCAL/Calib/Data");
+ if (entry) fCalibData = (AliEMCALCalibData*) entry->GetObject();
+ }
+
+ if(!fCalibData)
+ AliFatal("Calibration parameters not found in CDB!");
+
+}
+
+//____________________________________________________________________________
+void AliEMCALClusterizer::GetCaloCalibPedestal()
+{
+ // Set calibration parameters:
+ // if calibration database exists, they are read from database,
+ // otherwise, they are taken from digitizer.
+ //
+ // It is a user responsilibity to open CDB before reconstruction,
+ // for example:
+ // AliCDBStorage* storage = AliCDBManager::Instance()->GetStorage("local://CalibDB");
+
+ //Check if calibration is stored in data base
+
+ if(!fCaloPed)
+ {
+ AliCDBEntry *entry = (AliCDBEntry*)
+ AliCDBManager::Instance()->Get("EMCAL/Calib/Pedestals");
+ if (entry) fCaloPed = (AliCaloCalibPedestal*) entry->GetObject();
+ }
+
+ if(!fCaloPed)
+ AliFatal("Pedestal info not found in CDB!");
+
+}
+
+//____________________________________________________________________________
+void AliEMCALClusterizer::Init()
+{
+ // Make all memory allocations which can not be done in default constructor.
+ // Attach the Clusterizer task to the list of EMCAL tasks
+
+ AliRunLoader *rl = AliRunLoader::Instance();
+ if (rl->GetAliRun() && rl->GetAliRun()->GetDetector("EMCAL"))
+ fGeom = dynamic_cast<AliEMCAL*>(rl->GetAliRun()->GetDetector("EMCAL"))->GetGeometry();
+ else
+ fGeom = AliEMCALGeometry::GetInstance(AliEMCALGeometry::GetDefaultGeometryName());
+
+ AliDebug(1,Form("geom %p",fGeom));
+
+ if(!gMinuit)
+ gMinuit = new TMinuit(100) ;
+
+}
+
+//____________________________________________________________________________
+void AliEMCALClusterizer::InitParameters()
+{
+ // Initializes the parameters for the Clusterizer
+ fNumberOfECAClusters = 0 ;
+ fCalibData = 0 ;
+ fCaloPed = 0 ;
+
+ const AliEMCALRecParam* recParam = AliEMCALReconstructor::GetRecParam();
+ if(!recParam) {
+ AliFatal("Reconstruction parameters for EMCAL not set!");
+ } else {
+ fECAClusteringThreshold = recParam->GetClusteringThreshold();
+ fECAW0 = recParam->GetW0();
+ fMinECut = recParam->GetMinECut();
+ fToUnfold = recParam->GetUnfold();
+ if(fToUnfold) AliWarning("Cluster Unfolding ON. Implementing only for eta=0 case!!!");
+ fECALocMaxCut = recParam->GetLocMaxCut();
+ fTimeCut = recParam->GetTimeCut();
+ fTimeMin = recParam->GetTimeMin();
+ fTimeMax = recParam->GetTimeMax();
+
+ AliDebug(1,Form("Reconstruction parameters: fECAClusteringThreshold=%.3f GeV, fECAW=%.3f, fMinECut=%.3f GeV, fToUnfold=%d, fECALocMaxCut=%.3f GeV, fTimeCut=%e s,fTimeMin=%e s,fTimeMax=%e s",
+ fECAClusteringThreshold,fECAW0,fMinECut,fToUnfold,fECALocMaxCut,fTimeCut, fTimeMin, fTimeMax));
+ }
+
+}
+
+
+//____________________________________________________________________________
+void AliEMCALClusterizer::Print(Option_t * /*option*/)const
+{
+ // Print clusterizer parameters
+
+ TString message("\n") ;
+
+ if( strcmp(GetName(), "") !=0 ){
+
+ // Print parameters
+
+ TString taskName(Version()) ;
+
+ printf("--------------- ");
+ printf("%s",taskName.Data()) ;
+ printf(" ");
+ printf("Clusterizing digits: ");
+ printf("\n ECA Local Maximum cut = %f", fECALocMaxCut);
+ printf("\n ECA Logarithmic weight = %f", fECAW0);
+ if(fToUnfold)
+ printf("\nUnfolding on\n");
+ else
+ printf("\nUnfolding off\n");
+
+ printf("------------------------------------------------------------------");
+ }
+ else
+ printf("AliEMCALClusterizer not initialized ") ;
+}
+
+//____________________________________________________________________________
+void AliEMCALClusterizer::PrintRecPoints(Option_t * option)
+{
+ // Prints list of RecPoints produced at the current pass of AliEMCALClusterizer
+ if(strstr(option,"deb")) {
+ printf("PrintRecPoints: Clusterization result:") ;
+
+ printf(" Found %d ECA Rec Points\n ",
+ fRecPoints->GetEntriesFast()) ;
+ }
+
+ if(strstr(option,"all")) {
+ if(strstr(option,"deb")) {
+ printf("\n-----------------------------------------------------------------------\n") ;
+ printf("Clusters in ECAL section\n") ;
+ printf("Index Ene(GeV) Multi Module GX GY GZ lX lY lZ Dispersion Lambda 1 Lambda 2 # of prim Primaries list\n") ;
+ }
+ Int_t index;
+ for (index = 0 ; index < fRecPoints->GetEntries() ; index++) {
+ AliEMCALRecPoint * rp = dynamic_cast<AliEMCALRecPoint * >(fRecPoints->At(index)) ;
+ TVector3 globalpos;
+ //rp->GetGlobalPosition(globalpos);
+ TVector3 localpos;
+ rp->GetLocalPosition(localpos);
+ Float_t lambda[2];
+ rp->GetElipsAxis(lambda);
+
+ Int_t nprimaries=0;
+ Int_t * primaries = rp->GetPrimaries(nprimaries);
+
+ if(strstr(option,"deb"))
+ printf("\n%6d %8.4f %3d %4.1f %4.1f %4.1f %4.1f %4.1f %4.1f %4.1f %4f %4f %2d : ",
+ rp->GetIndexInList(), rp->GetEnergy(), rp->GetMultiplicity(),
+ globalpos.X(), globalpos.Y(), globalpos.Z(), localpos.X(), localpos.Y(), localpos.Z(),
+ rp->GetDispersion(), lambda[0], lambda[1], nprimaries) ;
+ if(strstr(option,"deb")){
+ for (Int_t iprimary=0; iprimary<nprimaries; iprimary++) {
+ printf("%d ", primaries[iprimary] ) ;
+ }
+ }
+ }
+
+ if(strstr(option,"deb"))
+ printf("\n-----------------------------------------------------------------------\n");
+ }
+}
+
+//___________________________________________________________________
+void AliEMCALClusterizer::PrintRecoInfo()
+{
+ //print reco version
+ printf(" AliEMCALClusterizer::PrintRecoInfo() : version %s \n", Version() );
+
+}
+
//____________________________________________________________________________
void AliEMCALClusterizer::SetInput(TTree *digitsTree)
{
Int_t bufsize = 32000;
fTreeR->Branch("EMCALECARP", "TObjArray", &fRecPoints, bufsize, split);
}
+
+
+
// --- Standard library ---
// --- AliRoot header files ---
+class AliEMCALGeometry ;
+class AliEMCALCalibData ;
+class AliCaloCalibPedestal ;
class AliEMCALClusterizer : public TObject {
AliEMCALClusterizer() ; // default ctor
virtual ~AliEMCALClusterizer() ; // dtorEM
+ AliEMCALClusterizer(AliEMCALGeometry* geometry);
+ AliEMCALClusterizer(AliEMCALGeometry* geometry, AliEMCALCalibData * calib, AliCaloCalibPedestal * pedestal);
virtual void Digits2Clusters(Option_t *option) = 0;
- virtual Float_t GetTimeCut() const = 0;
-
- virtual void SetECAClusteringThreshold(Float_t) = 0;
- virtual void SetECALocalMaxCut(Float_t) = 0;
- virtual void SetECALogWeight(Float_t) = 0;
- virtual void SetTimeCut(Float_t) = 0;
- virtual void SetUnfolding(Bool_t) = 0;
-
- virtual const char * Version() const {Warning("Version", "Not Defined") ; return 0 ; }
+ virtual Float_t Calibrate(const Float_t amp, const Float_t time, const Int_t cellId) ; // Tranforms Amp to energy
+ virtual void Init() ;
+ virtual void InitParameters() ; //{ AliInfo("Overload this method."); }
+
+ //Get/Set reconstruction parameters
+ virtual void GetCalibrationParameters(void) ;
+ virtual void GetCaloCalibPedestal(void) ;
+ virtual void SetCalibrationParameters(AliEMCALCalibData * calib) { fCalibData = calib ; }
+ virtual void SetCaloCalibPedestal(AliCaloCalibPedestal * caloped) { fCaloPed = caloped ; }
+
+ virtual Float_t GetTimeMin() const { return fTimeMin ; }
+ virtual Float_t GetTimeMax() const { return fTimeMax ; }
+ virtual Float_t GetTimeCut() const { return fTimeCut ; }
+ virtual void GetNumberOfClustersFound(int numb )const { numb = fNumberOfECAClusters ;}
+ virtual Float_t GetECAClusteringThreshold() const { return fECAClusteringThreshold;}
+ virtual Float_t GetECALocalMaxCut() const { return fECALocMaxCut;}
+ virtual Float_t GetECALogWeight() const { return fECAW0;}
+ virtual Float_t GetMinECut() const { return fMinECut;}
+
+ virtual void SetTimeMin(Float_t t) { fTimeMin = t ;}
+ virtual void SetTimeMax(Float_t t) { fTimeMax = t ;}
+ virtual void SetTimeCut(Float_t t) { fTimeCut = t ;}
+ virtual void SetECAClusteringThreshold(Float_t th) { fECAClusteringThreshold = th ; }
+ virtual void SetMinECut(Float_t mine) { fMinECut = mine; }
+ virtual void SetECALocalMaxCut(Float_t cut) { fECALocMaxCut = cut ; }
+ virtual void SetECALogWeight(Float_t w) { fECAW0 = w ; }
+ virtual void SetUnfolding(Bool_t toUnfold = kTRUE ) {fToUnfold = toUnfold ;}
virtual void SetInput(TTree *digitsTree);
virtual void SetOutput(TTree *clustersTree);
- virtual void InitParameters() { AliInfo("Overload this method."); }
+
+ virtual void Print(Option_t * option)const ;
+ virtual void PrintRecPoints(Option_t * option);
+ virtual void PrintRecoInfo(); //*MENU*
+
+ virtual const char * Version() const {Warning("Version", "Not Defined") ; return 0 ; }
protected:
TClonesArray *fDigitsArr; // Array with EMCAL digits
TTree *fTreeR; // Tree with output clusters
TObjArray *fRecPoints; // Array with EMCAL clusters
-
+
+ AliEMCALGeometry * fGeom; //! pointer to geometry for utilities
+ AliEMCALCalibData * fCalibData ; //! Calibration database if aval
+ AliCaloCalibPedestal * fCaloPed ; //! Tower status map if aval
+
+ Float_t fADCchannelECA ; // width of one ADC channel for EC section (GeV)
+ Float_t fADCpedestalECA ; // pedestal of ADC for EC section (GeV)
+
+ Float_t fTimeMin ; // Minimum time of physical signal in a cell/digit
+ Float_t fTimeMax ; // Maximum time of physical signal in a cell/digit
+ Float_t fTimeCut ; // Maximum time difference between the digits inside EMC cluster
+
+ Bool_t fDefaultInit; //! Says if the task was created by defaut ctor (only parameters are initialized)
+ Bool_t fToUnfold ; // To perform unfolding
+ Int_t fNumberOfECAClusters ; // number of clusters found in EC section
+
+ Float_t fECAClusteringThreshold ; // minimum energy to seed a EC digit in a cluster
+ Float_t fECALocMaxCut ; // minimum energy difference to distinguish local maxima in a cluster
+ Float_t fECAW0 ; // logarithmic weight for the cluster center of gravity calculation
+ Float_t fMinECut; // Minimum energy for a digit to be a member of a cluster
+
private:
AliEMCALClusterizer(const AliEMCALClusterizer &); //copy ctor
AliEMCALClusterizer & operator = (const AliEMCALClusterizer &);
-
- ClassDef(AliEMCALClusterizer,1) // Clusterization algorithm class
+
+
+
+ ClassDef(AliEMCALClusterizer,2) // Clusterization algorithm class
} ;
#endif // AliEMCALCLUSTERIZER_H
--- /dev/null
+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * *
+ * Author: The ALICE Off-line Project. *
+ * Contributors are mentioned in the code where appropriate. *
+ * *
+ * Permission to use, copy, modify and distribute this software and its *
+ * documentation strictly for non-commercial purposes is hereby granted *
+ * without fee, provided that the above copyright notice appears in all *
+ * copies and that both the copyright notice and this permission notice *
+ * appear in the supporting documentation. The authors make no claims *
+ * about the suitability of this software for any purpose. It is *
+ * provided "as is" without express or implied warranty. *
+ **************************************************************************/
+
+// This class derives from AliEMCALClustrerizer but also keeps the API of AliEMCALClusterizerv1
+// Algorithm:
+// 1. peek the most energetic cell
+// 2. assign it as a center of the cluster and add cells surrounding it: 3x3, 5x5...
+// 3. remove the cells contributing to the cluster
+// 4. start from 1 for the remaining clusters
+// 5. cluster splitting (not implemented yet) - use the shape analysis to resolve the energy sharing
+// - for high energy clusters check the surrounding of the 3x3 clusters for extra energy
+// (merge 3x3 clusters and resolve the internal energy sharing - case for 2 clusters merged)
+// Use Case:
+// root [0] AliEMCALClusterizerNxN * cl = new AliEMCALClusterizerNxN("galice.root")
+// Warning in <TDatabasePDG::TDatabasePDG>: object already instantiated
+// //reads gAlice from header file "..."
+// root [1] cl->ExecuteTask()
+// //finds RecPoints in all events stored in galice.root
+// root [2] cl->SetDigitsBranch("digits2")
+// //sets another title for Digitis (input) branch
+// root [3] cl->SetRecPointsBranch("recp2")
+// //sets another title four output branches
+// root [4] cl->SetTowerLocalMaxCut(0.03)
+// //set clusterization parameters
+// root [5] cl->ExecuteTask("deb all time")
+// //once more finds RecPoints options are
+// // deb - print number of found rec points
+// // deb all - print number of found RecPoints and some their characteristics
+// // time - print benchmarking results
+
+// --- ROOT system ---
+#include <TMath.h>
+#include <TMinuit.h>
+#include <TTree.h>
+#include <TBenchmark.h>
+#include <TBrowser.h>
+#include <TROOT.h>
+#include <TClonesArray.h>
+
+// --- Standard library ---
+#include <cassert>
+
+// --- AliRoot header files ---
+#include "AliLog.h"
+#include "AliEMCALClusterizerNxN.h"
+#include "AliEMCALRecPoint.h"
+#include "AliEMCALDigit.h"
+#include "AliEMCALGeometry.h"
+#include "AliCaloCalibPedestal.h"
+#include "AliEMCALCalibData.h"
+#include "AliESDCaloCluster.h"
+
+ClassImp(AliEMCALClusterizerNxN)
+
+Bool_t AliEMCALClusterizerNxN::fgkIsInputCalibrated = kFALSE;
+
+//____________________________________________________________________________
+AliEMCALClusterizerNxN::AliEMCALClusterizerNxN()
+ : AliEMCALClusterizer()
+{
+ // ctor with the indication of the file where header Tree and digits Tree are stored
+}
+
+//____________________________________________________________________________
+AliEMCALClusterizerNxN::AliEMCALClusterizerNxN(AliEMCALGeometry* geometry)
+ : AliEMCALClusterizer(geometry)
+{
+ // ctor with the indication of the file where header Tree and digits Tree are stored
+ // use this contructor to avoid usage of Init() which uses runloader
+ // change needed by HLT - MP
+
+}
+
+//____________________________________________________________________________
+AliEMCALClusterizerNxN::AliEMCALClusterizerNxN(AliEMCALGeometry* geometry, AliEMCALCalibData * calib, AliCaloCalibPedestal * caloped)
+: AliEMCALClusterizer(geometry, calib, caloped)
+{
+ // ctor, geometry and calibration are initialized elsewhere.
+
+}
+
+
+//____________________________________________________________________________
+ AliEMCALClusterizerNxN::~AliEMCALClusterizerNxN()
+{
+ // dtor
+}
+
+
+//____________________________________________________________________________
+void AliEMCALClusterizerNxN::Digits2Clusters(Option_t * option)
+{
+ // Steering method to perform clusterization for the current event
+ // in AliEMCALLoader
+
+ if(strstr(option,"tim"))
+ gBenchmark->Start("EMCALClusterizer");
+
+ if(strstr(option,"print"))
+ Print("") ;
+
+ //Get calibration parameters from file or digitizer default values.
+ GetCalibrationParameters() ;
+
+ //Get dead channel map from file or digitizer default values.
+ GetCaloCalibPedestal() ;
+
+ fNumberOfECAClusters = 0;
+
+ MakeClusters() ; //only the real clusters
+
+ if(fToUnfold)
+ MakeUnfolding() ;
+
+ Int_t index ;
+
+ //Evaluate position, dispersion and other RecPoint properties for EC section
+ for(index = 0; index < fRecPoints->GetEntries(); index++)
+ {
+ dynamic_cast<AliEMCALRecPoint *>(fRecPoints->At(index))->EvalAll(fECAW0,fDigitsArr) ;
+ AliDebug(5, Form("MAX INDEX %d ", dynamic_cast<AliEMCALRecPoint *>(fRecPoints->At(index))->GetMaximalEnergyIndex()));
+ //For each rec.point set the distance to the nearest bad crystal
+ dynamic_cast<AliEMCALRecPoint *>(fRecPoints->At(index))->EvalDistanceToBadChannels(fCaloPed);
+ }
+
+ fRecPoints->Sort() ;
+
+ for(index = 0; index < fRecPoints->GetEntries(); index++)
+ {
+ (dynamic_cast<AliEMCALRecPoint *>(fRecPoints->At(index)))->SetIndexInList(index) ;
+ (dynamic_cast<AliEMCALRecPoint *>(fRecPoints->At(index)))->Print();
+ }
+
+ fTreeR->Fill();
+
+ if(strstr(option,"deb") || strstr(option,"all"))
+ PrintRecPoints(option) ;
+
+ AliDebug(1,Form("EMCAL Clusterizer found %d Rec Points",fRecPoints->GetEntriesFast()));
+
+ fRecPoints->Delete();
+
+ if(strstr(option,"tim")){
+ gBenchmark->Stop("EMCALClusterizer");
+ printf("Exec took %f seconds for Clusterizing",
+ gBenchmark->GetCpuTime("EMCALClusterizer"));
+ }
+}
+
+//____________________________________________________________________________
+Int_t AliEMCALClusterizerNxN::AreNeighbours(AliEMCALDigit * d1, AliEMCALDigit * d2, Bool_t & shared) const
+{
+ // Gives the neighbourness of two digits = 0 are not neighbour ; continue searching
+ // = 1 are neighbour
+ // = 2 is in different SM; continue searching
+ // In case it is in different SM, but same phi rack, check if neigbours at eta=0
+ // neighbours are defined as digits having at least a common side
+ // The order of d1 and d2 is important: first (d1) should be a digit already in a cluster
+ // which is compared to a digit (d2) not yet in a cluster
+
+ static Int_t nSupMod1=0, nModule1=0, nIphi1=0, nIeta1=0, iphi1=0, ieta1=0;
+ static Int_t nSupMod2=0, nModule2=0, nIphi2=0, nIeta2=0, iphi2=0, ieta2=0;
+ static Int_t rowdiff=0, coldiff=0;
+
+ shared = kFALSE;
+
+ fGeom->GetCellIndex(d1->GetId(), nSupMod1,nModule1,nIphi1,nIeta1);
+ fGeom->GetCellIndex(d2->GetId(), nSupMod2,nModule2,nIphi2,nIeta2);
+ fGeom->GetCellPhiEtaIndexInSModule(nSupMod1,nModule1,nIphi1,nIeta1, iphi1,ieta1);
+ fGeom->GetCellPhiEtaIndexInSModule(nSupMod2,nModule2,nIphi2,nIeta2, iphi2,ieta2);
+
+ //If different SM, check if they are in the same phi, then consider cells close to eta=0 as neighbours; May 2010
+ if(nSupMod1 != nSupMod2 )
+ {
+ //Check if the 2 SM are in the same PHI position (0,1), (2,3), ...
+ Float_t smPhi1 = fGeom->GetEMCGeometry()->GetPhiCenterOfSM(nSupMod1);
+ Float_t smPhi2 = fGeom->GetEMCGeometry()->GetPhiCenterOfSM(nSupMod2);
+
+ if(!TMath::AreEqualAbs(smPhi1, smPhi2, 1e-3)) return 2; //Not phi rack equal, not neighbours
+
+ // In case of a shared cluster, index of SM in C side, columns start at 48 and ends at 48*2
+ // C Side impair SM, nSupMod%2=1; A side pair SM nSupMod%2=0
+ if(nSupMod1%2) ieta1+=AliEMCALGeoParams::fgkEMCALCols;
+ else ieta2+=AliEMCALGeoParams::fgkEMCALCols;
+
+ shared = kTRUE; // maybe a shared cluster, we know this later, set it for the moment.
+
+ }//Different SM, same phi
+
+ rowdiff = TMath::Abs(iphi1 - iphi2);
+ coldiff = TMath::Abs(ieta1 - ieta2) ;
+
+ // neighbours +-1 in col and row
+ if ( TMath::Abs(coldiff) < 2 && TMath::Abs(rowdiff) < 2)
+ {
+
+ AliDebug(9, Form("AliEMCALClusterizerNxN::AreNeighbours(): id1=%d, (row %d, col %d) ; id2=%d, (row %d, col %d), shared %d \n",
+ d1->GetId(), iphi1,ieta1, d2->GetId(), iphi2,ieta2, shared));
+
+ return 1;
+ }//Neighbours
+ else
+ {
+ AliDebug(9, Form("NOT AliEMCALClusterizerNxN::AreNeighbours(): id1=%d, (row %d, col %d) ; id2=%d, (row %d, col %d), shared %d \n",
+ d1->GetId(), iphi1,ieta1, d2->GetId(), iphi2,ieta2, shared));
+ shared = kFALSE;
+ return 2 ;
+ }//Not neighbours
+}
+
+//____________________________________________________________________________
+void AliEMCALClusterizerNxN::MakeClusters()
+{
+ // Steering method to construct the clusters stored in a list of Reconstructed Points
+ // A cluster is defined as a list of neighbour digits
+ // Mar 03, 2007 by PAI
+
+ if (fGeom==0) AliFatal("Did not get geometry from EMCALLoader");
+
+ fRecPoints->Clear();
+
+ // Set up TObjArray with pointers to digits to work on
+ //TObjArray *digitsC = new TObjArray();
+ TObjArray digitsC;
+ TIter nextdigit(fDigitsArr);
+ AliEMCALDigit *digit = 0;
+ while ( (digit = dynamic_cast<AliEMCALDigit*>(nextdigit())) ) {
+ digitsC.AddLast(digit);
+ }
+
+ TIter nextdigitC(&digitsC);
+
+ AliDebug(1,Form("MakeClusters: Number of digits %d -> (e %f)\n",
+ fDigitsArr->GetEntries(),fMinECut));
+
+ Bool_t bDone = kFALSE;
+ while ( bDone != kTRUE )
+ {
+ //first sort the digits:
+ Int_t iMaxEnergyDigit = -1;
+ Float_t dMaxEnergyDigit = -1;
+ AliEMCALDigit *pMaxEnergyDigit = 0;
+ nextdigitC.Reset();
+ while ( (digit = dynamic_cast<AliEMCALDigit *>(nextdigitC())) )
+ { // scan over the list of digitsC
+ Float_t dEnergyCalibrated = Calibrate(digit->GetAmplitude(), digit->GetTime(),digit->GetId());
+ //AliDebug(5, Form("-> Digit ENERGY: %1.5f", dEnergyCalibrated));
+
+ //if(fGeom->CheckAbsCellId(digit->GetId()) && dEnergyCalibrated > fECAClusteringThreshold )
+ if(fGeom->CheckAbsCellId(digit->GetId()) && dEnergyCalibrated > 0.0) // no threshold!
+ {
+ if (dEnergyCalibrated > dMaxEnergyDigit)
+ {
+ dMaxEnergyDigit = dEnergyCalibrated;
+ iMaxEnergyDigit = digit->GetId();
+ pMaxEnergyDigit = digit;
+ }
+ }
+ }
+
+ if (iMaxEnergyDigit < 0 || digitsC.GetEntries() <= 0)
+ {
+ bDone = kTRUE;
+ continue;
+ }
+
+ AliDebug (2, Form("Max digit found: %1.2f AbsId: %d", dMaxEnergyDigit, iMaxEnergyDigit));
+ AliDebug(5, Form("Max Digit ENERGY: %1.5f", dMaxEnergyDigit));
+
+ // keep the candidate digits in a list
+ TList clusterDigitList;
+ clusterDigitList.SetOwner(kFALSE);
+ clusterDigitList.AddLast(pMaxEnergyDigit);
+
+ Double_t clusterCandidateEnergy = dMaxEnergyDigit;
+
+ // now loop over the resto of the digits and cluster into NxN cluster
+ // we do not actually cluster yet: we keep them in the list clusterDigitList
+ nextdigitC.Reset();
+ while ( (digit = dynamic_cast<AliEMCALDigit *>(nextdigitC())) )
+ { // scan over the list of digitsC
+ if (digit == pMaxEnergyDigit) continue;
+ Float_t dEnergyCalibrated = Calibrate(digit->GetAmplitude(), digit->GetTime(),digit->GetId());
+ AliDebug(5, Form("-> Digit ENERGY: %1.5f", dEnergyCalibrated));
+ //if(fGeom->CheckAbsCellId(digit->GetId()) && dEnergyCalibrated > fECAClusteringThreshold )
+ if(fGeom->CheckAbsCellId(digit->GetId()) && dEnergyCalibrated > 0.0 )
+ {
+ Float_t time = pMaxEnergyDigit->GetTime(); //Time or TimeR?
+ if(TMath::Abs(time - digit->GetTime()) > fTimeCut ) continue; //Time or TimeR?
+ Bool_t shared = kFALSE; //cluster shared by 2 SuperModules?
+ if (AreNeighbours(pMaxEnergyDigit, digit, shared) == 1) // call (digit,digitN) in THAT order !!!!!
+ {
+ clusterDigitList.AddLast(digit) ;
+ clusterCandidateEnergy += dEnergyCalibrated;
+ }
+ }
+ }// loop over the next digits
+
+ // start a cluster here only if a cluster energy is larger than clustering threshold
+ //if (clusterCandidateEnergy > 0.1)
+ AliDebug(5, Form("Clusterization threshold is %f MeV", fECAClusteringThreshold));
+ if (clusterCandidateEnergy > fECAClusteringThreshold)
+ {
+ if(fNumberOfECAClusters >= fRecPoints->GetSize()) fRecPoints->Expand(2*fNumberOfECAClusters+1) ;
+
+ AliEMCALRecPoint *recPoint = new AliEMCALRecPoint("") ;
+ fRecPoints->AddAt(recPoint, fNumberOfECAClusters) ;
+ recPoint = dynamic_cast<AliEMCALRecPoint *>(fRecPoints->At(fNumberOfECAClusters)) ;
+ fNumberOfECAClusters++ ;
+ recPoint->SetClusterType(AliESDCaloCluster::kEMCALClusterv1);
+
+ AliDebug(9, Form("Number of cells per cluster (max is 9!): %d", clusterDigitList.GetEntries()));
+ for (Int_t idig = 0; idig < clusterDigitList.GetEntries(); idig++)
+ {
+
+ digit = (AliEMCALDigit*)clusterDigitList.At(idig);
+ Float_t dEnergyCalibrated = Calibrate(digit->GetAmplitude(), digit->GetTime(),digit->GetId());
+ AliDebug(5, Form(" Adding digit %d", digit->GetId()));
+ // note: this way the sharing info is lost!
+ recPoint->AddDigit(*digit, dEnergyCalibrated, kFALSE) ; //Time or TimeR?
+ digitsC.Remove(digit);
+ }
+ }
+ else
+ {
+ // we do not want to start clustering in the same spot!
+ // but in this case we may NOT reuse this seed for another cluster!
+ // need a better bookeeping?
+ digitsC.Remove(pMaxEnergyDigit);
+ }
+
+ AliDebug (2, Form("Number of digits left: %d", digitsC.GetEntries()));
+ } // while ! done
+
+ //delete digitsC ; //nope we use an object
+
+ AliDebug(1,Form("total no of clusters %d from %d digits",fNumberOfECAClusters,fDigitsArr->GetEntriesFast()));
+}
+
+//____________________________________________________________________________
+void AliEMCALClusterizerNxN::MakeUnfolding()
+{
+ // Unfolds clusters using the shape of an ElectroMagnetic shower
+ // Performs unfolding of all clusters
+
+ if(fNumberOfECAClusters > 0){
+ if (fGeom==0)
+ AliFatal("Did not get geometry from EMCALLoader") ;
+ Int_t nModulesToUnfold = fGeom->GetNCells();
+
+ Int_t numberofNotUnfolded = fNumberOfECAClusters ;
+ Int_t index ;
+ for(index = 0 ; index < numberofNotUnfolded ; index++){
+
+ AliEMCALRecPoint * recPoint = dynamic_cast<AliEMCALRecPoint *>( fRecPoints->At(index) ) ;
+
+ TVector3 gpos;
+ Int_t absId;
+ recPoint->GetGlobalPosition(gpos);
+ fGeom->GetAbsCellIdFromEtaPhi(gpos.Eta(),gpos.Phi(),absId);
+ if(absId > nModulesToUnfold)
+ break ;
+
+ Int_t nMultipl = recPoint->GetMultiplicity() ;
+ AliEMCALDigit ** maxAt = new AliEMCALDigit*[nMultipl] ;
+ Float_t * maxAtEnergy = new Float_t[nMultipl] ;
+ Int_t nMax = recPoint->GetNumberOfLocalMax(maxAt, maxAtEnergy,fECALocMaxCut,fDigitsArr) ;
+
+ if( nMax > 1 ) { // if cluster is very flat (no pronounced maximum) then nMax = 0
+ //UnfoldCluster(recPoint, nMax, maxAt, maxAtEnergy) ;
+ fRecPoints->Remove(recPoint);
+ fRecPoints->Compress() ;
+ index-- ;
+ fNumberOfECAClusters-- ;
+ numberofNotUnfolded-- ;
+ }
+ else{
+ recPoint->SetNExMax(1) ; //Only one local maximum
+ }
+
+ delete[] maxAt ;
+ delete[] maxAtEnergy ;
+ }
+ }
+ // End of Unfolding of clusters
+}
+
+//____________________________________________________________________________
+Double_t AliEMCALClusterizerNxN::ShowerShape(Double_t x, Double_t y)
+{
+ // Shape of the shower
+ // If you change this function, change also the gradient evaluation in ChiSquare()
+ //
+ Double_t r = sqrt(x*x+y*y);
+ Double_t r133 = TMath::Power(r, 1.33) ;
+ Double_t r669 = TMath::Power(r, 6.69) ;
+ Double_t shape = TMath::Exp( -r133 * (1. / (1.57 + 0.0860 * r133) - 0.55 / (1 + 0.000563 * r669) ) ) ;
+ return shape ;
+}
+
+//____________________________________________________________________________
+void AliEMCALClusterizerNxN::UnfoldCluster(AliEMCALRecPoint * /*iniTower*/,
+ Int_t /*nMax*/,
+ AliEMCALDigit ** /*maxAt*/,
+ Float_t * /*maxAtEnergy*/)
+{
+ //
+ // Performs the unfolding of a cluster with nMax overlapping showers
+ //
+ AliWarning("Not implemented. To be.");
+}
+
+//___________________________________________________________________
+void AliEMCALClusterizerNxN::SetInputCalibrated(Bool_t val)
+{
+ //
+ // input is calibrated - the case when we run already on ESD
+ //
+ AliEMCALClusterizerNxN::fgkIsInputCalibrated = val;
+}
--- /dev/null
+#ifndef ALIEMCALCLUSTERIZERNXN_H
+#define ALIEMCALCLUSTERIZERNXN_H
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice */
+
+/* $Id: AliEMCALClusterizerNxN.h 41181 2010-05-12 13:58:06Z gconesab $ */
+
+//_________________________________________________________________________
+// This class derives from AliEMCALClustrerizer but also keeps the API of AliEMCALClusterizerv1
+// Algorithm:
+// 1. peek the most energetic cell
+// 2. assign it as a center of the cluster and add cells surrounding it: 3x3, 5x5...
+// 3. remove the cells contributing to the cluster
+// 4. start from 1 for the remaining clusters
+// 5. cluster splitting (not implemented yet) - use the shape analysis to resolve the energy sharing
+// - for high energy clusters check the surrounding of the 3x3 clusters for extra energy
+// (merge 3x3 clusters and resolve the internal energy sharing - case for 2 clusters merged)
+
+#include "AliEMCALClusterizer.h"
+class AliEMCALRecPoint ;
+class AliEMCALDigit ;
+
+class AliEMCALClusterizerNxN : public AliEMCALClusterizer {
+
+public:
+
+ AliEMCALClusterizerNxN() ;
+ AliEMCALClusterizerNxN(AliEMCALGeometry* geometry);
+ AliEMCALClusterizerNxN(AliEMCALGeometry* geometry, AliEMCALCalibData * calib, AliCaloCalibPedestal * pedestal);
+
+ virtual ~AliEMCALClusterizerNxN() ;
+
+ virtual Int_t AreNeighbours(AliEMCALDigit * d1, AliEMCALDigit * d2, Bool_t & shared)const ;
+ // Checks if digits are in neighbour cells
+
+ virtual void Digits2Clusters(Option_t *option); // Does the job
+
+
+ static Double_t ShowerShape(Double_t x, Double_t y) ; // Shape of EM shower used in unfolding;
+ //class member function (not object member function)
+ static void SetInputCalibrated(Bool_t val);
+
+ virtual const char * Version() const { return "clu-NxN" ; }
+
+
+protected:
+
+ virtual void MakeClusters();
+
+ static Bool_t fgkIsInputCalibrated; // to enable reclusterization from ESD cells
+
+private:
+ AliEMCALClusterizerNxN(const AliEMCALClusterizerNxN &); //copy ctor
+ AliEMCALClusterizerNxN & operator = (const AliEMCALClusterizerNxN &);
+
+
+ virtual void MakeUnfolding();
+ void UnfoldCluster(AliEMCALRecPoint * iniEmc, Int_t Nmax,
+ AliEMCALDigit ** maxAt,
+ Float_t * maxAtEnergy ); //Unfolds cluster using TMinuit package
+
+private:
+
+ ClassDef(AliEMCALClusterizerNxN,1) // Clusterizer implementation version 1
+
+};
+
+#endif // AliEMCALCLUSTERIZERNXN_H
/* $Id$ */
//-- Author: Yves Schutz (SUBATECH) & Dmitri Peressounko (SUBATECH & Kurchatov Institute)
-// August 2002 Yves Schutz: clone PHOS as closely as possible and intoduction
-// of new IO (à la PHOS)
-// Mar 2007, Aleksei Pavlinov - new algoritmh of pseudo clusters
+//-- Gustavo Conesa (LPSC-Grenoble), move common clusterizer functionalities to mother class
//////////////////////////////////////////////////////////////////////////////
// Clusterization class. Performs clusterization (collects neighbouring active cells) and
// unfolds the clusters having several local maxima.
// Results are stored in TreeR#, branches EMCALTowerRP (EMC recPoints),
// EMCALPreShoRP (CPV RecPoints) and AliEMCALClusterizer (Clusterizer with all
// parameters including input digits branch title, thresholds etc.)
-// This TTask is normally called from Reconstructioner, but can as well be used in
-// standalone mode.
-// Use Case:
-// root [0] AliEMCALClusterizerv1 * cl = new AliEMCALClusterizerv1("galice.root")
-// Warning in <TDatabasePDG::TDatabasePDG>: object already instantiated
-// //reads gAlice from header file "..."
-// root [1] cl->ExecuteTask()
-// //finds RecPoints in all events stored in galice.root
-// root [2] cl->SetDigitsBranch("digits2")
-// //sets another title for Digitis (input) branch
-// root [3] cl->SetRecPointsBranch("recp2")
-// //sets another title four output branches
-// root [4] cl->SetTowerLocalMaxCut(0.03)
-// //set clusterization parameters
-// root [5] cl->ExecuteTask("deb all time")
-// //once more finds RecPoints options are
-// // deb - print number of found rec points
-// // deb all - print number of found RecPoints and some their characteristics
-// // time - print benchmarking results
+//
// --- ROOT system ---
-#include <cassert>
-class TROOT;
-#include <TH1.h>
#include <TFile.h>
-class TFolder;
#include <TMath.h>
#include <TMinuit.h>
#include <TTree.h>
-class TSystem;
#include <TBenchmark.h>
#include <TBrowser.h>
#include <TROOT.h>
+#include <TList.h>
+#include <TClonesArray.h>
// --- Standard library ---
-
+#include <cassert>
// --- AliRoot header files ---
-#include "AliRunLoader.h"
-#include "AliRun.h"
-#include "AliESD.h"
+#include "AliLog.h"
#include "AliEMCALClusterizerv1.h"
#include "AliEMCALRecPoint.h"
#include "AliEMCALDigit.h"
-#include "AliEMCALDigitizer.h"
-#include "AliEMCAL.h"
#include "AliEMCALGeometry.h"
-#include "AliEMCALRecParam.h"
-#include "AliEMCALReconstructor.h"
-#include "AliCDBManager.h"
#include "AliCaloCalibPedestal.h"
#include "AliEMCALCalibData.h"
-class AliCDBStorage;
-#include "AliCDBEntry.h"
+#include "AliESDCaloCluster.h"
ClassImp(AliEMCALClusterizerv1)
//____________________________________________________________________________
-AliEMCALClusterizerv1::AliEMCALClusterizerv1()
- : AliEMCALClusterizer(),
- fGeom(0),
- fDefaultInit(kFALSE),
- fToUnfold(kFALSE),
- fNumberOfECAClusters(0),fCalibData(0),fCaloPed(0),
- fADCchannelECA(0.),fADCpedestalECA(0.),fECAClusteringThreshold(0.),fECALocMaxCut(0.),
- fECAW0(0.),fTimeCut(1.),fTimeMin(-1.),fTimeMax(1.),fMinECut(0.)
+AliEMCALClusterizerv1::AliEMCALClusterizerv1(): AliEMCALClusterizer()
{
// ctor with the indication of the file where header Tree and digits Tree are stored
//____________________________________________________________________________
AliEMCALClusterizerv1::AliEMCALClusterizerv1(AliEMCALGeometry* geometry)
- : AliEMCALClusterizer(),
- fGeom(geometry),
- fDefaultInit(kFALSE),
- fToUnfold(kFALSE),
- fNumberOfECAClusters(0),fCalibData(0), fCaloPed(0),
- fADCchannelECA(0.),fADCpedestalECA(0.),fECAClusteringThreshold(0.),fECALocMaxCut(0.),
- fECAW0(0.),fTimeCut(1.),fTimeMin(-1.),fTimeMax(1.),fMinECut(0.)
+ : AliEMCALClusterizer(geometry)
{
// ctor with the indication of the file where header Tree and digits Tree are stored
// use this contructor to avoid usage of Init() which uses runloader
// change needed by HLT - MP
- // Note for the future: the use on runloader should be avoided or optional at least
- // another way is to make Init virtual and protected at least such that the deriving classes can overload
- // Init() ;
- //
-
- if (!fGeom)
- {
- AliFatal("Geometry not initialized.");
- }
-
}
//____________________________________________________________________________
AliEMCALClusterizerv1::AliEMCALClusterizerv1(AliEMCALGeometry* geometry, AliEMCALCalibData * calib, AliCaloCalibPedestal * caloped)
-: AliEMCALClusterizer(),
-fGeom(geometry),
-fDefaultInit(kFALSE),
-fToUnfold(kFALSE),
-fNumberOfECAClusters(0),fCalibData(calib), fCaloPed(caloped),
-fADCchannelECA(0.),fADCpedestalECA(0.),fECAClusteringThreshold(0.),fECALocMaxCut(0.),
-fECAW0(0.),fTimeCut(1.),fTimeMin(-1.),fTimeMax(1.),fMinECut(0.)
+: AliEMCALClusterizer(geometry, calib, caloped)
{
// ctor, geometry and calibration are initialized elsewhere.
-
- if (!fGeom)
- AliFatal("Geometry not initialized.");
-
+
}
// dtor
}
-//____________________________________________________________________________
-Float_t AliEMCALClusterizerv1::Calibrate(const Float_t amp, const Float_t time, const Int_t absId)
-{
-
- // Convert digitized amplitude into energy.
- // Calibration parameters are taken from calibration data base for raw data,
- // or from digitizer parameters for simulated data.
-
- if(fCalibData){
-
- if (fGeom==0)
- AliFatal("Did not get geometry from EMCALLoader") ;
-
- Int_t iSupMod = -1;
- Int_t nModule = -1;
- Int_t nIphi = -1;
- Int_t nIeta = -1;
- Int_t iphi = -1;
- Int_t ieta = -1;
-
- Bool_t bCell = fGeom->GetCellIndex(absId, iSupMod, nModule, nIphi, nIeta) ;
- if(!bCell) {
- fGeom->PrintGeometry();
- Error("Calibrate()"," Wrong cell id number : %i", absId);
- assert(0);
- }
-
- fGeom->GetCellPhiEtaIndexInSModule(iSupMod,nModule,nIphi, nIeta,iphi,ieta);
-
- // Check if channel is bad (dead or hot), in this case return 0.
- // Gustavo: 15-12-09 In case of RAW data this selection is already done, but not in simulation.
- // for the moment keep it here but remember to do the selection at the sdigitizer level
- // and remove it from here
- Int_t channelStatus = (Int_t)(fCaloPed->GetDeadMap(iSupMod))->GetBinContent(ieta,iphi);
- if(channelStatus == AliCaloCalibPedestal::kHot || channelStatus == AliCaloCalibPedestal::kDead) {
- AliDebug(2,Form("Tower from SM %d, ieta %d, iphi %d is BAD : status %d !!!",iSupMod,ieta,iphi, channelStatus));
- return 0;
- }
- //Check if time is too large or too small, indication of a noisy channel, remove in this case
- if(time > fTimeMax || time < fTimeMin) return 0;
-
- fADCchannelECA = fCalibData->GetADCchannel (iSupMod,ieta,iphi);
- fADCpedestalECA = fCalibData->GetADCpedestal(iSupMod,ieta,iphi);
-
- return -fADCpedestalECA + amp * fADCchannelECA ;
-
- }
- else //Return energy with default parameters if calibration is not available
- return -fADCpedestalECA + amp * fADCchannelECA ;
-
-}
-
//____________________________________________________________________________
void AliEMCALClusterizerv1::Digits2Clusters(Option_t * option)
{
return kFALSE ;
}
for(index = 0; index < nPar; index++){
- Double_t err ;
- Double_t val ;
+ Double_t err = 0. ;
+ Double_t val = 0. ;
gMinuit->GetParameter(index, val, err) ; // Returns value and error of parameter index
fitparameters[index] = val ;
}
}
-//____________________________________________________________________________
-void AliEMCALClusterizerv1::GetCalibrationParameters()
-{
- // Set calibration parameters:
- // if calibration database exists, they are read from database,
- // otherwise, they are taken from digitizer.
- //
- // It is a user responsilibity to open CDB before reconstruction,
- // for example:
- // AliCDBStorage* storage = AliCDBManager::Instance()->GetStorage("local://CalibDB");
-
- //Check if calibration is stored in data base
-
- if(!fCalibData)
- {
- AliCDBEntry *entry = (AliCDBEntry*)
- AliCDBManager::Instance()->Get("EMCAL/Calib/Data");
- if (entry) fCalibData = (AliEMCALCalibData*) entry->GetObject();
- }
-
- if(!fCalibData)
- AliFatal("Calibration parameters not found in CDB!");
-
-}
-
-//____________________________________________________________________________
-void AliEMCALClusterizerv1::GetCaloCalibPedestal()
-{
- // Set calibration parameters:
- // if calibration database exists, they are read from database,
- // otherwise, they are taken from digitizer.
- //
- // It is a user responsilibity to open CDB before reconstruction,
- // for example:
- // AliCDBStorage* storage = AliCDBManager::Instance()->GetStorage("local://CalibDB");
-
- //Check if calibration is stored in data base
-
- if(!fCaloPed)
- {
- AliCDBEntry *entry = (AliCDBEntry*)
- AliCDBManager::Instance()->Get("EMCAL/Calib/Pedestals");
- if (entry) fCaloPed = (AliCaloCalibPedestal*) entry->GetObject();
- }
-
- if(!fCaloPed)
- AliFatal("Pedestal info not found in CDB!");
-
-}
-
-
-//____________________________________________________________________________
-void AliEMCALClusterizerv1::Init()
-{
- // Make all memory allocations which can not be done in default constructor.
- // Attach the Clusterizer task to the list of EMCAL tasks
-
- AliRunLoader *rl = AliRunLoader::Instance();
- if (rl->GetAliRun() && rl->GetAliRun()->GetDetector("EMCAL"))
- fGeom = dynamic_cast<AliEMCAL*>(rl->GetAliRun()->GetDetector("EMCAL"))->GetGeometry();
- else
- fGeom = AliEMCALGeometry::GetInstance(AliEMCALGeometry::GetDefaultGeometryName());
-
- AliDebug(1,Form("geom %p",fGeom));
-
- if(!gMinuit)
- gMinuit = new TMinuit(100) ;
-
-}
-
-//____________________________________________________________________________
-void AliEMCALClusterizerv1::InitParameters()
-{
- // Initializes the parameters for the Clusterizer
- fNumberOfECAClusters = 0;
-
- fCalibData = 0 ;
- fCaloPed = 0 ;
-
- const AliEMCALRecParam* recParam = AliEMCALReconstructor::GetRecParam();
- if(!recParam) {
- AliFatal("Reconstruction parameters for EMCAL not set!");
- } else {
- fECAClusteringThreshold = recParam->GetClusteringThreshold();
- fECAW0 = recParam->GetW0();
- fMinECut = recParam->GetMinECut();
- fToUnfold = recParam->GetUnfold();
- if(fToUnfold) AliWarning("Cluster Unfolding ON. Implementing only for eta=0 case!!!");
- fECALocMaxCut = recParam->GetLocMaxCut();
- fTimeCut = recParam->GetTimeCut();
- fTimeMin = recParam->GetTimeMin();
- fTimeMax = recParam->GetTimeMax();
-
- AliDebug(1,Form("Reconstruction parameters: fECAClusteringThreshold=%.3f GeV, fECAW=%.3f, fMinECut=%.3f GeV, fToUnfold=%d, fECALocMaxCut=%.3f GeV, fTimeCut=%e s,fTimeMin=%e s,fTimeMax=%e s",
- fECAClusteringThreshold,fECAW0,fMinECut,fToUnfold,fECALocMaxCut,fTimeCut, fTimeMin, fTimeMax));
- }
-
-}
-
//____________________________________________________________________________
Int_t AliEMCALClusterizerv1::AreNeighbours(AliEMCALDigit * d1, AliEMCALDigit * d2, Bool_t & shared) const
{
static Int_t nSupMod1=0, nModule1=0, nIphi1=0, nIeta1=0, iphi1=0, ieta1=0;
static Int_t nSupMod2=0, nModule2=0, nIphi2=0, nIeta2=0, iphi2=0, ieta2=0;
- static Int_t rowdiff, coldiff;
shared = kFALSE;
}//Different SM, same phi
- rowdiff = TMath::Abs(iphi1 - iphi2);
- coldiff = TMath::Abs(ieta1 - ieta2) ;
+ Int_t rowdiff = TMath::Abs(iphi1 - iphi2);
+ Int_t coldiff = TMath::Abs(ieta1 - ieta2) ;
// neighbours with at least common side; May 11, 2007
if ((coldiff==0 && TMath::Abs(rowdiff)==1) || (rowdiff==0 && TMath::Abs(coldiff)==1)) {
AliEMCALRecPoint * recPoint = dynamic_cast<AliEMCALRecPoint *>( fRecPoints->At(index) ) ;
TVector3 gpos;
- Int_t absId;
+ Int_t absId = -1;
recPoint->GetGlobalPosition(gpos);
fGeom->GetAbsCellIdFromEtaPhi(gpos.Eta(),gpos.Phi(),absId);
if(absId > nModulesToUnfold)
AliEMCALDigit * digit = 0 ;
Int_t * digitsList = iniTower->GetDigitsList() ;
-
- Int_t iparam ;
+
+ Int_t iparam = 0 ;
Int_t iDigit ;
for(iDigit = 0 ; iDigit < nDigits ; iDigit ++){
digit = dynamic_cast<AliEMCALDigit*>( fDigitsArr->At(digitsList[iDigit] ) ) ;
fGeom->RelPosCellInSModule(digit->GetId(), yDigit, xDigit, zDigit);
efit[iDigit] = 0;
- iparam = 0 ;
while(iparam < nPar ){
xpar = fitparameters[iparam] ;
zpar = fitparameters[iparam+1] ;
// to its contribution to efit
Float_t * energiesList = iniTower->GetEnergiesList() ;
- Float_t ratio ;
+ Float_t ratio = 0 ;
iparam = 0 ;
while(iparam < nPar ){
fNumberOfECAClusters++ ;
recPoint->SetNExMax((Int_t)nPar/3) ;
- Float_t eDigit ;
+ Float_t eDigit = 0. ;
for(iDigit = 0 ; iDigit < nDigits ; iDigit ++){
digit = dynamic_cast<AliEMCALDigit*>( fDigitsArr->At( digitsList[iDigit] ) ) ;
fGeom->RelPosCellInSModule(digit->GetId(), yDigit, xDigit, zDigit);
for(iparam = 0 ; iparam < nPar ; iparam++)
Grad[iparam] = 0 ; // Will evaluate gradient
- Double_t efit ;
+ Double_t efit = 0. ;
AliEMCALDigit * digit ;
Int_t iDigit ;
if(iflag == 2){ // calculate gradient
Int_t iParam = 0 ;
- efit = 0 ;
+ efit = 0. ;
while(iParam < nPar ){
Double_t dx = (xDigit - x[iParam]) ;
iParam++ ;
// Here we assume, that sigma = sqrt(E)
}
}
-//____________________________________________________________________________
-void AliEMCALClusterizerv1::Print(Option_t * /*option*/)const
-{
- // Print clusterizer parameters
-
- TString message("\n") ;
-
- if( strcmp(GetName(), "") !=0 ){
-
- // Print parameters
-
- TString taskName(Version()) ;
-
- printf("--------------- ");
- printf("%s",taskName.Data()) ;
- printf(" ");
- printf("Clusterizing digits: ");
- printf("\n ECA Local Maximum cut = %f", fECALocMaxCut);
- printf("\n ECA Logarithmic weight = %f", fECAW0);
- if(fToUnfold)
- printf("\nUnfolding on\n");
- else
- printf("\nUnfolding off\n");
-
- printf("------------------------------------------------------------------");
- }
- else
- printf("AliEMCALClusterizerv1 not initialized ") ;
-}
-
-//____________________________________________________________________________
-void AliEMCALClusterizerv1::PrintRecPoints(Option_t * option)
-{
- // Prints list of RecPoints produced at the current pass of AliEMCALClusterizer
- if(strstr(option,"deb")) {
- printf("PrintRecPoints: Clusterization result:") ;
-
- printf(" Found %d ECA Rec Points\n ",
- fRecPoints->GetEntriesFast()) ;
- }
-
- if(strstr(option,"all")) {
- if(strstr(option,"deb")) {
- printf("\n-----------------------------------------------------------------------\n") ;
- printf("Clusters in ECAL section\n") ;
- printf("Index Ene(GeV) Multi Module GX GY GZ lX lY lZ Dispersion Lambda 1 Lambda 2 # of prim Primaries list\n") ;
- }
- Int_t index =0;
-
- for (index = 0 ; index < fRecPoints->GetEntries() ; index++) {
- AliEMCALRecPoint * rp = dynamic_cast<AliEMCALRecPoint * >(fRecPoints->At(index)) ;
- TVector3 globalpos;
- //rp->GetGlobalPosition(globalpos);
- TVector3 localpos;
- rp->GetLocalPosition(localpos);
- Float_t lambda[2];
- rp->GetElipsAxis(lambda);
- Int_t * primaries;
- Int_t nprimaries;
- primaries = rp->GetPrimaries(nprimaries);
- if(strstr(option,"deb"))
- printf("\n%6d %8.4f %3d %4.1f %4.1f %4.1f %4.1f %4.1f %4.1f %4.1f %4f %4f %2d : ",
- rp->GetIndexInList(), rp->GetEnergy(), rp->GetMultiplicity(),
- globalpos.X(), globalpos.Y(), globalpos.Z(), localpos.X(), localpos.Y(), localpos.Z(),
- rp->GetDispersion(), lambda[0], lambda[1], nprimaries) ;
- if(strstr(option,"deb")){
- for (Int_t iprimary=0; iprimary<nprimaries; iprimary++) {
- printf("%d ", primaries[iprimary] ) ;
- }
- }
- }
-
- if(strstr(option,"deb"))
- printf("\n-----------------------------------------------------------------------\n");
- }
-}
-
-//___________________________________________________________________
-void AliEMCALClusterizerv1::PrintRecoInfo()
-{
- printf(" AliEMCALClusterizerv1::PrintRecoInfo() : version %s \n", Version() );
-
-}
// Implementation version 1 of the clusterization algorithm
// Performs clusterization (collects neighbouring active cells) and
// unfolding of the clusters with several local maxima.
-// results are stored in TreeR#, branches PHOSEmcRP (EMC recPoints),
-// PHOSCpvRP (CPV RecPoints) and AliPHOSClusterizer
+// results are stored in TreeR
//
//*-- Author: Yves Schutz (SUBATECH)
-// Modif:
-// August 2002 Yves Schutz: clone PHOS as closely as possible and intoduction
-// of new IO (à la PHOS)
+//-- Gustavo Conesa (LPSC-Grenoble), move common clusterizer functionalities to mother class
+
// --- ROOT system ---
// --- Standard library ---
// --- AliRoot header files ---
-
#include "AliEMCALClusterizer.h"
-class TH1F;
class AliEMCALRecPoint ;
class AliEMCALDigit ;
-class AliEMCALDigitizer ;
-class AliEMCALGeometry ;
-class AliEMCALCalibData ;
-class AliCaloCalibPedestal ;
class AliEMCALClusterizerv1 : public AliEMCALClusterizer {
virtual ~AliEMCALClusterizerv1() ;
virtual Int_t AreNeighbours(AliEMCALDigit * d1, AliEMCALDigit * d2, Bool_t & shared)const ;
- // Checks if digits are in neighbour cells
-
- virtual Float_t Calibrate(const Float_t amp, const Float_t time, const Int_t cellId) ; // Tranforms Amp to energy
-
- virtual void GetNumberOfClustersFound(int numb )const { numb = fNumberOfECAClusters ;}
- virtual Float_t GetECAClusteringThreshold() const { return fECAClusteringThreshold;}
- virtual Float_t GetECALocalMaxCut() const { return fECALocMaxCut;}
- virtual Float_t GetECALogWeight() const { return fECAW0;}
- virtual Float_t GetMinECut() const { return fMinECut;}
-
- virtual Float_t GetTimeCut() const { return fTimeCut ; }
- virtual Float_t GetTimeMin() const { return fTimeMin ; }
- virtual Float_t GetTimeMax() const { return fTimeMax ; }
-
+ // Checks if digits are in neighbour cells
virtual void Digits2Clusters(Option_t *option); // Does the job
- virtual void Print(Option_t * option)const ;
-
- virtual void SetECAClusteringThreshold(Float_t cluth) { fECAClusteringThreshold = cluth ; }
- virtual void SetMinECut(Float_t mine) { fMinECut = mine; }
- virtual void SetECALocalMaxCut(Float_t cut) { fECALocMaxCut = cut ; }
- virtual void SetECALogWeight(Float_t w) { fECAW0 = w ; }
- virtual void SetTimeCut(Float_t t) { fTimeCut = t ;}
- virtual void SetTimeMin(Float_t t) { fTimeMin = t ;}
- virtual void SetTimeMax(Float_t t) { fTimeMax = t ;}
-
- virtual void SetUnfolding(Bool_t toUnfold = kTRUE ) {fToUnfold = toUnfold ;}
static Double_t ShowerShape(Double_t x, Double_t y) ; // Shape of EM shower used in unfolding;
//class member function (not object member function)
static void UnfoldingChiSquare(Int_t & nPar, Double_t * Grad, Double_t & fret, Double_t * x, Int_t iflag) ;
// Chi^2 of the fit. Should be static to be passes to MINUIT
virtual const char * Version() const { return "clu-v1" ; }
- void PrintRecoInfo(); //*MENU*
- void SetCalibrationParameters(AliEMCALCalibData * calib) { fCalibData = calib ; }
- void SetCaloCalibPedestal(AliCaloCalibPedestal * caloped) { fCaloPed = caloped ; }
-
protected:
virtual void MakeClusters();
AliEMCALClusterizerv1(const AliEMCALClusterizerv1 &); //copy ctor
AliEMCALClusterizerv1 & operator = (const AliEMCALClusterizerv1 &);
- void GetCalibrationParameters(void) ;
- void GetCaloCalibPedestal(void) ;
Bool_t FindFit(AliEMCALRecPoint * emcRP, AliEMCALDigit ** MaxAt, const Float_t * maxAtEnergy,
Int_t NPar, Float_t * FitParametres) const; //Used in UnfoldClusters, calls TMinuit
- void Init() ;
- void InitParameters();
virtual void MakeUnfolding();
void UnfoldCluster(AliEMCALRecPoint * iniEmc, Int_t Nmax,
- AliEMCALDigit ** maxAt,
- Float_t * maxAtEnergy ); //Unfolds cluster using TMinuit package
- void PrintRecPoints(Option_t * option) ;
+ AliEMCALDigit ** maxAt,
+ Float_t * maxAtEnergy ); //Unfolds cluster using TMinuit package
private:
- AliEMCALGeometry* fGeom; //! pointer to geometry for utilities
-
- Bool_t fDefaultInit; //! Says if the task was created by defaut ctor (only parameters are initialized)
- Bool_t fToUnfold ; // To perform unfolding
- Int_t fNumberOfECAClusters ; // number of clusters found in EC section
-
- //Calibration parameters... to be replaced by database
-
- AliEMCALCalibData * fCalibData ; //! Calibration database if aval
- AliCaloCalibPedestal * fCaloPed ; //! Tower status map if aval
- Float_t fADCchannelECA ; // width of one ADC channel for EC section (GeV)
- Float_t fADCpedestalECA ; // pedestal of ADC for EC section (GeV)
-
- Float_t fECAClusteringThreshold ; // minimum energy to seed a EC digit in a cluster
- Float_t fECALocMaxCut ; // minimum energy difference to distinguish local maxima in a cluster
- Float_t fECAW0 ; // logarithmic weight for the cluster center of gravity calculation
- Float_t fTimeCut ; // Maximum time difference between the digits inside EMC cluster
- Float_t fTimeMin ; // Minimum time of physical signal in a cell/digiy
- Float_t fTimeMax ; // Maximum time of physical signal in a cell/digit
- Float_t fMinECut; // Minimum energy for a digit to be a member of a cluster
- ClassDef(AliEMCALClusterizerv1,9) // Clusterizer implementation version 1
+ ClassDef(AliEMCALClusterizerv1,10) // Clusterizer implementation version 1
};
ClassImp(AliEMCALRecParam)
- TObjArray* AliEMCALRecParam::fgkMaps =0; //ALTRO mappings
+TObjArray* AliEMCALRecParam::fgkMaps =0; //ALTRO mappings
AliEMCALRecParam::AliEMCALRecParam() :
AliDetectorRecoParam(),
fTimeCut(1.),// high value, accept all
fTimeMin(-1.),// small value, accept all
fTimeMax(1.),// high value, accept all//clustering
+ fClusterizerFlag(AliEMCALRecParam::kClusterizerv1),
fTrkCutX(6.0),
fTrkCutY(6.0),
fTrkCutZ(6.0),
// Pb Pb
// as a first step, all array elements are initialized to 0.0
- Int_t i, j;
+ Int_t i=0, j=0;
for (i = 0; i < 6; i++) {
for (j = 0; j < 6; j++) {
fGamma[i][j] = fPiZero[i][j] = fHadron[i][j] = 0.;
fTimeCut(rp.fTimeCut),
fTimeMin(rp.fTimeMin),
fTimeMax(rp.fTimeMax),//clustering
+ fClusterizerFlag(rp.fClusterizerFlag),
fTrkCutX(rp.fTrkCutX),
fTrkCutY(rp.fTrkCutY),
fTrkCutZ(rp.fTrkCutZ),
//copy constructor
//PID values
- Int_t i, j;
+ Int_t i=0, j=0;
for (i = 0; i < 6; i++) {
for (j = 0; j < 6; j++) {
fGamma[i][j] = rp.fGamma[i][j];
fTimeCut = rp.fTimeCut;
fTimeMax = rp.fTimeMax;
fTimeMin = rp.fTimeMin;//clustering
+ fClusterizerFlag = rp.fClusterizerFlag;
fTrkCutX = rp.fTrkCutX;
fTrkCutY = rp.fTrkCutY;
fTrkCutZ = rp.fTrkCutZ;
fFitLEDEvents = rp.fFitLEDEvents;//raw signal
//PID values
- Int_t i, j;
+ Int_t i=0, j=0;
for (i = 0; i < 6; i++) {
for (j = 0; j < 6; j++) {
fGamma[i][j] = rp.fGamma[i][j];
//PID parameters for pp implemented
// as a first step, all array elements are initialized to 0.0
- Int_t i, j;
+ Int_t i=0, j=0;
for (i = 0; i < 6; i++) {
for (j = 0; j < 6; j++) {
params->SetGamma(i,j,0.);
// if nothing is specified print all, if "reco" just clusterization/track matching
// if "pid", just PID parameters, if "raw", just raw utils parameters.
if(!strcmp("",opt) || !strcmp("reco",opt)){
+ AliInfo(Form("Clusterizer selected: %d", fClusterizerFlag));
AliInfo(Form("Clusterization parameters :\n fClusteringThreshold=%.3f,\n fW0=%.3f,\n fMinECut=%.3f,\n fUnfold=%d,\n fLocMaxCut=%.3f,\n fTimeCut=%2.1f ns\n fTimeMin=%2.1f ns\n fTimeMax=%2.1f ns\n",
fClusteringThreshold,fW0,fMinECut,fUnfold,fLocMaxCut,fTimeCut*1.e9,fTimeMin*1e9,fTimeMax*1e9));
class AliEMCALRecParam : public AliDetectorRecoParam
{
public:
+
+ enum AliEMCALClusterizerFlag
+ {
+ kClusterizerv1 = 0,
+ kClusterizerNxN = 1
+ };
AliEMCALRecParam() ;
AliEMCALRecParam(const AliEMCALRecParam& recParam);
Bool_t UseFALTRO() const {return fUseFALTRO; }
Bool_t FitLEDEvents() const {return fFitLEDEvents; }
- virtual void Print(Option_t * option="") const ;
+ virtual void Print(Option_t * option="") const ;
static AliEMCALRecParam* GetDefaultParameters();
static AliEMCALRecParam* GetLowFluxParam();
static AliEMCALRecParam* GetHighFluxParam();
static AliEMCALRecParam* GetCalibParam();
static AliEMCALRecParam* GetCosmicParam();
-
+
static const TObjArray* GetMappings();
+ void SetClusterizerFlag(Short_t val) { fClusterizerFlag = val; }
+ Short_t GetClusterizerFlag() const { return fClusterizerFlag; }
+
private:
//Clustering
Float_t fClusteringThreshold ; // Minimum energy to seed a EC digit in a cluster
Float_t fTimeCut ; // Maximum time of digits with respect to EMC cluster max.
Float_t fTimeMin ; // Minimum time of digits
Float_t fTimeMax ; // Maximum time of digits
+ Short_t fClusterizerFlag ; // Choice of the clusterizer; Default selection (v1) is zero
//PID (Guenole)
Double_t fGamma[6][6]; // Parameter to Compute PID for photons
static TObjArray* fgkMaps; // ALTRO mappings for RCU0..RCUX
- ClassDef(AliEMCALRecParam,12) // Reconstruction parameters
+ ClassDef(AliEMCALRecParam,13) // Reconstruction parameters
} ;
#include "AliEMCALRawUtils.h"
#include "AliEMCALDigit.h"
#include "AliEMCALClusterizerv1.h"
+#include "AliEMCALClusterizerNxN.h"
#include "AliEMCALRecPoint.h"
#include "AliEMCALPID.h"
#include "AliEMCALTrigger.h"
if (entry) fPedestalData = (AliCaloCalibPedestal*) entry->GetObject();
}
- if(!fPedestalData)
- AliFatal("Dead map not found in CDB!");
-
-
- //Init the clusterizer with geometry and calibration pointers, avoid doing it twice.
- fgClusterizer = new AliEMCALClusterizerv1(fGeom, fCalibData,fPedestalData);
+ if(!fPedestalData)
+ AliFatal("Dead map not found in CDB!");
+ InitClusterizer();
+
if(!fGeom) AliFatal(Form("Could not get geometry!"));
AliEMCALTriggerDCSConfigDB* dcsConfigDB = AliEMCALTriggerDCSConfigDB::Instance();
// fList = AliEMCALHistoUtilities::GetTriggersListOfHists(kTRUE);
// }
+//____________________________________________________________________________
+void AliEMCALReconstructor::InitClusterizer()
+{
+ //Init the clusterizer with geometry and calibration pointers, avoid doing it twice.
+
+ AliEMCALRecParam *recParam = NULL;
+ AliCDBEntry *entry = (AliCDBEntry*)
+ AliCDBManager::Instance()->Get("EMCAL/Calib/RecoParam");
+ //Get The reco param for the default event specie
+ if (entry)
+ recParam = (AliEMCALRecParam*)((TObjArray *) entry->GetObject())->At(0);
+
+ if(!recParam)
+ AliFatal("RecoParam not found in CDB!");
+
+ if (recParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerv1)
+ {
+ fgClusterizer = new AliEMCALClusterizerv1(fGeom, fCalibData,fPedestalData);
+ }
+ else
+ {
+ fgClusterizer = new AliEMCALClusterizerNxN(fGeom, fCalibData,fPedestalData);
+ }
+
+}
+
//____________________________________________________________________________
void AliEMCALReconstructor::Reconstruct(TTree* digitsTree, TTree* clustersTree) const
{
AliCodeTimerAuto("",0)
ReadDigitsArrayFromTree(digitsTree);
+
fgClusterizer->InitParameters();
fgClusterizer->SetOutput(clustersTree);
for (Int_t idig = 0 ; idig < nDigits ; idig++) {
const AliEMCALDigit * dig = (const AliEMCALDigit*)digits->At(idig);
if(dig->GetAmplitude() > 0 ){
- energy = (static_cast<AliEMCALClusterizerv1*> (fgClusterizer))->Calibrate(dig->GetAmplitude(),dig->GetTime(),dig->GetId()); //TimeR or Time?
+ energy = fgClusterizer->Calibrate(dig->GetAmplitude(),dig->GetTime(),dig->GetId()); //TimeR or Time?
if(energy > 0){ //Digits tagged as bad (dead, hot, not alive) are set to 0 in calibrate, remove them
emcCells.SetCell(idignew,dig->GetId(),energy, dig->GetTime());
idignew++;
ec->SetPosition(xyz);
ec->SetE(clust->GetEnergy());
- //Distance to the nearest bad crystal
- ec->SetDistanceToBadChannel(clust->GetDistanceToBadTower());
+ //Distance to the nearest bad crystal
+ ec->SetDistanceToBadChannel(clust->GetDistanceToBadTower());
ec->SetNCells(newCellMult);
//Change type of list from short to ushort
UShort_t *newAbsIdList = new UShort_t[newCellMult];
- Double_t *newFracList = new Double_t[newCellMult];
+ Double_t *newFracList = new Double_t[newCellMult];
for(Int_t i = 0; i < newCellMult ; i++) {
newAbsIdList[i]=absIdList[i];
newFracList[i]=fracList[i];
virtual ~AliEMCALReconstructor() ; //dtor
virtual void Init() {;}
-
+ virtual void InitClusterizer();
+
Bool_t Debug() const { return fDebug ; }
using AliReconstructor::FillESD;
AliTracker* CreateTracker () const
{return new AliEMCALTracker;}
using AliReconstructor::Reconstruct;
- virtual void Reconstruct(TTree* digitsTree, TTree* clustersTree) const;
+ virtual void Reconstruct(TTree* digitsTree, TTree* clustersTree) const ;
virtual Bool_t HasDigitConversion() const {return kTRUE;};
virtual void ConvertDigits(AliRawReader* rawReader, TTree* digitsTree) const;
#pragma link C++ class AliEMCALReconstructor+;
#pragma link C++ class AliEMCALClusterizer+;
#pragma link C++ class AliEMCALClusterizerv1+;
+#pragma link C++ class AliEMCALClusterizerNxN+;
#pragma link C++ class AliEMCALTrack+;
#pragma link C++ class AliEMCALTracker+;
#pragma link C++ class AliEMCALPID+;
AliEMCALReconstructor.cxx \
AliEMCALClusterizer.cxx \
AliEMCALClusterizerv1.cxx \
+AliEMCALClusterizerNxN.cxx \
AliEMCALTrack.cxx \
AliEMCALTracker.cxx \
AliEMCALPID.cxx \
# have to tune
EINCLUDE:=PYTHIA6 RAW EVGEN THijing VZERO
-EXPORT:=AliEMCALTrack.h
-
ifeq (win32gcc,$(ALICE_TARGET))
PACKSOFLAGS:= $(SOFLAGS) -L$(ALICE_ROOT)/lib/tgt_$(ALICE_TARGET) \
-lEMCALUtils -lEMCALbase -lEMCALsim -lSTEER -lESD -lCDB -lSTEERBase \
//AliEMCALRecParam *recParamDB = AliEMCALRecParam::GetDefaultParameters();
AliEMCALRecParam* params = AliEMCALRecParam::GetDefaultParameters();
+ //params->SetClusterizerFlag(AliEMCALRecParam::kClusterizerNxN);
+ params->SetClusterizerFlag(AliEMCALRecParam::kClusterizerv1);
//Clusterization
params->SetClusteringThreshold(0.5);
//AliEMCALRecParam *recParamDB = AliEMCALRecParam::GetLowFluxParam();
AliEMCALRecParam* params = AliEMCALRecParam::GetDefaultParameters();
+ //params->SetClusterizerFlag(AliEMCALRecParam::kClusterizerNxN);
+ params->SetClusterizerFlag(AliEMCALRecParam::kClusterizerv1);
+
params->SetClusteringThreshold(0.1); // 100 MeV
params->SetMinECut(0.01); //10 MeV
params->SetTimeCut(250e-9);//250 ns
// Author: Gustavo Conesa (INFN-LNF)
-void PrintEMCALRecParam(char * file = "Run0_999999999_v0_s0.root"){
+void PrintEMCALRecParam(char * file = "$ALICE_ROOT/OCDB/EMCAL/Calib/RecoParam/Run0_999999999_v0_s0.root"){
TFile * f = new TFile(file,"READ");
cout<<"================================================"<<endl;
-
-rparam->Print("");
-
+//rparam->Print("reco");//Print only clusterizer parameters
+//rparam->Print("pid");//Print only pid parameters
+//rparam->Print("raw");//Print only raw digitization parameters
+rparam->Print("");// Print all
}
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff