New TrackReconstructor for Kalman tracking (Philippe Pillot)
[u/mrichter/AliRoot.git] / MUON / AliMUONReconstructor.cxx
CommitLineData
cac2eb58 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 **************************************************************************/
cac2eb58 15/* $Id$ */
16
9265505b 17//-----------------------------
18// Class AliMUONReconstructor //
19//----------------------------- //
20// Class for the
21// MUON track reconstruction
30178c30 22
cf464691 23#include "AliMUONReconstructor.h"
b2d7df0b 24
25#include "AliESD.h"
26#include "AliESDMuonTrack.h"
94f6fba9 27#include "AliLog.h"
b2d7df0b 28#include "AliMUON.h"
96fdfe9a 29#include "AliMUONConstants.h"
b2d7df0b 30#include "AliMUONCalibrationData.h"
346fb008 31#include "AliMUONClusterFinderAZ.h"
b2d7df0b 32#include "AliMUONClusterReconstructor.h"
33#include "AliMUONData.h"
34#include "AliMUONDigitCalibrator.h"
cc87ebcd 35#include "AliMUONEventRecoCombi.h"
2cf44ef3 36#include "AliMUONDigitMaker.h"
cac2eb58 37#include "AliMUONTrack.h"
38#include "AliMUONTrackParam.h"
b2d7df0b 39#include "AliMUONTrackReconstructor.h"
cac2eb58 40#include "AliMUONTriggerTrack.h"
96fdfe9a 41#include "AliMUONTriggerCircuitNew.h"
42#include "AliMUONTriggerCrateStore.h"
9265505b 43
44#include "AliMpSegmentation.h"
96fdfe9a 45
cf464691 46#include "AliRawReader.h"
b2d7df0b 47#include "AliRun.h"
48#include "AliRunLoader.h"
49#include "TTask.h"
94f6fba9 50#include "TStopwatch.h"
29f1b13a 51
9265505b 52/// \cond CLASSIMP
cac2eb58 53ClassImp(AliMUONReconstructor)
9265505b 54/// \endcond
b2d7df0b 55
cac2eb58 56//_____________________________________________________________________________
57AliMUONReconstructor::AliMUONReconstructor()
aa6ecf89 58 : AliReconstructor(),
59 fRunLoader(0x0),
60 fDigitMaker(new AliMUONDigitMaker()),
96fdfe9a 61 fCalibrationData(0x0),
62 fCrateManager(new AliMUONTriggerCrateStore()),
9c4b1ee7 63 fTriggerCircuit(new TClonesArray("AliMUONTriggerCircuitNew", 234)),
9c4b1ee7 64 fTransformer(new AliMUONGeometryTransformer(kTRUE))
65
8789635b 66{
d19b6003 67/// Default constructor
68
94f6fba9 69 AliDebug(1,"");
96fdfe9a 70 // Crate manager
71 fCrateManager->ReadFromFile();
72
73 // set to digit maker
74 fDigitMaker->SetCrateManager(fCrateManager);
75
9c4b1ee7 76 // transformater
77 fTransformer->ReadGeometryData("volpath.dat", "geometry.root");
78
96fdfe9a 79 // trigger circuit
80 for (Int_t i = 0; i < AliMUONConstants::NTriggerCircuit(); i++) {
81 AliMUONTriggerCircuitNew* c = new AliMUONTriggerCircuitNew();
9c4b1ee7 82 c->SetTransformer(fTransformer);
96fdfe9a 83 c->Init(i,*fCrateManager);
84 TClonesArray& circuit = *fTriggerCircuit;
85 new(circuit[circuit.GetEntriesFast()])AliMUONTriggerCircuitNew(*c);
86 delete c;
9c4b1ee7 87 }
88
89
8789635b 90}
94f6fba9 91
8789635b 92//_____________________________________________________________________________
93AliMUONReconstructor::~AliMUONReconstructor()
cac2eb58 94{
d19b6003 95/// Destructor
96
94f6fba9 97 AliDebug(1,"");
b2d7df0b 98 delete fCalibrationData;
aa6ecf89 99 delete fDigitMaker;
96fdfe9a 100 delete fCrateManager;
101 delete fTriggerCircuit;
9c4b1ee7 102 delete fTransformer;
cac2eb58 103}
b2d7df0b 104
105//_____________________________________________________________________________
106TTask*
107AliMUONReconstructor::GetCalibrationTask(AliMUONData* data) const
108{
d19b6003 109/// Create the calibration task(s).
b2d7df0b 110
111 const AliRun* run = fRunLoader->GetAliRun();
9f9729c3 112
113 AliInfo("Calibration will occur.");
114 Int_t runNumber = run->GetRunNumber();
115 fCalibrationData = new AliMUONCalibrationData(runNumber);
116 if ( !fCalibrationData->IsValid() )
b2d7df0b 117 {
118 AliError("Could not retrieve calibrations !");
119 delete fCalibrationData;
120 fCalibrationData = 0x0;
121 return 0x0;
122 }
9f9729c3 123 TTask* calibration = new TTask("MUONCalibrator","MUON Digit calibrator");
124 calibration->Add(new AliMUONDigitCalibrator(data,fCalibrationData));
125 //FIXME: calibration->Add(something about dead channels should go here).
126 return calibration;
127
b2d7df0b 128}
129
130//_____________________________________________________________________________
131void
132AliMUONReconstructor::Init(AliRunLoader* runLoader)
133{
d19b6003 134/// Initialize
135
b2d7df0b 136 fRunLoader = runLoader;
137}
138
cac2eb58 139//_____________________________________________________________________________
140void AliMUONReconstructor::Reconstruct(AliRunLoader* runLoader) const
141{
d19b6003 142/// Reconstruct
143/// \todo add more
144
cac2eb58 145 AliLoader* loader = runLoader->GetLoader("MUONLoader");
146 Int_t nEvents = runLoader->GetNumberOfEvents();
147
3bc8b580 148 AliMUONData* data = new AliMUONData(loader,"MUON","MUON");
149
150// passing loader as argument.
9c4b1ee7 151 AliMUONTrackReconstructor* recoEvent = new AliMUONTrackReconstructor(data);
96fdfe9a 152 recoEvent->SetTriggerCircuit(fTriggerCircuit);
3bc8b580 153
154 if (strstr(GetOption(),"Original"))
9c4b1ee7 155 recoEvent->SetTrackMethod(kOriginal); // Original tracking
3bc8b580 156 else if (strstr(GetOption(),"Combi"))
9c4b1ee7 157 recoEvent->SetTrackMethod(kCombi); // Combined cluster / track
3bc8b580 158 else
9c4b1ee7 159 recoEvent->SetTrackMethod(kKalman); // Kalman
3bc8b580 160
94f6fba9 161 AliMUONClusterReconstructor* recoCluster = new AliMUONClusterReconstructor(data);
3bc8b580 162
7e4a628d 163 AliMUONClusterFinderVS *recModel = recoCluster->GetRecoModel();
3bc8b580 164
8e0ae46c 165 if (!strstr(GetOption(),"VS")) {
1af223d7 166 recModel = (AliMUONClusterFinderVS*) new AliMUONClusterFinderAZ();
9ffe3ef4 167 recoCluster->SetRecoModel(recModel);
9ffe3ef4 168 }
7e4a628d 169 recModel->SetGhostChi2Cut(10);
cac2eb58 170
171 loader->LoadDigits("READ");
172 loader->LoadRecPoints("RECREATE");
173 loader->LoadTracks("RECREATE");
cc87ebcd 174
b2d7df0b 175 TTask* calibration = GetCalibrationTask(data);
176
9c4b1ee7 177 Int_t chBeg = (recoEvent->GetTrackMethod() == kCombi ? 6 : 0);
178
cc87ebcd 179 // Loop over events
cac2eb58 180 for(Int_t ievent = 0; ievent < nEvents; ievent++) {
94f6fba9 181
182 AliDebug(1,Form("Event %d",ievent));
183
cac2eb58 184 runLoader->GetEvent(ievent);
185
8210cb40 186 //----------------------- digit2cluster & Trigger2Trigger -------------------
cac2eb58 187 if (!loader->TreeR()) loader->MakeRecPointsContainer();
cc87ebcd 188
cac2eb58 189 // tracking branch
9c4b1ee7 190 if (recoEvent->GetTrackMethod() != kCombi) {
3bc8b580 191 data->MakeBranch("RC");
192 data->SetTreeAddress("D,RC");
cc87ebcd 193 } else {
3bc8b580 194 data->SetTreeAddress("D");
195 data->SetTreeAddress("RCC");
cc87ebcd 196 }
addbf88a 197 // Important for avoiding a memory leak when reading digits ( to be investigated more in detail)
cc87ebcd 198 // In any case the reading of GLT is needed for the Trigger2Tigger method below
3bc8b580 199 data->SetTreeAddress("GLT");
cc87ebcd 200
3bc8b580 201 data->GetDigits();
b2d7df0b 202
203 if ( calibration )
204 {
205 calibration->ExecuteTask();
206 }
207
cc87ebcd 208 recoCluster->Digits2Clusters(chBeg);
b2d7df0b 209
9c4b1ee7 210 if (recoEvent->GetTrackMethod() == kCombi) {
cc87ebcd 211 // Combined cluster / track finder
3bc8b580 212 AliMUONEventRecoCombi::Instance()->FillEvent(data, (AliMUONClusterFinderAZ*)recModel);
cc87ebcd 213 ((AliMUONClusterFinderAZ*) recModel)->SetReco(2);
214 }
3bc8b580 215 else data->Fill("RC");
cac2eb58 216
217 // trigger branch
3bc8b580 218 data->MakeBranch("TC");
219 data->SetTreeAddress("TC");
cc87ebcd 220 recoCluster->Trigger2Trigger();
3bc8b580 221 data->Fill("TC");
cac2eb58 222
cc87ebcd 223 //AZ loader->WriteRecPoints("OVERWRITE");
cac2eb58 224
225 //---------------------------- Track & TriggerTrack ---------------------
226 if (!loader->TreeT()) loader->MakeTracksContainer();
227
228 // trigger branch
3bc8b580 229 data->MakeBranch("RL"); //trigger track
230 data->SetTreeAddress("RL");
cac2eb58 231 recoEvent->EventReconstructTrigger();
3bc8b580 232 data->Fill("RL");
cac2eb58 233
234 // tracking branch
3bc8b580 235 data->MakeBranch("RT"); //track
236 data->SetTreeAddress("RT");
cac2eb58 237 recoEvent->EventReconstruct();
3bc8b580 238 data->Fill("RT");
cac2eb58 239
cc87ebcd 240 loader->WriteTracks("OVERWRITE");
241
9c4b1ee7 242 if (recoEvent->GetTrackMethod() == kCombi) {
cc87ebcd 243 // Combined cluster / track
244 ((AliMUONClusterFinderAZ*) recModel)->SetReco(1);
3bc8b580 245 data->MakeBranch("RC");
246 data->SetTreeAddress("RC");
247 AliMUONEventRecoCombi::Instance()->FillRecP(data, recoEvent);
248 data->Fill("RC");
cc87ebcd 249 }
250 loader->WriteRecPoints("OVERWRITE");
251
cac2eb58 252 //--------------------------- Resetting branches -----------------------
3bc8b580 253 data->ResetDigits();
254 data->ResetRawClusters();
255 data->ResetTrigger();
3bc8b580 256 data->ResetRecTracks();
257 data->ResetRecTriggerTracks();
29f1b13a 258
cac2eb58 259 }
260 loader->UnloadDigits();
261 loader->UnloadRecPoints();
262 loader->UnloadTracks();
263
264 delete recoCluster;
265 delete recoEvent;
3bc8b580 266 delete data;
b2d7df0b 267 delete calibration;
cac2eb58 268}
cf464691 269
270//_____________________________________________________________________________
271void AliMUONReconstructor::Reconstruct(AliRunLoader* runLoader, AliRawReader* rawReader) const
272{
d19b6003 273/// Recontruct
274/// \todo add more
275
94f6fba9 276 // AliLoader
cf464691 277 AliLoader* loader = runLoader->GetLoader("MUONLoader");
94f6fba9 278 AliMUONData data(loader,"MUON","MUON");
cf464691 279
94f6fba9 280 // passing loader as argument.
9c4b1ee7 281 AliMUONTrackReconstructor recoEvent(&data);
96fdfe9a 282 recoEvent.SetTriggerCircuit(fTriggerCircuit);
cf464691 283
aa6ecf89 284 fDigitMaker->SetMUONData(&data);
1197ff51 285
94f6fba9 286 AliMUONClusterReconstructor recoCluster(&data);
30eabfb8 287
288 if (strstr(GetOption(),"Original"))
94f6fba9 289 {
9c4b1ee7 290 recoEvent.SetTrackMethod(kOriginal); // Original tracking
94f6fba9 291 }
30eabfb8 292 else
94f6fba9 293 {
9c4b1ee7 294 recoEvent.SetTrackMethod(kKalman);
94f6fba9 295 }
296
297 AliMUONClusterFinderVS *recModel = recoCluster.GetRecoModel();
298 if (!strstr(GetOption(),"VS"))
299 {
30eabfb8 300 recModel = (AliMUONClusterFinderVS*) new AliMUONClusterFinderAZ();
94f6fba9 301 recoCluster.SetRecoModel(recModel);
30eabfb8 302 }
7e4a628d 303 recModel->SetGhostChi2Cut(10);
cf464691 304
94f6fba9 305 TTask* calibration = GetCalibrationTask(&data);
b2d7df0b 306
cf464691 307 loader->LoadRecPoints("RECREATE");
308 loader->LoadTracks("RECREATE");
94f6fba9 309 loader->LoadDigits("READ");
310
cc87ebcd 311 // Loop over events
8210cb40 312 Int_t iEvent = 0;
94f6fba9 313
314 TStopwatch totalTimer;
315 TStopwatch rawTimer;
316 TStopwatch calibTimer;
317 TStopwatch clusterTimer;
318 TStopwatch trackingTimer;
319
320 rawTimer.Start(kTRUE); rawTimer.Stop();
321 calibTimer.Start(kTRUE); calibTimer.Stop();
322 clusterTimer.Start(kTRUE); clusterTimer.Stop();
323 trackingTimer.Start(kTRUE); trackingTimer.Stop();
324
325 totalTimer.Start(kTRUE);
326
327 while (rawReader->NextEvent())
328 {
329 AliDebug(1,Form("Event %d",iEvent));
330
cf464691 331 runLoader->GetEvent(iEvent++);
ce6a659c 332
333 //----------------------- raw2digits & raw2trigger-------------------
94f6fba9 334 if (!loader->TreeD())
335 {
336 AliDebug(1,Form("Making Digit Container for event %d",iEvent));
337 loader->MakeDigitsContainer();
338 }
339
340 data.SetTreeAddress("D,GLT");
341 rawTimer.Start(kFALSE);
aa6ecf89 342 fDigitMaker->Raw2Digits(rawReader);
94f6fba9 343 rawTimer.Stop();
344
b2d7df0b 345 if ( calibration )
346 {
94f6fba9 347 calibTimer.Start(kFALSE);
b2d7df0b 348 calibration->ExecuteTask();
94f6fba9 349 calibTimer.Stop();
b2d7df0b 350 }
94f6fba9 351
8210cb40 352 //----------------------- digit2cluster & Trigger2Trigger -------------------
94f6fba9 353 clusterTimer.Start(kFALSE);
354
cf464691 355 if (!loader->TreeR()) loader->MakeRecPointsContainer();
cc87ebcd 356
cf464691 357 // tracking branch
94f6fba9 358 data.MakeBranch("RC");
359 data.SetTreeAddress("RC");
360 recoCluster.Digits2Clusters();
361 data.Fill("RC");
cf464691 362
363 // trigger branch
94f6fba9 364 data.MakeBranch("TC");
365 data.SetTreeAddress("TC");
b59ce4c6 366 recoCluster.Trigger2Trigger();
94f6fba9 367 data.Fill("TC");
368
cf464691 369 loader->WriteRecPoints("OVERWRITE");
370
94f6fba9 371 clusterTimer.Stop();
372
cf464691 373 //---------------------------- Track & TriggerTrack ---------------------
94f6fba9 374 trackingTimer.Start(kFALSE);
cf464691 375 if (!loader->TreeT()) loader->MakeTracksContainer();
376
377 // trigger branch
94f6fba9 378 data.MakeBranch("RL"); //trigger track
379 data.SetTreeAddress("RL");
380 recoEvent.EventReconstructTrigger();
381 data.Fill("RL");
cf464691 382
383 // tracking branch
94f6fba9 384 data.MakeBranch("RT"); //track
385 data.SetTreeAddress("RT");
386 recoEvent.EventReconstruct();
387 data.Fill("RT");
cf464691 388
cc87ebcd 389 loader->WriteTracks("OVERWRITE");
94f6fba9 390 trackingTimer.Stop();
391
cf464691 392 //--------------------------- Resetting branches -----------------------
94f6fba9 393 data.ResetDigits();
394 data.ResetRawClusters();
395 data.ResetTrigger();
94f6fba9 396 data.ResetRecTracks();
397 data.ResetRecTriggerTracks();
cc87ebcd 398
cf464691 399 }
94f6fba9 400
401 totalTimer.Stop();
402
cf464691 403 loader->UnloadRecPoints();
404 loader->UnloadTracks();
ce6a659c 405 loader->UnloadDigits();
cf464691 406
94f6fba9 407 AliInfo(Form("Execution time for converting RAW data to digits in MUON : R:%.2fs C:%.2fs",
408 rawTimer.RealTime(),rawTimer.CpuTime()));
409 AliInfo(Form("Execution time for calibrating MUON : R:%.2fs C:%.2fs",
410 calibTimer.RealTime(),calibTimer.CpuTime()));
411 AliInfo(Form("Execution time for clusterizing MUON : R:%.2fs C:%.2fs",
412 clusterTimer.RealTime(),clusterTimer.CpuTime()));
413 AliInfo(Form("Execution time for tracking MUON : R:%.2fs C:%.2fs",
414 trackingTimer.RealTime(),trackingTimer.CpuTime()));
415 AliInfo(Form("Total Execution time for Reconstruct(from raw) MUON : R:%.2fs C:%.2fs",
416 totalTimer.RealTime(),totalTimer.CpuTime()));
cf464691 417}
418
cac2eb58 419//_____________________________________________________________________________
420void AliMUONReconstructor::FillESD(AliRunLoader* runLoader, AliESD* esd) const
421{
d19b6003 422/// Fill ESD
423/// \todo add more
424
cac2eb58 425 TClonesArray* recTracksArray = 0;
426 TClonesArray* recTrigTracksArray = 0;
cc87ebcd 427
cac2eb58 428 AliLoader* loader = runLoader->GetLoader("MUONLoader");
429 loader->LoadTracks("READ");
cac2eb58 430 AliMUONData* muonData = new AliMUONData(loader,"MUON","MUON");
431
cc87ebcd 432 // declaration
a2da7817 433 Int_t iEvent;// nPart;
434 Int_t nTrackHits;// nPrimary;
cac2eb58 435 Double_t fitFmin;
cac2eb58 436
437 Double_t bendingSlope, nonBendingSlope, inverseBendingMomentum;
438 Double_t xRec, yRec, zRec, chi2MatchTrigger;
439 Bool_t matchTrigger;
440
b002ac70 441 // setting pointer for tracks, triggertracks & trackparam at vertex
cac2eb58 442 AliMUONTrack* recTrack = 0;
443 AliMUONTrackParam* trackParam = 0;
444 AliMUONTriggerTrack* recTriggerTrack = 0;
2cf44ef3 445
cc87ebcd 446 iEvent = runLoader->GetEventNumber();
cac2eb58 447 runLoader->GetEvent(iEvent);
448
94f6fba9 449 // Get vertex
450 Double_t vertex[3] = {0};
451 const AliESDVertex *esdVert = esd->GetVertex();
43939bd8 452 if (esdVert->GetNContributors()) {
453 esdVert->GetXYZ(vertex);
454 printf("find vertex\n");
455 }
cac2eb58 456 // setting ESD MUON class
30178c30 457 AliESDMuonTrack* theESDTrack = new AliESDMuonTrack() ;
cac2eb58 458
459 //-------------------- trigger tracks-------------
460 Long_t trigPat = 0;
461 muonData->SetTreeAddress("RL");
462 muonData->GetRecTriggerTracks();
463 recTrigTracksArray = muonData->RecTriggerTracks();
464
465 // ready global trigger pattern from first track
cc87ebcd 466 if (recTrigTracksArray)
cac2eb58 467 recTriggerTrack = (AliMUONTriggerTrack*) recTrigTracksArray->First();
468 if (recTriggerTrack) trigPat = recTriggerTrack->GetGTPattern();
469
470 //printf(">>> Event %d Number of Recconstructed tracks %d \n",iEvent, nrectracks);
cc87ebcd 471
cac2eb58 472 // -------------------- tracks-------------
473 muonData->SetTreeAddress("RT");
474 muonData->GetRecTracks();
475 recTracksArray = muonData->RecTracks();
cc87ebcd 476
cac2eb58 477 Int_t nRecTracks = 0;
478 if (recTracksArray)
479 nRecTracks = (Int_t) recTracksArray->GetEntriesFast(); //
cc87ebcd 480
cac2eb58 481 // loop over tracks
482 for (Int_t iRecTracks = 0; iRecTracks < nRecTracks; iRecTracks++) {
483
484 // reading info from tracks
485 recTrack = (AliMUONTrack*) recTracksArray->At(iRecTracks);
486
889a0215 487 trackParam = (AliMUONTrackParam*) (recTrack->GetTrackParamAtHit())->First();
211c52eb 488
43939bd8 489 if (esdVert->GetNContributors())
211c52eb 490 trackParam->ExtrapToVertex(vertex[0],vertex[1],vertex[2]);
cac2eb58 491
492 bendingSlope = trackParam->GetBendingSlope();
493 nonBendingSlope = trackParam->GetNonBendingSlope();
211c52eb 494 inverseBendingMomentum = trackParam->GetInverseBendingMomentum();
cac2eb58 495 xRec = trackParam->GetNonBendingCoor();
496 yRec = trackParam->GetBendingCoor();
497 zRec = trackParam->GetZ();
498
499 nTrackHits = recTrack->GetNTrackHits();
500 fitFmin = recTrack->GetFitFMin();
501 matchTrigger = recTrack->GetMatchTrigger();
502 chi2MatchTrigger = recTrack->GetChi2MatchTrigger();
503
504 // setting data member of ESD MUON
30178c30 505 theESDTrack->SetInverseBendingMomentum(inverseBendingMomentum);
506 theESDTrack->SetThetaX(TMath::ATan(nonBendingSlope));
507 theESDTrack->SetThetaY(TMath::ATan(bendingSlope));
889a0215 508 theESDTrack->SetZ(zRec);
509 theESDTrack->SetBendingCoor(yRec); // calculate vertex at ESD or Tracking level ?
510 theESDTrack->SetNonBendingCoor(xRec);
30178c30 511 theESDTrack->SetChi2(fitFmin);
512 theESDTrack->SetNHit(nTrackHits);
513 theESDTrack->SetMatchTrigger(matchTrigger);
514 theESDTrack->SetChi2MatchTrigger(chi2MatchTrigger);
cac2eb58 515
cc87ebcd 516 // storing ESD MUON Track into ESD Event
517 if (nRecTracks != 0)
30178c30 518 esd->AddMuonTrack(theESDTrack);
cac2eb58 519 } // end loop tracks
520
cac2eb58 521 // reset muondata
522 muonData->ResetRecTracks();
523 muonData->ResetRecTriggerTracks();
524
cc87ebcd 525 //} // end loop on event
526 loader->UnloadTracks();
2cf44ef3 527
30178c30 528 delete theESDTrack;
cac2eb58 529 delete muonData;
a2da7817 530}//_____________________________________________________________________________
531void AliMUONReconstructor::FillESD(AliRunLoader* runLoader, AliRawReader* /*rawReader*/, AliESD* esd) const
532{
d19b6003 533/// Fill ESD
534/// \todo add more
535
a2da7817 536 // don't need rawReader ???
537 FillESD(runLoader, esd);
cac2eb58 538}