Line data Source code
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 :
16 : // $Id$
17 :
18 : #include "AliMUONSimpleClusterServer.h"
19 :
20 : #include "AliCodeTimer.h"
21 : #include "AliLog.h"
22 : #include "AliMUONCluster.h"
23 : #include "AliMUONConstants.h"
24 : #include "AliMUONGeometryTransformer.h"
25 : #include "AliMUONPad.h"
26 : #include "AliMUONTriggerTrackToTrackerClusters.h"
27 : #include "AliMUONVCluster.h"
28 : #include "AliMUONVClusterFinder.h"
29 : #include "AliMUONVClusterStore.h"
30 : #include "AliMUONVDigitStore.h"
31 : #include "AliMUONRecoParam.h"
32 : #include "AliMpArea.h"
33 : #include "AliMpDEIterator.h"
34 : #include "AliMpDEManager.h"
35 : #include "AliMpExMap.h"
36 : #include "AliMpExMapIterator.h"
37 : #include "AliMpPad.h"
38 : #include "AliMpSegmentation.h"
39 : #include "AliMpVSegmentation.h"
40 : #include <Riostream.h>
41 : #include <TObjArray.h>
42 : #include <TString.h>
43 : #include <float.h>
44 :
45 : /// \class AliMUONSimpleClusterServer
46 : ///
47 : /// Implementation of AliMUONVClusterServer interface
48 : ///
49 : ///
50 : /// \author Laurent Aphecetche, Subatech
51 :
52 : using std::endl;
53 : using std::cout;
54 : /// \cond CLASSIMP
55 18 : ClassImp(AliMUONSimpleClusterServer)
56 : /// \endcond
57 :
58 : namespace
59 : {
60 : TString AsString(const AliMpArea& area)
61 : {
62 0 : return Form("(X,Y)=(%7.3f,%7.3f) (DX,DY)=(%7.3f,%7.3f)",
63 0 : area.GetPositionX(),
64 0 : area.GetPositionY(),
65 0 : area.GetDimensionX(), /// TBCL was Y !!!
66 0 : area.GetDimensionY());
67 : }
68 : }
69 :
70 : //_____________________________________________________________________________
71 : AliMUONSimpleClusterServer::AliMUONSimpleClusterServer(AliMUONVClusterFinder* clusterFinder,
72 : const AliMUONGeometryTransformer& transformer)
73 2 : : AliMUONVClusterServer(),
74 2 : fDigitStore(0x0),
75 2 : fClusterFinder(clusterFinder),
76 2 : fkTransformer(transformer),
77 2 : fPads(),
78 2 : fTriggerTrackStore(0x0),
79 2 : fBypass(0x0)
80 10 : {
81 : /// Ctor
82 : /// Note that we take ownership of the clusterFinder
83 :
84 6 : fPads[0] = new AliMpExMap;
85 6 : fPads[1] = new AliMpExMap;
86 :
87 4 : fPadsIterator[0] = fPads[0]->CreateIterator();
88 4 : fPadsIterator[1] = fPads[1]->CreateIterator();
89 4 : }
90 :
91 : //_____________________________________________________________________________
92 : AliMUONSimpleClusterServer::~AliMUONSimpleClusterServer()
93 12 : {
94 : /// Dtor
95 4 : delete fClusterFinder;
96 4 : delete fPads[0];
97 4 : delete fPads[1];
98 4 : delete fPadsIterator[0];
99 4 : delete fPadsIterator[1];
100 2 : delete fBypass;
101 6 : }
102 :
103 : //_____________________________________________________________________________
104 : Int_t
105 : AliMUONSimpleClusterServer::Clusterize(Int_t chamberId,
106 : AliMUONVClusterStore& clusterStore,
107 : const AliMpArea& area,
108 : const AliMUONRecoParam* recoParam)
109 : {
110 : /// Area is in absolute coordinate. If not valid, means clusterize all
111 : /// the chamber.
112 : ///
113 : /// We first find out the list of DE that have a non-zero overlap with area,
114 : /// and then use the clusterfinder to find clusters in those areas (and DE).
115 :
116 160 : AliCodeTimerAuto(Form("Chamber %d",chamberId),0);
117 :
118 80 : if ( fTriggerTrackStore && chamberId >= 6 )
119 : {
120 0 : return fBypass->GenerateClusters(chamberId,clusterStore);
121 : }
122 :
123 80 : if (!recoParam) {
124 0 : AliError("Reconstruction parameters are missing: unable to clusterize");
125 0 : return 0;
126 : }
127 :
128 80 : AliMpDEIterator it;
129 :
130 80 : it.First(chamberId);
131 :
132 : Int_t nofAddedClusters(0);
133 80 : Int_t fNCluster = clusterStore.GetSize();
134 :
135 400 : AliDebug(1,Form("chamberId = %2d NofClusters before = %d searchArea=%s",
136 : chamberId,fNCluster,AsString(area).Data()));
137 :
138 3904 : while ( !it.IsDone() )
139 : {
140 1248 : Int_t detElemId = it.CurrentDEId();
141 :
142 1248 : TObjArray* pads[2] =
143 3744 : {
144 2496 : static_cast<TObjArray*>(fPads[0]->GetValue(detElemId)),
145 2496 : static_cast<TObjArray*>(fPads[1]->GetValue(detElemId))
146 : };
147 :
148 1920 : if ( ( pads[0] && pads[0]->GetLast()>=0 ) ||
149 1360 : ( pads[1] && pads[1]->GetLast()>=0 ) )
150 : {
151 144 : AliMpArea deArea; // area in DE-local-coordinates
152 : Bool_t ok(kTRUE);
153 :
154 144 : if ( area.IsValid() )
155 : {
156 0 : ok = Overlap(detElemId,area,deArea);
157 0 : }
158 :
159 144 : if ( ok )
160 : {
161 144 : const AliMpVSegmentation* seg[2] =
162 576 : { AliMpSegmentation::Instance()->GetMpSegmentation(detElemId,AliMp::kCath0),
163 288 : AliMpSegmentation::Instance()->GetMpSegmentation(detElemId,AliMp::kCath1)
164 : };
165 :
166 288 : fClusterFinder->SetChargeHints(recoParam->LowestPadCharge(),
167 144 : recoParam->LowestClusterCharge());
168 :
169 288 : if ( fClusterFinder->NeedSegmentation() )
170 : {
171 288 : fClusterFinder->Prepare(detElemId,pads,deArea,seg);
172 : }
173 : else
174 : {
175 0 : fClusterFinder->Prepare(detElemId,pads,deArea);
176 : }
177 :
178 720 : AliDebug(1,Form("Clusterizing DE %04d with %3d pads (cath0) and %3d pads (cath1)",
179 : detElemId,
180 : (pads[0] ? pads[0]->GetLast()+1 : 0),
181 : (pads[1] ? pads[1]->GetLast()+1 : 0)));
182 :
183 : AliMUONCluster* cluster;
184 :
185 924 : while ( ( cluster = fClusterFinder->NextCluster() ) )
186 : {
187 : // add new cluster to the store with information to build its ID
188 : // increment the number of clusters into the store
189 164 : AliMUONVCluster* rawCluster = clusterStore.Add(chamberId, detElemId, fNCluster++);
190 :
191 164 : ++nofAddedClusters;
192 :
193 : // fill array of Id of digits attached to this cluster
194 164 : Int_t nPad = cluster->Multiplicity();
195 164 : if (nPad < 1) AliWarning("no pad attached to the cluster");
196 :
197 3188 : for (Int_t iPad=0; iPad<nPad; iPad++)
198 : {
199 1430 : AliMUONPad *pad = cluster->Pad(iPad);
200 :
201 : // skip virtual pads
202 1536 : if (!pad->IsReal()) continue;
203 :
204 2648 : rawCluster->AddDigitId(pad->GetUniqueID());
205 1324 : }
206 :
207 : // fill charge and other cluster informations
208 328 : rawCluster->SetCharge(cluster->Charge());
209 164 : rawCluster->SetChi2(cluster->Chi2());
210 :
211 164 : Double_t xg, yg, zg;
212 492 : fkTransformer.Local2Global(detElemId,
213 492 : cluster->Position().X(), cluster->Position().Y(),
214 : 0, xg, yg, zg);
215 164 : rawCluster->SetXYZ(xg, yg, zg);
216 164 : rawCluster->SetErrXY(recoParam->GetDefaultNonBendingReso(chamberId),recoParam->GetDefaultBendingReso(chamberId));
217 :
218 : // Set MC label
219 492 : if (fDigitStore && fDigitStore->HasMCInformation())
220 : {
221 164 : rawCluster->SetMCLabel(FindMCLabel(*cluster, detElemId, seg));
222 : }
223 :
224 820 : AliDebug(1,Form("Adding RawCluster detElemId %4d mult %2d charge %e (xl,yl,zl)=(%e,%e,%e) (xg,yg,zg)=(%e,%e,%e) label %d",
225 : detElemId,rawCluster->GetNDigits(),rawCluster->GetCharge(),
226 : cluster->Position().X(),cluster->Position().Y(),0.0,
227 : xg,yg,zg,rawCluster->GetMCLabel()));
228 164 : }
229 144 : }
230 144 : }
231 1248 : it.Next();
232 1248 : }
233 :
234 400 : AliDebug(1,Form("chamberId = %2d NofClusters after = %d",chamberId,fNCluster));
235 :
236 : return nofAddedClusters;
237 160 : }
238 :
239 :
240 : //_____________________________________________________________________________
241 : void
242 : AliMUONSimpleClusterServer::Global2Local(Int_t detElemId, const AliMpArea& globalArea,
243 : AliMpArea& localArea) const
244 : {
245 : /// Convert a global area in local area for a given DE
246 :
247 0 : Double_t xl,yl,zl;
248 :
249 0 : Int_t chamberId = AliMpDEManager::GetChamberId(detElemId);
250 0 : if ( chamberId < 0 ) {
251 0 : AliErrorStream() << "Cannot get chamberId from detElemId=" << detElemId << endl;
252 0 : return;
253 : }
254 0 : Double_t zg = AliMUONConstants::DefaultChamberZ(chamberId);
255 :
256 0 : fkTransformer.Global2Local(detElemId,
257 0 : globalArea.GetPositionX(),globalArea.GetPositionY(),zg,
258 : xl,yl,zl);
259 :
260 0 : localArea = AliMpArea(xl,yl, globalArea.GetDimensionX(), globalArea.GetDimensionY());
261 0 : }
262 :
263 : //_____________________________________________________________________________
264 : Bool_t
265 : AliMUONSimpleClusterServer::Overlap(Int_t detElemId,
266 : const AliMpArea& area,
267 : AliMpArea& deArea) const
268 : {
269 : /// Check whether (global) area overlaps with the given DE.
270 : /// If it is, set deArea to the overlap region and convert it
271 : /// in the local coordinate system of that DE.
272 :
273 : Bool_t overlap(kFALSE);
274 :
275 0 : AliMpArea* globalDEArea = fkTransformer.GetDEArea(detElemId);
276 :
277 0 : if (!globalDEArea) return kFALSE;
278 :
279 0 : AliMpArea overlapArea;
280 :
281 0 : if ( area.Overlap(*globalDEArea) )
282 : {
283 0 : overlapArea = area.Intersect(*globalDEArea);
284 0 : Global2Local(detElemId,overlapArea,deArea);
285 : overlap = kTRUE;
286 0 : }
287 : else
288 : {
289 0 : deArea = AliMpArea();
290 : }
291 :
292 0 : AliDebug(1,Form("DE %04d area %s globalDEArea %s overlapArea %s deArea %s overlap=%d",
293 : detElemId,
294 : AsString(area).Data(),
295 : AsString(*globalDEArea).Data(),
296 : AsString(overlapArea).Data(),
297 : AsString(deArea).Data(),
298 : overlap));
299 :
300 0 : return overlap;
301 0 : }
302 :
303 : //_____________________________________________________________________________
304 : TObjArray*
305 : AliMUONSimpleClusterServer::PadArray(Int_t detElemId, Int_t cathode) const
306 : {
307 : /// Return array for given cathode of given DE
308 :
309 2648 : return static_cast<TObjArray*>(fPads[cathode]->GetValue(detElemId));
310 : }
311 :
312 : //_____________________________________________________________________________
313 : Bool_t
314 : AliMUONSimpleClusterServer::UseTriggerTrackStore(AliMUONVTriggerTrackStore* trackStore)
315 : {
316 : /// Tells us to use trigger track store, and thus to bypass St45 clusters
317 0 : fTriggerTrackStore = trackStore; // not owner
318 0 : delete fBypass;
319 0 : fBypass = new AliMUONTriggerTrackToTrackerClusters(fkTransformer,fTriggerTrackStore);
320 0 : return kTRUE;
321 0 : }
322 :
323 : //_____________________________________________________________________________
324 : void
325 : AliMUONSimpleClusterServer::UseDigits(TIter& next, AliMUONVDigitStore* digitStore)
326 : {
327 : /// Convert digitStore into two arrays of AliMUONPads
328 :
329 16 : fDigitStore = digitStore;
330 :
331 : // Clear pads arrays in the maps
332 48 : for ( Int_t i=0; i<2; i++ ) {
333 16 : fPadsIterator[i]->Reset();
334 16 : Int_t key; TObject* obj;
335 672 : while ( ( obj = fPadsIterator[i]->Next(key) ) ) {
336 : //cout << "clearing array for detElemId " << key << " ";
337 320 : obj->Clear();
338 : }
339 16 : }
340 :
341 : AliMUONVDigit* d;
342 2980 : while ( ( d = static_cast<AliMUONVDigit*>(next()) ) )
343 : {
344 1640 : if ( ! (d->Charge() > 0.) ) continue; // skip void digits.
345 1548 : if ( ! d->IsTracker() ) continue; // skip trigger digits
346 1324 : Int_t ix = d->PadX();
347 1324 : Int_t iy = d->PadY();
348 1324 : Int_t cathode = d->Cathode();
349 1324 : Int_t detElemId = d->DetElemId();
350 2648 : const AliMpVSegmentation* seg = AliMpSegmentation::Instance()->
351 1324 : GetMpSegmentation(detElemId,AliMp::GetCathodType(cathode));
352 1324 : AliMpPad pad = seg->PadByIndices(ix,iy);
353 :
354 1324 : TObjArray* padArray = PadArray(detElemId,cathode);
355 1324 : if (!padArray)
356 : {
357 448 : padArray = new TObjArray(100);
358 224 : padArray->SetOwner(kTRUE);
359 224 : fPads[cathode]->Add(detElemId,padArray);
360 : }
361 :
362 3972 : AliMUONPad* mpad = new AliMUONPad(detElemId,cathode,
363 1324 : ix,iy,pad.GetPositionX(),pad.GetPositionY(),
364 1324 : pad.GetDimensionX(),pad.GetDimensionY(),
365 2648 : d->Charge());
366 2648 : if ( d->IsSaturated() ) mpad->SetSaturated(kTRUE);
367 2648 : mpad->SetUniqueID(d->GetUniqueID());
368 1324 : padArray->Add(mpad);
369 1324 : }
370 8 : }
371 :
372 : //_____________________________________________________________________________
373 : Int_t
374 : AliMUONSimpleClusterServer::FindMCLabel(const AliMUONCluster& cluster, Int_t detElemId, const AliMpVSegmentation* seg[2]) const
375 : {
376 : /// Find the label of the most contributing MC track (-1 in case of failure)
377 : /// The data member fDigitStore must be set
378 :
379 : // --- get the digit (if any) located at the cluster position on both cathods ---
380 164 : Int_t nTracks[2] = {0, 0};
381 82 : AliMUONVDigit* digit[2] = {0x0, 0x0};
382 492 : for (Int_t iCath = 0; iCath < 2; iCath++) {
383 164 : AliMpPad pad
384 492 : = seg[AliMp::GetCathodType(iCath)]->PadByPosition(cluster.Position().X(), cluster.Position().Y(),kFALSE);
385 164 : if (pad.IsValid()) {
386 656 : digit[iCath] = fDigitStore->FindObject(detElemId, pad.GetManuId(), pad.GetManuChannel(), iCath);
387 492 : if (digit[iCath]) nTracks[iCath] = digit[iCath]->Ntracks();
388 : }
389 164 : }
390 :
391 82 : if (nTracks[0] + nTracks[1] == 0) return -1;
392 :
393 : // --- build the list of contributing tracks and of the associated charge ---
394 82 : Int_t* trackId = new Int_t[nTracks[0] + nTracks[1]];
395 82 : Float_t* trackCharge = new Float_t[nTracks[0] + nTracks[1]];
396 : Int_t nTracksTot = 0;
397 :
398 : // fill with contributing tracks on first cathod
399 330 : for (Int_t iTrack1 = 0; iTrack1 < nTracks[0]; iTrack1++) {
400 83 : trackId[iTrack1] = digit[0]->Track(iTrack1);
401 83 : trackCharge[iTrack1] = digit[0]->TrackCharge(iTrack1);
402 : }
403 : nTracksTot = nTracks[0];
404 :
405 : // complement with contributing tracks on second cathod
406 330 : for (Int_t iTrack2 = 0; iTrack2 < nTracks[1]; iTrack2++) {
407 83 : Int_t trackId2 = digit[1]->Track(iTrack2);
408 : // check if track exist
409 : Bool_t trackExist = kFALSE;
410 168 : for (Int_t iTrack1 = 0; iTrack1 < nTracks[0]; iTrack1++) {
411 84 : if (trackId2 == trackId[iTrack1]) {
412 : // complement existing track
413 83 : trackCharge[iTrack1] += digit[1]->TrackCharge(iTrack2);
414 : trackExist = kTRUE;
415 83 : break;
416 : }
417 : }
418 : // add the new track
419 83 : if (!trackExist) {
420 0 : trackId[nTracksTot] = trackId2;
421 0 : trackCharge[nTracksTot] = digit[1]->TrackCharge(iTrack2);
422 0 : nTracksTot++;
423 0 : }
424 : }
425 :
426 : // --- Find the most contributing track ---
427 : Int_t mainTrackId = -1;
428 : Float_t maxCharge = 0.;
429 330 : for (Int_t iTrack = 0; iTrack < nTracksTot; iTrack++) {
430 83 : if (trackCharge[iTrack] > maxCharge) {
431 83 : mainTrackId = trackId[iTrack];
432 : maxCharge = trackCharge[iTrack];
433 83 : }
434 : }
435 :
436 164 : delete[] trackId;
437 164 : delete[] trackCharge;
438 :
439 : return mainTrackId;
440 82 : }
441 :
442 : //_____________________________________________________________________________
443 : void
444 : AliMUONSimpleClusterServer::Print(Option_t*) const
445 : {
446 : /// Printout for debug only
447 :
448 0 : AliMpDEIterator it;
449 :
450 0 : it.First();
451 :
452 0 : while ( !it.IsDone() )
453 : {
454 0 : Int_t detElemId = it.CurrentDEId();
455 :
456 : // printout the number of pads / de, and number of used pads / de
457 :
458 0 : if ( ( PadArray(detElemId,0) && PadArray(detElemId,0)->GetLast() >= 0 ) ||
459 0 : ( PadArray(detElemId,1) && PadArray(detElemId,1)->GetLast() >= 0 ) )
460 : {
461 0 : cout << Form("---- DE %04d",detElemId) << endl;
462 :
463 0 : for ( Int_t cathode = 0; cathode < 2; ++cathode )
464 : {
465 0 : cout << Form(" -- Cathode %1d",cathode) << endl;
466 :
467 0 : TObjArray* padArray = PadArray(detElemId,cathode);
468 :
469 0 : if (!padArray)
470 : {
471 0 : cout << "no pad array" << endl;
472 : }
473 0 : else if ( padArray->GetLast() < 0 )
474 : {
475 0 : cout << "no pads" << endl;
476 : }
477 : else
478 : {
479 0 : TIter next(padArray);
480 : AliMUONPad* pad;
481 0 : while ( ( pad = static_cast<AliMUONPad*>(next()) ) )
482 : {
483 0 : pad->Print("full");
484 : }
485 0 : }
486 : }
487 0 : }
488 0 : it.Next();
489 : }
490 0 : }
491 :
492 :
|
