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