Adding a reminder for coders
[u/mrichter/AliRoot.git] / EMCAL / AliEMCALAfterBurnerUF.cxx
CommitLineData
fc645679 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15// After-burner for the EMCAL cluster unfolding algorithm
16//
15cb4c51 17// Input: TObjArray *clusArray -- array of AliVClusters;
18// AliVCaloCells *cellsEMCAL -- EMCAL cells.
fc645679 19//
20// Output is appended to clusArray, the original (unfolded or not) clusters
21// are deleted or moved to another position in clusArray.
22//
23// If you want to use particular geometry, you must initialize it _before_
24// creating AliEMCALAfterBurnerUF instance. Add this or similar line to the
25// initialization section:
26//
27// AliEMCALGeometry::GetInstance("EMCAL_FIRSTYEARV1");
28//
29// gGeoManager must be initialized for this code to work! Do it yourself or
30// provide geometry.root file in the current directory so that
31// AliEMCALAfterBurnerUF will take it by itself.
32// How to use:
33//
34// // add this lines to the initialization section of your analysis
35// AliEMCALAfterBurnerUF *abuf = new AliEMCALAfterBurnerUF();
36// TObjArray *clusArray = new TObjArray(100);
37//
38//
004d0978 39// AliVEvent *event = InputEvent();
40// AliVCaloCells *cellsEMCAL = event->GetEMCALCells();
fc645679 41//
42// for (Int_t i = 0; i < event->GetNumberOfCaloClusters(); i++)
43// {
004d0978 44// AliVCluster *clus = event->GetCaloCluster(i);
fc645679 45//
46// clusArray->Add(clus->Clone()); // NOTE _CLONE_ in this line
47// }
48//
49// abuf->UnfoldClusters(clusArray, cellsEMCAL);
50//
51// // do an analysis with clusArray
52// // ....
53//
54// // prevent memory leak
55// clusArray->Delete();
56//
57//----
58// Author: Olga Driga (SUBATECH)
59
60// --- ROOT system ---
61#include <TObjArray.h>
62#include <TClonesArray.h>
63#include <TGeoManager.h>
64
65// --- Standard library --
66
67// --- AliRoot header files ---
68#include <AliEMCALAfterBurnerUF.h>
69#include <AliEMCALGeometry.h>
70#include <AliEMCALUnfolding.h>
71#include <AliAODCaloCluster.h>
15cb4c51 72#include <AliVCaloCells.h>
fc645679 73#include <AliEMCALRecPoint.h>
74#include <AliEMCALDigit.h>
75
76ClassImp(AliEMCALAfterBurnerUF)
77
78//------------------------------------------------------------------------
79AliEMCALAfterBurnerUF::AliEMCALAfterBurnerUF():
80 fGeom(NULL),
81 fLogWeight(4.5), // correct?
82 fECALocMaxCut(0.03), // value suggested by Adam
83 fRecPoints(NULL),
84 fDigitsArr(NULL),
85 fClusterUnfolding(NULL)
86{
87 // Use this constructor, if unsure
88
89 Init();
90}
91
92//------------------------------------------------------------------------
93AliEMCALAfterBurnerUF::AliEMCALAfterBurnerUF(Float_t logWeight, Float_t ECALocMaxCut):
94 fGeom(NULL),
95 fLogWeight(logWeight),
96 fECALocMaxCut(ECALocMaxCut),
97 fRecPoints(NULL),
98 fDigitsArr(NULL),
99 fClusterUnfolding(NULL)
100{
101 // This constructor allows to set parameters
102 // Recommended values:
103 // Float_t logWeight = 4.5, ECALocMaxCut = 0.03
104
105 Init();
106}
107
108//------------------------------------------------------------------------
109void AliEMCALAfterBurnerUF::Init()
110{
111 // After-burner initialization
112 // Imports geometry.root (if required), creates unfolding class instance
113 //
114 // TODO: geometry.root does not allow to use the method AliEMCALRecPoint::EvalAll();
115 // for this to work, the OCDB geometry must be imported instead
116
117 if (!gGeoManager)
4b7fc3d0 118 Warning("AliEMCALAfterBurnerUF::Init","GeoManager not found, please import the geometry.root file or pass to the geometry the misalignment matrices");
119// TGeoManager::Import("geometry.root");
fc645679 120
121 // required for global cluster position recalculation
122 if (!gGeoManager)
0d0d6b98 123 Info("AliEMCALAfterBurnerUF::Init", "gGeoManager was not set, be careful");
fc645679 124
125 // initialize geometry, if not yet initialized
126 if (!AliEMCALGeometry::GetInstance()) {
0d0d6b98 127 Warning("AliEMCALAfterBurnerUF::Init", "AliEMCALGeometry is not yet initialized. Initializing with EMCAL_COMPLETEV1");
128 AliEMCALGeometry::GetInstance("EMCAL_COMPLETEV1");
fc645679 129 }
130
131 // AliEMCALRecPoint is using exactly this call
132 fGeom = AliEMCALGeometry::GetInstance();
133 if (!fGeom)
134 Fatal("AliEMCALAfterBurnerUF::AliEMCALAfterBurnerUF", "could not get EMCAL geometry");
135
136 fClusterUnfolding = new AliEMCALUnfolding(fGeom);
137 fClusterUnfolding->SetECALocalMaxCut(fECALocMaxCut);
138
139 // clusters --> recPoints, cells --> digits and back ;)
140 fRecPoints = new TObjArray(100);
141 fDigitsArr = new TClonesArray("AliEMCALDigit",1152);
142}
143
144//------------------------------------------------------------------------
145AliEMCALAfterBurnerUF::~AliEMCALAfterBurnerUF()
146{
147 if (fClusterUnfolding) delete fClusterUnfolding;
148
0d0d6b98 149 if (fRecPoints) {
150 fRecPoints->Delete();
151 delete fRecPoints;
152 }
4b7fc3d0 153 if (fDigitsArr) {
154 fDigitsArr->Clear("C");
155 delete fDigitsArr;
156 }
fc645679 157}
158
159//------------------------------------------------------------------------
0d0d6b98 160void AliEMCALAfterBurnerUF::Clear()
161{
162 //Clean the arrays
163
164 if (fRecPoints) fRecPoints->Delete(); // do not Clear(), it leaks, why?
165 if (fDigitsArr) fDigitsArr->Clear("C");
166
167}
168//------------------------------------------------------------------------
fc645679 169void AliEMCALAfterBurnerUF::RecPoints2Clusters(TObjArray *clusArray)
170{
171 // Restore clusters from recPoints
172 // Cluster energy, global position, cells and their amplitude fractions are restored
173 //
174 // TODO: restore time and other parameters
175
176 for(Int_t i = 0; i < fRecPoints->GetEntriesFast(); i++)
177 {
178 AliEMCALRecPoint *recPoint = (AliEMCALRecPoint*) fRecPoints->At(i);
179
03b08eaf 180 const Int_t ncells = recPoint->GetMultiplicity();
fc645679 181 Int_t ncells_true = 0;
182
183 // cells and their amplitude fractions
184 UShort_t absIds[ncells]; // NOTE: unfolding must not give recPoints with no cells!
185 Double32_t ratios[ncells];
186
187 for (Int_t c = 0; c < ncells; c++) {
188 AliEMCALDigit *digit = (AliEMCALDigit*) fDigitsArr->At(recPoint->GetDigitsList()[c]);
189
190 absIds[ncells_true] = digit->GetId();
191 ratios[ncells_true] = recPoint->GetEnergiesList()[c]/digit->GetAmplitude();
192
193 if (ratios[ncells_true] > 0.001) ncells_true++;
194 }
195
196 if (ncells_true < 1) {
197 Warning("AliEMCALAfterBurnerUF::RecPoints2Clusters", "skipping cluster with no cells");
198 continue;
199 }
200
201 TVector3 gpos;
202 Float_t g[3];
203
204 // calculate new cluster position
205 recPoint->EvalGlobalPosition(fLogWeight, fDigitsArr);
206 recPoint->GetGlobalPosition(gpos);
207 gpos.GetXYZ(g);
208
209 // create a new cluster
210 AliAODCaloCluster *clus = new AliAODCaloCluster();
15cb4c51 211 clus->SetType(AliVCluster::kEMCALClusterv1);
fc645679 212 clus->SetE(recPoint->GetEnergy());
213 clus->SetPosition(g);
214 clus->SetNCells(ncells_true);
215 clus->SetCellsAbsId(absIds);
216 clus->SetCellsAmplitudeFraction(ratios);
217 // TODO: time not stored
218 // TODO: some other properties not stored
219
220 clusArray->Add(clus);
221 } // recPoints loop
222
223}
224
225//------------------------------------------------------------------------
004d0978 226void AliEMCALAfterBurnerUF::UnfoldClusters(TObjArray *clusArray, AliVCaloCells *cellsEMCAL)
fc645679 227{
228 // Unfolds clusters.
229 //
230 // Input: TObjArray of clusters, EMCAL cells.
231 // Output is added to the same array, original clusters are _deleted_ or moved to another position.
232
233 Int_t ndigits = 0;
234
235 Int_t nclus = clusArray->GetEntriesFast();
236
237 /* Fill recPoints with digits
238 */
239 for (Int_t i = 0; i < nclus; i++)
240 {
004d0978 241 AliVCluster *clus = (AliVCluster*) clusArray->At(i);
fc645679 242 if (!clus->IsEMCAL()) continue;
243
244 // new recPoint
245 AliEMCALRecPoint *recPoint = new AliEMCALRecPoint("");
15cb4c51 246 recPoint->SetClusterType(AliVCluster::kEMCALClusterv1);
fc645679 247 fRecPoints->Add(recPoint);
248
249 // fill digits
250 for (Int_t c = 0; c < clus->GetNCells(); c++) {
251 Int_t absId = clus->GetCellAbsId(c);
252 Double_t amp = cellsEMCAL->GetCellAmplitude(absId);
253 Double_t time = cellsEMCAL->GetCellTime(absId);
254
255 // NOTE: it is easy to include cells recalibration here:
256 // amp *= factor;
257
258 AliEMCALDigit *digit = (AliEMCALDigit*) fDigitsArr->New(ndigits);
259
260 digit->SetId(absId);
261 digit->SetAmplitude(amp);
262 digit->SetTime(time);
263 digit->SetTimeR(time);
264 digit->SetIndexInList(ndigits);
265
266 recPoint->AddDigit(*digit, amp, kFALSE);
267
268 ndigits++;
269 }
270
271 // this cluster will be substituted with the result of unfolding
272 clusArray->RemoveAt(i);
273 delete clus;
274 } // cluster loop
275
276
277 /* Peform unfolding
278 */
279 fClusterUnfolding->SetInput(fRecPoints->GetEntriesFast(), fRecPoints, fDigitsArr);
280 fClusterUnfolding->MakeUnfolding();
281
282 /* Restore clusters from recPoints.
283 */
284 RecPoints2Clusters(clusArray);
285
286 // clean up
0d0d6b98 287 fRecPoints->Delete(); // do not Clear(), it leaks, why?
4b7fc3d0 288 fDigitsArr->Clear("C");
fc645679 289
290 clusArray->Compress();
291
292}