* 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] AliEMCALClusterizerFixedWindow * cl = new AliEMCALClusterizerFixedWindow("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
+// This class derives from AliEMCALClustrerizer
-#include "AliEMCALClusterizerFixedWindow.h"
-
-// --- ROOT system ---
+// --- Root ---
#include <TMath.h>
#include <TMinuit.h>
#include <TTree.h>
#include <TClonesArray.h>
#include <TH1I.h>
-// --- Standard library ---
-#include <cassert>
-//#include <iostream>
-//#include <fstream>
-
-// --- AliRoot header files ---
+// --- AliRoot ---
#include "AliLog.h"
#include "AliEMCALRecPoint.h"
#include "AliEMCALDigit.h"
#include "AliEMCALCalibData.h"
#include "AliESDCaloCluster.h"
#include "AliEMCALUnfolding.h"
-#include "AliEMCALFixedWindowClusterInfo.h"
+
+#include "AliEMCALClusterizerFixedWindow.h"
ClassImp(AliEMCALClusterizerFixedWindow)
-//____________________________________________________________________________
-AliEMCALClusterizerFixedWindow::AliEMCALClusterizerFixedWindow()
-: AliEMCALClusterizer(),
-nPhi(4),
-nEta(4),
-shiftPhi(2),
-shiftEta(2),
-fTRUshift(0),
-clusters_array(0),
-fClustersInfo(new AliEMCALFixedWindowClusterInfo("clustersInfo"))
+//__________________________________________________________________________________________
+AliEMCALClusterizerFixedWindow::AliEMCALClusterizerFixedWindow() :
+ AliEMCALClusterizer(),
+ fNphi(4),
+ fNeta(4),
+ fShiftPhi(2),
+ fShiftEta(2),
+ fTRUshift(0),
+ fNEtaDigitsSupMod(0),
+ fNPhiDigitsSupMod(0),
+ fNTRUPhi(0),
+ fNTRUEta(0),
+ fNEtaDigits(0),
+ fNPhiDigits(0),
+ fMaxShiftPhi(0),
+ fMaxShiftEta(0),
+ fNDigitsCluster(0),
+ fNClusEtaNoShift(0),
+ fNClusPhiNoShift(0),
+ fNClusters(0),
+ fNTotalClus(0),
+ fClustersArray(0),
+ fInitialized(0)
{
- // ctor with the indication of the file where header Tree and digits Tree are stored
-
+ // Constructor
}
-//____________________________________________________________________________
-AliEMCALClusterizerFixedWindow::AliEMCALClusterizerFixedWindow(AliEMCALGeometry* geometry)
-: AliEMCALClusterizer(geometry),
-nPhi(4),
-nEta(4),
-shiftPhi(2),
-shiftEta(2),
-fTRUshift(0),
-clusters_array(0),
-fClustersInfo(new AliEMCALFixedWindowClusterInfo("clustersInfo"))
+//__________________________________________________________________________________________
+AliEMCALClusterizerFixedWindow::AliEMCALClusterizerFixedWindow(AliEMCALGeometry* geometry) :
+ AliEMCALClusterizer(geometry),
+ fNphi(4),
+ fNeta(4),
+ fShiftPhi(2),
+ fShiftEta(2),
+ fTRUshift(0),
+ fNEtaDigitsSupMod(0),
+ fNPhiDigitsSupMod(0),
+ fNTRUPhi(0),
+ fNTRUEta(0),
+ fNEtaDigits(0),
+ fNPhiDigits(0),
+ fMaxShiftPhi(0),
+ fMaxShiftEta(0),
+ fNDigitsCluster(0),
+ fNClusEtaNoShift(0),
+ fNClusPhiNoShift(0),
+ fNClusters(0),
+ fNTotalClus(0),
+ fClustersArray(0),
+ fInitialized(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
+ // Constructor
}
-//____________________________________________________________________________
-AliEMCALClusterizerFixedWindow::AliEMCALClusterizerFixedWindow(AliEMCALGeometry* geometry, AliEMCALCalibData * calib, AliCaloCalibPedestal * caloped)
-: AliEMCALClusterizer(geometry, calib, caloped),
-nPhi(4),
-nEta(4),
-shiftPhi(2),
-shiftEta(2),
-fTRUshift(0),
-clusters_array(0),
-fClustersInfo(new AliEMCALFixedWindowClusterInfo("clustersInfo"))
+//__________________________________________________________________________________________
+AliEMCALClusterizerFixedWindow::AliEMCALClusterizerFixedWindow(AliEMCALGeometry* geometry, AliEMCALCalibData * calib, AliCaloCalibPedestal * caloped) :
+ AliEMCALClusterizer(geometry, calib, caloped),
+ fNphi(4),
+ fNeta(4),
+ fShiftPhi(2),
+ fShiftEta(2),
+ fTRUshift(0),
+ fNEtaDigitsSupMod(0),
+ fNPhiDigitsSupMod(0),
+ fNTRUPhi(0),
+ fNTRUEta(0),
+ fNEtaDigits(0),
+ fNPhiDigits(0),
+ fMaxShiftPhi(0),
+ fMaxShiftEta(0),
+ fNDigitsCluster(0),
+ fNClusEtaNoShift(0),
+ fNClusPhiNoShift(0),
+ fNClusters(0),
+ fNTotalClus(0),
+ fClustersArray(0),
+ fInitialized(0)
{
- // ctor, geometry and calibration are initialized elsewhere.
+ // Constructor
}
-//____________________________________________________________________________
+//__________________________________________________________________________________________
AliEMCALClusterizerFixedWindow::~AliEMCALClusterizerFixedWindow()
{
- // dtor
- delete fClustersInfo;
+ // Destructor
+
+ if (fClustersArray) {
+ for (Int_t i = 0; i < fNTotalClus; i++) {
+ if (fClustersArray[i]) {
+ delete[] fClustersArray[i];
+ fClustersArray[i] = 0;
+ }
+ delete[] fClustersArray;
+ fClustersArray = 0;
+ }
+ }
+
+}
+
+//__________________________________________________________________________________________
+void AliEMCALClusterizerFixedWindow::SetNphi (Int_t n)
+{
+ // Set fNphi; if clusterizer already initialized gives a warning and does nothing
+
+ if (fInitialized != 0)
+ AliWarning("Clusterizer already initialized. Unable to change the parameters.");
+ else
+ fNphi = n;
+}
+
+//__________________________________________________________________________________________
+void AliEMCALClusterizerFixedWindow::SetNeta (Int_t n)
+{
+ // Set fNeta; if clusterizer already initialized gives a warning and does nothing
+
+ if (fInitialized != 0)
+ AliWarning("Clusterizer already initialized. Unable to change the parameters.");
+ else
+ fNeta = n;
+}
+
+//__________________________________________________________________________________________
+void AliEMCALClusterizerFixedWindow::SetShiftPhi (Int_t s)
+{
+ // Set fShiftPhi; if clusterizer already initialized gives a warning and does nothing
+
+ if (fInitialized != 0)
+ AliWarning("Clusterizer already initialized. Unable to change the parameters.");
+ else
+ fShiftPhi = s;
+}
+
+//__________________________________________________________________________________________
+void AliEMCALClusterizerFixedWindow::SetShiftEta (Int_t s)
+{
+ // Set fShiftEta; if clusterizer already initialized gives a warning and does nothing
+
+ if (fInitialized != 0)
+ AliWarning("Clusterizer already initialized. Unable to change the parameters.");
+ else
+ fShiftEta = s;
+}
+
+//__________________________________________________________________________________________
+void AliEMCALClusterizerFixedWindow::SetTRUshift(Bool_t b)
+{
+ // Set fTRUshift; if clusterizer already initialized gives a warning and does nothing
+
+ if (fInitialized != 0)
+ AliWarning("Clusterizer already initialized. Unable to change the parameters.");
+ else
+ fTRUshift = b;
}
-//____________________________________________________________________________
+
+//__________________________________________________________________________________________
void AliEMCALClusterizerFixedWindow::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("");
+ // Steering method to perform clusterization for the current event
+
+ static Float_t cputime = 0;
+ static Float_t realtime = 0;
+
+ if (strstr(option,"tim"))
+ gBenchmark->Start("EMCALClusterizer");
- //Get calibration parameters from file or digitizer default values.
- GetCalibrationParameters();
+ if (strstr(option,"print"))
+ Print("");
- //Get dead channel map from file or digitizer default values.
- GetCaloCalibPedestal();
+ //Get calibration parameters from file or digitizer default values.
+ GetCalibrationParameters();
- MakeClusters(); //only the real clusters
+ //Get dead channel map from file or digitizer default values.
+ GetCaloCalibPedestal();
- if (fToUnfold) {
- fClusterUnfolding->SetInput(fNumberOfECAClusters,fRecPoints,fDigitsArr);
- fClusterUnfolding->MakeUnfolding();
- }
+ MakeClusters(); //only the real clusters
- //Evaluate position, dispersion and other RecPoint properties for EC section
- for (Int_t index = 0; index < fRecPoints->GetEntries(); index++) {
- AliEMCALRecPoint * rp = dynamic_cast<AliEMCALRecPoint *>(fRecPoints->At(index));
- if (rp) {
- rp->EvalAll(fECAW0,fDigitsArr,fJustClusters);
- AliDebug(5, Form("MAX INDEX %d ", rp->GetMaximalEnergyIndex()));
- //For each rec.point set the distance to the nearest bad crystal
- if (fCaloPed)
- rp->EvalDistanceToBadChannels(fCaloPed);
- }
- }
+ if (fToUnfold) {
+ fClusterUnfolding->SetInput(fNumberOfECAClusters,fRecPoints,fDigitsArr);
+ fClusterUnfolding->MakeUnfolding();
+ }
- //fRecPoints->Sort();
+ //Evaluate position, dispersion and other RecPoint properties for EC section
+ for (Int_t index = 0; index < fRecPoints->GetEntries(); index++) {
+ AliEMCALRecPoint * rp = dynamic_cast<AliEMCALRecPoint *>(fRecPoints->At(index));
+ if (rp) {
+ rp->EvalAll(fECAW0,fDigitsArr,fJustClusters);
+ AliDebug(5, Form("MAX INDEX %d ", rp->GetMaximalEnergyIndex()));
+ //For each rec.point set the distance to the nearest bad crystal
+ if (fCaloPed)
+ rp->EvalDistanceToBadChannels(fCaloPed);
+ }
+ }
+
+ fRecPoints->Sort();
- for (Int_t index = 0; index < fRecPoints->GetEntries(); index++) {
-
AliEMCALRecPoint *rp = dynamic_cast<AliEMCALRecPoint *>(fRecPoints->At(index));
- if (rp) {
- rp->SetIndexInList(index);
- }
- else AliFatal("RecPoint NULL!!");
- }
+ for (Int_t index = 0; index < fRecPoints->GetEntries(); index++) {
+ AliEMCALRecPoint *rp = dynamic_cast<AliEMCALRecPoint *>(fRecPoints->At(index));
+ if (rp) {
+ rp->SetIndexInList(index);
+ }
+ else AliFatal("RecPoint NULL!!");
+ }
- if (fTreeR)
- fTreeR->Fill();
+ if (fTreeR)
+ fTreeR->Fill();
- if (strstr(option,"deb") || strstr(option,"all"))
- PrintRecPoints(option);
+ if (strstr(option,"deb") || strstr(option,"all"))
+ PrintRecPoints(option);
- AliDebug(1,Form("EMCAL Clusterizer found %d Rec Points",fRecPoints->GetEntriesFast()));
+ AliDebug(1,Form("EMCAL Clusterizer found %d Rec Points",fRecPoints->GetEntriesFast()));
- if (strstr(option,"tim")) {
- gBenchmark->Stop("EMCALClusterizer");
- printf("Exec took %f seconds for Clusterizing",
- gBenchmark->GetCpuTime("EMCALClusterizer"));
- }
+ if (strstr(option,"tim")) {
+ gBenchmark->Stop("EMCALClusterizer");
+ Printf("Exec took %f CPU time (%f real time) for clusterizing",
+ gBenchmark->GetCpuTime("EMCALClusterizer")-cputime,gBenchmark->GetRealTime("EMCALClusterizer")-realtime);
+ cputime = gBenchmark->GetCpuTime("EMCALClusterizer");
+ realtime = gBenchmark->GetRealTime("EMCALClusterizer");
+ }
}
-//____________________________________________________________________________
-void AliEMCALClusterizerFixedWindow::MakeClusters()
+//__________________________________________________________________________________________
+void AliEMCALClusterizerFixedWindow::ExecOnce()
{
- // Make clusters
-
- if (fGeom == 0)
- AliFatal("Did not get geometry from EMCALLoader");
-
- fNumberOfECAClusters = 0;
- fRecPoints->Delete();
+ // Initialize clusterizer.
- if (fClustersInfo->GetLastElementId() > 0)
- fClustersInfo->Clear();
-
- Int_t nSupMod=0, nModule=0, nIphi=0, nIeta=0, iphi=0, ieta=0;
+ fInitialized = -1;
+
+ if (!fGeom) {
+ AliError("Did not get geometry!");
+ return;
+ }
- // Defining geometry and clusterization parameter
- Int_t nEtaDigitsSupMod = fGeom->GetNEta() * fGeom->GetNETAdiv(); // always 48?;
- Int_t nPhiDigitsSupMod = fGeom->GetNPhi() * fGeom->GetNPHIdiv(); // always 24?;
+ // Defining geometry and clusterization parameter
+ fNEtaDigitsSupMod = fGeom->GetNEta() * fGeom->GetNETAdiv(); // always 48?;
+ fNPhiDigitsSupMod = fGeom->GetNPhi() * fGeom->GetNPHIdiv(); // always 24?;
- Int_t nTRUPhi = 1;
- Int_t nTRUEta = 1;
+ fNTRUPhi = 1;
+ fNTRUEta = 1;
- Int_t nEtaDigits = nEtaDigitsSupMod * fGeom->GetNumberOfSuperModules() / fGeom->GetNPhiSuperModule();
- Int_t nPhiDigits = nPhiDigitsSupMod * fGeom->GetNPhiSuperModule();
+ fNEtaDigits = fNEtaDigitsSupMod * fGeom->GetNumberOfSuperModules() / fGeom->GetNPhiSuperModule();
+ fNPhiDigits = fNPhiDigitsSupMod * fGeom->GetNPhiSuperModule();
- if (fTRUshift)
- {
- nTRUPhi = fGeom->GetNPhiSuperModule() * 3;
- nTRUEta = fGeom->GetNumberOfSuperModules() / fGeom->GetNPhiSuperModule();
- nEtaDigits /= nTRUEta;
- nPhiDigits /= nTRUPhi;
+ if (fTRUshift){
+ fNTRUPhi = fGeom->GetNPhiSuperModule() * 3;
+ fNTRUEta = fGeom->GetNumberOfSuperModules() / fGeom->GetNPhiSuperModule();
+ fNEtaDigits /= fNTRUEta;
+ fNPhiDigits /= fNTRUPhi;
}
// Check if clusterizer parameter are compatible with calorimeter geometry
- if (nEtaDigits < nEta)
- {
- AliFatal(Form("Error: nEta = %d is greater than nEtaDigits = %d.",nEta,nEtaDigits));
- return;
- }
- if (nPhiDigits < nPhi)
- {
- AliFatal(Form("Error: nPhi = %d is greater than nPhiDigits = %d.",nPhi,nPhiDigits));
- return;
- }
- if (nEtaDigits % shiftEta != 0)
- {
- AliFatal(Form("Error: shiftEta = %d is such that clusters cannot slide the whole calorimeter (nEtaDigits = %d).",shiftEta,nEtaDigits));
- return;
- }
- if (nPhiDigits % shiftPhi != 0)
- {
- AliFatal(Form("Error: shiftPhi = %d is such that clusters cannot slide the whole calorimeter (nPhiDigits = %d).",shiftPhi,nPhiDigits));
- return;
- }
- if (nEta % shiftEta != 0)
- {
- AliFatal(Form("Error: shiftEta = %d is not divisor of nEta = %d.",shiftEta,nEta));
- return;
- }
- if (nPhi % shiftPhi != 0)
- {
- AliFatal(Form("Error: shiftPhi = %d is not divisor of nPhi = %d).",shiftPhi,nPhi));
- return;
- }
-
- Int_t maxiShiftPhi = nPhi / shiftPhi;
- Int_t maxiShiftEta = nEta / shiftEta;
-
- Int_t nDigitsCluster = nPhi * nEta;
-
- Int_t nClusEtaNoShift = nEtaDigits / nEta;
- Int_t nClusPhiNoShift = nPhiDigits / nPhi;
-
- Int_t nClusters = nClusEtaNoShift * nClusPhiNoShift * nTRUEta * nTRUPhi;
+ if (fNEtaDigits < fNeta){
+ AliError(Form("Error: fNeta = %d is greater than nEtaDigits = %d.",fNeta,fNEtaDigits));
+ return;
+ }
+ if (fNPhiDigits < fNphi){
+ AliError(Form("Error: fNphi = %d is greater than nPhiDigits = %d.",fNphi,fNPhiDigits));
+ return;
+ }
+ if (fNEtaDigits % fShiftEta != 0){
+ AliError(Form("Error: fShiftEta = %d is such that clusters cannot slide the whole calorimeter (nEtaDigits = %d).",fShiftEta,fNEtaDigits));
+ return;
+ }
+ if (fNPhiDigits % fShiftPhi != 0){
+ AliError(Form("Error: fShiftPhi = %d is such that clusters cannot slide the whole calorimeter (nPhiDigits = %d).",fShiftPhi,fNPhiDigits));
+ return;
+ }
+ if (fNeta % fShiftEta != 0){
+ AliError(Form("Error: fShiftEta = %d is not divisor of fNeta = %d.",fShiftEta,fNeta));
+ return;
+ }
+ if (fNphi % fShiftPhi != 0){
+ AliError(Form("Error: fShiftPhi = %d is not divisor of fNphi = %d).",fShiftPhi,fNphi));
+ return;
+ }
- Int_t nTotalClus = nClusters * maxiShiftEta * maxiShiftPhi;
+ fMaxShiftPhi = fNphi / fShiftPhi;
+ fMaxShiftEta = fNeta / fShiftEta;
- if (!clusters_array)
- {
- clusters_array = new AliEMCALDigit**[nTotalClus];
- for (Int_t i = 0; i < nTotalClus; i++)
- {
- clusters_array[i] = NULL;
+ fNClusEtaNoShift = fNEtaDigits / fNeta;
+ fNClusPhiNoShift = fNPhiDigits / fNphi;
+
+ fNClusters = fNClusEtaNoShift * fNClusPhiNoShift * fNTRUEta * fNTRUPhi;
+
+ fNTotalClus = fNClusters * fMaxShiftEta * fMaxShiftPhi;
+
+ fNDigitsCluster = fNphi * fNeta;
+
+ if (fClustersArray) {
+ for (Int_t i = 0; i < fNTotalClus; i++) {
+ if (fClustersArray[i]) {
+ delete[] fClustersArray[i];
+ fClustersArray[i] = 0;
+ }
+ delete[] fClustersArray;
+ fClustersArray = 0;
+ }
+ }
+
+ fClustersArray = new AliEMCALDigit**[fNTotalClus];
+ for (Int_t i = 0; i < fNTotalClus; i++) {
+ fClustersArray[i] = new AliEMCALDigit*[fNDigitsCluster];
+ for (Int_t j = 0; j < fNDigitsCluster; j++) {
+ fClustersArray[i][j] = 0;
}
}
+
+ AliDebug(1,Form("****ExecOnce*****\n"
+ "fNphi = %d, fNeta = %d, fShiftPhi = %d, fShiftEta = %d, fTRUshift = %d\n"
+ "fNEtaDigitsSupMod = %d, fNPhiDigitsSupMod = %d, fNTRUPhi = %d, fNTRUEta = %d, fNEtaDigits = %d, fNPhiDigits = %d\n"
+ "fMaxShiftPhi = %d, fMaxShiftEta = %d, fNDigitsCluster = %d, fNClusEtaNoShift = %d, fNClusPhiNoShift = %d\n"
+ "fNClusters = %d, fNTotalClus = %d\n",
+ fNphi,fNeta,fShiftPhi,fShiftEta,fTRUshift,
+ fNEtaDigitsSupMod,fNPhiDigitsSupMod,fNTRUPhi,fNTRUEta,fNEtaDigits,fNPhiDigits,
+ fMaxShiftPhi,fMaxShiftEta,fNDigitsCluster,fNClusEtaNoShift,fNClusPhiNoShift,
+ fNClusters,fNTotalClus));
+
+ fInitialized = 1;
+}
- AliEMCALDigit *digit = 0;
+//__________________________________________________________________________________________
+void AliEMCALClusterizerFixedWindow::MakeClusters()
+{
+ // Make clusters.
+
+ fNumberOfECAClusters = 0;
+ fRecPoints->Delete();
+
+ if (fInitialized == 0)
+ ExecOnce();
+
+ if (fInitialized == -1) {
+ AliError(Form("%s: error initializing the clusterizer. No clusterization will be performed.",GetName()));
+ return;
+ }
+
+ // Set up TObjArray with pointers to digits to work on calibrated digits
+ TObjArray *digitsC = new TObjArray();
+ AliEMCALDigit *digit;
+ Float_t dEnergyCalibrated = 0.0, ehs = 0.0, time = 0.0;
+ TIter nextdigit(fDigitsArr);
+ while ((digit = dynamic_cast<AliEMCALDigit*>(nextdigit()))) { // calibrate and clean up digits
+ dEnergyCalibrated = digit->GetAmplitude();
+ time = digit->GetTime();
+ Calibrate(dEnergyCalibrated, time, digit->GetId());
+ digit->SetCalibAmp(dEnergyCalibrated);
+ digit->SetTime(time);
+ if (dEnergyCalibrated < fMinECut || time > fTimeMax || time < fTimeMin) {
+ continue;
+ }
+ else if (!fGeom->CheckAbsCellId(digit->GetId())) {
+ continue;
+ }
+ else {
+ ehs += dEnergyCalibrated;
+ digitsC->AddLast(digit);
+ }
+ }
- for (Int_t ishiftPhi = 0; ishiftPhi < maxiShiftPhi; ishiftPhi++)
- {
- Int_t nClusPhi = (nPhiDigits - shiftPhi * ishiftPhi) / nPhi;
+ AliDebug(1,Form("MakeClusters: Number of digits %d -> (e %f), ehs %f\n",
+ fDigitsArr->GetEntries(),fMinECut,ehs));
+
+ Int_t nSupMod=0, nModule=0, nIphi=0, nIeta=0;
+ Int_t iphi=0, ieta=0; // cell eta-phi indexes in SM
+
+ for (Int_t ishiftPhi = 0; ishiftPhi < fMaxShiftPhi; ishiftPhi++){
+ Int_t nClusPhi = (fNPhiDigits - fShiftPhi * ishiftPhi) / fNphi;
- for (Int_t ishiftEta = 0; ishiftEta < maxiShiftEta; ishiftEta++)
- {
-
- Int_t nClusEta = (nEtaDigits - shiftEta * ishiftEta) / nEta;
-
- Int_t iTotalClus = nClusters * (ishiftPhi * maxiShiftEta + ishiftEta);
+ for (Int_t ishiftEta = 0; ishiftEta < fMaxShiftEta; ishiftEta++) {
- TIter nextdigit(fDigitsArr);
+ Int_t nClusEta = (fNEtaDigits - fShiftEta * ishiftEta) / fNeta;
- nextdigit.Reset();
+ Int_t iTotalClus = fNClusters * (ishiftPhi * fMaxShiftEta + ishiftEta);
- while (digit = static_cast<AliEMCALDigit*>(nextdigit()))
- {
+ TIter nextdigitC(digitsC);
+ while ((digit = dynamic_cast<AliEMCALDigit*>(nextdigitC()))) { // scan over the list of digitsC
+
fGeom->GetCellIndex (digit->GetId(), nSupMod, nModule, nIphi, nIeta);
fGeom->GetCellPhiEtaIndexInSModule (nSupMod, nModule, nIphi, nIeta, iphi, ieta);
- Int_t iphi_eff = iphi - shiftPhi * ishiftPhi + nPhiDigitsSupMod * (nSupMod / 2); // N supermodules along phi
+ Int_t iphi_eff = iphi - fShiftPhi * ishiftPhi + fNPhiDigitsSupMod * (nSupMod / 2); // N supermodules along phi
- Int_t iTRUphi = iphi_eff / nPhiDigits;
+ Int_t iTRUphi = iphi_eff / fNPhiDigits;
- iphi_eff -= iTRUphi * nPhiDigits;
+ iphi_eff -= iTRUphi * fNPhiDigits;
- Int_t iClusPhi = iphi_eff / nPhi;
+ Int_t iClusPhi = iphi_eff / fNphi;
if (iphi_eff < 0 || iClusPhi >= nClusPhi)
continue;
- Int_t ieta_eff = ieta - shiftEta * ishiftEta + nEtaDigitsSupMod * (nSupMod % 2); // 2 supermodules along eta
+ Int_t ieta_eff = ieta - fShiftEta * ishiftEta + fNEtaDigitsSupMod * (nSupMod % 2); // 2 supermodules along eta
- Int_t iTRUeta = ieta_eff / nEtaDigits;
+ Int_t iTRUeta = ieta_eff / fNEtaDigits;
- ieta_eff -= iTRUeta * nEtaDigits;
+ ieta_eff -= iTRUeta * fNEtaDigits;
- Int_t iClusEta = ieta_eff / nEta;
+ Int_t iClusEta = ieta_eff / fNeta;
if (ieta_eff < 0 || iClusEta >= nClusEta)
continue;
- iphi_eff += iTRUphi * nPhiDigits;
- iClusPhi = iphi_eff / nPhi;
+ iphi_eff += iTRUphi * fNPhiDigits;
+ iClusPhi = iphi_eff / fNphi;
- ieta_eff += iTRUeta * nEtaDigits;
- iClusEta = ieta_eff / nEta;
+ ieta_eff += iTRUeta * fNEtaDigits;
+ iClusEta = ieta_eff / fNeta;
- Int_t iCluster = iClusPhi + iClusEta * nClusPhiNoShift * nTRUPhi;
- Int_t iDigit = iphi_eff % nPhi + (ieta_eff % nEta) * nPhi;
+ Int_t iCluster = iClusPhi + iClusEta * fNClusPhiNoShift * fNTRUPhi;
+ Int_t iDigit = iphi_eff % fNphi + (ieta_eff % fNeta) * fNphi;
-
- if (iCluster >= nClusters)
- {
- AliFatal(Form("ERROR: iCluster out of range! iCluster = %d, nClusters = %d", iCluster, nClusters));
+ if (iCluster >= fNClusters){
+ AliError(Form("iCluster out of range! iCluster = %d, fNClusters = %d (should never happen...)", iCluster, fNClusters));
return;
}
iCluster += iTotalClus;
-
- if (clusters_array[iCluster] == NULL)
- {
- fNumberOfECAClusters++;
- clusters_array[iCluster] = new AliEMCALDigit*[nDigitsCluster];
- for (Int_t i = 0; i < nDigitsCluster; i++)
- {
- clusters_array[iCluster][i] = NULL;
- }
-
- fClustersInfo->Add(iCluster, -1, iClusEta, iClusPhi);
+
+ if (iCluster >= fNTotalClus){
+ AliError(Form("iCluster out of range! iCluster = %d, fNTotalClus = %d (should never happen...)", iCluster, fNTotalClus));
+ return;
}
-
- if (clusters_array[iCluster][iDigit] != NULL)
- {
- AliFatal("ERROR: digit already added!");
+
+ if (iDigit >= fNDigitsCluster){
+ AliError(Form("iDigit out of range! iDigit = %d, fNDigitsCluster = %d (should never happen...)", iDigit, fNDigitsCluster));
+ return;
+ }
+
+ if (fClustersArray[iCluster][iDigit] != 0){
+ AliError("Digit already added! (should never happen...)");
return;
}
- clusters_array[iCluster][iDigit] = digit;
+ fClustersArray[iCluster][iDigit] = digit;
} // loop on digit
} // loop on eta shift
- } // loop on phi shift
-
- Int_t iRecPoint = 0;
- for (Int_t iCluster = 0; iCluster < nTotalClus; iCluster++)
- {
-
- if (clusters_array[iCluster] == NULL) continue;
-
- (*fRecPoints)[iRecPoint] = new AliEMCALRecPoint("");
- AliEMCALRecPoint *recPoint = dynamic_cast<AliEMCALRecPoint*> (fRecPoints->At(iRecPoint));
-
- if (recPoint)
- {
- if (fClustersInfo->ContainsIndex(iRecPoint))
- AliFatal(Form("ERROR: index present already, %d", iRecPoint));
-
- fClustersInfo->SetIndexFromId(iCluster, iRecPoint);
-
- iRecPoint++;
- recPoint->SetClusterType(AliVCluster::kEMCALClusterv1);
- // note: this way the sharing info is lost!
- for (Int_t iDigit = 0; iDigit < nDigitsCluster; iDigit++)
- {
- if (clusters_array[iCluster][iDigit] == NULL) continue;
- digit = clusters_array[iCluster][iDigit];
- Float_t dEnergyCalibrated = digit->GetAmplitude();
- Float_t time = digit->GetTime();
- Calibrate(dEnergyCalibrated,time,digit->GetId());
- digit->SetCalibAmp(dEnergyCalibrated);
- recPoint->AddDigit(*digit, dEnergyCalibrated, kFALSE); //Time or TimeR?
- clusters_array[iCluster][iDigit] = NULL;
- }
- }
-
- delete[] clusters_array[iCluster];
- clusters_array[iCluster] = NULL;
- }
+ } // loop on phi shift
- AliDebug(1,Form("MakeClusters: Number of digits %d -> (e %f)\n",
- fDigitsArr->GetEntries(),fMinECut));
-
- AliDebug(1,Form("total no of clusters %d from %d digits",fNumberOfECAClusters,fDigitsArr->GetEntriesFast()));
+ AliEMCALRecPoint *recPoint = 0;
+ Bool_t recPointOk = kFALSE;
+ for (Int_t iCluster = 0; iCluster < fNTotalClus; iCluster++) {
+
+ if (!recPoint) {
+ if(fNumberOfECAClusters >= fRecPoints->GetSize()) fRecPoints->Expand(fNumberOfECAClusters*2+1);
+
+ (*fRecPoints)[fNumberOfECAClusters] = new AliEMCALRecPoint("");
+ recPoint = static_cast<AliEMCALRecPoint*>(fRecPoints->At(fNumberOfECAClusters));
+ }
+
+ if (recPoint) {
+ recPoint->SetClusterType(AliVCluster::kEMCALClusterv1);
+ recPoint->SetUniqueID(iCluster);
+ fNumberOfECAClusters++;
+
+ for (Int_t iDigit = 0; iDigit < fNDigitsCluster; iDigit++) {
+ if (fClustersArray[iCluster][iDigit] == 0) continue;
+ digit = fClustersArray[iCluster][iDigit];
+ recPoint->AddDigit(*digit, digit->GetCalibAmp(), kFALSE); //Time or TimeR?
+ fClustersArray[iCluster][iDigit] = 0;
+ recPointOk = kTRUE;
+ }
+
+ if (recPointOk) { // unset the pointer so that a new rec point will be allocated in the next iteration
+ recPoint = 0;
+ recPointOk = kFALSE;
+ }
+ }
+ else {
+ AliError("Error allocating rec points!");
+ break;
+ }
+ }
+
+ if (!recPointOk) {
+ fRecPoints->RemoveLast();
+ fNumberOfECAClusters--;
+ }
+
+ delete digitsC;
+ AliDebug(1, Form("MakeClusters: Number of digits %d -> (e %f)\n", fDigitsArr->GetEntries(),fMinECut));
+ AliDebug(1, Form("total no of clusters %d from %d digits", fNumberOfECAClusters, fDigitsArr->GetEntriesFast()));
}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff