]>
Commit | Line | Data |
---|---|---|
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 | ||
76 | ClassImp(AliEMCALAfterBurnerUF) | |
77 | ||
78 | //------------------------------------------------------------------------ | |
79 | AliEMCALAfterBurnerUF::AliEMCALAfterBurnerUF(): | |
80 | fGeom(NULL), | |
81 | fLogWeight(4.5), // correct? | |
82 | fECALocMaxCut(0.03), // value suggested by Adam | |
71332f0e | 83 | fMinECut(0.01), |
fc645679 | 84 | fRecPoints(NULL), |
85 | fDigitsArr(NULL), | |
86 | fClusterUnfolding(NULL) | |
87 | { | |
88 | // Use this constructor, if unsure | |
89 | ||
90 | Init(); | |
91 | } | |
92 | ||
93 | //------------------------------------------------------------------------ | |
71332f0e | 94 | AliEMCALAfterBurnerUF::AliEMCALAfterBurnerUF(Float_t logWeight, Float_t ecaLocMaxCut, Float_t minECut): |
fc645679 | 95 | fGeom(NULL), |
96 | fLogWeight(logWeight), | |
71332f0e | 97 | fECALocMaxCut(ecaLocMaxCut), |
98 | fMinECut(minECut), | |
fc645679 | 99 | fRecPoints(NULL), |
100 | fDigitsArr(NULL), | |
101 | fClusterUnfolding(NULL) | |
102 | { | |
103 | // This constructor allows to set parameters | |
104 | // Recommended values: | |
105 | // Float_t logWeight = 4.5, ECALocMaxCut = 0.03 | |
106 | ||
107 | Init(); | |
108 | } | |
109 | ||
110 | //------------------------------------------------------------------------ | |
111 | void AliEMCALAfterBurnerUF::Init() | |
112 | { | |
113 | // After-burner initialization | |
114 | // Imports geometry.root (if required), creates unfolding class instance | |
115 | // | |
116 | // TODO: geometry.root does not allow to use the method AliEMCALRecPoint::EvalAll(); | |
117 | // for this to work, the OCDB geometry must be imported instead | |
118 | ||
119 | if (!gGeoManager) | |
4b7fc3d0 | 120 | Warning("AliEMCALAfterBurnerUF::Init","GeoManager not found, please import the geometry.root file or pass to the geometry the misalignment matrices"); |
121 | // TGeoManager::Import("geometry.root"); | |
fc645679 | 122 | |
123 | // required for global cluster position recalculation | |
124 | if (!gGeoManager) | |
0d0d6b98 | 125 | Info("AliEMCALAfterBurnerUF::Init", "gGeoManager was not set, be careful"); |
fc645679 | 126 | |
127 | // initialize geometry, if not yet initialized | |
128 | if (!AliEMCALGeometry::GetInstance()) { | |
0d0d6b98 | 129 | Warning("AliEMCALAfterBurnerUF::Init", "AliEMCALGeometry is not yet initialized. Initializing with EMCAL_COMPLETEV1"); |
130 | AliEMCALGeometry::GetInstance("EMCAL_COMPLETEV1"); | |
fc645679 | 131 | } |
132 | ||
133 | // AliEMCALRecPoint is using exactly this call | |
134 | fGeom = AliEMCALGeometry::GetInstance(); | |
135 | if (!fGeom) | |
136 | Fatal("AliEMCALAfterBurnerUF::AliEMCALAfterBurnerUF", "could not get EMCAL geometry"); | |
137 | ||
138 | fClusterUnfolding = new AliEMCALUnfolding(fGeom); | |
139 | fClusterUnfolding->SetECALocalMaxCut(fECALocMaxCut); | |
71332f0e | 140 | fClusterUnfolding->SetThreshold(fMinECut); |
141 | ||
fc645679 | 142 | // clusters --> recPoints, cells --> digits and back ;) |
143 | fRecPoints = new TObjArray(100); | |
144 | fDigitsArr = new TClonesArray("AliEMCALDigit",1152); | |
145 | } | |
146 | ||
147 | //------------------------------------------------------------------------ | |
148 | AliEMCALAfterBurnerUF::~AliEMCALAfterBurnerUF() | |
149 | { | |
1186cd2b | 150 | // |
151 | // destructor | |
152 | // | |
153 | ||
fc645679 | 154 | if (fClusterUnfolding) delete fClusterUnfolding; |
155 | ||
0d0d6b98 | 156 | if (fRecPoints) { |
157 | fRecPoints->Delete(); | |
158 | delete fRecPoints; | |
159 | } | |
4b7fc3d0 | 160 | if (fDigitsArr) { |
161 | fDigitsArr->Clear("C"); | |
162 | delete fDigitsArr; | |
163 | } | |
fc645679 | 164 | } |
165 | ||
0d0d6b98 | 166 | //------------------------------------------------------------------------ |
167 | void AliEMCALAfterBurnerUF::Clear() | |
168 | { | |
169 | //Clean the arrays | |
170 | ||
171 | if (fRecPoints) fRecPoints->Delete(); // do not Clear(), it leaks, why? | |
172 | if (fDigitsArr) fDigitsArr->Clear("C"); | |
173 | ||
174 | } | |
fc645679 | 175 | //------------------------------------------------------------------------ |
176 | void AliEMCALAfterBurnerUF::RecPoints2Clusters(TObjArray *clusArray) | |
177 | { | |
178 | // Restore clusters from recPoints | |
179 | // Cluster energy, global position, cells and their amplitude fractions are restored | |
180 | // | |
181 | // TODO: restore time and other parameters | |
182 | ||
183 | for(Int_t i = 0; i < fRecPoints->GetEntriesFast(); i++) | |
184 | { | |
185 | AliEMCALRecPoint *recPoint = (AliEMCALRecPoint*) fRecPoints->At(i); | |
186 | ||
03b08eaf | 187 | const Int_t ncells = recPoint->GetMultiplicity(); |
1186cd2b | 188 | Int_t ncellsTrue = 0; |
fc645679 | 189 | |
190 | // cells and their amplitude fractions | |
191 | UShort_t absIds[ncells]; // NOTE: unfolding must not give recPoints with no cells! | |
192 | Double32_t ratios[ncells]; | |
193 | ||
194 | for (Int_t c = 0; c < ncells; c++) { | |
195 | AliEMCALDigit *digit = (AliEMCALDigit*) fDigitsArr->At(recPoint->GetDigitsList()[c]); | |
196 | ||
1186cd2b | 197 | absIds[ncellsTrue] = digit->GetId(); |
198 | ratios[ncellsTrue] = recPoint->GetEnergiesList()[c]/digit->GetAmplitude(); | |
fc645679 | 199 | |
1186cd2b | 200 | if (ratios[ncellsTrue] > 0.001) ncellsTrue++; |
fc645679 | 201 | } |
202 | ||
1186cd2b | 203 | if (ncellsTrue < 1) { |
fc645679 | 204 | Warning("AliEMCALAfterBurnerUF::RecPoints2Clusters", "skipping cluster with no cells"); |
205 | continue; | |
206 | } | |
207 | ||
208 | TVector3 gpos; | |
209 | Float_t g[3]; | |
210 | ||
211 | // calculate new cluster position | |
212 | recPoint->EvalGlobalPosition(fLogWeight, fDigitsArr); | |
213 | recPoint->GetGlobalPosition(gpos); | |
214 | gpos.GetXYZ(g); | |
215 | ||
216 | // create a new cluster | |
217 | AliAODCaloCluster *clus = new AliAODCaloCluster(); | |
15cb4c51 | 218 | clus->SetType(AliVCluster::kEMCALClusterv1); |
fc645679 | 219 | clus->SetE(recPoint->GetEnergy()); |
220 | clus->SetPosition(g); | |
1186cd2b | 221 | clus->SetNCells(ncellsTrue); |
fc645679 | 222 | clus->SetCellsAbsId(absIds); |
223 | clus->SetCellsAmplitudeFraction(ratios); | |
224 | // TODO: time not stored | |
225 | // TODO: some other properties not stored | |
226 | ||
227 | clusArray->Add(clus); | |
228 | } // recPoints loop | |
229 | ||
230 | } | |
231 | ||
232 | //------------------------------------------------------------------------ | |
004d0978 | 233 | void AliEMCALAfterBurnerUF::UnfoldClusters(TObjArray *clusArray, AliVCaloCells *cellsEMCAL) |
fc645679 | 234 | { |
235 | // Unfolds clusters. | |
236 | // | |
237 | // Input: TObjArray of clusters, EMCAL cells. | |
238 | // Output is added to the same array, original clusters are _deleted_ or moved to another position. | |
239 | ||
240 | Int_t ndigits = 0; | |
241 | ||
242 | Int_t nclus = clusArray->GetEntriesFast(); | |
243 | ||
244 | /* Fill recPoints with digits | |
245 | */ | |
246 | for (Int_t i = 0; i < nclus; i++) | |
247 | { | |
004d0978 | 248 | AliVCluster *clus = (AliVCluster*) clusArray->At(i); |
fc645679 | 249 | if (!clus->IsEMCAL()) continue; |
250 | ||
251 | // new recPoint | |
252 | AliEMCALRecPoint *recPoint = new AliEMCALRecPoint(""); | |
15cb4c51 | 253 | recPoint->SetClusterType(AliVCluster::kEMCALClusterv1); |
fc645679 | 254 | fRecPoints->Add(recPoint); |
255 | ||
256 | // fill digits | |
257 | for (Int_t c = 0; c < clus->GetNCells(); c++) { | |
258 | Int_t absId = clus->GetCellAbsId(c); | |
259 | Double_t amp = cellsEMCAL->GetCellAmplitude(absId); | |
260 | Double_t time = cellsEMCAL->GetCellTime(absId); | |
261 | ||
262 | // NOTE: it is easy to include cells recalibration here: | |
263 | // amp *= factor; | |
264 | ||
265 | AliEMCALDigit *digit = (AliEMCALDigit*) fDigitsArr->New(ndigits); | |
266 | ||
267 | digit->SetId(absId); | |
268 | digit->SetAmplitude(amp); | |
269 | digit->SetTime(time); | |
270 | digit->SetTimeR(time); | |
271 | digit->SetIndexInList(ndigits); | |
272 | ||
273 | recPoint->AddDigit(*digit, amp, kFALSE); | |
274 | ||
275 | ndigits++; | |
276 | } | |
277 | ||
278 | // this cluster will be substituted with the result of unfolding | |
279 | clusArray->RemoveAt(i); | |
280 | delete clus; | |
281 | } // cluster loop | |
282 | ||
283 | ||
284 | /* Peform unfolding | |
285 | */ | |
286 | fClusterUnfolding->SetInput(fRecPoints->GetEntriesFast(), fRecPoints, fDigitsArr); | |
287 | fClusterUnfolding->MakeUnfolding(); | |
288 | ||
289 | /* Restore clusters from recPoints. | |
290 | */ | |
291 | RecPoints2Clusters(clusArray); | |
292 | ||
293 | // clean up | |
0d0d6b98 | 294 | fRecPoints->Delete(); // do not Clear(), it leaks, why? |
4b7fc3d0 | 295 | fDigitsArr->Clear("C"); |
fc645679 | 296 | |
297 | clusArray->Compress(); | |
298 | ||
299 | } |