]>
Commit | Line | Data |
---|---|---|
5502f2ea | 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 | /* $Id$ */ | |
16 | //_________________________________________________________________________ | |
17 | // Implementation version 1 of algorithm class to construct EMCAL track segments | |
18 | // Track segment for EMCAL is list of | |
19 | // ECAL RecPoint + (possibly) PRE RecPoint + (possibly) HCAL RecPoint | |
20 | // To find TrackSegments we do the following: | |
21 | // for each ECAL RecPoint we look for PRE and HC RecPoints with same direction within fSame. | |
22 | // If there is such a PRE or ECAL RecPoint, | |
23 | // we make a "Link": indexes of ECAL and PRE, HCAL RecPoints and their scalar product. | |
24 | // Then we sort "Links", starting from the | |
25 | // least "Link" pointing to the unassigned RecPoints assigning them to a new TrackSegment. | |
26 | // If there is no PRE, HCAL RecPoint we make a TrackSegment | |
27 | // consisting from ECAL alone. There is no TrackSegments without ECAL RecPoint. | |
28 | //// In principle this class should be called from AliEMCALReconstructioner, but | |
29 | // one can use it as well in standalone mode. | |
30 | // | |
31 | //*-- Author: Dmitri Peressounko (RRC Ki & SUBATECH) & Yves Schutz (SUBATECH) | |
32 | // | |
33 | ||
34 | // --- ROOT system --- | |
5502f2ea | 35 | #include "TTree.h" |
5502f2ea | 36 | #include "TBenchmark.h" |
37 | ||
38 | // --- Standard library --- | |
39 | ||
40 | // --- AliRoot header files --- | |
41 | ||
42 | #include "AliEMCALTrackSegmentMakerv1.h" | |
5502f2ea | 43 | #include "AliEMCALTrackSegment.h" |
44 | #include "AliEMCALLink.h" | |
45 | #include "AliEMCALGetter.h" | |
5502f2ea | 46 | |
47 | ClassImp( AliEMCALTrackSegmentMakerv1) | |
48 | ||
49 | ||
50 | //____________________________________________________________________________ | |
51 | AliEMCALTrackSegmentMakerv1::AliEMCALTrackSegmentMakerv1() : AliEMCALTrackSegmentMaker() | |
52 | { | |
53 | // default ctor (to be used mainly by Streamer) | |
54 | ||
55 | InitParameters() ; | |
5502f2ea | 56 | fDefaultInit = kTRUE ; |
57 | } | |
58 | ||
59 | //____________________________________________________________________________ | |
88cb7938 | 60 | AliEMCALTrackSegmentMakerv1::AliEMCALTrackSegmentMakerv1(const TString alirunFileName, const TString eventFolderName) |
61 | :AliEMCALTrackSegmentMaker(alirunFileName, eventFolderName) | |
5502f2ea | 62 | { |
63 | // ctor | |
64 | ||
65 | InitParameters() ; | |
5502f2ea | 66 | Init() ; |
5502f2ea | 67 | fDefaultInit = kFALSE ; |
68 | ||
69 | } | |
70 | ||
71 | //____________________________________________________________________________ | |
72 | AliEMCALTrackSegmentMakerv1::~AliEMCALTrackSegmentMakerv1() | |
73 | { | |
74 | // dtor | |
75 | // fDefaultInit = kTRUE if TrackSegmentMaker created by default ctor (to get just the parameters) | |
fdebddeb | 76 | |
5502f2ea | 77 | } |
78 | ||
5502f2ea | 79 | //____________________________________________________________________________ |
80 | const TString AliEMCALTrackSegmentMakerv1::BranchName() const | |
81 | { | |
88cb7938 | 82 | return GetName() ; |
83 | ||
5502f2ea | 84 | } |
85 | ||
86 | //____________________________________________________________________________ | |
70a93198 | 87 | Float_t AliEMCALTrackSegmentMakerv1::HowClose(AliEMCALRecPoint * ec, AliEMCALRecPoint * rp, Bool_t &toofar)const |
5502f2ea | 88 | { |
89 | // Calculates the distance between the EMC RecPoint and the PPSD RecPoint | |
90 | // Clusters are sorted in "rows" and "columns" of width 1 cm | |
91 | ||
92 | Float_t r = -1. ; | |
93 | Float_t delta = 10. ; // large enough to be ineffective ??! | |
94 | ||
95 | ||
70a93198 | 96 | TVector3 vecEC; |
97 | TVector3 vecRP; | |
98 | ec->GetGlobalPosition(vecEC); | |
99 | rp->GetGlobalPosition(vecRP); | |
5502f2ea | 100 | |
101 | Float_t pro = TMath::Abs(1 - (vecEC * vecRP / ( vecEC.Mag() * vecRP.Mag() ))) ; | |
102 | ||
103 | if(pro <= delta) { | |
104 | r = pro ; | |
105 | toofar = kFALSE ; | |
106 | } | |
107 | else | |
108 | toofar = kTRUE ; | |
109 | ||
110 | if (gDebug == 2 ) | |
fdebddeb | 111 | printf("HowClose: ec = %d, rp = %d pro = %f, toofar=%d", ec->GetIndexInList(), rp->GetIndexInList(), pro, toofar ) ; |
5502f2ea | 112 | |
113 | return r ; | |
114 | } | |
115 | ||
116 | //____________________________________________________________________________ | |
117 | void AliEMCALTrackSegmentMakerv1::Init() | |
118 | { | |
119 | // Make all memory allocations that are not possible in default constructor | |
120 | ||
88cb7938 | 121 | AliEMCALGetter* gime = AliEMCALGetter::Instance(GetTitle(), fEventFolderName.Data()); |
fdebddeb | 122 | |
88cb7938 | 123 | if ( !gime->TrackSegmentMaker() ) { |
124 | gime->PostTrackSegmentMaker(this); | |
125 | } | |
5502f2ea | 126 | } |
127 | ||
128 | //____________________________________________________________________________ | |
129 | void AliEMCALTrackSegmentMakerv1::InitParameters() | |
130 | { | |
88cb7938 | 131 | fClose = 10e-3 ; |
88cb7938 | 132 | fTrackSegmentsInRun = 0 ; |
5502f2ea | 133 | } |
134 | ||
135 | ||
136 | //____________________________________________________________________________ | |
137 | void AliEMCALTrackSegmentMakerv1::MakeLinks()const | |
138 | { | |
139 | // Finds distances (links) between all PRE, EC and HC clusters, | |
140 | // which are not further apart from each other than fDangle | |
141 | // and sort them in accordance with this distance | |
142 | ||
fdebddeb | 143 | /* AliEMCALGetter * gime = AliEMCALGetter::Instance() ; |
88cb7938 | 144 | TObjArray * aECARecPoints = gime->ECARecPoints() ; |
fdebddeb | 145 | // TObjArray * aPRERecPoints = gime->PRERecPoints() ; |
146 | //TObjArray * aHCARecPoints = gime->HCARecPoints() ; | |
5502f2ea | 147 | |
148 | fPRELinkArray->Clear() ; | |
88cb7938 | 149 | fHCALinkArray->Clear() ; |
5502f2ea | 150 | |
70a93198 | 151 | AliEMCALRecPoint * pre ; |
152 | AliEMCALRecPoint * eca ; | |
153 | AliEMCALRecPoint * hca ; | |
5502f2ea | 154 | |
155 | Int_t iPRELink = 0 ; | |
88cb7938 | 156 | Int_t iHCALink = 0 ; |
5502f2ea | 157 | |
88cb7938 | 158 | Int_t iECARP; |
159 | for(iECARP = 0; iECARP < aECARecPoints->GetEntriesFast(); iECARP++ ) { | |
70a93198 | 160 | eca = dynamic_cast<AliEMCALRecPoint *>(aECARecPoints->At(iECARP)) ; |
5502f2ea | 161 | Bool_t toofar = kTRUE ; |
162 | Int_t iPRERP = 0 ; | |
163 | for(iPRERP = 0; iPRERP < aPRERecPoints->GetEntriesFast(); iPRERP++ ) { | |
70a93198 | 164 | pre = dynamic_cast<AliEMCALRecPoint *>(aPRERecPoints->At(iPRERP)) ; |
88cb7938 | 165 | Float_t prod = HowClose(eca, pre, toofar) ; |
5502f2ea | 166 | if(toofar) |
167 | break ; | |
168 | if(prod < fClose) { | |
88cb7938 | 169 | new ((*fPRELinkArray)[iPRELink++]) AliEMCALLink(prod, iECARP, iPRERP, 0) ; |
5502f2ea | 170 | } |
171 | } | |
172 | toofar = kTRUE ; | |
88cb7938 | 173 | Int_t iHCARP = 0 ; |
174 | for(iHCARP = 0; iHCARP < aHCARecPoints->GetEntriesFast(); iHCARP++ ) { | |
70a93198 | 175 | hca = dynamic_cast<AliEMCALRecPoint *>(aHCARecPoints->At(iHCARP)) ; |
88cb7938 | 176 | Float_t prod = HowClose(eca, hca, toofar) ; |
5502f2ea | 177 | if(toofar) |
178 | break ; | |
179 | if(prod < fClose) { | |
88cb7938 | 180 | new ((*fHCALinkArray)[iHCALink++]) AliEMCALLink(prod, iECARP, iHCARP, 1) ; |
5502f2ea | 181 | } |
182 | } | |
183 | } | |
184 | ||
185 | fPRELinkArray->Sort() ; //first links with largest scalar product | |
88cb7938 | 186 | fHCALinkArray->Sort() ; //first links with largest scalar product |
fdebddeb | 187 | */ |
5502f2ea | 188 | } |
189 | ||
190 | //____________________________________________________________________________ | |
191 | void AliEMCALTrackSegmentMakerv1::MakePairs() | |
192 | { | |
193 | // Using the previously made list of "links", we found the best link - i.e. | |
194 | // link with the largest scalar product (closest to one) to still | |
195 | // unassigned RecParticles. We assign these RecPoints to TrackSegment and | |
196 | // remove them from the list of "unassigned". | |
197 | ||
fdebddeb | 198 | /*AliEMCALGetter * gime = AliEMCALGetter::Instance() ; |
199 | TObjArray * aECARecPoints = gime->ECARecPoints() ; | |
5502f2ea | 200 | TObjArray * aPRERecPoints = gime->PRERecPoints() ; |
88cb7938 | 201 | TObjArray * aHCARecPoints = gime->HCARecPoints() ; |
202 | TClonesArray * trackSegments = gime->TrackSegments() ; | |
5502f2ea | 203 | |
204 | //Make arrays to mark clusters already chosen | |
88cb7938 | 205 | Int_t * ecaExist = 0; |
206 | Int_t nECA = aECARecPoints->GetEntriesFast() ; | |
207 | if (nECA) | |
208 | ecaExist = new Int_t[nECA] ; | |
5502f2ea | 209 | |
210 | Int_t index; | |
88cb7938 | 211 | for(index = 0; index < nECA; index ++) |
212 | ecaExist[index] = 1 ; | |
5502f2ea | 213 | |
214 | Bool_t * preExist = 0; | |
215 | Int_t nPRE = aPRERecPoints->GetEntriesFast() ; | |
216 | if(nPRE) | |
217 | preExist = new Bool_t[nPRE] ; | |
218 | for(index = 0; index < nPRE; index ++) | |
219 | preExist[index] = kTRUE ; | |
220 | ||
88cb7938 | 221 | Bool_t * hcaExist = 0; |
222 | Int_t nHCA = aHCARecPoints->GetEntriesFast() ; | |
223 | if(nHCA) | |
224 | hcaExist = new Bool_t[nHCA] ; | |
225 | for(index = 0; index < nHCA; index ++) | |
226 | hcaExist[index] = kTRUE ; | |
5502f2ea | 227 | |
70a93198 | 228 | AliEMCALRecPoint * null = 0 ; |
5502f2ea | 229 | // Finds the smallest links and makes pairs of PRE and ECAL clusters with largest scalar product |
230 | ||
231 | TIter nextPRE(fPRELinkArray) ; | |
232 | AliEMCALLink * linkPRE ; | |
233 | ||
234 | while ( (linkPRE = static_cast<AliEMCALLink *>(nextPRE()) ) ){ | |
235 | ||
88cb7938 | 236 | if(ecaExist[linkPRE->GetECA()] != -1){ //without PRE yet |
5502f2ea | 237 | |
238 | if(preExist[linkPRE->GetOther()]){ // PRE still exist | |
239 | ||
240 | new ((* trackSegments)[fNTrackSegments]) | |
70a93198 | 241 | AliEMCALTrackSegment(dynamic_cast<AliEMCALRecPoint *>(aECARecPoints->At(linkPRE->GetECA())) , |
242 | dynamic_cast<AliEMCALRecPoint *>(aPRERecPoints->At(linkPRE->GetOther())), null) ; | |
5502f2ea | 243 | (dynamic_cast<AliEMCALTrackSegment *>(trackSegments->At(fNTrackSegments)))->SetIndexInList(fNTrackSegments); |
244 | fNTrackSegments++ ; | |
245 | if (gDebug == 2 ) | |
fdebddeb | 246 | printf("MakePairs: ECAL section with PRE section") ; |
88cb7938 | 247 | ecaExist[linkPRE->GetECA()] = -1 ; //Mark ecal that pre was found |
5502f2ea | 248 | //mark PRE recpoint as already used |
249 | preExist[linkPRE->GetOther()] = kFALSE ; | |
250 | } //if PRE still exist | |
251 | } | |
252 | } | |
253 | ||
254 | // Finds the smallest links and makes pairs of HCAL and ECAL clusters with largest scalar product | |
255 | ||
88cb7938 | 256 | TIter nextHCA(fHCALinkArray) ; |
257 | AliEMCALLink * linkHCA ; | |
5502f2ea | 258 | |
88cb7938 | 259 | while ( (linkHCA = static_cast<AliEMCALLink *>(nextHCA()) ) ){ |
5502f2ea | 260 | |
88cb7938 | 261 | if(ecaExist[linkHCA->GetECA()] != -2){ //without HCAL yet |
5502f2ea | 262 | |
88cb7938 | 263 | if(hcaExist[linkHCA->GetOther()]){ // HCAL still exist |
5502f2ea | 264 | // search among the already existing track segments |
265 | Int_t ii ; | |
266 | Bool_t found = kFALSE ; | |
267 | AliEMCALTrackSegment * ts = 0 ; | |
268 | for ( ii = 0 ; ii < fNTrackSegments ; ii++ ) { | |
269 | ts = dynamic_cast<AliEMCALTrackSegment *>(trackSegments->At(ii)) ; | |
88cb7938 | 270 | if ( ts->GetECAIndex() == linkHCA->GetECA() ) { |
5502f2ea | 271 | found = kTRUE ; |
272 | break ; | |
273 | } | |
274 | } | |
275 | if (found){ | |
70a93198 | 276 | ts->SetHCARecPoint( dynamic_cast<AliEMCALRecPoint *>(aHCARecPoints->At(linkHCA->GetOther())) ) ; |
5502f2ea | 277 | if (gDebug == 2 ) |
fdebddeb | 278 | printf("MakePairs: ECAL section with PRE and HCAL sections") ; |
5502f2ea | 279 | } |
280 | if (!found) { | |
281 | new ((* trackSegments)[fNTrackSegments]) | |
70a93198 | 282 | AliEMCALTrackSegment(dynamic_cast<AliEMCALRecPoint *>(aECARecPoints->At(linkHCA->GetECA())), null, |
283 | dynamic_cast<AliEMCALRecPoint *>(aHCARecPoints->At(linkHCA->GetOther()))) ; | |
5502f2ea | 284 | (dynamic_cast<AliEMCALTrackSegment *>(trackSegments->At(fNTrackSegments)))->SetIndexInList(fNTrackSegments); |
285 | fNTrackSegments++ ; | |
286 | if (gDebug == 2 ) | |
fdebddeb | 287 | printf("MakePairs: ECAL section with HCAL section") ; |
5502f2ea | 288 | } |
88cb7938 | 289 | ecaExist[linkHCA->GetECA()] = -2 ; //Mark ecal that hcal was found |
5502f2ea | 290 | //mark HCAL recpoint as already used |
88cb7938 | 291 | hcaExist[linkHCA->GetOther()] = kFALSE ; |
5502f2ea | 292 | } //if HCAL still exist |
293 | } | |
294 | } | |
295 | ||
296 | ||
297 | //look through ECAL recPoints left without PRE/HCAL | |
88cb7938 | 298 | if(ecaExist){ //if there is ecal rec point |
299 | Int_t iECARP ; | |
300 | for(iECARP = 0; iECARP < nECA ; iECARP++ ){ | |
301 | if(ecaExist[iECARP] > 0 ){ | |
5502f2ea | 302 | new ((*trackSegments)[fNTrackSegments]) |
70a93198 | 303 | AliEMCALTrackSegment(dynamic_cast<AliEMCALRecPoint *>(aECARecPoints->At(iECARP)), null, null) ; |
5502f2ea | 304 | (dynamic_cast<AliEMCALTrackSegment *>(trackSegments->At(fNTrackSegments)))->SetIndexInList(fNTrackSegments); |
305 | fNTrackSegments++; | |
306 | if( gDebug == 2 ) | |
fdebddeb | 307 | printf("MakePairs: ECAL section alone") ; |
5502f2ea | 308 | } |
309 | } | |
310 | } | |
88cb7938 | 311 | delete [] ecaExist ; |
5502f2ea | 312 | delete [] preExist ; |
88cb7938 | 313 | delete [] hcaExist ; |
fdebddeb | 314 | */ |
5502f2ea | 315 | } |
316 | ||
317 | //____________________________________________________________________________ | |
318 | void AliEMCALTrackSegmentMakerv1::Exec(Option_t * option) | |
319 | { | |
320 | // STEERing method | |
321 | ||
5502f2ea | 322 | |
323 | if(strstr(option,"tim")) | |
324 | gBenchmark->Start("EMCALTSMaker"); | |
325 | ||
326 | if(strstr(option,"print")) { | |
327 | Print("") ; | |
328 | return ; | |
329 | } | |
330 | ||
88cb7938 | 331 | AliEMCALGetter * gime = AliEMCALGetter::Instance() ; |
5502f2ea | 332 | |
88cb7938 | 333 | Int_t nevents = gime->MaxEvent() ; |
5502f2ea | 334 | Int_t ievent ; |
335 | ||
336 | for(ievent = 0; ievent < nevents; ievent++){ | |
5502f2ea | 337 | gime->Event(ievent,"R") ; |
5502f2ea | 338 | //Make some initializations |
339 | fNTrackSegments = 0 ; | |
5502f2ea | 340 | |
88cb7938 | 341 | gime->TrackSegments()->Clear() ; |
5502f2ea | 342 | |
88cb7938 | 343 | MakeLinks() ; |
5502f2ea | 344 | MakePairs() ; |
345 | ||
9e5d2067 | 346 | WriteTrackSegments() ; |
5502f2ea | 347 | |
348 | if(strstr(option,"deb")) | |
349 | PrintTrackSegments(option) ; | |
350 | ||
351 | //increment the total number of track segments per run | |
88cb7938 | 352 | fTrackSegmentsInRun += gime->TrackSegments()->GetEntriesFast() ; |
5502f2ea | 353 | |
354 | } | |
355 | ||
356 | if(strstr(option,"tim")){ | |
357 | gBenchmark->Stop("EMCALTSMaker"); | |
fdebddeb | 358 | printf("Exec: took %f seconds for making TS %f seconds per event", |
5502f2ea | 359 | gBenchmark->GetCpuTime("EMCALTSMaker"), gBenchmark->GetCpuTime("EMCALTSMaker")/nevents) ; |
360 | } | |
88cb7938 | 361 | Unload(); |
362 | } | |
363 | ||
364 | //____________________________________________________________________________ | |
365 | void AliEMCALTrackSegmentMakerv1::Unload() | |
366 | { | |
fdebddeb | 367 | // Unloads the RecPoints and Tracks |
88cb7938 | 368 | AliEMCALGetter * gime = AliEMCALGetter::Instance() ; |
369 | gime->EmcalLoader()->UnloadRecPoints() ; | |
370 | gime->EmcalLoader()->UnloadTracks() ; | |
5502f2ea | 371 | } |
372 | ||
373 | //____________________________________________________________________________ | |
9e5d2067 | 374 | void AliEMCALTrackSegmentMakerv1::Print(Option_t * /*option*/)const |
5502f2ea | 375 | { |
376 | // Print TrackSegmentMaker parameters | |
377 | ||
fdebddeb | 378 | printf("Print: TrackSegmentMakerv1 parameters:") ; |
5502f2ea | 379 | if( strcmp(GetName(), "") != 0 ) { |
380 | printf("Making Track segments with parameters:\n") ; | |
381 | printf(" Allowed spred on the scalar product of two recpoints with same direction: %f\n", fClose) ; | |
382 | printf("============================================\n") ; | |
383 | } | |
384 | else | |
385 | printf("AliEMCALTrackSegmentMakerv1 not initialized ") ; | |
386 | } | |
387 | ||
388 | //____________________________________________________________________________ | |
9e5d2067 | 389 | void AliEMCALTrackSegmentMakerv1::WriteTrackSegments() |
5502f2ea | 390 | { |
391 | // Writes found TrackSegments to TreeR. Creates branches | |
392 | // "EMCALTS" and "AliEMCALTrackSegmentMaker" with the same title. | |
393 | // In the former branch found TrackSegments are stored, while | |
394 | // in the latter all parameters, with which TS were made. | |
395 | // ROOT does not allow overwriting existing branches, therefore | |
396 | // first we check, if branches with the same title already exist. | |
397 | // If yes - exits without writing. | |
398 | ||
88cb7938 | 399 | AliEMCALGetter *gime = AliEMCALGetter::Instance() ; |
5502f2ea | 400 | |
401 | TClonesArray * trackSegments = gime->TrackSegments() ; | |
402 | trackSegments->Expand(trackSegments->GetEntriesFast()) ; | |
88cb7938 | 403 | |
404 | TTree * treeT = gime->TreeT(); | |
5502f2ea | 405 | |
406 | //First TS | |
407 | Int_t bufferSize = 32000 ; | |
88cb7938 | 408 | TBranch * tsBranch = treeT->Branch("EMCALTS",&trackSegments,bufferSize); |
409 | tsBranch->Fill() ; | |
410 | ||
411 | gime->WriteTracks("OVERWRITE"); | |
412 | gime->WriteTrackSegmentMaker("OVERWRITE"); | |
5502f2ea | 413 | } |
414 | ||
415 | ||
416 | //____________________________________________________________________________ | |
417 | void AliEMCALTrackSegmentMakerv1::PrintTrackSegments(Option_t * option) | |
418 | { | |
419 | // option deb - prints # of found TrackSegments | |
420 | // option deb all - prints as well indexed of found RecParticles assigned to the TS | |
5502f2ea | 421 | |
88cb7938 | 422 | TClonesArray * trackSegments = AliEMCALGetter::Instance()->TrackSegments() ; |
5502f2ea | 423 | |
424 | ||
fdebddeb | 425 | printf("PrintTrackSegments: Results from TrackSegmentMaker:") ; |
5502f2ea | 426 | printf("nevent: %d\n", gAlice->GetEvNumber()) ; |
427 | printf(" Found %d TrackSegments\n", trackSegments->GetEntriesFast() ); | |
428 | ||
429 | if(strstr(option,"all")) { // printing found TS | |
fdebddeb | 430 | printf("TrackSegment# ECAL RP# \n") ; |
5502f2ea | 431 | Int_t index; |
432 | for (index = 0 ; index < fNTrackSegments ; index++) { | |
433 | AliEMCALTrackSegment * ts = (AliEMCALTrackSegment * )trackSegments->At(index) ; | |
fdebddeb | 434 | printf(" %d %d \n", |
435 | ts->GetIndexInList(), ts->GetECAIndex()); | |
5502f2ea | 436 | } |
437 | } | |
438 | } |