Bug fix in the constructor (thanks to A.Marin)
[u/mrichter/AliRoot.git] / TRD / AliTRDReconstructor.cxx
CommitLineData
121a60bd 1/**************************************************************************
91834b8d 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**************************************************************************/
121a60bd 15
16/* $Id$ */
17
18///////////////////////////////////////////////////////////////////////////////
19// //
4fad09c9 20// Class for TRD reconstruction //
954135b1 21// //
22// For the special options which can be used during reconstruction and their //
23// default values pls. see function SetOption(). //
121a60bd 24// //
25///////////////////////////////////////////////////////////////////////////////
26
3a039a31 27#include <TObjString.h>
28#include <TObjArray.h>
d39401ad 29#include <TTreeStream.h>
30#include <TDirectory.h>
cb8b99ee 31#include <TRef.h>
121a60bd 32
3bc9d03e 33#include "AliRawReader.h"
cb8b99ee 34#include "AliLog.h"
3bc9d03e 35
36#include "AliTRDReconstructor.h"
3fe61b77 37#include "AliTRDclusterizer.h"
a8fbb161 38#include "AliTRDrawData.h"
8df1f8f5 39#include "AliTRDrawStream.h"
a8fbb161 40#include "AliTRDdigitsManager.h"
e4f2f73d 41#include "AliTRDtrackerV1.h"
cb8b99ee 42#include "AliESDEvent.h"
43#include "AliESDTrdTrack.h"
44#include "AliESDTrdTracklet.h"
45#include "AliTRDtrackletWord.h"
29f95561 46
97b6dfac 47#define SETFLG(n,f) ((n) |= f)
48#define CLRFLG(n,f) ((n) &= ~f)
49
121a60bd 50ClassImp(AliTRDReconstructor)
51
828c6f80 52TClonesArray *AliTRDReconstructor::fgClusters = NULL;
53TClonesArray *AliTRDReconstructor::fgTracklets = NULL;
cb8b99ee 54TClonesArray *AliTRDReconstructor::fgTracks = NULL;
828c6f80 55Char_t const * AliTRDReconstructor::fgSteerNames[kNsteer] = {
b72f4eaf 56 "DigitsConversion "
b72f4eaf 57 ,"Write Clusters "
58 ,"Write Online Tracklets "
b72f4eaf 59 ,"Stand Alone Tracking "
8fc736d7 60 ,"HLT Mode "
61 ,"Process Online Trklts "
62 ,"Debug Streaming "
63 ,"Cl. Radial Correction "
b72f4eaf 64};
828c6f80 65Char_t const * AliTRDReconstructor::fgSteerFlags[kNsteer] = {
b72f4eaf 66 "dc"// digits conversion [false]
b72f4eaf 67 ,"cw"// write clusters [true]
68 ,"tw"// write online tracklets [false]
b72f4eaf 69 ,"sa"// track seeding (stand alone tracking) [true]
b72f4eaf 70 ,"hlt"// HLT reconstruction [false]
a5b99acd 71 ,"tp"// also use online tracklets for reconstruction [false]
a2fbb6ec 72 ,"deb"// Write debug stream [false]
b30d8c09 73 ,"cc" // Cluster radial correction during reconstruction [false]
b72f4eaf 74};
828c6f80 75Char_t const * AliTRDReconstructor::fgTaskNames[AliTRDrecoParam::kTRDreconstructionTasks] = {
a2fbb6ec 76 "Clusterizer"
b72f4eaf 77 ,"Tracker"
78 ,"PID"
79};
828c6f80 80Char_t const * AliTRDReconstructor::fgTaskFlags[AliTRDrecoParam::kTRDreconstructionTasks] = {
a2fbb6ec 81 "cl"
b72f4eaf 82 ,"tr"
83 ,"pd"
84};
c6f7c6cb 85Int_t AliTRDReconstructor::fgNTimeBins = -1;
952051c5 86const Float_t AliTRDReconstructor::fgkMinClustersInTrack = 0.5; //
87const Float_t AliTRDReconstructor::fgkLabelFraction = 0.8; //
88const Double_t AliTRDReconstructor::fgkMaxChi2 = 12.0; //
89const Double_t AliTRDReconstructor::fgkMaxSnp = 0.95; // Maximum local sine of the azimuthal angle
90const Double_t AliTRDReconstructor::fgkMaxStep = 2.0; // Maximal step size in propagation
91const Double_t AliTRDReconstructor::fgkEpsilon = 1.e-5; // Precision of radial coordinate
b72f4eaf 92
eb38ed55 93//_____________________________________________________________________________
3a039a31 94AliTRDReconstructor::AliTRDReconstructor()
95 :AliReconstructor()
97b6dfac 96 ,fSteerParam(0)
5127281e 97 ,fClusterizer(NULL)
3a039a31 98{
97b6dfac 99 // setting default "ON" steering parameters
89971837 100 // owner of debug streamers
101 SETFLG(fSteerParam, kOwner);
97b6dfac 102 // write clusters [cw]
103 SETFLG(fSteerParam, kWriteClusters);
104 // track seeding (stand alone tracking) [sa]
9dcc64cc 105 //SETFLG(fSteerParam, kSeeding);
8fc736d7 106 // Cluster radial correction during reconstruction [cc]
b30d8c09 107 //SETFLG(fSteerParam, kClRadialCorr);
a2fbb6ec 108 memset(fDebugStream, 0, sizeof(TTreeSRedirector *) * AliTRDrecoParam::kTRDreconstructionTasks);
eb38ed55 109}
110
3a039a31 111//_____________________________________________________________________________
48f8adf3 112AliTRDReconstructor::~AliTRDReconstructor()
113{
d39401ad 114 //
115 // Destructor
116 //
117
48f8adf3 118 if(fgClusters) {
6d54fcbf 119 fgClusters->Delete();
120 delete fgClusters;
121 fgClusters = NULL;
48f8adf3 122 }
a5b99acd 123 if(fgTracklets) {
6d54fcbf 124 fgTracklets->Delete();
125 delete fgTracklets;
126 fgTracklets = NULL;
a5b99acd 127 }
cb8b99ee 128 if(fgTracks) {
129 fgTracks->Delete();
130 delete fgTracks;
131 fgTracks = NULL;
132 }
f29f13a6 133 if(fSteerParam&kOwner){
a2fbb6ec 134 for(Int_t itask = 0; itask < AliTRDrecoParam::kTRDreconstructionTasks; itask++)
29f95561 135 if(fDebugStream[itask]) delete fDebugStream[itask];
136 }
5127281e 137 if(fClusterizer){
138 delete fClusterizer;
139 fClusterizer = NULL;
140 }
18c05eb3 141}
48f8adf3 142
91834b8d 143
48f8adf3 144//_____________________________________________________________________________
91834b8d 145void AliTRDReconstructor::Init(){
18c05eb3 146 //
147 // Init Options
148 //
18c05eb3 149 SetOption(GetOption());
a2fbb6ec 150 Options(fSteerParam);
151
5127281e 152 if(!fClusterizer){
153 fClusterizer = new AliTRDclusterizer(fgTaskNames[AliTRDrecoParam::kClusterizer], fgTaskNames[AliTRDrecoParam::kClusterizer]);
154 fClusterizer->SetReconstructor(this);
155 }
156
a2fbb6ec 157 // Make Debug Streams when Debug Streaming
158 if(IsDebugStreaming()){
159 for(Int_t task = 0; task < AliTRDrecoParam::kTRDreconstructionTasks; task++){
160 TDirectory *savedir = gDirectory;
161 fDebugStream[task] = new TTreeSRedirector(Form("TRD.Debug%s.root", fgTaskNames[task]));
162 savedir->cd();
163 SETFLG(fSteerParam, kOwner);
164 }
165 }
acd241e9 166}
0397ce53 167
121a60bd 168//_____________________________________________________________________________
a8fbb161 169void AliTRDReconstructor::ConvertDigits(AliRawReader *rawReader
91834b8d 170 , TTree *digitsTree) const
121a60bd 171{
3bc9d03e 172 //
25ca55ce 173 // Convert raw data digits into digit objects in a root tree
3bc9d03e 174 //
121a60bd 175
ad5dd6ec 176 //AliInfo("Convert raw data digits into digit objects [RawReader -> Digit TTree]");
91834b8d 177
a8fbb161 178 AliTRDrawData rawData;
3a0f6479 179 rawReader->Reset();
180 rawReader->Select("TRD");
a8fbb161 181 AliTRDdigitsManager *manager = rawData.Raw2Digits(rawReader);
182 manager->MakeBranch(digitsTree);
183 manager->WriteDigits();
c0ab62ff 184 delete manager;
0ee00e25 185
cb8b99ee 186 // take over ownership of online tracklets
187 fgTracklets = rawData.TrackletsArray();
188 rawData.SetTrackletsOwner(0x0);
189
190 // take over GTU tracks
191 fgTracks = rawData.TracksArray();
192 rawData.SetTracksOwner(0x0);
121a60bd 193}
194
195//_____________________________________________________________________________
25ca55ce 196void AliTRDReconstructor::Reconstruct(AliRawReader *rawReader
197 , TTree *clusterTree) const
198{
199 //
200 // Reconstruct clusters
201 //
202
91834b8d 203 //AliInfo("Reconstruct TRD clusters from RAW data [RawReader -> Cluster TTree]");
204
205
3a0f6479 206 rawReader->Reset();
207 rawReader->Select("TRD");
208
5127281e 209 if(!fClusterizer){
210 AliFatal("Clusterizer not available!");
211 return;
212 }
213
da0a4b87 214 fClusterizer->ResetRecPoints();
215
5127281e 216 fClusterizer->OpenOutput(clusterTree);
5127281e 217 fClusterizer->SetUseLabels(kFALSE);
218 fClusterizer->Raw2ClustersChamber(rawReader);
66f6bfd9 219
3e1f6592 220 fgNTimeBins = fClusterizer->GetNTimeBins();
221
cb8b99ee 222 // take over ownership of online tracklets
223 fgTracklets = fClusterizer->TrackletsArray();
224 fClusterizer->SetTrackletsOwner(kFALSE);
225
226 // take over GTU tracks
227 fgTracks = fClusterizer->TracksArray();
228 fClusterizer->SetTracksOwner(kFALSE);
229
66f6bfd9 230 if(IsWritingClusters()) return;
d6e36038 231
66f6bfd9 232 // take over ownership of clusters
5127281e 233 fgClusters = fClusterizer->RecPoints();
234 fClusterizer->SetClustersOwner(kFALSE);
25ca55ce 235}
236
237//_____________________________________________________________________________
a8fbb161 238void AliTRDReconstructor::Reconstruct(TTree *digitsTree
239 , TTree *clusterTree) const
25ca55ce 240{
241 //
242 // Reconstruct clusters
243 //
e4f2f73d 244
91834b8d 245 //AliInfo("Reconstruct TRD clusters from Digits [Digit TTree -> Cluster TTree]");
a2fbb6ec 246
c665eaa7 247 AliTRDclusterizer clusterer(fgTaskNames[AliTRDrecoParam::kClusterizer], fgTaskNames[AliTRDrecoParam::kClusterizer]);
248 clusterer.SetReconstructor(this);
249 clusterer.OpenOutput(clusterTree);
250 clusterer.ReadDigits(digitsTree);
251 clusterer.MakeClusters();
d6e36038 252
cb8b99ee 253 // read tracklets and tracks if not done during reading of raw data
254 if (!fgTracklets) {
255 clusterer.ReadTracklets();
256 fgTracklets = clusterer.TrackletsArray();
257 clusterer.SetTrackletsOwner(kFALSE);
258 }
259 if (!fgTracks) {
260 clusterer.ReadTracks();
261 fgTracks = clusterer.TracksArray();
262 clusterer.SetTracksOwner(kFALSE);
263 }
264
3e1f6592 265 fgNTimeBins = clusterer.GetNTimeBins();
266
66f6bfd9 267 if(IsWritingClusters()) return;
268
269 // take over ownership of clusters
c665eaa7 270 fgClusters = clusterer.RecPoints();
271 clusterer.SetClustersOwner(kFALSE);
25ca55ce 272}
273
274//_____________________________________________________________________________
44dbae42 275AliTracker *AliTRDReconstructor::CreateTracker() const
121a60bd 276{
3bc9d03e 277 //
278 // Create a TRD tracker
279 //
121a60bd 280
e4f2f73d 281 //return new AliTRDtracker(NULL);
3a039a31 282 AliTRDtrackerV1 *tracker = new AliTRDtrackerV1();
283 tracker->SetReconstructor(this);
284 return tracker;
25ca55ce 285
25ca55ce 286}
287
288//_____________________________________________________________________________
289void AliTRDReconstructor::FillESD(TTree* /*digitsTree*/
91834b8d 290 , TTree* /*clusterTree*/
cb8b99ee 291 , AliESDEvent* esd) const
25ca55ce 292{
293 //
af26ce80 294 // Fill ESD
25ca55ce 295 //
296
cb8b99ee 297 // ----- filling tracklets -----
298 Int_t trackletIndex[1080] = { 0 };
299 AliDebug(1, Form("Filling tracklets from %p (%i)",
300 fgTracklets, fgTracklets ? fgTracklets->GetEntriesFast() : 0));
301 if (fgTracklets) {
302 Int_t nTracklets = fgTracklets->GetEntriesFast();
303
304 Int_t lastHC = -1;
305 for (Int_t iTracklet = 0; iTracklet < nTracklets; iTracklet++) {
306 AliTRDtrackletBase *trkl = (AliTRDtrackletBase*) ((*fgTracklets)[iTracklet]);
307 Int_t hc = trkl->GetHCId();
308 // hc += trkl->GetY() > 0 ? 1 : 0;
309 if ((hc < 0) || (hc >= 1080)) {
310 AliError(Form("HC for tracklet: 0x%08x out of range: %i", trkl->GetTrackletWord(), trkl->GetHCId()));
311 continue;
312}
313 AliDebug(5, Form("hc: %4i : 0x%08x z: %2i", hc, trkl->GetTrackletWord(), trkl->GetZbin()));
314 if (hc != lastHC) {
315 AliDebug(2, Form("set tracklet index for HC %i to %i", hc, iTracklet));
316 trackletIndex[hc] = iTracklet + 1;
317 lastHC = hc;
318 }
319 }
320
321 for (Int_t iDet = 0; iDet < 540; iDet++) {
322 Int_t trklIndexA = trackletIndex[2*iDet + 0] - 1;
323 Int_t trklIndexB = trackletIndex[2*iDet + 1] - 1;
324 Int_t trklIndex = esd->GetNumberOfTrdTracklets();
325 AliTRDtrackletBase *trklA = trklIndexA > -1 ? (AliTRDtrackletBase*) ((*fgTracklets)[trklIndexA]) : 0x0;
326 AliTRDtrackletBase *trklB = trklIndexB > -1 ? (AliTRDtrackletBase*) ((*fgTracklets)[trklIndexB]) : 0x0;
327 while (trklA != 0x0 || trklB != 0x0) {
328 AliDebug(5, Form("det %i - A: %i/%i -> %p, B: %i/%i -> %p",
329 iDet, trklIndexA, nTracklets, trklA, trklIndexB, nTracklets, trklB));
330 if (trklA == 0x0) {
331 AliESDTrdTracklet esdTracklet(trklB->GetTrackletWord(), trklB->GetHCId());
332 esd->AddTrdTracklet(&esdTracklet);
333 trklIndexB++;
334 trklB = trklIndexB < nTracklets ? (AliTRDtrackletBase*) ((*fgTracklets)[trklIndexB]) : 0x0;
335 if (trklB && trklB->GetHCId() != 2*iDet + 1)
336 trklB = 0x0;
337 }
338 else if (trklB == 0x0) {
339 AliESDTrdTracklet esdTracklet(trklA->GetTrackletWord(), trklA->GetHCId());
340 esd->AddTrdTracklet(&esdTracklet);
341 trklIndexA++;
342 trklA = trklIndexA < nTracklets ? (AliTRDtrackletBase*) ((*fgTracklets)[trklIndexA]) : 0x0;
343 if (trklA && trklA->GetHCId() != 2*iDet)
344 trklA = 0x0;
345 }
346 else {
347 if ((trklA->GetZbin() < trklB->GetZbin()) ||
348 ((trklA->GetZbin() == trklB->GetZbin()) && (trklA->GetYbin() < trklB->GetYbin()))) {
349 AliESDTrdTracklet esdTracklet(trklA->GetTrackletWord(), trklA->GetHCId());
350 esd->AddTrdTracklet(&esdTracklet);
351 trklIndexA++;
352 trklA = trklIndexA < nTracklets ? (AliTRDtrackletBase*) ((*fgTracklets)[trklIndexA]) : 0x0;
353 if (trklA && trklA->GetHCId() != 2*iDet)
354 trklA = 0x0;
355 }
356 else {
357 AliESDTrdTracklet esdTracklet(trklB->GetTrackletWord(), trklB->GetHCId());
358 esd->AddTrdTracklet(&esdTracklet);
359 trklIndexB++;
360 trklB = trklIndexB < nTracklets ? (AliTRDtrackletBase*) ((*fgTracklets)[trklIndexB]) : 0x0;
361 if (trklB && trklB->GetHCId() != 2*iDet + 1)
362 trklB = 0x0;
363 }
364 }
365 }
366 // updating tracklet indices as in ESD
367 if (esd->GetNumberOfTrdTracklets() != trklIndex) {
368 trackletIndex[2*iDet + 0] = trackletIndex[2*iDet + 1] = trklIndex;
369 }
370 else
371 trackletIndex[2*iDet + 0] = trackletIndex[2*iDet + 1] = -1;
372 }
373 }
374
375 // ----- filling GTU tracks -----
376 AliDebug(1, Form("Now filling ESD with GTU tracks from %p (%i)",
377 fgTracks, fgTracks ? fgTracks->GetEntriesFast() : 0));
378 if (fgTracks) {
379 for (Int_t iTrack = 0; iTrack < fgTracks->GetEntriesFast(); iTrack++) {
380 AliESDTrdTrack *trdTrack = (AliESDTrdTrack*) ((*fgTracks)[iTrack]);
381
382 UInt_t mask = trdTrack->GetLayerMask();
383 UInt_t stack = trdTrack->GetStack();
384
385 for (Int_t iLayer = 0; iLayer < 6; iLayer++) {
386 if (mask & (1 << iLayer)) {
387
388 Int_t det = trdTrack->GetSector()*30 + stack*6 + iLayer;
389 Int_t idx = trdTrack->GetTrackletIndex(iLayer);
390
391 if ((det < 0) || (det > 539)) {
392 AliError(Form("Invalid detector no. from track: %i", 2*det));
393 continue;
394 }
395 if (trackletIndex[2*det] >= 0) {
396 if ((trackletIndex[2*det] + idx > -1) &&
397 (trackletIndex[2*det] + idx < esd->GetNumberOfTrdTracklets())) {
398 AliESDTrdTracklet *trkl = esd->GetTrdTracklet(trackletIndex[2*det] + idx);
399 if (trkl) {
400 AliDebug(5, Form("adding tracklet: 0x%08x", trkl->GetTrackletWord()));
401 trdTrack->AddTrackletReference(trkl, iLayer);
402 }
403 }
404 }
405 }
406 }
407 // only add the track when it's complete (including tracklet references)
408 esd->AddTrdTrack(trdTrack);
409 }
410 }
411
412 // clearing variables for next event
413 fgTracklets = 0x0;
414 fgTracks = 0x0;
25ca55ce 415}
3a039a31 416
828c6f80 417//_____________________________________________________________________________
3a039a31 418void AliTRDReconstructor::SetOption(Option_t *opt)
419{
d39401ad 420 //
421 // Read option string into the steer param.
422 //
954135b1 423 // The following string options are available during reconstruction.
424 // In square brackets the default values are given.
425 // "dc" : digits conversion [false]
426 // "cw" : write clusters [true]
427 // "tw" : write online tracklets [false]
428 // "sa" : track seeding (stand alone tracking) [true]
429 // "hlt" : HLT reconstruction [false]
430 // "tp" : also use online tracklets for reconstruction [false]
431 // "deb" : Write debug stream [false]
432 // "cc" : Cluster radial correction during reconstruction [false]
433 //
434 // To check the actual options used during reconstruction include the following line in your rec.C script
435 // AliLog::SetClassDebugLevel("AliTRDReconstructor", 1);
97b6dfac 436
437 AliReconstructor::SetOption(opt);
3a039a31 438
439 TString s(opt);
440 TObjArray *opar = s.Tokenize(",");
441 for(Int_t ipar=0; ipar<opar->GetEntriesFast(); ipar++){
d39401ad 442 Bool_t processed = kFALSE;
3a039a31 443 TString sopt(((TObjString*)(*opar)[ipar])->String());
b72f4eaf 444 for(Int_t iopt=0; iopt<kNsteer; iopt++){
445 if(!sopt.Contains(fgSteerFlags[iopt])) continue;
446 SETFLG(fSteerParam, BIT(iopt));
447 if(sopt.Contains("!")) CLRFLG(fSteerParam, BIT(iopt));
d39401ad 448 processed = kTRUE;
b72f4eaf 449 break;
450 }
d39401ad 451 if(processed) continue;
b72f4eaf 452
453 AliWarning(Form("Unknown option flag %s.", sopt.Data()));
3a039a31 454 }
455}
456
acd241e9 457//_____________________________________________________________________________
a2fbb6ec 458void AliTRDReconstructor::Options(UInt_t steer)
b72f4eaf 459{
d39401ad 460 //
461 // Print the options
462 //
463
b72f4eaf 464 for(Int_t iopt=0; iopt<kNsteer; iopt++){
9e8fdff4 465 AliDebugGeneral("AliTRDReconstructor", 1, Form(" %s[%s]%s", fgSteerNames[iopt], fgSteerFlags[iopt], steer ?(((steer>>iopt)&1)?" : ON":" : OFF"):""));
b72f4eaf 466 }
b72f4eaf 467}
468