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 **************************************************************************/
18 ///////////////////////////////////////////////////////////////////////////////
20 // class for running the reconstruction //
22 // Clusters and tracks are created for all detectors and all events by //
25 // AliReconstruction rec; //
28 // The Run method returns kTRUE in case of successful execution. //
30 // If the input to the reconstruction are not simulated digits but raw data, //
31 // this can be specified by an argument of the Run method or by the method //
33 // rec.SetInput("..."); //
35 // The input formats and the corresponding argument are: //
36 // - DDL raw data files: directory name, ends with "/" //
37 // - raw data root file: root file name, extension ".root" //
38 // - raw data DATE file: DATE file name, any other non-empty string //
39 // - MC root files : empty string, default //
41 // By default all events are reconstructed. The reconstruction can be //
42 // limited to a range of events by giving the index of the first and the //
43 // last event as an argument to the Run method or by calling //
45 // rec.SetEventRange(..., ...); //
47 // The index -1 (default) can be used for the last event to indicate no //
48 // upper limit of the event range. //
50 // In case of raw-data reconstruction the user can modify the default //
51 // number of events per digits/clusters/tracks file. In case the option //
52 // is not used the number is set 1. In case the user provides 0, than //
53 // the number of events is equal to the number of events inside the //
54 // raw-data file (i.e. one digits/clusters/tracks file): //
56 // rec.SetNumberOfEventsPerFile(...); //
59 // The name of the galice file can be changed from the default //
60 // "galice.root" by passing it as argument to the AliReconstruction //
61 // constructor or by //
63 // rec.SetGAliceFile("..."); //
65 // The local reconstruction can be switched on or off for individual //
68 // rec.SetRunLocalReconstruction("..."); //
70 // The argument is a (case sensitive) string with the names of the //
71 // detectors separated by a space. The special string "ALL" selects all //
72 // available detectors. This is the default. //
74 // The reconstruction of the primary vertex position can be switched off by //
76 // rec.SetRunVertexFinder(kFALSE); //
78 // The tracking and the creation of ESD tracks can be switched on for //
79 // selected detectors by //
81 // rec.SetRunTracking("..."); //
83 // Uniform/nonuniform field tracking switches (default: uniform field) //
85 // rec.SetUniformFieldTracking(); ( rec.SetUniformFieldTracking(kFALSE); ) //
87 // The filling of additional ESD information can be steered by //
89 // rec.SetFillESD("..."); //
91 // Again, for both methods the string specifies the list of detectors. //
92 // The default is "ALL". //
94 // The call of the shortcut method //
96 // rec.SetRunReconstruction("..."); //
98 // is equivalent to calling SetRunLocalReconstruction, SetRunTracking and //
99 // SetFillESD with the same detector selecting string as argument. //
101 // The reconstruction requires digits or raw data as input. For the creation //
102 // of digits and raw data have a look at the class AliSimulation. //
104 // The input data of a detector can be replaced by the corresponding HLT //
105 // data by calling (usual detector string) //
106 // SetUseHLTData("..."); //
109 ///////////////////////////////////////////////////////////////////////////////
116 #include <TPluginManager.h>
117 #include <TGeoManager.h>
118 #include <TLorentzVector.h>
121 #include <TObjArray.h>
125 #include <TProofOutputFile.h>
126 #include <TParameter.h>
128 #include "AliReconstruction.h"
129 #include "AliCodeTimer.h"
130 #include "AliReconstructor.h"
132 #include "AliRunLoader.h"
134 #include "AliRawReaderFile.h"
135 #include "AliRawReaderDate.h"
136 #include "AliRawReaderRoot.h"
137 #include "AliRawEventHeaderBase.h"
138 #include "AliRawEvent.h"
139 #include "AliESDEvent.h"
140 #include "AliESDMuonTrack.h"
141 #include "AliESDfriend.h"
142 #include "AliESDVertex.h"
143 #include "AliESDcascade.h"
144 #include "AliESDkink.h"
145 #include "AliESDtrack.h"
146 #include "AliESDCaloCluster.h"
147 #include "AliESDCaloCells.h"
148 #include "AliMultiplicity.h"
149 #include "AliTracker.h"
150 #include "AliVertexer.h"
151 #include "AliVertexerTracks.h"
152 #include "AliV0vertexer.h"
153 #include "AliCascadeVertexer.h"
154 #include "AliHeader.h"
155 #include "AliGenEventHeader.h"
157 #include "AliESDpid.h"
158 #include "AliESDtrack.h"
159 #include "AliESDPmdTrack.h"
161 #include "AliESDTagCreator.h"
163 #include "AliGeomManager.h"
164 #include "AliTrackPointArray.h"
165 #include "AliCDBManager.h"
166 #include "AliCDBStorage.h"
167 #include "AliCDBEntry.h"
168 #include "AliAlignObj.h"
170 #include "AliCentralTrigger.h"
171 #include "AliTriggerConfiguration.h"
172 #include "AliTriggerClass.h"
173 #include "AliTriggerCluster.h"
174 #include "AliCTPRawStream.h"
176 #include "AliQADataMakerRec.h"
177 #include "AliGlobalQADataMaker.h"
179 #include "AliQADataMakerSteer.h"
181 #include "AliPlaneEff.h"
183 #include "AliSysInfo.h" // memory snapshots
184 #include "AliRawHLTManager.h"
186 #include "AliMagWrapCheb.h"
188 #include "AliDetectorRecoParam.h"
189 #include "AliGRPRecoParam.h"
190 #include "AliRunInfo.h"
191 #include "AliEventInfo.h"
195 ClassImp(AliReconstruction)
197 //_____________________________________________________________________________
198 const char* AliReconstruction::fgkDetectorName[AliReconstruction::fgkNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "HMPID", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "T0", "VZERO", "ACORDE", "HLT"};
200 //_____________________________________________________________________________
201 AliReconstruction::AliReconstruction(const char* gAliceFilename) :
203 fUniformField(kFALSE),
204 fForcedFieldMap(NULL),
205 fRunVertexFinder(kTRUE),
206 fRunVertexFinderTracks(kTRUE),
207 fRunHLTTracking(kFALSE),
208 fRunMuonTracking(kFALSE),
210 fRunCascadeFinder(kTRUE),
211 fStopOnError(kFALSE),
212 fWriteAlignmentData(kFALSE),
213 fWriteESDfriend(kFALSE),
214 fFillTriggerESD(kTRUE),
222 fRunLocalReconstruction("ALL"),
226 fUseTrackingErrorsForAlignment(""),
227 fGAliceFileName(gAliceFilename),
232 fNumberOfEventsPerFile(1),
234 fLoadAlignFromCDB(kTRUE),
235 fLoadAlignData("ALL"),
243 fParentRawReader(NULL),
248 fDiamondProfile(NULL),
249 fDiamondProfileTPC(NULL),
250 fMeanVertexConstraint(kTRUE),
254 fAlignObjArray(NULL),
257 fInitCDBCalled(kFALSE),
258 fSetRunNumberFromDataCalled(kFALSE),
264 fSameQACycle(kFALSE),
266 fRunPlaneEff(kFALSE),
275 fIsNewRunLoader(kFALSE),
279 // create reconstruction object with default parameters
282 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
283 fReconstructor[iDet] = NULL;
284 fLoader[iDet] = NULL;
285 fTracker[iDet] = NULL;
287 for (Int_t iDet = 0; iDet < AliQA::kNDET; iDet++) {
288 fQACycles[iDet] = 999999 ;
289 fQAWriteExpert[iDet] = kFALSE ;
295 //_____________________________________________________________________________
296 AliReconstruction::AliReconstruction(const AliReconstruction& rec) :
298 fUniformField(rec.fUniformField),
299 fForcedFieldMap(NULL),
300 fRunVertexFinder(rec.fRunVertexFinder),
301 fRunVertexFinderTracks(rec.fRunVertexFinderTracks),
302 fRunHLTTracking(rec.fRunHLTTracking),
303 fRunMuonTracking(rec.fRunMuonTracking),
304 fRunV0Finder(rec.fRunV0Finder),
305 fRunCascadeFinder(rec.fRunCascadeFinder),
306 fStopOnError(rec.fStopOnError),
307 fWriteAlignmentData(rec.fWriteAlignmentData),
308 fWriteESDfriend(rec.fWriteESDfriend),
309 fFillTriggerESD(rec.fFillTriggerESD),
311 fCleanESD(rec.fCleanESD),
312 fV0DCAmax(rec.fV0DCAmax),
313 fV0CsPmin(rec.fV0CsPmin),
317 fRunLocalReconstruction(rec.fRunLocalReconstruction),
318 fRunTracking(rec.fRunTracking),
319 fFillESD(rec.fFillESD),
320 fLoadCDB(rec.fLoadCDB),
321 fUseTrackingErrorsForAlignment(rec.fUseTrackingErrorsForAlignment),
322 fGAliceFileName(rec.fGAliceFileName),
323 fRawInput(rec.fRawInput),
324 fEquipIdMap(rec.fEquipIdMap),
325 fFirstEvent(rec.fFirstEvent),
326 fLastEvent(rec.fLastEvent),
327 fNumberOfEventsPerFile(rec.fNumberOfEventsPerFile),
329 fLoadAlignFromCDB(rec.fLoadAlignFromCDB),
330 fLoadAlignData(rec.fLoadAlignData),
331 fESDPar(rec.fESDPar),
332 fUseHLTData(rec.fUseHLTData),
338 fParentRawReader(NULL),
340 fRecoParam(rec.fRecoParam),
343 fDiamondProfile(rec.fDiamondProfile),
344 fDiamondProfileTPC(rec.fDiamondProfileTPC),
345 fMeanVertexConstraint(rec.fMeanVertexConstraint),
349 fAlignObjArray(rec.fAlignObjArray),
350 fCDBUri(rec.fCDBUri),
352 fInitCDBCalled(rec.fInitCDBCalled),
353 fSetRunNumberFromDataCalled(rec.fSetRunNumberFromDataCalled),
354 fQADetectors(rec.fQADetectors),
356 fQATasks(rec.fQATasks),
358 fRunGlobalQA(rec.fRunGlobalQA),
359 fSameQACycle(rec.fSameQACycle),
360 fRunPlaneEff(rec.fRunPlaneEff),
369 fIsNewRunLoader(rec.fIsNewRunLoader),
375 for (Int_t i = 0; i < rec.fOptions.GetEntriesFast(); i++) {
376 if (rec.fOptions[i]) fOptions.Add(rec.fOptions[i]->Clone());
378 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
379 fReconstructor[iDet] = NULL;
380 fLoader[iDet] = NULL;
381 fTracker[iDet] = NULL;
384 for (Int_t iDet = 0; iDet < AliQA::kNDET; iDet++) {
385 fQACycles[iDet] = rec.fQACycles[iDet];
386 fQAWriteExpert[iDet] = rec.fQAWriteExpert[iDet] ;
389 for (Int_t i = 0; i < rec.fSpecCDBUri.GetEntriesFast(); i++) {
390 if (rec.fSpecCDBUri[i]) fSpecCDBUri.Add(rec.fSpecCDBUri[i]->Clone());
394 //_____________________________________________________________________________
395 AliReconstruction& AliReconstruction::operator = (const AliReconstruction& rec)
397 // assignment operator
398 // Used in PROOF mode
399 // Be very careful while modifing it!
400 // Simple rules to follow:
401 // for persistent data members - use their assignment operators
402 // for non-persistent ones - do nothing or take the default values from constructor
403 // TSelector members should not be touched
404 if(&rec == this) return *this;
406 fUniformField = rec.fUniformField;
407 fForcedFieldMap = NULL;
408 fRunVertexFinder = rec.fRunVertexFinder;
409 fRunVertexFinderTracks = rec.fRunVertexFinderTracks;
410 fRunHLTTracking = rec.fRunHLTTracking;
411 fRunMuonTracking = rec.fRunMuonTracking;
412 fRunV0Finder = rec.fRunV0Finder;
413 fRunCascadeFinder = rec.fRunCascadeFinder;
414 fStopOnError = rec.fStopOnError;
415 fWriteAlignmentData = rec.fWriteAlignmentData;
416 fWriteESDfriend = rec.fWriteESDfriend;
417 fFillTriggerESD = rec.fFillTriggerESD;
419 fCleanESD = rec.fCleanESD;
420 fV0DCAmax = rec.fV0DCAmax;
421 fV0CsPmin = rec.fV0CsPmin;
425 fRunLocalReconstruction = rec.fRunLocalReconstruction;
426 fRunTracking = rec.fRunTracking;
427 fFillESD = rec.fFillESD;
428 fLoadCDB = rec.fLoadCDB;
429 fUseTrackingErrorsForAlignment = rec.fUseTrackingErrorsForAlignment;
430 fGAliceFileName = rec.fGAliceFileName;
431 fRawInput = rec.fRawInput;
432 fEquipIdMap = rec.fEquipIdMap;
433 fFirstEvent = rec.fFirstEvent;
434 fLastEvent = rec.fLastEvent;
435 fNumberOfEventsPerFile = rec.fNumberOfEventsPerFile;
437 for (Int_t i = 0; i < rec.fOptions.GetEntriesFast(); i++) {
438 if (rec.fOptions[i]) fOptions.Add(rec.fOptions[i]->Clone());
441 fLoadAlignFromCDB = rec.fLoadAlignFromCDB;
442 fLoadAlignData = rec.fLoadAlignData;
443 fESDPar = rec.fESDPar;
444 fUseHLTData = rec.fUseHLTData;
446 delete fRunInfo; fRunInfo = NULL;
447 if (rec.fRunInfo) fRunInfo = new AliRunInfo(*rec.fRunInfo);
449 fEventInfo = rec.fEventInfo;
453 fParentRawReader = NULL;
455 fRecoParam = rec.fRecoParam;
457 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
458 delete fReconstructor[iDet]; fReconstructor[iDet] = NULL;
459 delete fLoader[iDet]; fLoader[iDet] = NULL;
460 delete fTracker[iDet]; fTracker[iDet] = NULL;
463 for (Int_t iDet = 0; iDet < AliQA::kNDET; iDet++) {
464 fQACycles[iDet] = rec.fQACycles[iDet];
465 fQAWriteExpert[iDet] = rec.fQAWriteExpert[iDet] ;
469 delete fDiamondProfile; fDiamondProfile = NULL;
470 if (rec.fDiamondProfile) fDiamondProfile = new AliESDVertex(*rec.fDiamondProfile);
471 delete fDiamondProfileTPC; fDiamondProfileTPC = NULL;
472 if (rec.fDiamondProfileTPC) fDiamondProfileTPC = new AliESDVertex(*rec.fDiamondProfileTPC);
473 fMeanVertexConstraint = rec.fMeanVertexConstraint;
475 delete fGRPData; fGRPData = NULL;
476 if (rec.fGRPData) fGRPData = (TMap*)((rec.fGRPData)->Clone());
478 delete fAlignObjArray; fAlignObjArray = NULL;
481 fSpecCDBUri.Delete();
482 fInitCDBCalled = rec.fInitCDBCalled;
483 fSetRunNumberFromDataCalled = rec.fSetRunNumberFromDataCalled;
484 fQADetectors = rec.fQADetectors;
486 fQATasks = rec.fQATasks;
488 fRunGlobalQA = rec.fRunGlobalQA;
489 fSameQACycle = rec.fSameQACycle;
490 fRunPlaneEff = rec.fRunPlaneEff;
499 fIsNewRunLoader = rec.fIsNewRunLoader;
506 //_____________________________________________________________________________
507 AliReconstruction::~AliReconstruction()
512 delete fForcedFieldMap;
514 if (fAlignObjArray) {
515 fAlignObjArray->Delete();
516 delete fAlignObjArray;
518 fSpecCDBUri.Delete();
520 AliCodeTimer::Instance()->Print();
523 //_____________________________________________________________________________
524 void AliReconstruction::InitCDB()
526 // activate a default CDB storage
527 // First check if we have any CDB storage set, because it is used
528 // to retrieve the calibration and alignment constants
529 AliCodeTimerAuto("");
531 if (fInitCDBCalled) return;
532 fInitCDBCalled = kTRUE;
534 AliCDBManager* man = AliCDBManager::Instance();
535 if (man->IsDefaultStorageSet())
537 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
538 AliWarning("Default CDB storage has been already set !");
539 AliWarning(Form("Ignoring the default storage declared in AliReconstruction: %s",fCDBUri.Data()));
540 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
541 fCDBUri = man->GetDefaultStorage()->GetURI();
544 if (fCDBUri.Length() > 0)
546 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
547 AliDebug(2, Form("Default CDB storage is set to: %s", fCDBUri.Data()));
548 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
550 fCDBUri="local://$ALICE_ROOT";
551 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
552 AliWarning("Default CDB storage not yet set !!!!");
553 AliWarning(Form("Setting it now to: %s", fCDBUri.Data()));
554 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
557 man->SetDefaultStorage(fCDBUri);
560 // Now activate the detector specific CDB storage locations
561 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
562 TObject* obj = fSpecCDBUri[i];
564 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
565 AliDebug(2, Form("Specific CDB storage for %s is set to: %s",obj->GetName(),obj->GetTitle()));
566 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
567 man->SetSpecificStorage(obj->GetName(), obj->GetTitle());
569 AliSysInfo::AddStamp("InitCDB");
572 //_____________________________________________________________________________
573 void AliReconstruction::SetDefaultStorage(const char* uri) {
574 // Store the desired default CDB storage location
575 // Activate it later within the Run() method
581 //_____________________________________________________________________________
582 void AliReconstruction::SetSpecificStorage(const char* calibType, const char* uri) {
583 // Store a detector-specific CDB storage location
584 // Activate it later within the Run() method
586 AliCDBPath aPath(calibType);
587 if(!aPath.IsValid()){
588 // if calibType is not wildcard but it is a valid detector, add "/*" to make it a valid path
589 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
590 if(!strcmp(calibType, fgkDetectorName[iDet])) {
591 aPath.SetPath(Form("%s/*", calibType));
592 AliInfo(Form("Path for specific storage set to %s", aPath.GetPath().Data()));
596 if(!aPath.IsValid()){
597 AliError(Form("Not a valid path or detector: %s", calibType));
602 // // check that calibType refers to a "valid" detector name
603 // Bool_t isDetector = kFALSE;
604 // for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
605 // TString detName = fgkDetectorName[iDet];
606 // if(aPath.GetLevel0() == detName) {
607 // isDetector = kTRUE;
613 // AliError(Form("Not a valid detector: %s", aPath.GetLevel0().Data()));
617 TObject* obj = fSpecCDBUri.FindObject(aPath.GetPath().Data());
618 if (obj) fSpecCDBUri.Remove(obj);
619 fSpecCDBUri.Add(new TNamed(aPath.GetPath().Data(), uri));
623 //_____________________________________________________________________________
624 Bool_t AliReconstruction::SetRunNumberFromData()
626 // The method is called in Run() in order
627 // to set a correct run number.
628 // In case of raw data reconstruction the
629 // run number is taken from the raw data header
631 if (fSetRunNumberFromDataCalled) return kTRUE;
632 fSetRunNumberFromDataCalled = kTRUE;
634 AliCDBManager* man = AliCDBManager::Instance();
637 if(fRawReader->NextEvent()) {
638 if(man->GetRun() > 0) {
639 AliWarning("Run number is taken from raw-event header! Ignoring settings in AliCDBManager!");
641 man->SetRun(fRawReader->GetRunNumber());
642 fRawReader->RewindEvents();
645 if(man->GetRun() > 0) {
646 AliWarning("No raw-data events are found ! Using settings in AliCDBManager !");
649 AliWarning("Neither raw events nor settings in AliCDBManager are found !");
655 AliRunLoader *rl = AliRunLoader::Open(fGAliceFileName.Data());
657 AliError(Form("No run loader found in file %s", fGAliceFileName.Data()));
662 // read run number from gAlice
663 if(rl->GetHeader()) {
664 man->SetRun(rl->GetHeader()->GetRun());
669 AliError("Neither run-loader header nor RawReader objects are found !");
681 //_____________________________________________________________________________
682 void AliReconstruction::SetCDBLock() {
683 // Set CDB lock: from now on it is forbidden to reset the run number
684 // or the default storage or to activate any further storage!
686 AliCDBManager::Instance()->SetLock(1);
689 //_____________________________________________________________________________
690 Bool_t AliReconstruction::MisalignGeometry(const TString& detectors)
692 // Read the alignment objects from CDB.
693 // Each detector is supposed to have the
694 // alignment objects in DET/Align/Data CDB path.
695 // All the detector objects are then collected,
696 // sorted by geometry level (starting from ALIC) and
697 // then applied to the TGeo geometry.
698 // Finally an overlaps check is performed.
700 // Load alignment data from CDB and fill fAlignObjArray
701 if(fLoadAlignFromCDB){
703 TString detStr = detectors;
704 TString loadAlObjsListOfDets = "";
706 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
707 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
708 loadAlObjsListOfDets += fgkDetectorName[iDet];
709 loadAlObjsListOfDets += " ";
710 } // end loop over detectors
711 loadAlObjsListOfDets.Prepend("GRP "); //add alignment objects for non-sensitive modules
712 AliGeomManager::ApplyAlignObjsFromCDB(loadAlObjsListOfDets.Data());
713 AliCDBManager::Instance()->UnloadFromCache("*/Align/*");
715 // Check if the array with alignment objects was
716 // provided by the user. If yes, apply the objects
717 // to the present TGeo geometry
718 if (fAlignObjArray) {
719 if (gGeoManager && gGeoManager->IsClosed()) {
720 if (AliGeomManager::ApplyAlignObjsToGeom(*fAlignObjArray) == kFALSE) {
721 AliError("The misalignment of one or more volumes failed!"
722 "Compare the list of simulated detectors and the list of detector alignment data!");
727 AliError("Can't apply the misalignment! gGeoManager doesn't exist or it is still opened!");
733 if (fAlignObjArray) {
734 fAlignObjArray->Delete();
735 delete fAlignObjArray; fAlignObjArray=NULL;
741 //_____________________________________________________________________________
742 void AliReconstruction::SetGAliceFile(const char* fileName)
744 // set the name of the galice file
746 fGAliceFileName = fileName;
749 //_____________________________________________________________________________
750 void AliReconstruction::SetInput(const char* input)
752 // In case the input string starts with 'mem://', we run in an online mode
753 // and AliRawReaderDateOnline object is created. In all other cases a raw-data
754 // file is assumed. One can give as an input:
755 // mem://: - events taken from DAQ monitoring libs online
757 // mem://<filename> - emulation of the above mode (via DATE monitoring libs)
758 if (input) fRawInput = input;
761 //_____________________________________________________________________________
762 void AliReconstruction::SetOption(const char* detector, const char* option)
764 // set options for the reconstruction of a detector
766 TObject* obj = fOptions.FindObject(detector);
767 if (obj) fOptions.Remove(obj);
768 fOptions.Add(new TNamed(detector, option));
771 //_____________________________________________________________________________
772 void AliReconstruction::SetRecoParam(const char* detector, AliDetectorRecoParam *par)
774 // Set custom reconstruction parameters for a given detector
775 // Single set of parameters for all the events
777 // First check if the reco-params are global
778 if(!strcmp(detector, "GRP")) {
780 fRecoParam.AddDetRecoParam(fgkNDetectors,par);
784 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
785 if(!strcmp(detector, fgkDetectorName[iDet])) {
787 fRecoParam.AddDetRecoParam(iDet,par);
794 //_____________________________________________________________________________
795 Bool_t AliReconstruction::SetFieldMap(Float_t l3Current, Float_t diCurrent, Float_t factor, const char *path) {
796 //------------------------------------------------
797 // The magnetic field map, defined externally...
798 // L3 current 30000 A -> 0.5 T
799 // L3 current 12000 A -> 0.2 T
800 // dipole current 6000 A
801 // The polarities must be the same
802 //------------------------------------------------
803 const Float_t l3NominalCurrent1=30000.; // (A)
804 const Float_t l3NominalCurrent2=12000.; // (A)
805 const Float_t diNominalCurrent =6000. ; // (A)
807 const Float_t tolerance=0.03; // relative current tolerance
808 const Float_t zero=77.; // "zero" current (A)
811 Bool_t dipoleON=kFALSE;
813 TString s=(factor < 0) ? "L3: -" : "L3: +";
815 l3Current = TMath::Abs(l3Current);
816 if (TMath::Abs(l3Current-l3NominalCurrent1)/l3NominalCurrent1 < tolerance) {
817 map=AliMagWrapCheb::k5kG;
820 if (TMath::Abs(l3Current-l3NominalCurrent2)/l3NominalCurrent2 < tolerance) {
821 map=AliMagWrapCheb::k2kG;
824 if (l3Current < zero) {
825 map=AliMagWrapCheb::k2kG;
827 factor=0.; // in fact, this is a global factor...
828 fUniformField=kTRUE; // track with the uniform (zero) B field
830 AliError(Form("Wrong L3 current (%f A)!",l3Current));
834 diCurrent = TMath::Abs(diCurrent);
835 if (TMath::Abs(diCurrent-diNominalCurrent)/diNominalCurrent < tolerance) {
836 // 3% current tolerance...
840 if (diCurrent < zero) { // some small current..
844 AliError(Form("Wrong dipole current (%f A)!",diCurrent));
848 delete fForcedFieldMap;
850 new AliMagWrapCheb("B field map ",s,2,factor,10.,map,dipoleON,path);
852 fForcedFieldMap->Print();
854 AliTracker::SetFieldMap(fForcedFieldMap,fUniformField);
860 Bool_t AliReconstruction::InitGRP() {
861 //------------------------------------
862 // Initialization of the GRP entry
863 //------------------------------------
864 AliCDBEntry* entry = AliCDBManager::Instance()->Get("GRP/GRP/Data");
867 fGRPData = dynamic_cast<TMap*>(entry->GetObject());
869 AliCDBManager::Instance()->UnloadFromCache("GRP/GRP/Data");
873 AliError("No GRP entry found in OCDB!");
877 TObjString *lhcState=
878 dynamic_cast<TObjString*>(fGRPData->GetValue("fLHCState"));
880 AliError("GRP/GRP/Data entry: missing value for the LHC state ! Using UNKNOWN");
883 TObjString *beamType=
884 dynamic_cast<TObjString*>(fGRPData->GetValue("fAliceBeamType"));
886 AliError("GRP/GRP/Data entry: missing value for the beam type ! Using UNKNOWN");
889 TObjString *beamEnergyStr=
890 dynamic_cast<TObjString*>(fGRPData->GetValue("fAliceBeamEnergy"));
891 if (!beamEnergyStr) {
892 AliError("GRP/GRP/Data entry: missing value for the beam energy ! Using 0");
896 dynamic_cast<TObjString*>(fGRPData->GetValue("fRunType"));
898 AliError("GRP/GRP/Data entry: missing value for the run type ! Using UNKNOWN");
901 TObjString *activeDetectors=
902 dynamic_cast<TObjString*>(fGRPData->GetValue("fDetectorMask"));
903 if (!activeDetectors) {
904 AliError("GRP/GRP/Data entry: missing value for the detector mask ! Using 1074790399");
907 fRunInfo = new AliRunInfo(lhcState ? lhcState->GetString().Data() : "UNKNOWN",
908 beamType ? beamType->GetString().Data() : "UNKNOWN",
909 beamEnergyStr ? beamEnergyStr->GetString().Atof() : 0,
910 runType ? runType->GetString().Data() : "UNKNOWN",
911 activeDetectors ? activeDetectors->GetString().Atoi() : 1074790399);
913 // Process the list of active detectors
914 if (activeDetectors && activeDetectors->GetString().IsDigit()) {
915 UInt_t detMask = activeDetectors->GetString().Atoi();
916 fLoadCDB.Form("%s %s %s %s",
917 fRunLocalReconstruction.Data(),
920 fQADetectors.Data());
921 fRunLocalReconstruction = MatchDetectorList(fRunLocalReconstruction,detMask);
922 fRunTracking = MatchDetectorList(fRunTracking,detMask);
923 fFillESD = MatchDetectorList(fFillESD,detMask);
924 fQADetectors = MatchDetectorList(fQADetectors,detMask);
925 fLoadCDB = MatchDetectorList(fLoadCDB,detMask);
926 fLoadAlignData = MatchDetectorList(fLoadAlignData,detMask);
929 AliInfo("===================================================================================");
930 AliInfo(Form("Running local reconstruction for detectors: %s",fRunLocalReconstruction.Data()));
931 AliInfo(Form("Running tracking for detectors: %s",fRunTracking.Data()));
932 AliInfo(Form("Filling ESD for detectors: %s",fFillESD.Data()));
933 AliInfo(Form("Quality assurance is active for detectors: %s",fQADetectors.Data()));
934 AliInfo(Form("CDB and reconstruction parameters are loaded for detectors: %s",fLoadCDB.Data()));
935 AliInfo("===================================================================================");
937 //*** Dealing with the magnetic field map
938 if (AliTracker::GetFieldMap()) {
939 AliInfo("Running with the externally set B field !");
941 // Construct the field map out of the information retrieved from GRP.
946 TObjString *l3Current=
947 dynamic_cast<TObjString*>(fGRPData->GetValue("fL3Current"));
949 AliError("GRP/GRP/Data entry: missing value for the L3 current !");
952 TObjString *l3Polarity=
953 dynamic_cast<TObjString*>(fGRPData->GetValue("fL3Polarity"));
955 AliError("GRP/GRP/Data entry: missing value for the L3 polarity !");
960 TObjString *diCurrent=
961 dynamic_cast<TObjString*>(fGRPData->GetValue("fDipoleCurrent"));
963 AliError("GRP/GRP/Data entry: missing value for the dipole current !");
966 TObjString *diPolarity=
967 dynamic_cast<TObjString*>(fGRPData->GetValue("fDipolePolarity"));
969 AliError("GRP/GRP/Data entry: missing value for the dipole polarity !");
974 Float_t l3Cur=TMath::Abs(atof(l3Current->GetName()));
975 Float_t diCur=TMath::Abs(atof(diCurrent->GetName()));
976 Float_t l3Pol=atof(l3Polarity->GetName());
978 if (l3Pol != 0.) factor=-1.;
981 if (!SetFieldMap(l3Cur, diCur, factor)) {
982 AliFatal("Failed to creat a B field map ! Exiting...");
984 AliInfo("Running with the B field constructed out of GRP !");
987 AliFatal("B field is neither set nor constructed from GRP ! Exitig...");
992 //*** Get the diamond profile from OCDB
993 entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertex");
995 if (fMeanVertexConstraint)
996 fDiamondProfile = dynamic_cast<AliESDVertex*> (entry->GetObject());
998 AliError("No diamond profile found in OCDB!");
1001 entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertexTPC");
1003 if (fMeanVertexConstraint)
1004 fDiamondProfileTPC = dynamic_cast<AliESDVertex*> (entry->GetObject());
1006 AliError("No diamond profile found in OCDB!");
1012 //_____________________________________________________________________________
1013 Bool_t AliReconstruction::LoadCDB()
1015 AliCodeTimerAuto("");
1017 AliCDBManager::Instance()->Get("GRP/CTP/Config");
1019 TString detStr = fLoadCDB;
1020 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1021 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1022 AliCDBManager::Instance()->GetAll(Form("%s/Calib/*",fgkDetectorName[iDet]));
1027 //_____________________________________________________________________________
1028 Bool_t AliReconstruction::Run(const char* input)
1031 AliCodeTimerAuto("");
1034 if (GetAbort() != TSelector::kContinue) return kFALSE;
1036 TChain *chain = NULL;
1037 if (fRawReader && (chain = fRawReader->GetChain())) {
1040 gProof->AddInput(this);
1042 outputFile.SetProtocol("root",kTRUE);
1043 outputFile.SetHost(gSystem->HostName());
1044 outputFile.SetFile(Form("%s/AliESDs.root",gSystem->pwd()));
1045 AliInfo(Form("Output file with ESDs is %s",outputFile.GetUrl()));
1046 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE",outputFile.GetUrl()));
1048 chain->Process("AliReconstruction");
1051 chain->Process(this);
1056 if (GetAbort() != TSelector::kContinue) return kFALSE;
1058 if (GetAbort() != TSelector::kContinue) return kFALSE;
1059 //******* The loop over events
1061 while ((iEvent < fRunLoader->GetNumberOfEvents()) ||
1062 (fRawReader && fRawReader->NextEvent())) {
1063 if (!ProcessEvent(iEvent)) {
1064 Abort("ProcessEvent",TSelector::kAbortFile);
1070 if (GetAbort() != TSelector::kContinue) return kFALSE;
1072 if (GetAbort() != TSelector::kContinue) return kFALSE;
1078 //_____________________________________________________________________________
1079 void AliReconstruction::InitRawReader(const char* input)
1081 AliCodeTimerAuto("");
1083 // Init raw-reader and
1084 // set the input in case of raw data
1085 if (input) fRawInput = input;
1086 fRawReader = AliRawReader::Create(fRawInput.Data());
1088 AliInfo("Reconstruction will run over digits");
1090 if (!fEquipIdMap.IsNull() && fRawReader)
1091 fRawReader->LoadEquipmentIdsMap(fEquipIdMap);
1093 if (!fUseHLTData.IsNull()) {
1094 // create the RawReaderHLT which performs redirection of HLT input data for
1095 // the specified detectors
1096 AliRawReader* pRawReader=AliRawHLTManager::CreateRawReaderHLT(fRawReader, fUseHLTData.Data());
1098 fParentRawReader=fRawReader;
1099 fRawReader=pRawReader;
1101 AliError(Form("can not create Raw Reader for HLT input %s", fUseHLTData.Data()));
1104 AliSysInfo::AddStamp("CreateRawReader");
1107 //_____________________________________________________________________________
1108 void AliReconstruction::InitRun(const char* input)
1110 // Initialization of raw-reader,
1111 // run number, CDB etc.
1112 AliCodeTimerAuto("");
1113 AliSysInfo::AddStamp("Start");
1115 // Initialize raw-reader if any
1116 InitRawReader(input);
1118 // Initialize the CDB storage
1121 // Set run number in CDBManager (if it is not already set by the user)
1122 if (!SetRunNumberFromData()) {
1123 Abort("SetRunNumberFromData", TSelector::kAbortProcess);
1127 // Set CDB lock: from now on it is forbidden to reset the run number
1128 // or the default storage or to activate any further storage!
1133 //_____________________________________________________________________________
1134 void AliReconstruction::Begin(TTree *)
1136 // Initialize AlReconstruction before
1137 // going into the event loop
1138 // Should follow the TSelector convention
1139 // i.e. initialize only the object on the client side
1140 AliCodeTimerAuto("");
1142 AliReconstruction *reco = NULL;
1144 if ((reco = (AliReconstruction*)fInput->FindObject("AliReconstruction"))) {
1147 AliSysInfo::AddStamp("ReadInputInBegin");
1150 // Import ideal TGeo geometry and apply misalignment
1152 TString geom(gSystem->DirName(fGAliceFileName));
1153 geom += "/geometry.root";
1154 AliGeomManager::LoadGeometry(geom.Data());
1156 Abort("LoadGeometry", TSelector::kAbortProcess);
1159 AliSysInfo::AddStamp("LoadGeom");
1160 TString detsToCheck=fRunLocalReconstruction;
1161 if(!AliGeomManager::CheckSymNamesLUT(detsToCheck.Data())) {
1162 Abort("CheckSymNamesLUT", TSelector::kAbortProcess);
1165 AliSysInfo::AddStamp("CheckGeom");
1169 Abort("InitGRP", TSelector::kAbortProcess);
1172 AliSysInfo::AddStamp("InitGRP");
1174 if (!MisalignGeometry(fLoadAlignData)) {
1175 Abort("MisalignGeometry", TSelector::kAbortProcess);
1178 AliCDBManager::Instance()->UnloadFromCache("GRP/Geometry/Data");
1179 AliSysInfo::AddStamp("MisalignGeom");
1182 Abort("LoadCDB", TSelector::kAbortProcess);
1185 AliSysInfo::AddStamp("LoadCDB");
1187 // Read the reconstruction parameters from OCDB
1188 if (!InitRecoParams()) {
1189 AliWarning("Not all detectors have correct RecoParam objects initialized");
1191 AliSysInfo::AddStamp("InitRecoParams");
1194 if (reco) *reco = *this;
1195 fInput->Add(gGeoManager);
1197 fInput->Add(const_cast<TMap*>(AliCDBManager::Instance()->GetEntryCache()));
1198 fInput->Add(new TParameter<Int_t>("RunNumber",AliCDBManager::Instance()->GetRun()));
1199 AliMagF *magFieldMap = (AliMagF*)AliTracker::GetFieldMap();
1200 magFieldMap->SetName("MagneticFieldMap");
1201 fInput->Add(magFieldMap);
1206 //_____________________________________________________________________________
1207 void AliReconstruction::SlaveBegin(TTree*)
1209 // Initialization related to run-loader,
1210 // vertexer, trackers, recontructors
1211 // In proof mode it is executed on the slave
1212 AliCodeTimerAuto("");
1214 TProofOutputFile *outProofFile = NULL;
1216 if (AliReconstruction *reco = (AliReconstruction*)fInput->FindObject("AliReconstruction")) {
1219 if (TGeoManager *tgeo = (TGeoManager*)fInput->FindObject("Geometry")) {
1221 AliGeomManager::SetGeometry(tgeo);
1223 if (TMap *entryCache = (TMap*)fInput->FindObject("CDBEntryCache")) {
1224 Int_t runNumber = -1;
1225 if (TProof::GetParameter(fInput,"RunNumber",runNumber) == 0) {
1226 AliCDBManager *man = AliCDBManager::Instance(entryCache,runNumber);
1227 man->SetCacheFlag(kTRUE);
1228 man->SetLock(kTRUE);
1232 if (AliMagF *map = (AliMagF*)fInput->FindObject("MagneticFieldMap")) {
1233 AliTracker::SetFieldMap(map,fUniformField);
1235 if (TNamed *outputFileName = (TNamed *) fInput->FindObject("PROOF_OUTPUTFILE")) {
1236 outProofFile = new TProofOutputFile(gSystem->BaseName(TUrl(outputFileName->GetTitle()).GetFile()));
1237 outProofFile->SetOutputFileName(outputFileName->GetTitle());
1238 fOutput->Add(outProofFile);
1240 AliSysInfo::AddStamp("ReadInputInSlaveBegin");
1243 // get the run loader
1244 if (!InitRunLoader()) {
1245 Abort("InitRunLoader", TSelector::kAbortProcess);
1248 AliSysInfo::AddStamp("LoadLoader");
1250 ftVertexer = new AliVertexerTracks(AliTracker::GetBz());
1251 if(fDiamondProfile && fMeanVertexConstraint) ftVertexer->SetVtxStart(fDiamondProfile);
1254 if (fRunVertexFinder && !CreateVertexer()) {
1255 Abort("CreateVertexer", TSelector::kAbortProcess);
1258 AliSysInfo::AddStamp("CreateVertexer");
1261 if (!fRunTracking.IsNull() && !CreateTrackers(fRunTracking)) {
1262 Abort("CreateTrackers", TSelector::kAbortProcess);
1265 AliSysInfo::AddStamp("CreateTrackers");
1267 // create the ESD output file and tree
1268 if (!outProofFile) {
1269 ffile = TFile::Open("AliESDs.root", "RECREATE");
1270 ffile->SetCompressionLevel(2);
1271 if (!ffile->IsOpen()) {
1272 Abort("OpenESDFile", TSelector::kAbortProcess);
1277 if (!(ffile = outProofFile->OpenFile("RECREATE"))) {
1278 Abort(Form("Problems opening output PROOF file: %s/%s",
1279 outProofFile->GetDir(), outProofFile->GetFileName()),
1280 TSelector::kAbortProcess);
1285 ftree = new TTree("esdTree", "Tree with ESD objects");
1286 fesd = new AliESDEvent();
1287 fesd->CreateStdContent();
1288 fesd->WriteToTree(ftree);
1290 fhlttree = new TTree("HLTesdTree", "Tree with HLT ESD objects");
1291 fhltesd = new AliESDEvent();
1292 fhltesd->CreateStdContent();
1293 fhltesd->WriteToTree(fhlttree);
1296 if (fWriteESDfriend) {
1297 fesdf = new AliESDfriend();
1298 TBranch *br=ftree->Branch("ESDfriend.","AliESDfriend", &fesdf);
1299 br->SetFile("AliESDfriends.root");
1300 fesd->AddObject(fesdf);
1303 ProcInfo_t ProcInfo;
1304 gSystem->GetProcInfo(&ProcInfo);
1305 AliInfo(Form("Current memory usage %d %d", ProcInfo.fMemResident, ProcInfo.fMemVirtual));
1308 //Initialize the QA and start of cycle
1310 fQASteer = new AliQADataMakerSteer("rec") ;
1311 fQASteer->SetActiveDetectors(fQADetectors) ;
1312 for (Int_t det = 0 ; det < AliQA::kNDET ; det++) {
1313 fQASteer->SetCycleLength(AliQA::DETECTORINDEX_t(det), fQACycles[det]) ;
1314 if (fQAWriteExpert[det])
1315 fQASteer->SetWriteExpert(AliQA::DETECTORINDEX_t(det)) ;
1318 if (!fRawReader && fQATasks.Contains(AliQA::kRAWS))
1319 fQATasks.ReplaceAll(Form("%d",AliQA::kRAWS), "") ;
1320 fQASteer->SetTasks(fQATasks) ;
1321 fQASteer->InitQADataMaker(AliCDBManager::Instance()->GetRun(), fRecoParam) ;
1325 Bool_t sameCycle = kFALSE ;
1326 if (!fQASteer) fQASteer = new AliQADataMakerSteer("rec") ;
1327 AliQADataMaker *qadm = fQASteer->GetQADataMaker(AliQA::kGLOBAL);
1328 AliInfo(Form("Initializing the global QA data maker"));
1329 if (fQATasks.Contains(Form("%d", AliQA::kRECPOINTS))) {
1330 qadm->StartOfCycle(AliQA::kRECPOINTS, AliCDBManager::Instance()->GetRun(), sameCycle) ;
1331 TObjArray *arr=qadm->Init(AliQA::kRECPOINTS);
1332 AliTracker::SetResidualsArray(arr);
1335 if (fQATasks.Contains(Form("%d", AliQA::kESDS))) {
1336 qadm->StartOfCycle(AliQA::kESDS, AliCDBManager::Instance()->GetRun(), sameCycle) ;
1337 qadm->Init(AliQA::kESDS);
1341 //Initialize the Plane Efficiency framework
1342 if (fRunPlaneEff && !InitPlaneEff()) {
1343 Abort("InitPlaneEff", TSelector::kAbortProcess);
1347 if (strcmp(gProgName,"alieve") == 0)
1348 fRunAliEVE = InitAliEVE();
1353 //_____________________________________________________________________________
1354 Bool_t AliReconstruction::Process(Long64_t entry)
1356 // run the reconstruction over a single entry
1357 // from the chain with raw data
1358 AliCodeTimerAuto("");
1360 TTree *currTree = fChain->GetTree();
1361 AliRawEvent *event = new AliRawEvent;
1362 currTree->SetBranchAddress("rawevent",&event);
1363 currTree->GetEntry(entry);
1364 fRawReader = new AliRawReaderRoot(event);
1365 fStatus = ProcessEvent(fRunLoader->GetNumberOfEvents());
1373 //_____________________________________________________________________________
1374 void AliReconstruction::Init(TTree *tree)
1377 AliError("The input tree is not found!");
1383 //_____________________________________________________________________________
1384 Bool_t AliReconstruction::ProcessEvent(Int_t iEvent)
1386 // run the reconstruction over a single event
1387 // The event loop is steered in Run method
1389 AliCodeTimerAuto("");
1391 if (iEvent >= fRunLoader->GetNumberOfEvents()) {
1392 fRunLoader->SetEventNumber(iEvent);
1393 fRunLoader->GetHeader()->Reset(fRawReader->GetRunNumber(),
1395 fRunLoader->TreeE()->Fill();
1396 if (fRawReader && fRawReader->UseAutoSaveESD())
1397 fRunLoader->TreeE()->AutoSave("SaveSelf");
1400 if ((iEvent < fFirstEvent) || ((fLastEvent >= 0) && (iEvent > fLastEvent))) {
1404 AliInfo(Form("processing event %d", iEvent));
1406 // Fill Event-info object
1408 fRecoParam.SetEventSpecie(fRunInfo,fEventInfo);
1410 // Set the reco-params
1412 TString detStr = fLoadCDB;
1413 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1414 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1415 AliReconstructor *reconstructor = GetReconstructor(iDet);
1416 if (reconstructor && fRecoParam.GetDetRecoParamArray(iDet)) {
1417 const AliDetectorRecoParam *par = fRecoParam.GetDetRecoParam(iDet);
1418 reconstructor->SetRecoParam(par);
1423 fRunLoader->GetEvent(iEvent);
1427 fQASteer->RunOneEvent(fRawReader) ;
1429 // local single event reconstruction
1430 if (!fRunLocalReconstruction.IsNull()) {
1431 TString detectors=fRunLocalReconstruction;
1432 // run HLT event reconstruction first
1433 // ;-( IsSelected changes the string
1434 if (IsSelected("HLT", detectors) &&
1435 !RunLocalEventReconstruction("HLT")) {
1436 if (fStopOnError) {CleanUp(); return kFALSE;}
1438 detectors=fRunLocalReconstruction;
1439 detectors.ReplaceAll("HLT", "");
1440 if (!RunLocalEventReconstruction(detectors)) {
1441 if (fStopOnError) {CleanUp(); return kFALSE;}
1445 fesd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
1446 fhltesd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
1447 fesd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
1448 fhltesd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
1450 // Set magnetic field from the tracker
1451 fesd->SetMagneticField(AliTracker::GetBz());
1452 fhltesd->SetMagneticField(AliTracker::GetBz());
1454 // Set most probable pt, for B=0 tracking
1455 const AliGRPRecoParam *grpRecoParam = dynamic_cast<const AliGRPRecoParam*>(fRecoParam.GetDetRecoParam(fgkNDetectors));
1456 if (grpRecoParam) AliExternalTrackParam::SetMostProbablePt(grpRecoParam->GetMostProbablePt());
1458 // Fill raw-data error log into the ESD
1459 if (fRawReader) FillRawDataErrorLog(iEvent,fesd);
1462 if (fRunVertexFinder) {
1463 if (!RunVertexFinder(fesd)) {
1464 if (fStopOnError) {CleanUp(); return kFALSE;}
1469 if (!fRunTracking.IsNull()) {
1470 if (fRunMuonTracking) {
1471 if (!RunMuonTracking(fesd)) {
1472 if (fStopOnError) {CleanUp(); return kFALSE;}
1478 if (!fRunTracking.IsNull()) {
1479 if (!RunTracking(fesd)) {
1480 if (fStopOnError) {CleanUp(); return kFALSE;}
1485 if (!fFillESD.IsNull()) {
1486 TString detectors=fFillESD;
1487 // run HLT first and on hltesd
1488 // ;-( IsSelected changes the string
1489 if (IsSelected("HLT", detectors) &&
1490 !FillESD(fhltesd, "HLT")) {
1491 if (fStopOnError) {CleanUp(); return kFALSE;}
1494 // Temporary fix to avoid problems with HLT that overwrites the offline ESDs
1495 if (detectors.Contains("ALL")) {
1497 for (Int_t idet=0; idet<fgkNDetectors; ++idet){
1498 detectors += fgkDetectorName[idet];
1502 detectors.ReplaceAll("HLT", "");
1503 if (!FillESD(fesd, detectors)) {
1504 if (fStopOnError) {CleanUp(); return kFALSE;}
1508 // fill Event header information from the RawEventHeader
1509 if (fRawReader){FillRawEventHeaderESD(fesd);}
1512 AliESDpid::MakePID(fesd);
1514 if (fFillTriggerESD) {
1515 if (!FillTriggerESD(fesd)) {
1516 if (fStopOnError) {CleanUp(); return kFALSE;}
1523 // Propagate track to the beam pipe (if not already done by ITS)
1525 const Int_t ntracks = fesd->GetNumberOfTracks();
1526 const Double_t kBz = fesd->GetMagneticField();
1527 const Double_t kRadius = 2.8; //something less than the beam pipe radius
1530 UShort_t *selectedIdx=new UShort_t[ntracks];
1532 for (Int_t itrack=0; itrack<ntracks; itrack++){
1533 const Double_t kMaxStep = 5; //max step over the material
1536 AliESDtrack *track = fesd->GetTrack(itrack);
1537 if (!track) continue;
1539 AliExternalTrackParam *tpcTrack =
1540 (AliExternalTrackParam *)track->GetTPCInnerParam();
1544 PropagateTrackTo(tpcTrack,kRadius,track->GetMass(),kMaxStep,kTRUE);
1547 Int_t n=trkArray.GetEntriesFast();
1548 selectedIdx[n]=track->GetID();
1549 trkArray.AddLast(tpcTrack);
1552 //Tracks refitted by ITS should already be at the SPD vertex
1553 if (track->IsOn(AliESDtrack::kITSrefit)) continue;
1556 PropagateTrackTo(track,kRadius,track->GetMass(),kMaxStep,kTRUE);
1557 track->RelateToVertex(fesd->GetPrimaryVertexSPD(), kBz, kVeryBig);
1562 // Improve the reconstructed primary vertex position using the tracks
1564 Bool_t runVertexFinderTracks = fRunVertexFinderTracks;
1565 if(fesd->GetPrimaryVertexSPD()) {
1566 TString vtitle = fesd->GetPrimaryVertexSPD()->GetTitle();
1567 if(vtitle.Contains("cosmics")) {
1568 runVertexFinderTracks=kFALSE;
1572 if (runVertexFinderTracks) {
1573 // Get the global reco-params. They are atposition 16 inside the array of detectors in fRecoParam
1574 const AliGRPRecoParam *grpRecoParam = dynamic_cast<const AliGRPRecoParam*>(fRecoParam.GetDetRecoParam(fgkNDetectors));
1576 // TPC + ITS primary vertex
1577 ftVertexer->SetITSMode();
1578 // get cuts for vertexer from AliGRPRecoParam
1580 Int_t nCutsVertexer = grpRecoParam->GetVertexerTracksNCuts();
1581 Double_t *cutsVertexer = new Double_t[nCutsVertexer];
1582 grpRecoParam->GetVertexerTracksCutsITS(cutsVertexer);
1583 ftVertexer->SetCuts(cutsVertexer);
1584 delete [] cutsVertexer; cutsVertexer = NULL;
1586 if(fDiamondProfile && fMeanVertexConstraint) {
1587 ftVertexer->SetVtxStart(fDiamondProfile);
1589 ftVertexer->SetConstraintOff();
1591 AliESDVertex *pvtx=ftVertexer->FindPrimaryVertex(fesd);
1593 if (pvtx->GetStatus()) {
1594 fesd->SetPrimaryVertex(pvtx);
1595 for (Int_t i=0; i<ntracks; i++) {
1596 AliESDtrack *t = fesd->GetTrack(i);
1597 t->RelateToVertex(pvtx, kBz, kVeryBig);
1602 // TPC-only primary vertex
1603 ftVertexer->SetTPCMode();
1604 // get cuts for vertexer from AliGRPRecoParam
1606 Int_t nCutsVertexer = grpRecoParam->GetVertexerTracksNCuts();
1607 Double_t *cutsVertexer = new Double_t[nCutsVertexer];
1608 grpRecoParam->GetVertexerTracksCutsTPC(cutsVertexer);
1609 ftVertexer->SetCuts(cutsVertexer);
1610 delete [] cutsVertexer; cutsVertexer = NULL;
1612 if(fDiamondProfileTPC && fMeanVertexConstraint) {
1613 ftVertexer->SetVtxStart(fDiamondProfileTPC);
1615 ftVertexer->SetConstraintOff();
1617 pvtx=ftVertexer->FindPrimaryVertex(&trkArray,selectedIdx);
1619 if (pvtx->GetStatus()) {
1620 fesd->SetPrimaryVertexTPC(pvtx);
1621 for (Int_t i=0; i<ntracks; i++) {
1622 AliESDtrack *t = fesd->GetTrack(i);
1623 t->RelateToVertexTPC(pvtx, kBz, kVeryBig);
1629 delete[] selectedIdx;
1631 if(fDiamondProfile) fesd->SetDiamond(fDiamondProfile);
1636 AliV0vertexer vtxer;
1637 vtxer.Tracks2V0vertices(fesd);
1639 if (fRunCascadeFinder) {
1641 AliCascadeVertexer cvtxer;
1642 cvtxer.V0sTracks2CascadeVertices(fesd);
1647 if (fCleanESD) CleanESD(fesd);
1650 fQASteer->RunOneEvent(fesd) ;
1653 AliQADataMaker *qadm = fQASteer->GetQADataMaker(AliQA::kGLOBAL);
1654 if (qadm && fQATasks.Contains(Form("%d", AliQA::kESDS)))
1655 qadm->Exec(AliQA::kESDS, fesd);
1658 if (fWriteESDfriend) {
1659 fesdf->~AliESDfriend();
1660 new (fesdf) AliESDfriend(); // Reset...
1661 fesd->GetESDfriend(fesdf);
1665 // Auto-save the ESD tree in case of prompt reco @P2
1666 if (fRawReader && fRawReader->UseAutoSaveESD()) {
1667 ftree->AutoSave("SaveSelf");
1668 TFile *friendfile = (TFile *)(gROOT->GetListOfFiles()->FindObject("AliESDfriends.root"));
1669 if (friendfile) friendfile->Save();
1676 if (fRunAliEVE) RunAliEVE();
1680 if (fWriteESDfriend) {
1681 fesdf->~AliESDfriend();
1682 new (fesdf) AliESDfriend(); // Reset...
1685 ProcInfo_t ProcInfo;
1686 gSystem->GetProcInfo(&ProcInfo);
1687 AliInfo(Form("Event %d -> Current memory usage %d %d",iEvent, ProcInfo.fMemResident, ProcInfo.fMemVirtual));
1690 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1691 if (fReconstructor[iDet])
1692 fReconstructor[iDet]->SetRecoParam(NULL);
1695 fQASteer->Increment() ;
1699 //_____________________________________________________________________________
1700 void AliReconstruction::SlaveTerminate()
1702 // Finalize the run on the slave side
1703 // Called after the exit
1704 // from the event loop
1705 AliCodeTimerAuto("");
1707 if (fIsNewRunLoader) { // galice.root didn't exist
1708 fRunLoader->WriteHeader("OVERWRITE");
1709 fRunLoader->CdGAFile();
1710 fRunLoader->Write(0, TObject::kOverwrite);
1713 ftree->GetUserInfo()->Add(fesd);
1714 fhlttree->GetUserInfo()->Add(fhltesd);
1716 const TMap *cdbMap = AliCDBManager::Instance()->GetStorageMap();
1717 const TList *cdbList = AliCDBManager::Instance()->GetRetrievedIds();
1719 TMap *cdbMapCopy = new TMap(cdbMap->GetEntries());
1720 cdbMapCopy->SetOwner(1);
1721 cdbMapCopy->SetName("cdbMap");
1722 TIter iter(cdbMap->GetTable());
1725 while((pair = dynamic_cast<TPair*> (iter.Next()))){
1726 TObjString* keyStr = dynamic_cast<TObjString*> (pair->Key());
1727 TObjString* valStr = dynamic_cast<TObjString*> (pair->Value());
1728 cdbMapCopy->Add(new TObjString(keyStr->GetName()), new TObjString(valStr->GetName()));
1731 TList *cdbListCopy = new TList();
1732 cdbListCopy->SetOwner(1);
1733 cdbListCopy->SetName("cdbList");
1735 TIter iter2(cdbList);
1738 while((id = dynamic_cast<AliCDBId*> (iter2.Next()))){
1739 cdbListCopy->Add(new TObjString(id->ToString().Data()));
1742 ftree->GetUserInfo()->Add(cdbMapCopy);
1743 ftree->GetUserInfo()->Add(cdbListCopy);
1746 if(fESDPar.Contains("ESD.par")){
1747 AliInfo("Attaching ESD.par to Tree");
1748 TNamed *fn = CopyFileToTNamed(fESDPar.Data(),"ESD.par");
1749 ftree->GetUserInfo()->Add(fn);
1755 if (fWriteESDfriend)
1756 ftree->SetBranchStatus("ESDfriend*",0);
1757 // we want to have only one tree version number
1758 ftree->Write(ftree->GetName(),TObject::kOverwrite);
1761 // Finish with Plane Efficiency evaluation: before of CleanUp !!!
1762 if (fRunPlaneEff && !FinishPlaneEff()) {
1763 AliWarning("Finish PlaneEff evaluation failed");
1766 // End of cycle for the in-loop
1768 fQASteer->EndOfCycle() ;
1771 AliQADataMaker *qadm = fQASteer->GetQADataMaker(AliQA::kGLOBAL);
1773 if (fQATasks.Contains(Form("%d", AliQA::kRECPOINTS)))
1774 qadm->EndOfCycle(AliQA::kRECPOINTS);
1775 if (fQATasks.Contains(Form("%d", AliQA::kESDS)))
1776 qadm->EndOfCycle(AliQA::kESDS);
1784 //_____________________________________________________________________________
1785 void AliReconstruction::Terminate()
1787 // Create tags for the events in the ESD tree (the ESD tree is always present)
1788 // In case of empty events the tags will contain dummy values
1789 AliCodeTimerAuto("");
1791 AliESDTagCreator *esdtagCreator = new AliESDTagCreator();
1792 esdtagCreator->CreateESDTags(fFirstEvent,fLastEvent,fGRPData);
1794 // Cleanup of CDB manager: cache and active storages!
1795 AliCDBManager::Instance()->ClearCache();
1798 //_____________________________________________________________________________
1799 Bool_t AliReconstruction::RunLocalEventReconstruction(const TString& detectors)
1801 // run the local reconstruction
1803 static Int_t eventNr=0;
1804 AliCodeTimerAuto("")
1806 TString detStr = detectors;
1807 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1808 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1809 AliReconstructor* reconstructor = GetReconstructor(iDet);
1810 if (!reconstructor) continue;
1811 AliLoader* loader = fLoader[iDet];
1812 // Matthias April 2008: temporary fix to run HLT reconstruction
1813 // although the HLT loader is missing
1814 if (strcmp(fgkDetectorName[iDet], "HLT")==0) {
1816 reconstructor->Reconstruct(fRawReader, NULL);
1819 reconstructor->Reconstruct(dummy, NULL);
1824 AliWarning(Form("No loader is defined for %s!",fgkDetectorName[iDet]));
1827 // conversion of digits
1828 if (fRawReader && reconstructor->HasDigitConversion()) {
1829 AliInfo(Form("converting raw data digits into root objects for %s",
1830 fgkDetectorName[iDet]));
1831 // AliCodeTimerAuto(Form("converting raw data digits into root objects for %s",
1832 // fgkDetectorName[iDet]));
1833 loader->LoadDigits("update");
1834 loader->CleanDigits();
1835 loader->MakeDigitsContainer();
1836 TTree* digitsTree = loader->TreeD();
1837 reconstructor->ConvertDigits(fRawReader, digitsTree);
1838 loader->WriteDigits("OVERWRITE");
1839 loader->UnloadDigits();
1841 // local reconstruction
1842 AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
1843 //AliCodeTimerAuto(Form("running reconstruction for %s", fgkDetectorName[iDet]));
1844 loader->LoadRecPoints("update");
1845 loader->CleanRecPoints();
1846 loader->MakeRecPointsContainer();
1847 TTree* clustersTree = loader->TreeR();
1848 if (fRawReader && !reconstructor->HasDigitConversion()) {
1849 reconstructor->Reconstruct(fRawReader, clustersTree);
1851 loader->LoadDigits("read");
1852 TTree* digitsTree = loader->TreeD();
1854 AliError(Form("Can't get the %s digits tree", fgkDetectorName[iDet]));
1855 if (fStopOnError) return kFALSE;
1857 reconstructor->Reconstruct(digitsTree, clustersTree);
1859 loader->UnloadDigits();
1862 TString detQAStr(fQADetectors) ;
1864 fQASteer->RunOneEventInOneDetector(iDet, clustersTree) ;
1866 loader->WriteRecPoints("OVERWRITE");
1867 loader->UnloadRecPoints();
1868 AliSysInfo::AddStamp(Form("LRec%s_%d",fgkDetectorName[iDet],eventNr), iDet,1,eventNr);
1870 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1871 AliError(Form("the following detectors were not found: %s",
1873 if (fStopOnError) return kFALSE;
1879 //_____________________________________________________________________________
1880 Bool_t AliReconstruction::RunVertexFinder(AliESDEvent*& esd)
1882 // run the barrel tracking
1884 AliCodeTimerAuto("")
1886 AliESDVertex* vertex = NULL;
1887 Double_t vtxPos[3] = {0, 0, 0};
1888 Double_t vtxErr[3] = {0.07, 0.07, 0.1};
1889 TArrayF mcVertex(3);
1890 if (fRunLoader->GetHeader() && fRunLoader->GetHeader()->GenEventHeader()) {
1891 fRunLoader->GetHeader()->GenEventHeader()->PrimaryVertex(mcVertex);
1892 for (Int_t i = 0; i < 3; i++) vtxPos[i] = mcVertex[i];
1896 AliInfo("running the ITS vertex finder");
1898 fLoader[0]->LoadRecPoints();
1899 TTree* cltree = fLoader[0]->TreeR();
1901 if(fDiamondProfile) fVertexer->SetVtxStart(fDiamondProfile);
1902 vertex = fVertexer->FindVertexForCurrentEvent(cltree);
1905 AliError("Can't get the ITS cluster tree");
1907 fLoader[0]->UnloadRecPoints();
1910 AliError("Can't get the ITS loader");
1913 AliWarning("Vertex not found");
1914 vertex = new AliESDVertex();
1915 vertex->SetName("default");
1918 vertex->SetName("reconstructed");
1922 AliInfo("getting the primary vertex from MC");
1923 vertex = new AliESDVertex(vtxPos, vtxErr);
1927 vertex->GetXYZ(vtxPos);
1928 vertex->GetSigmaXYZ(vtxErr);
1930 AliWarning("no vertex reconstructed");
1931 vertex = new AliESDVertex(vtxPos, vtxErr);
1933 esd->SetPrimaryVertexSPD(vertex);
1934 // if SPD multiplicity has been determined, it is stored in the ESD
1935 AliMultiplicity *mult = fVertexer->GetMultiplicity();
1936 if(mult)esd->SetMultiplicity(mult);
1938 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1939 if (fTracker[iDet]) fTracker[iDet]->SetVertex(vtxPos, vtxErr);
1946 //_____________________________________________________________________________
1947 Bool_t AliReconstruction::RunHLTTracking(AliESDEvent*& esd)
1949 // run the HLT barrel tracking
1951 AliCodeTimerAuto("")
1954 AliError("Missing runLoader!");
1958 AliInfo("running HLT tracking");
1960 // Get a pointer to the HLT reconstructor
1961 AliReconstructor *reconstructor = GetReconstructor(fgkNDetectors-1);
1962 if (!reconstructor) return kFALSE;
1965 for (Int_t iDet = 1; iDet >= 0; iDet--) {
1966 TString detName = fgkDetectorName[iDet];
1967 AliDebug(1, Form("%s HLT tracking", detName.Data()));
1968 reconstructor->SetOption(detName.Data());
1969 AliTracker *tracker = reconstructor->CreateTracker();
1971 AliWarning(Form("couldn't create a HLT tracker for %s", detName.Data()));
1972 if (fStopOnError) return kFALSE;
1976 Double_t vtxErr[3]={0.005,0.005,0.010};
1977 const AliESDVertex *vertex = esd->GetVertex();
1978 vertex->GetXYZ(vtxPos);
1979 tracker->SetVertex(vtxPos,vtxErr);
1981 fLoader[iDet]->LoadRecPoints("read");
1982 TTree* tree = fLoader[iDet]->TreeR();
1984 AliError(Form("Can't get the %s cluster tree", detName.Data()));
1987 tracker->LoadClusters(tree);
1989 if (tracker->Clusters2Tracks(esd) != 0) {
1990 AliError(Form("HLT %s Clusters2Tracks failed", fgkDetectorName[iDet]));
1994 tracker->UnloadClusters();
2002 //_____________________________________________________________________________
2003 Bool_t AliReconstruction::RunMuonTracking(AliESDEvent*& esd)
2005 // run the muon spectrometer tracking
2007 AliCodeTimerAuto("")
2010 AliError("Missing runLoader!");
2013 Int_t iDet = 7; // for MUON
2015 AliInfo("is running...");
2017 // Get a pointer to the MUON reconstructor
2018 AliReconstructor *reconstructor = GetReconstructor(iDet);
2019 if (!reconstructor) return kFALSE;
2022 TString detName = fgkDetectorName[iDet];
2023 AliDebug(1, Form("%s tracking", detName.Data()));
2024 AliTracker *tracker = reconstructor->CreateTracker();
2026 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
2031 fLoader[iDet]->LoadRecPoints("read");
2033 tracker->LoadClusters(fLoader[iDet]->TreeR());
2035 Int_t rv = tracker->Clusters2Tracks(esd);
2039 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
2043 fLoader[iDet]->UnloadRecPoints();
2045 tracker->UnloadClusters();
2053 //_____________________________________________________________________________
2054 Bool_t AliReconstruction::RunTracking(AliESDEvent*& esd)
2056 // run the barrel tracking
2057 static Int_t eventNr=0;
2058 AliCodeTimerAuto("")
2060 AliInfo("running tracking");
2062 //Fill the ESD with the T0 info (will be used by the TOF)
2063 if (fReconstructor[11] && fLoader[11]) {
2064 fLoader[11]->LoadRecPoints("READ");
2065 TTree *treeR = fLoader[11]->TreeR();
2066 GetReconstructor(11)->FillESD((TTree *)NULL,treeR,esd);
2069 // pass 1: TPC + ITS inwards
2070 for (Int_t iDet = 1; iDet >= 0; iDet--) {
2071 if (!fTracker[iDet]) continue;
2072 AliDebug(1, Form("%s tracking", fgkDetectorName[iDet]));
2075 fLoader[iDet]->LoadRecPoints("read");
2076 AliSysInfo::AddStamp(Form("RLoadCluster%s_%d",fgkDetectorName[iDet],eventNr),iDet,1, eventNr);
2077 TTree* tree = fLoader[iDet]->TreeR();
2079 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
2082 fTracker[iDet]->LoadClusters(tree);
2083 AliSysInfo::AddStamp(Form("TLoadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,2, eventNr);
2085 if (fTracker[iDet]->Clusters2Tracks(esd) != 0) {
2086 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
2089 // preliminary PID in TPC needed by the ITS tracker
2091 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
2092 AliESDpid::MakePID(esd);
2094 AliSysInfo::AddStamp(Form("Tracking0%s_%d",fgkDetectorName[iDet],eventNr), iDet,3,eventNr);
2097 // pass 2: ALL backwards
2099 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2100 if (!fTracker[iDet]) continue;
2101 AliDebug(1, Form("%s back propagation", fgkDetectorName[iDet]));
2104 if (iDet > 1) { // all except ITS, TPC
2106 fLoader[iDet]->LoadRecPoints("read");
2107 AliSysInfo::AddStamp(Form("RLoadCluster0%s_%d",fgkDetectorName[iDet],eventNr), iDet,1, eventNr);
2108 tree = fLoader[iDet]->TreeR();
2110 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
2113 fTracker[iDet]->LoadClusters(tree);
2114 AliSysInfo::AddStamp(Form("TLoadCluster0%s_%d",fgkDetectorName[iDet],eventNr), iDet,2, eventNr);
2118 if (iDet>1) // start filling residuals for the "outer" detectors
2119 if (fRunGlobalQA) AliTracker::SetFillResiduals(kTRUE);
2121 if (fTracker[iDet]->PropagateBack(esd) != 0) {
2122 AliError(Form("%s backward propagation failed", fgkDetectorName[iDet]));
2127 if (iDet > 3) { // all except ITS, TPC, TRD and TOF
2128 fTracker[iDet]->UnloadClusters();
2129 fLoader[iDet]->UnloadRecPoints();
2131 // updated PID in TPC needed by the ITS tracker -MI
2133 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
2134 AliESDpid::MakePID(esd);
2136 AliSysInfo::AddStamp(Form("Tracking1%s_%d",fgkDetectorName[iDet],eventNr), iDet,3, eventNr);
2138 //stop filling residuals for the "outer" detectors
2139 if (fRunGlobalQA) AliTracker::SetFillResiduals(kFALSE);
2141 // pass 3: TRD + TPC + ITS refit inwards
2143 for (Int_t iDet = 2; iDet >= 0; iDet--) {
2144 if (!fTracker[iDet]) continue;
2145 AliDebug(1, Form("%s inward refit", fgkDetectorName[iDet]));
2148 if (iDet<2) // start filling residuals for TPC and ITS
2149 if (fRunGlobalQA) AliTracker::SetFillResiduals(kTRUE);
2151 if (fTracker[iDet]->RefitInward(esd) != 0) {
2152 AliError(Form("%s inward refit failed", fgkDetectorName[iDet]));
2155 // run postprocessing
2156 if (fTracker[iDet]->PostProcess(esd) != 0) {
2157 AliError(Form("%s postprocessing failed", fgkDetectorName[iDet]));
2160 AliSysInfo::AddStamp(Form("Tracking2%s_%d",fgkDetectorName[iDet],eventNr), iDet,3, eventNr);
2163 // write space-points to the ESD in case alignment data output
2165 if (fWriteAlignmentData)
2166 WriteAlignmentData(esd);
2168 for (Int_t iDet = 3; iDet >= 0; iDet--) {
2169 if (!fTracker[iDet]) continue;
2171 fTracker[iDet]->UnloadClusters();
2172 AliSysInfo::AddStamp(Form("TUnloadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,4, eventNr);
2173 fLoader[iDet]->UnloadRecPoints();
2174 AliSysInfo::AddStamp(Form("RUnloadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,5, eventNr);
2176 // stop filling residuals for TPC and ITS
2177 if (fRunGlobalQA) AliTracker::SetFillResiduals(kFALSE);
2183 //_____________________________________________________________________________
2184 Bool_t AliReconstruction::CleanESD(AliESDEvent *esd){
2186 // Remove the data which are not needed for the physics analysis.
2189 Int_t nTracks=esd->GetNumberOfTracks();
2190 Int_t nV0s=esd->GetNumberOfV0s();
2192 (Form("Number of ESD tracks and V0s before cleaning: %d %d",nTracks,nV0s));
2194 Float_t cleanPars[]={fV0DCAmax,fV0CsPmin,fDmax,fZmax};
2195 Bool_t rc=esd->Clean(cleanPars);
2197 nTracks=esd->GetNumberOfTracks();
2198 nV0s=esd->GetNumberOfV0s();
2200 (Form("Number of ESD tracks and V0s after cleaning %d %d",nTracks,nV0s));
2205 //_____________________________________________________________________________
2206 Bool_t AliReconstruction::FillESD(AliESDEvent*& esd, const TString& detectors)
2208 // fill the event summary data
2210 AliCodeTimerAuto("")
2211 static Int_t eventNr=0;
2212 TString detStr = detectors;
2214 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2215 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
2216 AliReconstructor* reconstructor = GetReconstructor(iDet);
2217 if (!reconstructor) continue;
2218 AliDebug(1, Form("filling ESD for %s", fgkDetectorName[iDet]));
2219 TTree* clustersTree = NULL;
2220 if (fLoader[iDet]) {
2221 fLoader[iDet]->LoadRecPoints("read");
2222 clustersTree = fLoader[iDet]->TreeR();
2223 if (!clustersTree) {
2224 AliError(Form("Can't get the %s clusters tree",
2225 fgkDetectorName[iDet]));
2226 if (fStopOnError) return kFALSE;
2229 if (fRawReader && !reconstructor->HasDigitConversion()) {
2230 reconstructor->FillESD(fRawReader, clustersTree, esd);
2232 TTree* digitsTree = NULL;
2233 if (fLoader[iDet]) {
2234 fLoader[iDet]->LoadDigits("read");
2235 digitsTree = fLoader[iDet]->TreeD();
2237 AliError(Form("Can't get the %s digits tree",
2238 fgkDetectorName[iDet]));
2239 if (fStopOnError) return kFALSE;
2242 reconstructor->FillESD(digitsTree, clustersTree, esd);
2243 if (fLoader[iDet]) fLoader[iDet]->UnloadDigits();
2245 if (fLoader[iDet]) {
2246 fLoader[iDet]->UnloadRecPoints();
2250 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
2251 AliError(Form("the following detectors were not found: %s",
2253 if (fStopOnError) return kFALSE;
2255 AliSysInfo::AddStamp(Form("FillESD%d",eventNr), 0,1, eventNr);
2260 //_____________________________________________________________________________
2261 Bool_t AliReconstruction::FillTriggerESD(AliESDEvent*& esd)
2263 // Reads the trigger decision which is
2264 // stored in Trigger.root file and fills
2265 // the corresponding esd entries
2267 AliCodeTimerAuto("")
2269 AliInfo("Filling trigger information into the ESD");
2272 AliCTPRawStream input(fRawReader);
2273 if (!input.Next()) {
2274 AliWarning("No valid CTP (trigger) DDL raw data is found ! The trigger info is taken from the event header!");
2277 if (esd->GetTriggerMask() != input.GetClassMask())
2278 AliError(Form("Invalid trigger pattern found in CTP raw-data: %llx %llx",
2279 input.GetClassMask(),esd->GetTriggerMask()));
2280 if (esd->GetOrbitNumber() != input.GetOrbitID())
2281 AliError(Form("Invalid orbit id found in CTP raw-data: %x %x",
2282 input.GetOrbitID(),esd->GetOrbitNumber()));
2283 if (esd->GetBunchCrossNumber() != input.GetBCID())
2284 AliError(Form("Invalid bunch-crossing id found in CTP raw-data: %x %x",
2285 input.GetBCID(),esd->GetBunchCrossNumber()));
2288 // Here one has to add the filling of trigger inputs and
2289 // interaction records
2299 //_____________________________________________________________________________
2300 Bool_t AliReconstruction::FillRawEventHeaderESD(AliESDEvent*& esd)
2303 // Filling information from RawReader Header
2306 if (!fRawReader) return kFALSE;
2308 AliInfo("Filling information from RawReader Header");
2310 esd->SetBunchCrossNumber(fRawReader->GetBCID());
2311 esd->SetOrbitNumber(fRawReader->GetOrbitID());
2312 esd->SetPeriodNumber(fRawReader->GetPeriod());
2314 esd->SetTimeStamp(fRawReader->GetTimestamp());
2315 esd->SetEventType(fRawReader->GetType());
2321 //_____________________________________________________________________________
2322 Bool_t AliReconstruction::IsSelected(TString detName, TString& detectors) const
2324 // check whether detName is contained in detectors
2325 // if yes, it is removed from detectors
2327 // check if all detectors are selected
2328 if ((detectors.CompareTo("ALL") == 0) ||
2329 detectors.BeginsWith("ALL ") ||
2330 detectors.EndsWith(" ALL") ||
2331 detectors.Contains(" ALL ")) {
2336 // search for the given detector
2337 Bool_t result = kFALSE;
2338 if ((detectors.CompareTo(detName) == 0) ||
2339 detectors.BeginsWith(detName+" ") ||
2340 detectors.EndsWith(" "+detName) ||
2341 detectors.Contains(" "+detName+" ")) {
2342 detectors.ReplaceAll(detName, "");
2346 // clean up the detectors string
2347 while (detectors.Contains(" ")) detectors.ReplaceAll(" ", " ");
2348 while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
2349 while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
2354 //_____________________________________________________________________________
2355 Bool_t AliReconstruction::InitRunLoader()
2357 // get or create the run loader
2359 if (gAlice) delete gAlice;
2362 if (!gSystem->AccessPathName(fGAliceFileName.Data())) { // galice.root exists
2363 // load all base libraries to get the loader classes
2364 TString libs = gSystem->GetLibraries();
2365 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2366 TString detName = fgkDetectorName[iDet];
2367 if (detName == "HLT") continue;
2368 if (libs.Contains("lib" + detName + "base.so")) continue;
2369 gSystem->Load("lib" + detName + "base.so");
2371 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data());
2373 AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
2378 fRunLoader->CdGAFile();
2379 fRunLoader->LoadgAlice();
2381 //PH This is a temporary fix to give access to the kinematics
2382 //PH that is needed for the labels of ITS clusters
2383 fRunLoader->LoadHeader();
2384 fRunLoader->LoadKinematics();
2386 } else { // galice.root does not exist
2388 AliError(Form("the file %s does not exist", fGAliceFileName.Data()));
2390 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data(),
2391 AliConfig::GetDefaultEventFolderName(),
2394 AliError(Form("could not create run loader in file %s",
2395 fGAliceFileName.Data()));
2399 fIsNewRunLoader = kTRUE;
2400 fRunLoader->MakeTree("E");
2402 if (fNumberOfEventsPerFile > 0)
2403 fRunLoader->SetNumberOfEventsPerFile(fNumberOfEventsPerFile);
2405 fRunLoader->SetNumberOfEventsPerFile((UInt_t)-1);
2411 //_____________________________________________________________________________
2412 AliReconstructor* AliReconstruction::GetReconstructor(Int_t iDet)
2414 // get the reconstructor object and the loader for a detector
2416 if (fReconstructor[iDet]) {
2417 if (fRecoParam.GetDetRecoParamArray(iDet) && !AliReconstructor::GetRecoParam(iDet)) {
2418 const AliDetectorRecoParam *par = fRecoParam.GetDetRecoParam(iDet);
2419 fReconstructor[iDet]->SetRecoParam(par);
2421 return fReconstructor[iDet];
2424 // load the reconstructor object
2425 TPluginManager* pluginManager = gROOT->GetPluginManager();
2426 TString detName = fgkDetectorName[iDet];
2427 TString recName = "Ali" + detName + "Reconstructor";
2429 if (!fIsNewRunLoader && !fRunLoader->GetLoader(detName+"Loader") && (detName != "HLT")) return NULL;
2431 AliReconstructor* reconstructor = NULL;
2432 // first check if a plugin is defined for the reconstructor
2433 TPluginHandler* pluginHandler =
2434 pluginManager->FindHandler("AliReconstructor", detName);
2435 // if not, add a plugin for it
2436 if (!pluginHandler) {
2437 AliDebug(1, Form("defining plugin for %s", recName.Data()));
2438 TString libs = gSystem->GetLibraries();
2439 if (libs.Contains("lib" + detName + "base.so") ||
2440 (gSystem->Load("lib" + detName + "base.so") >= 0)) {
2441 pluginManager->AddHandler("AliReconstructor", detName,
2442 recName, detName + "rec", recName + "()");
2444 pluginManager->AddHandler("AliReconstructor", detName,
2445 recName, detName, recName + "()");
2447 pluginHandler = pluginManager->FindHandler("AliReconstructor", detName);
2449 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
2450 reconstructor = (AliReconstructor*) pluginHandler->ExecPlugin(0);
2452 if (reconstructor) {
2453 TObject* obj = fOptions.FindObject(detName.Data());
2454 if (obj) reconstructor->SetOption(obj->GetTitle());
2455 reconstructor->Init();
2456 fReconstructor[iDet] = reconstructor;
2459 // get or create the loader
2460 if (detName != "HLT") {
2461 fLoader[iDet] = fRunLoader->GetLoader(detName + "Loader");
2462 if (!fLoader[iDet]) {
2463 AliConfig::Instance()
2464 ->CreateDetectorFolders(fRunLoader->GetEventFolder(),
2466 // first check if a plugin is defined for the loader
2468 pluginManager->FindHandler("AliLoader", detName);
2469 // if not, add a plugin for it
2470 if (!pluginHandler) {
2471 TString loaderName = "Ali" + detName + "Loader";
2472 AliDebug(1, Form("defining plugin for %s", loaderName.Data()));
2473 pluginManager->AddHandler("AliLoader", detName,
2474 loaderName, detName + "base",
2475 loaderName + "(const char*, TFolder*)");
2476 pluginHandler = pluginManager->FindHandler("AliLoader", detName);
2478 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
2480 (AliLoader*) pluginHandler->ExecPlugin(2, detName.Data(),
2481 fRunLoader->GetEventFolder());
2483 if (!fLoader[iDet]) { // use default loader
2484 fLoader[iDet] = new AliLoader(detName, fRunLoader->GetEventFolder());
2486 if (!fLoader[iDet]) {
2487 AliWarning(Form("couldn't get loader for %s", detName.Data()));
2488 if (fStopOnError) return NULL;
2490 fRunLoader->AddLoader(fLoader[iDet]);
2491 fRunLoader->CdGAFile();
2492 if (gFile && !gFile->IsWritable()) gFile->ReOpen("UPDATE");
2493 fRunLoader->Write(0, TObject::kOverwrite);
2498 if (fRecoParam.GetDetRecoParamArray(iDet) && !AliReconstructor::GetRecoParam(iDet)) {
2499 const AliDetectorRecoParam *par = fRecoParam.GetDetRecoParam(iDet);
2500 reconstructor->SetRecoParam(par);
2502 return reconstructor;
2505 //_____________________________________________________________________________
2506 Bool_t AliReconstruction::CreateVertexer()
2508 // create the vertexer
2511 AliReconstructor* itsReconstructor = GetReconstructor(0);
2512 if (itsReconstructor) {
2513 fVertexer = itsReconstructor->CreateVertexer();
2516 AliWarning("couldn't create a vertexer for ITS");
2517 if (fStopOnError) return kFALSE;
2523 //_____________________________________________________________________________
2524 Bool_t AliReconstruction::CreateTrackers(const TString& detectors)
2526 // create the trackers
2528 TString detStr = detectors;
2529 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2530 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
2531 AliReconstructor* reconstructor = GetReconstructor(iDet);
2532 if (!reconstructor) continue;
2533 TString detName = fgkDetectorName[iDet];
2534 if (detName == "HLT") {
2535 fRunHLTTracking = kTRUE;
2538 if (detName == "MUON") {
2539 fRunMuonTracking = kTRUE;
2544 fTracker[iDet] = reconstructor->CreateTracker();
2545 if (!fTracker[iDet] && (iDet < 7)) {
2546 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
2547 if (fStopOnError) return kFALSE;
2549 AliSysInfo::AddStamp(Form("LTracker%s",fgkDetectorName[iDet]), iDet,0);
2555 //_____________________________________________________________________________
2556 void AliReconstruction::CleanUp()
2558 // delete trackers and the run loader and close and delete the file
2560 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2561 delete fReconstructor[iDet];
2562 fReconstructor[iDet] = NULL;
2563 fLoader[iDet] = NULL;
2564 delete fTracker[iDet];
2565 fTracker[iDet] = NULL;
2576 if(!(AliCDBManager::Instance()->GetCacheFlag())) {
2577 delete fDiamondProfile;
2578 fDiamondProfile = NULL;
2579 delete fDiamondProfileTPC;
2580 fDiamondProfileTPC = NULL;
2589 delete fParentRawReader;
2590 fParentRawReader=NULL;
2599 void AliReconstruction::WriteAlignmentData(AliESDEvent* esd)
2601 // Write space-points which are then used in the alignment procedures
2602 // For the moment only ITS, TPC, TRD and TOF
2604 Int_t ntracks = esd->GetNumberOfTracks();
2605 for (Int_t itrack = 0; itrack < ntracks; itrack++)
2607 AliESDtrack *track = esd->GetTrack(itrack);
2610 for (Int_t iDet = 3; iDet >= 0; iDet--) {// TOF, TRD, TPC, ITS clusters
2611 nsp += track->GetNcls(iDet);
2613 if (iDet==0) { // ITS "extra" clusters
2614 track->GetClusters(iDet,idx);
2615 for (Int_t i=6; i<12; i++) if(idx[i] >= 0) nsp++;
2620 AliTrackPointArray *sp = new AliTrackPointArray(nsp);
2621 track->SetTrackPointArray(sp);
2623 for (Int_t iDet = 3; iDet >= 0; iDet--) {
2624 AliTracker *tracker = fTracker[iDet];
2625 if (!tracker) continue;
2626 Int_t nspdet = track->GetClusters(iDet,idx);
2628 if (iDet==0) // ITS "extra" clusters
2629 for (Int_t i=6; i<12; i++) if(idx[i] >= 0) nspdet++;
2631 if (nspdet <= 0) continue;
2635 while (isp2 < nspdet) {
2636 Bool_t isvalid=kTRUE;
2638 Int_t index=idx[isp++];
2639 if (index < 0) continue;
2641 TString dets = fgkDetectorName[iDet];
2642 if ((fUseTrackingErrorsForAlignment.CompareTo(dets) == 0) ||
2643 fUseTrackingErrorsForAlignment.BeginsWith(dets+" ") ||
2644 fUseTrackingErrorsForAlignment.EndsWith(" "+dets) ||
2645 fUseTrackingErrorsForAlignment.Contains(" "+dets+" ")) {
2646 isvalid = tracker->GetTrackPointTrackingError(index,p,track);
2648 isvalid = tracker->GetTrackPoint(index,p);
2651 if (!isvalid) continue;
2652 sp->AddPoint(isptrack,&p); isptrack++;
2659 //_____________________________________________________________________________
2660 void AliReconstruction::FillRawDataErrorLog(Int_t iEvent, AliESDEvent* esd)
2662 // The method reads the raw-data error log
2663 // accumulated within the rawReader.
2664 // It extracts the raw-data errors related to
2665 // the current event and stores them into
2666 // a TClonesArray inside the esd object.
2668 if (!fRawReader) return;
2670 for(Int_t i = 0; i < fRawReader->GetNumberOfErrorLogs(); i++) {
2672 AliRawDataErrorLog *log = fRawReader->GetErrorLog(i);
2674 if (iEvent != log->GetEventNumber()) continue;
2676 esd->AddRawDataErrorLog(log);
2681 TNamed* AliReconstruction::CopyFileToTNamed(TString fPath,TString pName){
2682 // Dump a file content into a char in TNamed
2684 in.open(fPath.Data(),ios::in | ios::binary|ios::ate);
2685 Int_t kBytes = (Int_t)in.tellg();
2686 printf("Size: %d \n",kBytes);
2689 char* memblock = new char [kBytes];
2690 in.seekg (0, ios::beg);
2691 in.read (memblock, kBytes);
2693 TString fData(memblock,kBytes);
2694 fn = new TNamed(pName,fData);
2695 printf("fData Size: %d \n",fData.Sizeof());
2696 printf("pName Size: %d \n",pName.Sizeof());
2697 printf("fn Size: %d \n",fn->Sizeof());
2701 AliInfo(Form("Could not Open %s\n",fPath.Data()));
2707 void AliReconstruction::TNamedToFile(TTree* fTree, TString pName){
2708 // This is not really needed in AliReconstruction at the moment
2709 // but can serve as a template
2711 TList *fList = fTree->GetUserInfo();
2712 TNamed *fn = (TNamed*)fList->FindObject(pName.Data());
2713 printf("fn Size: %d \n",fn->Sizeof());
2715 TString fTmp(fn->GetName()); // to be 100% sure in principle pName also works
2716 const char* cdata = fn->GetTitle();
2717 printf("fTmp Size %d\n",fTmp.Sizeof());
2719 int size = fn->Sizeof()-fTmp.Sizeof()-sizeof(UChar_t)-sizeof(Int_t); // see dfinition of TString::SizeOf()...
2720 printf("calculated size %d\n",size);
2721 ofstream out(pName.Data(),ios::out | ios::binary);
2722 out.write(cdata,size);
2727 //_____________________________________________________________________________
2728 void AliReconstruction::CheckQA()
2730 // check the QA of SIM for this run and remove the detectors
2731 // with status Fatal
2733 TString newRunLocalReconstruction ;
2734 TString newRunTracking ;
2735 TString newFillESD ;
2737 for (Int_t iDet = 0; iDet < AliQA::kNDET; iDet++) {
2738 TString detName(AliQA::GetDetName(iDet)) ;
2739 AliQA * qa = AliQA::Instance(AliQA::DETECTORINDEX_t(iDet)) ;
2740 if ( qa->IsSet(AliQA::DETECTORINDEX_t(iDet), AliQA::kSIM, AliQA::kFATAL)) {
2741 AliInfo(Form("QA status for %s in Hits and/or SDIGITS and/or Digits was Fatal; No reconstruction performed", detName.Data())) ;
2743 if ( fRunLocalReconstruction.Contains(AliQA::GetDetName(iDet)) ||
2744 fRunLocalReconstruction.Contains("ALL") ) {
2745 newRunLocalReconstruction += detName ;
2746 newRunLocalReconstruction += " " ;
2748 if ( fRunTracking.Contains(AliQA::GetDetName(iDet)) ||
2749 fRunTracking.Contains("ALL") ) {
2750 newRunTracking += detName ;
2751 newRunTracking += " " ;
2753 if ( fFillESD.Contains(AliQA::GetDetName(iDet)) ||
2754 fFillESD.Contains("ALL") ) {
2755 newFillESD += detName ;
2760 fRunLocalReconstruction = newRunLocalReconstruction ;
2761 fRunTracking = newRunTracking ;
2762 fFillESD = newFillESD ;
2765 //_____________________________________________________________________________
2766 Int_t AliReconstruction::GetDetIndex(const char* detector)
2768 // return the detector index corresponding to detector
2770 for (index = 0; index < fgkNDetectors ; index++) {
2771 if ( strcmp(detector, fgkDetectorName[index]) == 0 )
2776 //_____________________________________________________________________________
2777 Bool_t AliReconstruction::FinishPlaneEff() {
2779 // Here execute all the necessary operationis, at the end of the tracking phase,
2780 // in case that evaluation of PlaneEfficiencies was required for some detector.
2781 // E.g., write into a DataBase file the PlaneEfficiency which have been evaluated.
2783 // This Preliminary version works only FOR ITS !!!!!
2784 // other detectors (TOF,TRD, etc. have to develop their specific codes)
2787 // Return: kTRUE if all operations have been done properly, kFALSE otherwise
2790 //for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2791 for (Int_t iDet = 0; iDet < 1; iDet++) { // for the time being only ITS
2792 //if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
2793 if(fTracker[iDet]) {
2794 AliPlaneEff *planeeff=fTracker[iDet]->GetPlaneEff();
2795 TString name=planeeff->GetName();
2797 TFile* pefile = TFile::Open(name, "RECREATE");
2798 ret=(Bool_t)planeeff->Write();
2800 if(planeeff->GetCreateHistos()) {
2801 TString hname=planeeff->GetName();
2802 hname+="Histo.root";
2803 ret*=planeeff->WriteHistosToFile(hname,"RECREATE");
2809 //_____________________________________________________________________________
2810 Bool_t AliReconstruction::InitPlaneEff() {
2812 // Here execute all the necessary operations, before of the tracking phase,
2813 // for the evaluation of PlaneEfficiencies, in case required for some detectors.
2814 // E.g., read from a DataBase file a first evaluation of the PlaneEfficiency
2815 // which should be updated/recalculated.
2817 // This Preliminary version will work only FOR ITS !!!!!
2818 // other detectors (TOF,TRD, etc. have to develop their specific codes)
2821 // Return: kTRUE if all operations have been done properly, kFALSE otherwise
2823 AliWarning(Form("Implementation of this method not yet done !! Method return kTRUE"));
2827 //_____________________________________________________________________________
2828 Bool_t AliReconstruction::InitAliEVE()
2830 // This method should be called only in case
2831 // AliReconstruction is run
2832 // within the alieve environment.
2833 // It will initialize AliEVE in a way
2834 // so that it can visualize event processed
2835 // by AliReconstruction.
2836 // The return flag shows whenever the
2837 // AliEVE initialization was successful or not.
2840 macroStr.Form("%s/EVE/macros/alieve_online.C",gSystem->ExpandPathName("$ALICE_ROOT"));
2841 AliInfo(Form("Loading AliEVE macro: %s",macroStr.Data()));
2842 if (gROOT->LoadMacro(macroStr.Data()) != 0) return kFALSE;
2844 gROOT->ProcessLine("if (!gAliEveEvent) {gAliEveEvent = new AliEveEventManager();gAliEveEvent->AddNewEventCommand(\"alieve_online_on_new_event()\");gEve->AddEvent(gAliEveEvent);};");
2845 gROOT->ProcessLine("alieve_online_init()");
2850 //_____________________________________________________________________________
2851 void AliReconstruction::RunAliEVE()
2853 // Runs AliEVE visualisation of
2854 // the current event.
2855 // Should be executed only after
2856 // successful initialization of AliEVE.
2858 AliInfo("Running AliEVE...");
2859 gROOT->ProcessLine(Form("gAliEveEvent->SetEvent((AliRunLoader*)%p,(AliRawReader*)%p,(AliESDEvent*)%p);",fRunLoader,fRawReader,fesd));
2863 //_____________________________________________________________________________
2864 Bool_t AliReconstruction::SetRunQA(TString detAndAction)
2866 // Allows to run QA for a selected set of detectors
2867 // and a selected set of tasks among RAWS, RECPOINTS and ESDS
2868 // all selected detectors run the same selected tasks
2870 if (!detAndAction.Contains(":")) {
2871 AliError( Form("%s is a wrong syntax, use \"DetectorList:ActionList\" \n", detAndAction.Data()) ) ;
2875 Int_t colon = detAndAction.Index(":") ;
2876 fQADetectors = detAndAction(0, colon) ;
2877 if (fQADetectors.Contains("ALL") )
2878 fQADetectors = fFillESD ;
2879 fQATasks = detAndAction(colon+1, detAndAction.Sizeof() ) ;
2880 if (fQATasks.Contains("ALL") ) {
2881 fQATasks = Form("%d %d %d", AliQA::kRAWS, AliQA::kRECPOINTS, AliQA::kESDS) ;
2883 fQATasks.ToUpper() ;
2885 if ( fQATasks.Contains("RAW") )
2886 tempo = Form("%d ", AliQA::kRAWS) ;
2887 if ( fQATasks.Contains("RECPOINT") )
2888 tempo += Form("%d ", AliQA::kRECPOINTS) ;
2889 if ( fQATasks.Contains("ESD") )
2890 tempo += Form("%d ", AliQA::kESDS) ;
2892 if (fQATasks.IsNull()) {
2893 AliInfo("No QA requested\n") ;
2898 TString tempo(fQATasks) ;
2899 tempo.ReplaceAll(Form("%d", AliQA::kRAWS), AliQA::GetTaskName(AliQA::kRAWS)) ;
2900 tempo.ReplaceAll(Form("%d", AliQA::kRECPOINTS), AliQA::GetTaskName(AliQA::kRECPOINTS)) ;
2901 tempo.ReplaceAll(Form("%d", AliQA::kESDS), AliQA::GetTaskName(AliQA::kESDS)) ;
2902 AliInfo( Form("QA will be done on \"%s\" for \"%s\"\n", fQADetectors.Data(), tempo.Data()) ) ;
2907 //_____________________________________________________________________________
2908 Bool_t AliReconstruction::InitRecoParams()
2910 // The method accesses OCDB and retrieves all
2911 // the available reco-param objects from there.
2913 Bool_t isOK = kTRUE;
2915 TString detStr = fLoadCDB;
2916 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2918 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
2920 if (fRecoParam.GetDetRecoParamArray(iDet)) {
2921 AliInfo(Form("Using custom reconstruction parameters for detector %s",fgkDetectorName[iDet]));
2925 AliInfo(Form("Loading reconstruction parameter objects for detector %s",fgkDetectorName[iDet]));
2927 AliCDBPath path(fgkDetectorName[iDet],"Calib","RecoParam");
2928 AliCDBEntry *entry=AliCDBManager::Instance()->Get(path.GetPath());
2930 AliWarning(Form("Couldn't find RecoParam entry in OCDB for detector %s",fgkDetectorName[iDet]));
2934 TObject *recoParamObj = entry->GetObject();
2935 if (dynamic_cast<TObjArray*>(recoParamObj)) {
2936 // The detector has a normal TobjArray of AliDetectorRecoParam objects
2937 // Registering them in AliRecoParam
2938 fRecoParam.AddDetRecoParamArray(iDet,dynamic_cast<TObjArray*>(recoParamObj));
2940 else if (dynamic_cast<AliDetectorRecoParam*>(recoParamObj)) {
2941 // The detector has only onse set of reco parameters
2942 // Registering it in AliRecoParam
2943 AliInfo(Form("Single set of reconstruction parameters found for detector %s",fgkDetectorName[iDet]));
2944 dynamic_cast<AliDetectorRecoParam*>(recoParamObj)->SetAsDefault();
2945 fRecoParam.AddDetRecoParam(iDet,dynamic_cast<AliDetectorRecoParam*>(recoParamObj));
2948 AliError(Form("No valid RecoParam object found in the OCDB for detector %s",fgkDetectorName[iDet]));
2952 AliCDBManager::Instance()->UnloadFromCache(path.GetPath());
2956 if (AliDebugLevel() > 0) fRecoParam.Print();
2961 //_____________________________________________________________________________
2962 Bool_t AliReconstruction::GetEventInfo()
2964 // Fill the event info object
2966 AliCodeTimerAuto("")
2968 AliCentralTrigger *aCTP = NULL;
2970 fEventInfo.SetEventType(fRawReader->GetType());
2972 ULong64_t mask = fRawReader->GetClassMask();
2973 fEventInfo.SetTriggerMask(mask);
2974 UInt_t clmask = fRawReader->GetDetectorPattern()[0];
2975 fEventInfo.SetTriggerCluster(AliDAQ::ListOfTriggeredDetectors(clmask));
2977 aCTP = new AliCentralTrigger();
2978 TString configstr("");
2979 if (!aCTP->LoadConfiguration(configstr)) { // Load CTP config from OCDB
2980 AliError("No trigger configuration found in OCDB! The trigger configuration information will not be used!");
2984 aCTP->SetClassMask(mask);
2985 aCTP->SetClusterMask(clmask);
2988 fEventInfo.SetEventType(AliRawEventHeaderBase::kPhysicsEvent);
2990 if (fRunLoader && (!fRunLoader->LoadTrigger())) {
2991 aCTP = fRunLoader->GetTrigger();
2992 fEventInfo.SetTriggerMask(aCTP->GetClassMask());
2993 fEventInfo.SetTriggerCluster(AliDAQ::ListOfTriggeredDetectors(aCTP->GetClusterMask()));
2996 AliWarning("No trigger can be loaded! The trigger information will not be used!");
3001 AliTriggerConfiguration *config = aCTP->GetConfiguration();
3003 AliError("No trigger configuration has been found! The trigger configuration information will not be used!");
3004 if (fRawReader) delete aCTP;
3008 UChar_t clustmask = 0;
3010 ULong64_t trmask = fEventInfo.GetTriggerMask();
3011 const TObjArray& classesArray = config->GetClasses();
3012 Int_t nclasses = classesArray.GetEntriesFast();
3013 for( Int_t iclass=0; iclass < nclasses; iclass++ ) {
3014 AliTriggerClass* trclass = (AliTriggerClass*)classesArray.At(iclass);
3016 Int_t trindex = (Int_t)TMath::Log2(trclass->GetMask());
3017 fesd->SetTriggerClass(trclass->GetName(),trindex);
3018 if (trmask & (1 << trindex)) {
3020 trclasses += trclass->GetName();
3022 clustmask |= trclass->GetCluster()->GetClusterMask();
3026 fEventInfo.SetTriggerClasses(trclasses);
3028 // Set the information in ESD
3029 fesd->SetTriggerMask(trmask);
3030 fesd->SetTriggerCluster(clustmask);
3032 if (!aCTP->CheckTriggeredDetectors()) {
3033 if (fRawReader) delete aCTP;
3037 if (fRawReader) delete aCTP;
3039 // We have to fill also the HLT decision here!!
3045 const char *AliReconstruction::MatchDetectorList(const char *detectorList, UInt_t detectorMask)
3047 // Match the detector list found in the rec.C or the default 'ALL'
3048 // to the list found in the GRP (stored there by the shuttle PP which
3049 // gets the information from ECS)
3050 static TString resultList;
3051 TString detList = detectorList;
3055 for(Int_t iDet = 0; iDet < (AliDAQ::kNDetectors-1); iDet++) {
3056 if ((detectorMask >> iDet) & 0x1) {
3057 TString det = AliDAQ::OfflineModuleName(iDet);
3058 if ((detList.CompareTo("ALL") == 0) ||
3059 detList.BeginsWith("ALL ") ||
3060 detList.EndsWith(" ALL") ||
3061 detList.Contains(" ALL ") ||
3062 (detList.CompareTo(det) == 0) ||
3063 detList.BeginsWith(det) ||
3064 detList.EndsWith(det) ||
3065 detList.Contains( " "+det+" " )) {
3066 if (!resultList.EndsWith(det + " ")) {
3075 if ((detectorMask >> AliDAQ::kHLTId) & 0x1) {
3076 TString hltDet = AliDAQ::OfflineModuleName(AliDAQ::kNDetectors-1);
3077 if ((detList.CompareTo("ALL") == 0) ||
3078 detList.BeginsWith("ALL ") ||
3079 detList.EndsWith(" ALL") ||
3080 detList.Contains(" ALL ") ||
3081 (detList.CompareTo(hltDet) == 0) ||
3082 detList.BeginsWith(hltDet) ||
3083 detList.EndsWith(hltDet) ||
3084 detList.Contains( " "+hltDet+" " )) {
3085 resultList += hltDet;
3089 return resultList.Data();
3093 //______________________________________________________________________________
3094 void AliReconstruction::Abort(const char *method, EAbort what)
3096 // Abort processing. If what = kAbortProcess, the Process() loop will be
3097 // aborted. If what = kAbortFile, the current file in a chain will be
3098 // aborted and the processing will continue with the next file, if there
3099 // is no next file then Process() will be aborted. Abort() can also be
3100 // called from Begin(), SlaveBegin(), Init() and Notify(). After abort
3101 // the SlaveTerminate() and Terminate() are always called. The abort flag
3102 // can be checked in these methods using GetAbort().
3104 // The method is overwritten in AliReconstruction for better handling of
3105 // reco specific errors
3107 if (!fStopOnError) return;
3111 TString whyMess = method;
3112 whyMess += " failed! Aborting...";
3114 AliError(whyMess.Data());
3117 TString mess = "Abort";
3118 if (fAbort == kAbortProcess)
3119 mess = "AbortProcess";
3120 else if (fAbort == kAbortFile)
3123 Info(mess, whyMess.Data());