1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
17 ///////////////////////////////////////////////////////////////////////////////
19 // class for MUON reconstruction //
21 ///////////////////////////////////////////////////////////////////////////////
23 #include "AliMUONReconstructor.h"
26 #include "AliESDMuonTrack.h"
27 #include "AliGenEventHeader.h"
28 #include "AliHeader.h"
31 #include "AliMUONCalibrationData.h"
32 #include "AliMUONClusterFinderAZ.h"
33 #include "AliMUONClusterFinderVS.h"
34 #include "AliMUONClusterReconstructor.h"
35 #include "AliMUONData.h"
36 #include "AliMUONDigitCalibrator.h"
37 #include "AliMUONEventRecoCombi.h"
38 #include "AliMUONRawReader.h"
39 #include "AliMUONTrack.h"
40 #include "AliMUONTrackParam.h"
41 #include "AliMUONTrackReconstructor.h"
42 #include "AliMUONTriggerTrack.h"
43 #include "AliRawReader.h"
45 #include "AliRunLoader.h"
48 #include <TParticle.h>
49 #include "TStopwatch.h"
51 ClassImp(AliMUONReconstructor)
53 //_____________________________________________________________________________
54 AliMUONReconstructor::AliMUONReconstructor()
55 : AliReconstructor(), fCalibrationData(0x0)
60 //_____________________________________________________________________________
61 AliMUONReconstructor::~AliMUONReconstructor()
64 delete fCalibrationData;
67 //_____________________________________________________________________________
69 AliMUONReconstructor::GetCalibrationTask(AliMUONData* data) const
72 // Create the calibration task(s).
75 const AliRun* run = fRunLoader->GetAliRun();
77 // Not really clean, but for the moment we must check whether the
78 // simulation has decalibrated the data or not...
79 const AliMUON* muon = static_cast<AliMUON*>(run->GetModule("MUON"));
80 if ( muon->DigitizerType().Contains("NewDigitizer") )
82 AliInfo("Calibration will occur.");
83 Int_t runNumber = run->GetRunNumber();
84 fCalibrationData = new AliMUONCalibrationData(runNumber);
85 if ( !fCalibrationData->IsValid() )
87 AliError("Could not retrieve calibrations !");
88 delete fCalibrationData;
89 fCalibrationData = 0x0;
92 TTask* calibration = new TTask("MUONCalibrator","MUON Digit calibrator");
93 calibration->Add(new AliMUONDigitCalibrator(data,fCalibrationData));
94 //FIXME: calibration->Add(something about dead channels should go here).
99 AliInfo("Detected the usage of old digitizer (w/o decalibration). "
100 "Will not calibrate then.");
105 //_____________________________________________________________________________
107 AliMUONReconstructor::Init(AliRunLoader* runLoader)
109 fRunLoader = runLoader;
112 //_____________________________________________________________________________
113 void AliMUONReconstructor::Reconstruct(AliRunLoader* runLoader) const
116 AliLoader* loader = runLoader->GetLoader("MUONLoader");
117 Int_t nEvents = runLoader->GetNumberOfEvents();
119 AliMUONData* data = new AliMUONData(loader,"MUON","MUON");
121 // passing loader as argument.
122 AliMUONTrackReconstructor* recoEvent = new AliMUONTrackReconstructor(loader, data);
124 if (strstr(GetOption(),"Original"))
125 recoEvent->SetTrackMethod(1); // Original tracking
126 else if (strstr(GetOption(),"Combi"))
127 recoEvent->SetTrackMethod(3); // Combined cluster / track
129 recoEvent->SetTrackMethod(2); // Kalman
131 AliMUONClusterReconstructor* recoCluster = new AliMUONClusterReconstructor(data);
133 AliMUONClusterFinderVS *recModel = recoCluster->GetRecoModel();
135 if (!strstr(GetOption(),"VS")) {
136 recModel = (AliMUONClusterFinderVS*) new AliMUONClusterFinderAZ();
137 recoCluster->SetRecoModel(recModel);
139 recModel->SetGhostChi2Cut(10);
141 loader->LoadDigits("READ");
142 loader->LoadRecPoints("RECREATE");
143 loader->LoadTracks("RECREATE");
145 TTask* calibration = GetCalibrationTask(data);
147 Int_t chBeg = recoEvent->GetTrackMethod() == 3 ? 6 : 0;
149 for(Int_t ievent = 0; ievent < nEvents; ievent++) {
151 AliDebug(1,Form("Event %d",ievent));
153 runLoader->GetEvent(ievent);
155 //----------------------- digit2cluster & Trigger2Trigger -------------------
156 if (!loader->TreeR()) loader->MakeRecPointsContainer();
159 if (recoEvent->GetTrackMethod() != 3) {
160 data->MakeBranch("RC");
161 data->SetTreeAddress("D,RC");
163 data->SetTreeAddress("D");
164 data->SetTreeAddress("RCC");
166 // Important for avoiding a memory leak when reading digits ( to be investigated more in detail)
167 // In any case the reading of GLT is needed for the Trigger2Tigger method below
168 data->SetTreeAddress("GLT");
174 calibration->ExecuteTask();
177 recoCluster->Digits2Clusters(chBeg);
179 if (recoEvent->GetTrackMethod() == 3) {
180 // Combined cluster / track finder
181 AliMUONEventRecoCombi::Instance()->FillEvent(data, (AliMUONClusterFinderAZ*)recModel);
182 ((AliMUONClusterFinderAZ*) recModel)->SetReco(2);
184 else data->Fill("RC");
187 data->MakeBranch("TC");
188 data->SetTreeAddress("TC");
189 recoCluster->Trigger2Trigger();
192 //AZ loader->WriteRecPoints("OVERWRITE");
194 //---------------------------- Track & TriggerTrack ---------------------
195 if (!loader->TreeT()) loader->MakeTracksContainer();
198 data->MakeBranch("RL"); //trigger track
199 data->SetTreeAddress("RL");
200 recoEvent->EventReconstructTrigger();
204 data->MakeBranch("RT"); //track
205 data->SetTreeAddress("RT");
206 recoEvent->EventReconstruct();
209 loader->WriteTracks("OVERWRITE");
211 if (recoEvent->GetTrackMethod() == 3) {
212 // Combined cluster / track
213 ((AliMUONClusterFinderAZ*) recModel)->SetReco(1);
214 data->MakeBranch("RC");
215 data->SetTreeAddress("RC");
216 AliMUONEventRecoCombi::Instance()->FillRecP(data, recoEvent);
219 loader->WriteRecPoints("OVERWRITE");
221 //--------------------------- Resetting branches -----------------------
223 data->ResetRawClusters();
224 data->ResetTrigger();
226 data->ResetRawClusters();
227 data->ResetTrigger();
228 data->ResetRecTracks();
229 data->ResetRecTriggerTracks();
232 loader->UnloadDigits();
233 loader->UnloadRecPoints();
234 loader->UnloadTracks();
242 //_____________________________________________________________________________
243 void AliMUONReconstructor::Reconstruct(AliRunLoader* runLoader, AliRawReader* rawReader) const
246 AliLoader* loader = runLoader->GetLoader("MUONLoader");
247 AliMUONData data(loader,"MUON","MUON");
249 // passing loader as argument.
250 AliMUONTrackReconstructor recoEvent(loader, &data);
252 AliMUONRawReader rawData(&data);
254 AliMUONClusterReconstructor recoCluster(&data);
256 if (strstr(GetOption(),"Original"))
258 recoEvent.SetTrackMethod(1); // Original tracking
262 recoEvent.SetTrackMethod(2); // Kalman
265 AliMUONClusterFinderVS *recModel = recoCluster.GetRecoModel();
266 if (!strstr(GetOption(),"VS"))
268 recModel = (AliMUONClusterFinderVS*) new AliMUONClusterFinderAZ();
269 recoCluster.SetRecoModel(recModel);
271 recModel->SetGhostChi2Cut(10);
273 TTask* calibration = GetCalibrationTask(&data);
275 loader->LoadRecPoints("RECREATE");
276 loader->LoadTracks("RECREATE");
277 loader->LoadDigits("READ");
282 TStopwatch totalTimer;
284 TStopwatch calibTimer;
285 TStopwatch clusterTimer;
286 TStopwatch trackingTimer;
288 rawTimer.Start(kTRUE); rawTimer.Stop();
289 calibTimer.Start(kTRUE); calibTimer.Stop();
290 clusterTimer.Start(kTRUE); clusterTimer.Stop();
291 trackingTimer.Start(kTRUE); trackingTimer.Stop();
293 totalTimer.Start(kTRUE);
295 while (rawReader->NextEvent())
297 AliDebug(1,Form("Event %d",iEvent));
299 runLoader->GetEvent(iEvent++);
301 //----------------------- raw2digits & raw2trigger-------------------
302 if (!loader->TreeD())
304 AliDebug(1,Form("Making Digit Container for event %d",iEvent));
305 loader->MakeDigitsContainer();
308 data.SetTreeAddress("D,GLT");
309 rawTimer.Start(kFALSE);
310 rawData.Raw2Digits(rawReader);
315 calibTimer.Start(kFALSE);
316 calibration->ExecuteTask();
320 //----------------------- digit2cluster & Trigger2Trigger -------------------
321 clusterTimer.Start(kFALSE);
323 if (!loader->TreeR()) loader->MakeRecPointsContainer();
326 data.MakeBranch("RC");
327 data.SetTreeAddress("RC");
328 recoCluster.Digits2Clusters();
332 data.MakeBranch("TC");
333 data.SetTreeAddress("TC");
334 // recoCluster.Trigger2Trigger();
337 loader->WriteRecPoints("OVERWRITE");
341 //---------------------------- Track & TriggerTrack ---------------------
342 trackingTimer.Start(kFALSE);
343 if (!loader->TreeT()) loader->MakeTracksContainer();
346 data.MakeBranch("RL"); //trigger track
347 data.SetTreeAddress("RL");
348 recoEvent.EventReconstructTrigger();
352 data.MakeBranch("RT"); //track
353 data.SetTreeAddress("RT");
354 recoEvent.EventReconstruct();
357 loader->WriteTracks("OVERWRITE");
358 trackingTimer.Stop();
360 //--------------------------- Resetting branches -----------------------
362 data.ResetRawClusters();
365 data.ResetRawClusters();
367 data.ResetRecTracks();
368 data.ResetRecTriggerTracks();
374 loader->UnloadRecPoints();
375 loader->UnloadTracks();
376 loader->UnloadDigits();
378 AliInfo(Form("Execution time for converting RAW data to digits in MUON : R:%.2fs C:%.2fs",
379 rawTimer.RealTime(),rawTimer.CpuTime()));
380 AliInfo(Form("Execution time for calibrating MUON : R:%.2fs C:%.2fs",
381 calibTimer.RealTime(),calibTimer.CpuTime()));
382 AliInfo(Form("Execution time for clusterizing MUON : R:%.2fs C:%.2fs",
383 clusterTimer.RealTime(),clusterTimer.CpuTime()));
384 AliInfo(Form("Execution time for tracking MUON : R:%.2fs C:%.2fs",
385 trackingTimer.RealTime(),trackingTimer.CpuTime()));
386 AliInfo(Form("Total Execution time for Reconstruct(from raw) MUON : R:%.2fs C:%.2fs",
387 totalTimer.RealTime(),totalTimer.CpuTime()));
390 //_____________________________________________________________________________
391 void AliMUONReconstructor::FillESD(AliRunLoader* runLoader, AliESD* esd) const
393 TClonesArray* recTracksArray = 0;
394 TClonesArray* recTrigTracksArray = 0;
396 AliLoader* loader = runLoader->GetLoader("MUONLoader");
397 loader->LoadTracks("READ");
398 AliMUONData* muonData = new AliMUONData(loader,"MUON","MUON");
401 Int_t iEvent;// nPart;
402 Int_t nTrackHits;// nPrimary;
405 Double_t bendingSlope, nonBendingSlope, inverseBendingMomentum;
406 Double_t xRec, yRec, zRec, chi2MatchTrigger;
409 // setting pointer for tracks, triggertracks & trackparam at vertex
410 AliMUONTrack* recTrack = 0;
411 AliMUONTrackParam* trackParam = 0;
412 AliMUONTriggerTrack* recTriggerTrack = 0;
413 // TParticle* particle = new TParticle();
414 // AliGenEventHeader* header = 0;
415 iEvent = runLoader->GetEventNumber();
416 runLoader->GetEvent(iEvent);
419 Double_t vertex[3] = {0};
420 const AliESDVertex *esdVert = esd->GetVertex();
421 if (esdVert) esdVert->GetXYZ(vertex);
424 // if ( (header = runLoader->GetHeader()->GenEventHeader()) ) {
425 // header->PrimaryVertex(vertex);
427 // runLoader->LoadKinematics("READ");
428 // runLoader->TreeK()->GetBranch("Particles")->SetAddress(&particle);
429 // nPart = (Int_t)runLoader->TreeK()->GetEntries();
430 // for(Int_t iPart = 0; iPart < nPart; iPart++) {
431 // runLoader->TreeK()->GetEvent(iPart);
432 // if (particle->GetFirstMother() == -1) {
433 // vertex[0] += particle->Vx();
434 // vertex[1] += particle->Vy();
435 // vertex[2] += particle->Vz();
439 // vertex[0] /= (double)nPrimary;
440 // vertex[1] /= (double)nPrimary;
441 // vertex[2] /= (double)nPrimary;
445 // setting ESD MUON class
446 AliESDMuonTrack* theESDTrack = new AliESDMuonTrack() ;
448 //-------------------- trigger tracks-------------
450 muonData->SetTreeAddress("RL");
451 muonData->GetRecTriggerTracks();
452 recTrigTracksArray = muonData->RecTriggerTracks();
454 // ready global trigger pattern from first track
455 if (recTrigTracksArray)
456 recTriggerTrack = (AliMUONTriggerTrack*) recTrigTracksArray->First();
457 if (recTriggerTrack) trigPat = recTriggerTrack->GetGTPattern();
459 //printf(">>> Event %d Number of Recconstructed tracks %d \n",iEvent, nrectracks);
461 // -------------------- tracks-------------
462 muonData->SetTreeAddress("RT");
463 muonData->GetRecTracks();
464 recTracksArray = muonData->RecTracks();
466 Int_t nRecTracks = 0;
468 nRecTracks = (Int_t) recTracksArray->GetEntriesFast(); //
471 for (Int_t iRecTracks = 0; iRecTracks < nRecTracks; iRecTracks++) {
473 // reading info from tracks
474 recTrack = (AliMUONTrack*) recTracksArray->At(iRecTracks);
476 trackParam = (AliMUONTrackParam*) (recTrack->GetTrackParamAtHit())->First();
477 trackParam->ExtrapToVertex(vertex[0],vertex[1],vertex[2]);
479 bendingSlope = trackParam->GetBendingSlope();
480 nonBendingSlope = trackParam->GetNonBendingSlope();
481 inverseBendingMomentum = trackParam->GetInverseBendingMomentum();
482 xRec = trackParam->GetNonBendingCoor();
483 yRec = trackParam->GetBendingCoor();
484 zRec = trackParam->GetZ();
486 nTrackHits = recTrack->GetNTrackHits();
487 fitFmin = recTrack->GetFitFMin();
488 matchTrigger = recTrack->GetMatchTrigger();
489 chi2MatchTrigger = recTrack->GetChi2MatchTrigger();
491 // setting data member of ESD MUON
492 theESDTrack->SetInverseBendingMomentum(inverseBendingMomentum);
493 theESDTrack->SetThetaX(TMath::ATan(nonBendingSlope));
494 theESDTrack->SetThetaY(TMath::ATan(bendingSlope));
495 theESDTrack->SetZ(zRec);
496 theESDTrack->SetBendingCoor(yRec); // calculate vertex at ESD or Tracking level ?
497 theESDTrack->SetNonBendingCoor(xRec);
498 theESDTrack->SetChi2(fitFmin);
499 theESDTrack->SetNHit(nTrackHits);
500 theESDTrack->SetMatchTrigger(matchTrigger);
501 theESDTrack->SetChi2MatchTrigger(chi2MatchTrigger);
503 // storing ESD MUON Track into ESD Event
505 esd->AddMuonTrack(theESDTrack);
508 // add global trigger pattern
510 esd->SetTrigger(trigPat);
513 muonData->ResetRecTracks();
514 muonData->ResetRecTriggerTracks();
516 //} // end loop on event
517 loader->UnloadTracks();
519 // runLoader->UnloadKinematics();
523 }//_____________________________________________________________________________
524 void AliMUONReconstructor::FillESD(AliRunLoader* runLoader, AliRawReader* /*rawReader*/, AliESD* esd) const
526 // don't need rawReader ???
527 FillESD(runLoader, esd);