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 // For debug purposes the method SetCheckPointLevel can be used. If the //
105 // argument is greater than 0, files with ESD events will be written after //
106 // selected steps of the reconstruction for each event: //
107 // level 1: after tracking and after filling of ESD (final) //
108 // level 2: in addition after each tracking step //
109 // level 3: in addition after the filling of ESD for each detector //
110 // If a final check point file exists for an event, this event will be //
111 // skipped in the reconstruction. The tracking and the filling of ESD for //
112 // a detector will be skipped as well, if the corresponding check point //
113 // file exists. The ESD event will then be loaded from the file instead. //
115 ///////////////////////////////////////////////////////////////////////////////
121 #include <TPluginManager.h>
122 #include <TGeoManager.h>
123 #include <TLorentzVector.h>
125 #include "AliReconstruction.h"
126 #include "AliCodeTimer.h"
127 #include "AliReconstructor.h"
129 #include "AliRunLoader.h"
131 #include "AliRawReaderFile.h"
132 #include "AliRawReaderDate.h"
133 #include "AliRawReaderRoot.h"
134 #include "AliRawEventHeaderBase.h"
135 #include "AliESDEvent.h"
136 #include "AliESDfriend.h"
137 #include "AliESDVertex.h"
138 #include "AliMultiplicity.h"
139 #include "AliTracker.h"
140 #include "AliVertexer.h"
141 #include "AliVertexerTracks.h"
142 #include "AliV0vertexer.h"
143 #include "AliCascadeVertexer.h"
144 #include "AliHeader.h"
145 #include "AliGenEventHeader.h"
147 #include "AliESDpid.h"
148 #include "AliESDtrack.h"
150 #include "AliTagCreator.h"
152 #include "AliGeomManager.h"
153 #include "AliTrackPointArray.h"
154 #include "AliCDBManager.h"
155 #include "AliCDBEntry.h"
156 #include "AliAlignObj.h"
158 #include "AliCentralTrigger.h"
159 #include "AliCTPRawStream.h"
161 #include "AliAODEvent.h"
162 #include "AliAODHeader.h"
163 #include "AliAODTrack.h"
164 #include "AliAODVertex.h"
165 #include "AliAODCluster.h"
168 ClassImp(AliReconstruction)
171 //_____________________________________________________________________________
172 const char* AliReconstruction::fgkDetectorName[AliReconstruction::fgkNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "HMPID", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "T0", "VZERO", "ACORDE", "HLT"};
174 //_____________________________________________________________________________
175 AliReconstruction::AliReconstruction(const char* gAliceFilename, const char* cdbUri,
176 const char* name, const char* title) :
179 fUniformField(kTRUE),
180 fRunVertexFinder(kTRUE),
181 fRunHLTTracking(kFALSE),
182 fRunMuonTracking(kFALSE),
183 fStopOnError(kFALSE),
184 fWriteAlignmentData(kFALSE),
185 fWriteESDfriend(kFALSE),
187 fFillTriggerESD(kTRUE),
189 fRunLocalReconstruction("ALL"),
192 fGAliceFileName(gAliceFilename),
197 fNumberOfEventsPerFile(1),
200 fLoadAlignFromCDB(kTRUE),
201 fLoadAlignData("ALL"),
208 fDiamondProfile(NULL),
210 fAlignObjArray(NULL),
214 // create reconstruction object with default parameters
216 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
217 fReconstructor[iDet] = NULL;
218 fLoader[iDet] = NULL;
219 fTracker[iDet] = NULL;
224 //_____________________________________________________________________________
225 AliReconstruction::AliReconstruction(const AliReconstruction& rec) :
228 fUniformField(rec.fUniformField),
229 fRunVertexFinder(rec.fRunVertexFinder),
230 fRunHLTTracking(rec.fRunHLTTracking),
231 fRunMuonTracking(rec.fRunMuonTracking),
232 fStopOnError(rec.fStopOnError),
233 fWriteAlignmentData(rec.fWriteAlignmentData),
234 fWriteESDfriend(rec.fWriteESDfriend),
235 fWriteAOD(rec.fWriteAOD),
236 fFillTriggerESD(rec.fFillTriggerESD),
238 fRunLocalReconstruction(rec.fRunLocalReconstruction),
239 fRunTracking(rec.fRunTracking),
240 fFillESD(rec.fFillESD),
241 fGAliceFileName(rec.fGAliceFileName),
243 fEquipIdMap(rec.fEquipIdMap),
244 fFirstEvent(rec.fFirstEvent),
245 fLastEvent(rec.fLastEvent),
246 fNumberOfEventsPerFile(rec.fNumberOfEventsPerFile),
249 fLoadAlignFromCDB(rec.fLoadAlignFromCDB),
250 fLoadAlignData(rec.fLoadAlignData),
251 fESDPar(rec.fESDPar),
257 fDiamondProfile(NULL),
259 fAlignObjArray(rec.fAlignObjArray),
260 fCDBUri(rec.fCDBUri),
265 for (Int_t i = 0; i < rec.fOptions.GetEntriesFast(); i++) {
266 if (rec.fOptions[i]) fOptions.Add(rec.fOptions[i]->Clone());
268 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
269 fReconstructor[iDet] = NULL;
270 fLoader[iDet] = NULL;
271 fTracker[iDet] = NULL;
273 for (Int_t i = 0; i < rec.fSpecCDBUri.GetEntriesFast(); i++) {
274 if (rec.fSpecCDBUri[i]) fSpecCDBUri.Add(rec.fSpecCDBUri[i]->Clone());
278 //_____________________________________________________________________________
279 AliReconstruction& AliReconstruction::operator = (const AliReconstruction& rec)
281 // assignment operator
283 this->~AliReconstruction();
284 new(this) AliReconstruction(rec);
288 //_____________________________________________________________________________
289 AliReconstruction::~AliReconstruction()
295 fSpecCDBUri.Delete();
297 AliCodeTimer::Instance()->Print();
300 //_____________________________________________________________________________
301 void AliReconstruction::InitCDBStorage()
303 // activate a default CDB storage
304 // First check if we have any CDB storage set, because it is used
305 // to retrieve the calibration and alignment constants
307 AliCDBManager* man = AliCDBManager::Instance();
308 if (man->IsDefaultStorageSet())
310 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
311 AliWarning("Default CDB storage has been already set !");
312 AliWarning(Form("Ignoring the default storage declared in AliReconstruction: %s",fCDBUri.Data()));
313 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
317 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
318 AliDebug(2, Form("Default CDB storage is set to: %s",fCDBUri.Data()));
319 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
320 man->SetDefaultStorage(fCDBUri);
323 // Now activate the detector specific CDB storage locations
324 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
325 TObject* obj = fSpecCDBUri[i];
327 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
328 AliDebug(2, Form("Specific CDB storage for %s is set to: %s",obj->GetName(),obj->GetTitle()));
329 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
330 man->SetSpecificStorage(obj->GetName(), obj->GetTitle());
335 //_____________________________________________________________________________
336 void AliReconstruction::SetDefaultStorage(const char* uri) {
337 // Store the desired default CDB storage location
338 // Activate it later within the Run() method
344 //_____________________________________________________________________________
345 void AliReconstruction::SetSpecificStorage(const char* calibType, const char* uri) {
346 // Store a detector-specific CDB storage location
347 // Activate it later within the Run() method
349 AliCDBPath aPath(calibType);
350 if(!aPath.IsValid()){
351 // if calibType is not wildcard but it is a valid detector, add "/*" to make it a valid path
352 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
353 if(!strcmp(calibType, fgkDetectorName[iDet])) {
354 aPath.SetPath(Form("%s/*", calibType));
355 AliInfo(Form("Path for specific storage set to %s", aPath.GetPath().Data()));
359 if(!aPath.IsValid()){
360 AliError(Form("Not a valid path or detector: %s", calibType));
365 // check that calibType refers to a "valid" detector name
366 Bool_t isDetector = kFALSE;
367 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
368 TString detName = fgkDetectorName[iDet];
369 if(aPath.GetLevel0() == detName) {
376 AliError(Form("Not a valid detector: %s", aPath.GetLevel0().Data()));
380 TObject* obj = fSpecCDBUri.FindObject(aPath.GetPath().Data());
381 if (obj) fSpecCDBUri.Remove(obj);
382 fSpecCDBUri.Add(new TNamed(aPath.GetPath().Data(), uri));
389 //_____________________________________________________________________________
390 Bool_t AliReconstruction::SetRunNumber()
392 // The method is called in Run() in order
393 // to set a correct run number.
394 // In case of raw data reconstruction the
395 // run number is taken from the raw data header
397 if(AliCDBManager::Instance()->GetRun() < 0) {
399 AliError("No run loader is found !");
402 // read run number from gAlice
403 if(fRunLoader->GetAliRun())
404 AliCDBManager::Instance()->SetRun(fRunLoader->GetAliRun()->GetRunNumber());
407 if(fRawReader->NextEvent()) {
408 AliCDBManager::Instance()->SetRun(fRawReader->GetRunNumber());
409 fRawReader->RewindEvents();
412 AliError("No raw-data events found !");
417 AliError("Neither gAlice nor RawReader objects are found !");
421 AliInfo(Form("CDB Run number: %d",AliCDBManager::Instance()->GetRun()));
426 //_____________________________________________________________________________
427 Bool_t AliReconstruction::MisalignGeometry(const TString& detectors)
429 // Read the alignment objects from CDB.
430 // Each detector is supposed to have the
431 // alignment objects in DET/Align/Data CDB path.
432 // All the detector objects are then collected,
433 // sorted by geometry level (starting from ALIC) and
434 // then applied to the TGeo geometry.
435 // Finally an overlaps check is performed.
437 // Load alignment data from CDB and fill fAlignObjArray
438 if(fLoadAlignFromCDB){
440 TString detStr = detectors;
441 TString loadAlObjsListOfDets = "";
443 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
444 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
445 loadAlObjsListOfDets += fgkDetectorName[iDet];
446 loadAlObjsListOfDets += " ";
447 } // end loop over detectors
448 AliGeomManager::ApplyAlignObjsFromCDB(loadAlObjsListOfDets.Data());
450 // Check if the array with alignment objects was
451 // provided by the user. If yes, apply the objects
452 // to the present TGeo geometry
453 if (fAlignObjArray) {
454 if (gGeoManager && gGeoManager->IsClosed()) {
455 if (AliGeomManager::ApplyAlignObjsToGeom(*fAlignObjArray) == kFALSE) {
456 AliError("The misalignment of one or more volumes failed!"
457 "Compare the list of simulated detectors and the list of detector alignment data!");
462 AliError("Can't apply the misalignment! gGeoManager doesn't exist or it is still opened!");
468 delete fAlignObjArray; fAlignObjArray=0;
473 //_____________________________________________________________________________
474 void AliReconstruction::SetGAliceFile(const char* fileName)
476 // set the name of the galice file
478 fGAliceFileName = fileName;
481 //_____________________________________________________________________________
482 void AliReconstruction::SetOption(const char* detector, const char* option)
484 // set options for the reconstruction of a detector
486 TObject* obj = fOptions.FindObject(detector);
487 if (obj) fOptions.Remove(obj);
488 fOptions.Add(new TNamed(detector, option));
492 //_____________________________________________________________________________
493 Bool_t AliReconstruction::Run(const char* input)
495 // run the reconstruction
500 if (!input) input = fInput.Data();
501 TString fileName(input);
502 if (fileName.EndsWith("/")) {
503 fRawReader = new AliRawReaderFile(fileName);
504 } else if (fileName.EndsWith(".root")) {
505 fRawReader = new AliRawReaderRoot(fileName);
506 } else if (!fileName.IsNull()) {
507 fRawReader = new AliRawReaderDate(fileName);
508 fRawReader->SelectEvents(7);
510 if (!fEquipIdMap.IsNull() && fRawReader)
511 fRawReader->LoadEquipmentIdsMap(fEquipIdMap);
514 // get the run loader
515 if (!InitRunLoader()) return kFALSE;
517 // Initialize the CDB storage
520 // Set run number in CDBManager (if it is not already set by the user)
521 if (!SetRunNumber()) if (fStopOnError) return kFALSE;
523 // Import ideal TGeo geometry and apply misalignment
525 TString geom(gSystem->DirName(fGAliceFileName));
526 geom += "/geometry.root";
527 AliGeomManager::LoadGeometry(geom.Data());
528 if (!gGeoManager) if (fStopOnError) return kFALSE;
531 AliCDBManager* man = AliCDBManager::Instance();
532 if (!MisalignGeometry(fLoadAlignData)) if (fStopOnError) return kFALSE;
534 // local reconstruction
535 if (!fRunLocalReconstruction.IsNull()) {
536 if (!RunLocalReconstruction(fRunLocalReconstruction)) {
537 if (fStopOnError) {CleanUp(); return kFALSE;}
540 // if (!fRunVertexFinder && fRunTracking.IsNull() &&
541 // fFillESD.IsNull()) return kTRUE;
544 if (fRunVertexFinder && !CreateVertexer()) {
552 if (!fRunTracking.IsNull() && !CreateTrackers(fRunTracking)) {
559 // get the possibly already existing ESD file and tree
560 AliESDEvent* esd = new AliESDEvent(); AliESDEvent* hltesd = new AliESDEvent();
561 TFile* fileOld = NULL;
562 TTree* treeOld = NULL; TTree *hlttreeOld = NULL;
563 if (!gSystem->AccessPathName("AliESDs.root")){
564 gSystem->CopyFile("AliESDs.root", "AliESDs.old.root", kTRUE);
565 fileOld = TFile::Open("AliESDs.old.root");
566 if (fileOld && fileOld->IsOpen()) {
567 treeOld = (TTree*) fileOld->Get("esdTree");
568 if (treeOld)esd->ReadFromTree(treeOld);
569 hlttreeOld = (TTree*) fileOld->Get("HLTesdTree");
570 if (hlttreeOld) hltesd->ReadFromTree(hlttreeOld);
574 // create the ESD output file and tree
575 TFile* file = TFile::Open("AliESDs.root", "RECREATE");
576 file->SetCompressionLevel(2);
577 if (!file->IsOpen()) {
578 AliError("opening AliESDs.root failed");
579 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
582 TTree* tree = new TTree("esdTree", "Tree with ESD objects");
583 esd = new AliESDEvent();
584 esd->CreateStdContent();
585 esd->WriteToTree(tree);
587 TTree* hlttree = new TTree("HLTesdTree", "Tree with HLT ESD objects");
588 hltesd = new AliESDEvent();
589 hltesd->CreateStdContent();
590 hltesd->WriteToTree(hlttree);
593 delete esd; delete hltesd;
594 esd = NULL; hltesd = NULL;
596 // create the branch with ESD additions
600 AliESDfriend *esdf = 0;
601 if (fWriteESDfriend) {
602 esdf = new AliESDfriend();
603 TBranch *br=tree->Branch("ESDfriend.","AliESDfriend", &esdf);
604 br->SetFile("AliESDfriends.root");
605 esd->AddObject(esdf);
609 // Get the diamond profile from OCDB
610 AliCDBEntry* entry = AliCDBManager::Instance()
611 ->Get("GRP/Calib/MeanVertex");
614 fDiamondProfile = dynamic_cast<AliESDVertex*> (entry->GetObject());
616 AliError("No diamond profile found in OCDB!");
619 AliVertexerTracks tVertexer(AliTracker::GetBz());
620 if(fDiamondProfile) tVertexer.SetVtxStart(fDiamondProfile);
623 if (fRawReader) fRawReader->RewindEvents();
625 for (Int_t iEvent = 0; iEvent < fRunLoader->GetNumberOfEvents(); iEvent++) {
626 if (fRawReader) fRawReader->NextEvent();
627 if ((iEvent < fFirstEvent) || ((fLastEvent >= 0) && (iEvent > fLastEvent))) {
628 // copy old ESD to the new one
630 esd->ReadFromTree(treeOld);
631 treeOld->GetEntry(iEvent);
635 esd->ReadFromTree(hlttreeOld);
636 hlttreeOld->GetEntry(iEvent);
642 AliInfo(Form("processing event %d", iEvent));
643 fRunLoader->GetEvent(iEvent);
646 sprintf(aFileName, "ESD_%d.%d_final.root",
647 fRunLoader->GetHeader()->GetRun(),
648 fRunLoader->GetHeader()->GetEventNrInRun());
649 if (!gSystem->AccessPathName(aFileName)) continue;
651 // local reconstruction
652 if (!fRunLocalReconstruction.IsNull()) {
653 if (!RunLocalEventReconstruction(fRunLocalReconstruction)) {
654 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
659 esd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
660 hltesd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
661 esd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
662 hltesd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
664 // Set magnetic field from the tracker
665 esd->SetMagneticField(AliTracker::GetBz());
666 hltesd->SetMagneticField(AliTracker::GetBz());
670 // Fill raw-data error log into the ESD
671 if (fRawReader) FillRawDataErrorLog(iEvent,esd);
674 if (fRunVertexFinder) {
675 if (!ReadESD(esd, "vertex")) {
676 if (!RunVertexFinder(esd)) {
677 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
679 if (fCheckPointLevel > 0) WriteESD(esd, "vertex");
684 if (!fRunTracking.IsNull()) {
685 if (fRunHLTTracking) {
686 hltesd->SetVertex(esd->GetVertex());
687 if (!RunHLTTracking(hltesd)) {
688 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
694 if (!fRunTracking.IsNull()) {
695 if (fRunMuonTracking) {
696 if (!RunMuonTracking(esd)) {
697 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
703 if (!fRunTracking.IsNull()) {
704 if (!ReadESD(esd, "tracking")) {
705 if (!RunTracking(esd)) {
706 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
708 if (fCheckPointLevel > 0) WriteESD(esd, "tracking");
713 if (!fFillESD.IsNull()) {
714 if (!FillESD(esd, fFillESD)) {
715 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
718 // fill Event header information from the RawEventHeader
719 if (fRawReader){FillRawEventHeaderESD(esd);}
722 AliESDpid::MakePID(esd);
723 if (fCheckPointLevel > 1) WriteESD(esd, "PID");
725 if (fFillTriggerESD) {
726 if (!ReadESD(esd, "trigger")) {
727 if (!FillTriggerESD(esd)) {
728 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
730 if (fCheckPointLevel > 1) WriteESD(esd, "trigger");
734 //Try to improve the reconstructed primary vertex position using the tracks
735 AliESDVertex *pvtx=0;
736 Bool_t dovertex=kTRUE;
737 TObject* obj = fOptions.FindObject("ITS");
739 TString optITS = obj->GetTitle();
740 if (optITS.Contains("cosmics") || optITS.Contains("COSMICS"))
743 if(dovertex) pvtx=tVertexer.FindPrimaryVertex(esd);
744 if(fDiamondProfile) esd->SetDiamond(fDiamondProfile);
747 if (pvtx->GetStatus()) {
748 // Store the improved primary vertex
749 esd->SetPrimaryVertex(pvtx);
750 // Propagate the tracks to the DCA to the improved primary vertex
751 Double_t somethingbig = 777.;
752 Double_t bz = esd->GetMagneticField();
753 Int_t nt=esd->GetNumberOfTracks();
755 AliESDtrack *t = esd->GetTrack(nt);
756 t->RelateToVertex(pvtx, bz, somethingbig);
763 vtxer.Tracks2V0vertices(esd);
766 AliCascadeVertexer cvtxer;
767 cvtxer.V0sTracks2CascadeVertices(esd);
771 if (fWriteESDfriend) {
772 new (esdf) AliESDfriend(); // Reset...
773 esd->GetESDfriend(esdf);
780 if (fCheckPointLevel > 0) WriteESD(esd, "final");
783 if (fWriteESDfriend) {
784 new (esdf) AliESDfriend(); // Reset...
787 // delete esdf; esdf = 0;
793 tree->GetUserInfo()->Add(esd);
794 hlttree->GetUserInfo()->Add(hltesd);
798 if(fESDPar.Contains("ESD.par")){
799 AliInfo("Attaching ESD.par to Tree");
800 TNamed *fn = CopyFileToTNamed(fESDPar.Data(),"ESD.par");
801 tree->GetUserInfo()->Add(fn);
807 tree->SetBranchStatus("ESDfriend*",0);
812 TFile *aodFile = TFile::Open("AliAOD.root", "RECREATE");
813 ESDFile2AODFile(file, aodFile);
818 CleanUp(file, fileOld);
820 // Create tags for the events in the ESD tree (the ESD tree is always present)
821 // In case of empty events the tags will contain dummy values
822 AliTagCreator *tagCreator = new AliTagCreator();
823 tagCreator->CreateESDTags(fFirstEvent,fLastEvent);
824 if (fWriteAOD) tagCreator->CreateAODTags(fFirstEvent,fLastEvent);
830 //_____________________________________________________________________________
831 Bool_t AliReconstruction::RunLocalReconstruction(const TString& detectors)
833 // run the local reconstruction
837 AliCDBManager* man = AliCDBManager::Instance();
838 Bool_t origCache = man->GetCacheFlag();
840 TString detStr = detectors;
841 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
842 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
843 AliReconstructor* reconstructor = GetReconstructor(iDet);
844 if (!reconstructor) continue;
845 if (reconstructor->HasLocalReconstruction()) continue;
847 AliCodeTimerStart(Form("running reconstruction for %s", fgkDetectorName[iDet]));
848 AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
850 AliCodeTimerStart(Form("Loading calibration data from OCDB for %s", fgkDetectorName[iDet]));
851 AliInfo(Form("Loading calibration data from OCDB for %s", fgkDetectorName[iDet]));
853 man->SetCacheFlag(kTRUE);
854 TString calibPath = Form("%s/Calib/*", fgkDetectorName[iDet]);
855 man->GetAll(calibPath); // entries are cached!
857 AliCodeTimerStop(Form("Loading calibration data from OCDB for %s", fgkDetectorName[iDet]));
860 fRawReader->RewindEvents();
861 reconstructor->Reconstruct(fRunLoader, fRawReader);
863 reconstructor->Reconstruct(fRunLoader);
866 AliCodeTimerStop(Form("running reconstruction for %s", fgkDetectorName[iDet]));
868 // unload calibration data
869 man->UnloadFromCache(calibPath);
873 man->SetCacheFlag(origCache);
875 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
876 AliError(Form("the following detectors were not found: %s",
878 if (fStopOnError) return kFALSE;
884 //_____________________________________________________________________________
885 Bool_t AliReconstruction::RunLocalEventReconstruction(const TString& detectors)
887 // run the local reconstruction
891 TString detStr = detectors;
892 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
893 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
894 AliReconstructor* reconstructor = GetReconstructor(iDet);
895 if (!reconstructor) continue;
896 AliLoader* loader = fLoader[iDet];
898 // conversion of digits
899 if (fRawReader && reconstructor->HasDigitConversion()) {
900 AliInfo(Form("converting raw data digits into root objects for %s",
901 fgkDetectorName[iDet]));
902 AliCodeTimerAuto(Form("converting raw data digits into root objects for %s",
903 fgkDetectorName[iDet]));
904 loader->LoadDigits("update");
905 loader->CleanDigits();
906 loader->MakeDigitsContainer();
907 TTree* digitsTree = loader->TreeD();
908 reconstructor->ConvertDigits(fRawReader, digitsTree);
909 loader->WriteDigits("OVERWRITE");
910 loader->UnloadDigits();
913 // local reconstruction
914 if (!reconstructor->HasLocalReconstruction()) continue;
915 AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
916 AliCodeTimerAuto(Form("running reconstruction for %s", fgkDetectorName[iDet]));
917 loader->LoadRecPoints("update");
918 loader->CleanRecPoints();
919 loader->MakeRecPointsContainer();
920 TTree* clustersTree = loader->TreeR();
921 if (fRawReader && !reconstructor->HasDigitConversion()) {
922 reconstructor->Reconstruct(fRawReader, clustersTree);
924 loader->LoadDigits("read");
925 TTree* digitsTree = loader->TreeD();
927 AliError(Form("Can't get the %s digits tree", fgkDetectorName[iDet]));
928 if (fStopOnError) return kFALSE;
930 reconstructor->Reconstruct(digitsTree, clustersTree);
932 loader->UnloadDigits();
934 loader->WriteRecPoints("OVERWRITE");
935 loader->UnloadRecPoints();
938 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
939 AliError(Form("the following detectors were not found: %s",
941 if (fStopOnError) return kFALSE;
947 //_____________________________________________________________________________
948 Bool_t AliReconstruction::RunVertexFinder(AliESDEvent*& esd)
950 // run the barrel tracking
954 AliESDVertex* vertex = NULL;
955 Double_t vtxPos[3] = {0, 0, 0};
956 Double_t vtxErr[3] = {0.07, 0.07, 0.1};
958 if (fRunLoader->GetHeader() && fRunLoader->GetHeader()->GenEventHeader()) {
959 fRunLoader->GetHeader()->GenEventHeader()->PrimaryVertex(mcVertex);
960 for (Int_t i = 0; i < 3; i++) vtxPos[i] = mcVertex[i];
964 if(fDiamondProfile) fVertexer->SetVtxStart(fDiamondProfile);
965 AliInfo("running the ITS vertex finder");
966 if (fLoader[0]) fLoader[0]->LoadRecPoints();
967 vertex = fVertexer->FindVertexForCurrentEvent(fRunLoader->GetEventNumber());
968 if (fLoader[0]) fLoader[0]->UnloadRecPoints();
970 AliWarning("Vertex not found");
971 vertex = new AliESDVertex();
972 vertex->SetName("default");
975 vertex->SetName("reconstructed");
979 AliInfo("getting the primary vertex from MC");
980 vertex = new AliESDVertex(vtxPos, vtxErr);
984 vertex->GetXYZ(vtxPos);
985 vertex->GetSigmaXYZ(vtxErr);
987 AliWarning("no vertex reconstructed");
988 vertex = new AliESDVertex(vtxPos, vtxErr);
990 esd->SetVertex(vertex);
991 // if SPD multiplicity has been determined, it is stored in the ESD
992 AliMultiplicity *mult = fVertexer->GetMultiplicity();
993 if(mult)esd->SetMultiplicity(mult);
995 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
996 if (fTracker[iDet]) fTracker[iDet]->SetVertex(vtxPos, vtxErr);
1003 //_____________________________________________________________________________
1004 Bool_t AliReconstruction::RunHLTTracking(AliESDEvent*& esd)
1006 // run the HLT barrel tracking
1008 AliCodeTimerAuto("")
1011 AliError("Missing runLoader!");
1015 AliInfo("running HLT tracking");
1017 // Get a pointer to the HLT reconstructor
1018 AliReconstructor *reconstructor = GetReconstructor(fgkNDetectors-1);
1019 if (!reconstructor) return kFALSE;
1022 for (Int_t iDet = 1; iDet >= 0; iDet--) {
1023 TString detName = fgkDetectorName[iDet];
1024 AliDebug(1, Form("%s HLT tracking", detName.Data()));
1025 reconstructor->SetOption(detName.Data());
1026 AliTracker *tracker = reconstructor->CreateTracker(fRunLoader);
1028 AliWarning(Form("couldn't create a HLT tracker for %s", detName.Data()));
1029 if (fStopOnError) return kFALSE;
1033 Double_t vtxErr[3]={0.005,0.005,0.010};
1034 const AliESDVertex *vertex = esd->GetVertex();
1035 vertex->GetXYZ(vtxPos);
1036 tracker->SetVertex(vtxPos,vtxErr);
1038 fLoader[iDet]->LoadRecPoints("read");
1039 TTree* tree = fLoader[iDet]->TreeR();
1041 AliError(Form("Can't get the %s cluster tree", detName.Data()));
1044 tracker->LoadClusters(tree);
1046 if (tracker->Clusters2Tracks(esd) != 0) {
1047 AliError(Form("HLT %s Clusters2Tracks failed", fgkDetectorName[iDet]));
1051 tracker->UnloadClusters();
1059 //_____________________________________________________________________________
1060 Bool_t AliReconstruction::RunMuonTracking(AliESDEvent*& esd)
1062 // run the muon spectrometer tracking
1064 AliCodeTimerAuto("")
1067 AliError("Missing runLoader!");
1070 Int_t iDet = 7; // for MUON
1072 AliInfo("is running...");
1074 // Get a pointer to the MUON reconstructor
1075 AliReconstructor *reconstructor = GetReconstructor(iDet);
1076 if (!reconstructor) return kFALSE;
1079 TString detName = fgkDetectorName[iDet];
1080 AliDebug(1, Form("%s tracking", detName.Data()));
1081 AliTracker *tracker = reconstructor->CreateTracker(fRunLoader);
1083 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
1088 fLoader[iDet]->LoadTracks("update");
1089 fLoader[iDet]->CleanTracks();
1090 fLoader[iDet]->MakeTracksContainer();
1093 fLoader[iDet]->LoadRecPoints("read");
1094 tracker->LoadClusters(fLoader[iDet]->TreeR());
1096 Int_t rv = tracker->Clusters2Tracks(esd);
1098 fLoader[iDet]->UnloadRecPoints();
1102 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
1106 tracker->UnloadClusters();
1108 fLoader[iDet]->UnloadRecPoints();
1110 fLoader[iDet]->WriteTracks("OVERWRITE");
1111 fLoader[iDet]->UnloadTracks();
1119 //_____________________________________________________________________________
1120 Bool_t AliReconstruction::RunTracking(AliESDEvent*& esd)
1122 // run the barrel tracking
1124 AliCodeTimerAuto("")
1126 AliInfo("running tracking");
1128 //Fill the ESD with the T0 info (will be used by the TOF)
1129 if (fReconstructor[11])
1130 GetReconstructor(11)->FillESD(fRunLoader, esd);
1132 // pass 1: TPC + ITS inwards
1133 for (Int_t iDet = 1; iDet >= 0; iDet--) {
1134 if (!fTracker[iDet]) continue;
1135 AliDebug(1, Form("%s tracking", fgkDetectorName[iDet]));
1138 fLoader[iDet]->LoadRecPoints("read");
1139 TTree* tree = fLoader[iDet]->TreeR();
1141 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
1144 fTracker[iDet]->LoadClusters(tree);
1147 if (fTracker[iDet]->Clusters2Tracks(esd) != 0) {
1148 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
1151 if (fCheckPointLevel > 1) {
1152 WriteESD(esd, Form("%s.tracking", fgkDetectorName[iDet]));
1154 // preliminary PID in TPC needed by the ITS tracker
1156 GetReconstructor(1)->FillESD(fRunLoader, esd);
1157 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
1158 AliESDpid::MakePID(esd);
1162 // pass 2: ALL backwards
1163 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1164 if (!fTracker[iDet]) continue;
1165 AliDebug(1, Form("%s back propagation", fgkDetectorName[iDet]));
1168 if (iDet > 1) { // all except ITS, TPC
1170 fLoader[iDet]->LoadRecPoints("read");
1171 tree = fLoader[iDet]->TreeR();
1173 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
1176 fTracker[iDet]->LoadClusters(tree);
1180 if (fTracker[iDet]->PropagateBack(esd) != 0) {
1181 AliError(Form("%s backward propagation failed", fgkDetectorName[iDet]));
1184 if (fCheckPointLevel > 1) {
1185 WriteESD(esd, Form("%s.back", fgkDetectorName[iDet]));
1189 if (iDet > 2) { // all except ITS, TPC, TRD
1190 fTracker[iDet]->UnloadClusters();
1191 fLoader[iDet]->UnloadRecPoints();
1193 // updated PID in TPC needed by the ITS tracker -MI
1195 GetReconstructor(1)->FillESD(fRunLoader, esd);
1196 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
1197 AliESDpid::MakePID(esd);
1201 // write space-points to the ESD in case alignment data output
1203 if (fWriteAlignmentData)
1204 WriteAlignmentData(esd);
1206 // pass 3: TRD + TPC + ITS refit inwards
1207 for (Int_t iDet = 2; iDet >= 0; iDet--) {
1208 if (!fTracker[iDet]) continue;
1209 AliDebug(1, Form("%s inward refit", fgkDetectorName[iDet]));
1212 if (fTracker[iDet]->RefitInward(esd) != 0) {
1213 AliError(Form("%s inward refit failed", fgkDetectorName[iDet]));
1216 if (fCheckPointLevel > 1) {
1217 WriteESD(esd, Form("%s.refit", fgkDetectorName[iDet]));
1221 fTracker[iDet]->UnloadClusters();
1222 fLoader[iDet]->UnloadRecPoints();
1225 // Propagate track to the vertex - if not done by ITS
1227 Int_t ntracks = esd->GetNumberOfTracks();
1228 for (Int_t itrack=0; itrack<ntracks; itrack++){
1229 const Double_t kRadius = 3; // beam pipe radius
1230 const Double_t kMaxStep = 5; // max step
1231 const Double_t kMaxD = 123456; // max distance to prim vertex
1232 Double_t fieldZ = AliTracker::GetBz(); //
1233 AliESDtrack * track = esd->GetTrack(itrack);
1234 if (!track) continue;
1235 if (track->IsOn(AliESDtrack::kITSrefit)) continue;
1236 track->PropagateTo(kRadius, fieldZ, track->GetMass(),kMaxStep,kTRUE);
1237 track->RelateToVertex(esd->GetVertex(),fieldZ, kMaxD);
1243 //_____________________________________________________________________________
1244 Bool_t AliReconstruction::FillESD(AliESDEvent*& esd, const TString& detectors)
1246 // fill the event summary data
1248 AliCodeTimerAuto("")
1250 TString detStr = detectors;
1251 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1252 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1253 AliReconstructor* reconstructor = GetReconstructor(iDet);
1254 if (!reconstructor) continue;
1256 if (!ReadESD(esd, fgkDetectorName[iDet])) {
1257 AliDebug(1, Form("filling ESD for %s", fgkDetectorName[iDet]));
1258 TTree* clustersTree = NULL;
1259 if (reconstructor->HasLocalReconstruction() && fLoader[iDet]) {
1260 fLoader[iDet]->LoadRecPoints("read");
1261 clustersTree = fLoader[iDet]->TreeR();
1262 if (!clustersTree) {
1263 AliError(Form("Can't get the %s clusters tree",
1264 fgkDetectorName[iDet]));
1265 if (fStopOnError) return kFALSE;
1268 if (fRawReader && !reconstructor->HasDigitConversion()) {
1269 reconstructor->FillESD(fRawReader, clustersTree, esd);
1271 TTree* digitsTree = NULL;
1272 if (fLoader[iDet]) {
1273 fLoader[iDet]->LoadDigits("read");
1274 digitsTree = fLoader[iDet]->TreeD();
1276 AliError(Form("Can't get the %s digits tree",
1277 fgkDetectorName[iDet]));
1278 if (fStopOnError) return kFALSE;
1281 reconstructor->FillESD(digitsTree, clustersTree, esd);
1282 if (fLoader[iDet]) fLoader[iDet]->UnloadDigits();
1284 if (reconstructor->HasLocalReconstruction() && fLoader[iDet]) {
1285 fLoader[iDet]->UnloadRecPoints();
1289 reconstructor->FillESD(fRunLoader, fRawReader, esd);
1291 reconstructor->FillESD(fRunLoader, esd);
1293 if (fCheckPointLevel > 2) WriteESD(esd, fgkDetectorName[iDet]);
1297 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1298 AliError(Form("the following detectors were not found: %s",
1300 if (fStopOnError) return kFALSE;
1306 //_____________________________________________________________________________
1307 Bool_t AliReconstruction::FillTriggerESD(AliESDEvent*& esd)
1309 // Reads the trigger decision which is
1310 // stored in Trigger.root file and fills
1311 // the corresponding esd entries
1313 AliCodeTimerAuto("")
1315 AliInfo("Filling trigger information into the ESD");
1318 AliCTPRawStream input(fRawReader);
1319 if (!input.Next()) {
1320 AliError("No valid CTP (trigger) DDL raw data is found ! The trigger information is not stored in the ESD !");
1323 esd->SetTriggerMask(input.GetClassMask());
1324 esd->SetTriggerCluster(input.GetClusterMask());
1327 AliRunLoader *runloader = AliRunLoader::GetRunLoader();
1329 if (!runloader->LoadTrigger()) {
1330 AliCentralTrigger *aCTP = runloader->GetTrigger();
1331 esd->SetTriggerMask(aCTP->GetClassMask());
1332 esd->SetTriggerCluster(aCTP->GetClusterMask());
1335 AliWarning("No trigger can be loaded! The trigger information is not stored in the ESD !");
1340 AliError("No run loader is available! The trigger information is not stored in the ESD !");
1352 //_____________________________________________________________________________
1353 Bool_t AliReconstruction::FillRawEventHeaderESD(AliESDEvent*& esd)
1356 // Filling information from RawReader Header
1359 AliInfo("Filling information from RawReader Header");
1360 esd->SetBunchCrossNumber(0);
1361 esd->SetOrbitNumber(0);
1362 esd->SetPeriodNumber(0);
1363 esd->SetTimeStamp(0);
1364 esd->SetEventType(0);
1365 const AliRawEventHeaderBase * eventHeader = fRawReader->GetEventHeader();
1368 const UInt_t *id = eventHeader->GetP("Id");
1369 esd->SetBunchCrossNumber((id)[1]&0x00000fff);
1370 esd->SetOrbitNumber((((id)[0]<<20)&0xf00000)|(((id)[1]>>12)&0xfffff));
1371 esd->SetPeriodNumber(((id)[0]>>4)&0x0fffffff);
1373 esd->SetTimeStamp((eventHeader->Get("Timestamp")));
1374 esd->SetEventType((eventHeader->Get("Type")));
1381 //_____________________________________________________________________________
1382 Bool_t AliReconstruction::IsSelected(TString detName, TString& detectors) const
1384 // check whether detName is contained in detectors
1385 // if yes, it is removed from detectors
1387 // check if all detectors are selected
1388 if ((detectors.CompareTo("ALL") == 0) ||
1389 detectors.BeginsWith("ALL ") ||
1390 detectors.EndsWith(" ALL") ||
1391 detectors.Contains(" ALL ")) {
1396 // search for the given detector
1397 Bool_t result = kFALSE;
1398 if ((detectors.CompareTo(detName) == 0) ||
1399 detectors.BeginsWith(detName+" ") ||
1400 detectors.EndsWith(" "+detName) ||
1401 detectors.Contains(" "+detName+" ")) {
1402 detectors.ReplaceAll(detName, "");
1406 // clean up the detectors string
1407 while (detectors.Contains(" ")) detectors.ReplaceAll(" ", " ");
1408 while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
1409 while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
1414 //_____________________________________________________________________________
1415 Bool_t AliReconstruction::InitRunLoader()
1417 // get or create the run loader
1419 if (gAlice) delete gAlice;
1422 if (!gSystem->AccessPathName(fGAliceFileName.Data())) { // galice.root exists
1423 // load all base libraries to get the loader classes
1424 TString libs = gSystem->GetLibraries();
1425 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1426 TString detName = fgkDetectorName[iDet];
1427 if (detName == "HLT") continue;
1428 if (libs.Contains("lib" + detName + "base.so")) continue;
1429 gSystem->Load("lib" + detName + "base.so");
1431 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data());
1433 AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
1437 fRunLoader->CdGAFile();
1438 if (gFile->GetKey(AliRunLoader::GetGAliceName())) {
1439 if (fRunLoader->LoadgAlice() == 0) {
1440 gAlice = fRunLoader->GetAliRun();
1441 AliTracker::SetFieldMap(gAlice->Field(),fUniformField);
1444 if (!gAlice && !fRawReader) {
1445 AliError(Form("no gAlice object found in file %s",
1446 fGAliceFileName.Data()));
1451 } else { // galice.root does not exist
1453 AliError(Form("the file %s does not exist", fGAliceFileName.Data()));
1457 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data(),
1458 AliConfig::GetDefaultEventFolderName(),
1461 AliError(Form("could not create run loader in file %s",
1462 fGAliceFileName.Data()));
1466 fRunLoader->MakeTree("E");
1468 while (fRawReader->NextEvent()) {
1469 fRunLoader->SetEventNumber(iEvent);
1470 fRunLoader->GetHeader()->Reset(fRawReader->GetRunNumber(),
1472 fRunLoader->MakeTree("H");
1473 fRunLoader->TreeE()->Fill();
1476 fRawReader->RewindEvents();
1477 if (fNumberOfEventsPerFile > 0)
1478 fRunLoader->SetNumberOfEventsPerFile(fNumberOfEventsPerFile);
1480 fRunLoader->SetNumberOfEventsPerFile(iEvent);
1481 fRunLoader->WriteHeader("OVERWRITE");
1482 fRunLoader->CdGAFile();
1483 fRunLoader->Write(0, TObject::kOverwrite);
1484 // AliTracker::SetFieldMap(???);
1490 //_____________________________________________________________________________
1491 AliReconstructor* AliReconstruction::GetReconstructor(Int_t iDet)
1493 // get the reconstructor object and the loader for a detector
1495 if (fReconstructor[iDet]) return fReconstructor[iDet];
1497 // load the reconstructor object
1498 TPluginManager* pluginManager = gROOT->GetPluginManager();
1499 TString detName = fgkDetectorName[iDet];
1500 TString recName = "Ali" + detName + "Reconstructor";
1501 if (gAlice && !gAlice->GetDetector(detName) && (detName != "HLT")) return NULL;
1503 if (detName == "HLT") {
1504 if (!gROOT->GetClass("AliLevel3")) {
1505 gSystem->Load("libAliHLTSrc.so");
1506 gSystem->Load("libAliHLTMisc.so");
1507 gSystem->Load("libAliHLTHough.so");
1508 gSystem->Load("libAliHLTComp.so");
1512 AliReconstructor* reconstructor = NULL;
1513 // first check if a plugin is defined for the reconstructor
1514 TPluginHandler* pluginHandler =
1515 pluginManager->FindHandler("AliReconstructor", detName);
1516 // if not, add a plugin for it
1517 if (!pluginHandler) {
1518 AliDebug(1, Form("defining plugin for %s", recName.Data()));
1519 TString libs = gSystem->GetLibraries();
1520 if (libs.Contains("lib" + detName + "base.so") ||
1521 (gSystem->Load("lib" + detName + "base.so") >= 0)) {
1522 pluginManager->AddHandler("AliReconstructor", detName,
1523 recName, detName + "rec", recName + "()");
1525 pluginManager->AddHandler("AliReconstructor", detName,
1526 recName, detName, recName + "()");
1528 pluginHandler = pluginManager->FindHandler("AliReconstructor", detName);
1530 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
1531 reconstructor = (AliReconstructor*) pluginHandler->ExecPlugin(0);
1533 if (reconstructor) {
1534 TObject* obj = fOptions.FindObject(detName.Data());
1535 if (obj) reconstructor->SetOption(obj->GetTitle());
1536 reconstructor->Init(fRunLoader);
1537 fReconstructor[iDet] = reconstructor;
1540 // get or create the loader
1541 if (detName != "HLT") {
1542 fLoader[iDet] = fRunLoader->GetLoader(detName + "Loader");
1543 if (!fLoader[iDet]) {
1544 AliConfig::Instance()
1545 ->CreateDetectorFolders(fRunLoader->GetEventFolder(),
1547 // first check if a plugin is defined for the loader
1549 pluginManager->FindHandler("AliLoader", detName);
1550 // if not, add a plugin for it
1551 if (!pluginHandler) {
1552 TString loaderName = "Ali" + detName + "Loader";
1553 AliDebug(1, Form("defining plugin for %s", loaderName.Data()));
1554 pluginManager->AddHandler("AliLoader", detName,
1555 loaderName, detName + "base",
1556 loaderName + "(const char*, TFolder*)");
1557 pluginHandler = pluginManager->FindHandler("AliLoader", detName);
1559 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
1561 (AliLoader*) pluginHandler->ExecPlugin(2, detName.Data(),
1562 fRunLoader->GetEventFolder());
1564 if (!fLoader[iDet]) { // use default loader
1565 fLoader[iDet] = new AliLoader(detName, fRunLoader->GetEventFolder());
1567 if (!fLoader[iDet]) {
1568 AliWarning(Form("couldn't get loader for %s", detName.Data()));
1569 if (fStopOnError) return NULL;
1571 fRunLoader->AddLoader(fLoader[iDet]);
1572 fRunLoader->CdGAFile();
1573 if (gFile && !gFile->IsWritable()) gFile->ReOpen("UPDATE");
1574 fRunLoader->Write(0, TObject::kOverwrite);
1579 return reconstructor;
1582 //_____________________________________________________________________________
1583 Bool_t AliReconstruction::CreateVertexer()
1585 // create the vertexer
1588 AliReconstructor* itsReconstructor = GetReconstructor(0);
1589 if (itsReconstructor) {
1590 fVertexer = itsReconstructor->CreateVertexer(fRunLoader);
1593 AliWarning("couldn't create a vertexer for ITS");
1594 if (fStopOnError) return kFALSE;
1600 //_____________________________________________________________________________
1601 Bool_t AliReconstruction::CreateTrackers(const TString& detectors)
1603 // create the trackers
1605 TString detStr = detectors;
1606 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1607 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1608 AliReconstructor* reconstructor = GetReconstructor(iDet);
1609 if (!reconstructor) continue;
1610 TString detName = fgkDetectorName[iDet];
1611 if (detName == "HLT") {
1612 fRunHLTTracking = kTRUE;
1615 if (detName == "MUON") {
1616 fRunMuonTracking = kTRUE;
1621 fTracker[iDet] = reconstructor->CreateTracker(fRunLoader);
1622 if (!fTracker[iDet] && (iDet < 7)) {
1623 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
1624 if (fStopOnError) return kFALSE;
1631 //_____________________________________________________________________________
1632 void AliReconstruction::CleanUp(TFile* file, TFile* fileOld)
1634 // delete trackers and the run loader and close and delete the file
1636 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1637 delete fReconstructor[iDet];
1638 fReconstructor[iDet] = NULL;
1639 fLoader[iDet] = NULL;
1640 delete fTracker[iDet];
1641 fTracker[iDet] = NULL;
1645 delete fDiamondProfile;
1646 fDiamondProfile = NULL;
1661 gSystem->Unlink("AliESDs.old.root");
1666 //_____________________________________________________________________________
1668 Bool_t AliReconstruction::ReadESD(AliESDEvent*& esd, const char* recStep) const
1670 // read the ESD event from a file
1672 if (!esd) return kFALSE;
1674 sprintf(fileName, "ESD_%d.%d_%s.root",
1675 esd->GetRunNumber(), esd->GetEventNumberInFile(), recStep);
1676 if (gSystem->AccessPathName(fileName)) return kFALSE;
1678 AliInfo(Form("reading ESD from file %s", fileName));
1679 AliDebug(1, Form("reading ESD from file %s", fileName));
1680 TFile* file = TFile::Open(fileName);
1681 if (!file || !file->IsOpen()) {
1682 AliError(Form("opening %s failed", fileName));
1689 esd = (AliESDEvent*) file->Get("ESD");
1698 //_____________________________________________________________________________
1699 void AliReconstruction::WriteESD(AliESDEvent* esd, const char* recStep) const
1701 // write the ESD event to a file
1705 sprintf(fileName, "ESD_%d.%d_%s.root",
1706 esd->GetRunNumber(), esd->GetEventNumberInFile(), recStep);
1708 AliDebug(1, Form("writing ESD to file %s", fileName));
1709 TFile* file = TFile::Open(fileName, "recreate");
1710 if (!file || !file->IsOpen()) {
1711 AliError(Form("opening %s failed", fileName));
1723 //_____________________________________________________________________________
1724 void AliReconstruction::ESDFile2AODFile(TFile* esdFile, TFile* aodFile)
1726 // write all files from the given esd file to an aod file
1728 // create an AliAOD object
1729 AliAODEvent *aod = new AliAODEvent();
1730 aod->CreateStdContent();
1736 TTree *aodTree = new TTree("AOD", "AliAOD tree");
1737 aodTree->Branch(aod->GetList());
1740 TTree *t = (TTree*) esdFile->Get("esdTree");
1741 AliESDEvent *esd = new AliESDEvent();
1742 esd->ReadFromTree(t);
1744 Int_t nEvents = t->GetEntries();
1746 // set arrays and pointers
1754 // loop over events and fill them
1755 for (Int_t iEvent = 0; iEvent < nEvents; ++iEvent) {
1756 t->GetEntry(iEvent);
1758 // Multiplicity information needed by the header (to be revised!)
1759 Int_t nTracks = esd->GetNumberOfTracks();
1760 Int_t nPosTracks = 0;
1761 for (Int_t iTrack=0; iTrack<nTracks; ++iTrack)
1762 if (esd->GetTrack(iTrack)->GetSign()> 0) nPosTracks++;
1764 // Update the header
1765 AliAODHeader* header = aod->GetHeader();
1767 header->SetRunNumber (esd->GetRunNumber() );
1768 header->SetBunchCrossNumber(esd->GetBunchCrossNumber());
1769 header->SetOrbitNumber (esd->GetOrbitNumber() );
1770 header->SetPeriodNumber (esd->GetPeriodNumber() );
1771 header->SetTriggerMask (esd->GetTriggerMask() );
1772 header->SetTriggerCluster (esd->GetTriggerCluster() );
1773 header->SetEventType (esd->GetEventType() );
1774 header->SetMagneticField (esd->GetMagneticField() );
1775 header->SetZDCN1Energy (esd->GetZDCN1Energy() );
1776 header->SetZDCP1Energy (esd->GetZDCP1Energy() );
1777 header->SetZDCN2Energy (esd->GetZDCN2Energy() );
1778 header->SetZDCP2Energy (esd->GetZDCP2Energy() );
1779 header->SetZDCEMEnergy (esd->GetZDCEMEnergy() );
1780 header->SetRefMultiplicity (nTracks);
1781 header->SetRefMultiplicityPos(nPosTracks);
1782 header->SetRefMultiplicityNeg(nTracks - nPosTracks);
1783 header->SetMuonMagFieldScale(-999.); // FIXME
1784 header->SetCentrality(-999.); // FIXME
1788 Int_t nV0s = esd->GetNumberOfV0s();
1789 Int_t nCascades = esd->GetNumberOfCascades();
1790 Int_t nKinks = esd->GetNumberOfKinks();
1791 Int_t nVertices = nV0s + nCascades + nKinks;
1793 aod->ResetStd(nTracks, nVertices);
1794 AliAODTrack *aodTrack;
1797 // Array to take into account the tracks already added to the AOD
1798 Bool_t * usedTrack = NULL;
1800 usedTrack = new Bool_t[nTracks];
1801 for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) usedTrack[iTrack]=kFALSE;
1803 // Array to take into account the V0s already added to the AOD
1804 Bool_t * usedV0 = NULL;
1806 usedV0 = new Bool_t[nV0s];
1807 for (Int_t iV0=0; iV0<nV0s; ++iV0) usedV0[iV0]=kFALSE;
1809 // Array to take into account the kinks already added to the AOD
1810 Bool_t * usedKink = NULL;
1812 usedKink = new Bool_t[nKinks];
1813 for (Int_t iKink=0; iKink<nKinks; ++iKink) usedKink[iKink]=kFALSE;
1816 // Access to the AOD container of vertices
1817 TClonesArray &vertices = *(aod->GetVertices());
1820 // Access to the AOD container of tracks
1821 TClonesArray &tracks = *(aod->GetTracks());
1824 // Add primary vertex. The primary tracks will be defined
1825 // after the loops on the composite objects (V0, cascades, kinks)
1826 const AliESDVertex *vtx = esd->GetPrimaryVertex();
1828 vtx->GetXYZ(pos); // position
1829 vtx->GetCovMatrix(covVtx); //covariance matrix
1831 AliAODVertex * primary = new(vertices[jVertices++])
1832 AliAODVertex(pos, covVtx, vtx->GetChi2toNDF(), NULL, AliAODVertex::kPrimary);
1834 // Create vertices starting from the most complex objects
1837 for (Int_t nCascade = 0; nCascade < nCascades; ++nCascade) {
1838 AliESDcascade *cascade = esd->GetCascade(nCascade);
1840 cascade->GetXYZ(pos[0], pos[1], pos[2]);
1841 cascade->GetPosCovXi(covVtx);
1843 // Add the cascade vertex
1844 AliAODVertex * vcascade = new(vertices[jVertices++]) AliAODVertex(pos,
1846 cascade->GetChi2Xi(), // = chi2/NDF since NDF = 2*2-3
1848 AliAODVertex::kCascade);
1850 primary->AddDaughter(vcascade);
1852 // Add the V0 from the cascade. The ESD class have to be optimized...
1853 // Now we have to search for the corresponding Vo in the list of V0s
1854 // using the indeces of the positive and negative tracks
1856 Int_t posFromV0 = cascade->GetPindex();
1857 Int_t negFromV0 = cascade->GetNindex();
1860 AliESDv0 * v0 = 0x0;
1863 for (Int_t iV0=0; iV0<nV0s; ++iV0) {
1865 v0 = esd->GetV0(iV0);
1866 Int_t posV0 = v0->GetPindex();
1867 Int_t negV0 = v0->GetNindex();
1869 if (posV0==posFromV0 && negV0==negFromV0) {
1875 AliAODVertex * vV0FromCascade = 0x0;
1877 if (indV0>-1 && !usedV0[indV0] ) {
1879 // the V0 exists in the array of V0s and is not used
1881 usedV0[indV0] = kTRUE;
1883 v0->GetXYZ(pos[0], pos[1], pos[2]);
1884 v0->GetPosCov(covVtx);
1886 vV0FromCascade = new(vertices[jVertices++]) AliAODVertex(pos,
1888 v0->GetChi2V0(), // = chi2/NDF since NDF = 2*2-3
1893 // the V0 doesn't exist in the array of V0s or was used
1894 cerr << "Error: event " << iEvent << " cascade " << nCascade
1895 << " The V0 " << indV0
1896 << " doesn't exist in the array of V0s or was used!" << endl;
1898 cascade->GetXYZ(pos[0], pos[1], pos[2]);
1899 cascade->GetPosCov(covVtx);
1901 vV0FromCascade = new(vertices[jVertices++]) AliAODVertex(pos,
1903 v0->GetChi2V0(), // = chi2/NDF since NDF = 2*2-3
1906 vcascade->AddDaughter(vV0FromCascade);
1909 // Add the positive tracks from the V0
1911 if (! usedTrack[posFromV0]) {
1913 usedTrack[posFromV0] = kTRUE;
1915 AliESDtrack *esdTrack = esd->GetTrack(posFromV0);
1916 esdTrack->GetPxPyPz(p);
1917 esdTrack->GetXYZ(pos);
1918 esdTrack->GetCovarianceXYZPxPyPz(covTr);
1919 esdTrack->GetESDpid(pid);
1921 vV0FromCascade->AddDaughter(aodTrack =
1922 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
1923 esdTrack->GetLabel(),
1929 (Short_t)esdTrack->GetSign(),
1930 esdTrack->GetITSClusterMap(),
1933 kTRUE, // check if this is right
1934 kFALSE, // check if this is right
1935 AliAODTrack::kSecondary)
1937 aodTrack->ConvertAliPIDtoAODPID();
1940 cerr << "Error: event " << iEvent << " cascade " << nCascade
1941 << " track " << posFromV0 << " has already been used!" << endl;
1944 // Add the negative tracks from the V0
1946 if (!usedTrack[negFromV0]) {
1948 usedTrack[negFromV0] = kTRUE;
1950 AliESDtrack *esdTrack = esd->GetTrack(negFromV0);
1951 esdTrack->GetPxPyPz(p);
1952 esdTrack->GetXYZ(pos);
1953 esdTrack->GetCovarianceXYZPxPyPz(covTr);
1954 esdTrack->GetESDpid(pid);
1956 vV0FromCascade->AddDaughter(aodTrack =
1957 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
1958 esdTrack->GetLabel(),
1964 (Short_t)esdTrack->GetSign(),
1965 esdTrack->GetITSClusterMap(),
1968 kTRUE, // check if this is right
1969 kFALSE, // check if this is right
1970 AliAODTrack::kSecondary)
1972 aodTrack->ConvertAliPIDtoAODPID();
1975 cerr << "Error: event " << iEvent << " cascade " << nCascade
1976 << " track " << negFromV0 << " has already been used!" << endl;
1979 // Add the bachelor track from the cascade
1981 Int_t bachelor = cascade->GetBindex();
1983 if(!usedTrack[bachelor]) {
1985 usedTrack[bachelor] = kTRUE;
1987 AliESDtrack *esdTrack = esd->GetTrack(bachelor);
1988 esdTrack->GetPxPyPz(p);
1989 esdTrack->GetXYZ(pos);
1990 esdTrack->GetCovarianceXYZPxPyPz(covTr);
1991 esdTrack->GetESDpid(pid);
1993 vcascade->AddDaughter(aodTrack =
1994 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
1995 esdTrack->GetLabel(),
2001 (Short_t)esdTrack->GetSign(),
2002 esdTrack->GetITSClusterMap(),
2005 kTRUE, // check if this is right
2006 kFALSE, // check if this is right
2007 AliAODTrack::kSecondary)
2009 aodTrack->ConvertAliPIDtoAODPID();
2012 cerr << "Error: event " << iEvent << " cascade " << nCascade
2013 << " track " << bachelor << " has already been used!" << endl;
2016 // Add the primary track of the cascade (if any)
2018 } // end of the loop on cascades
2022 for (Int_t nV0 = 0; nV0 < nV0s; ++nV0) {
2024 if (usedV0[nV0]) continue; // skip if aready added to the AOD
2026 AliESDv0 *v0 = esd->GetV0(nV0);
2028 v0->GetXYZ(pos[0], pos[1], pos[2]);
2029 v0->GetPosCov(covVtx);
2031 AliAODVertex * vV0 =
2032 new(vertices[jVertices++]) AliAODVertex(pos,
2034 v0->GetChi2V0(), // = chi2/NDF since NDF = 2*2-3
2037 primary->AddDaughter(vV0);
2039 Int_t posFromV0 = v0->GetPindex();
2040 Int_t negFromV0 = v0->GetNindex();
2042 // Add the positive tracks from the V0
2044 if (!usedTrack[posFromV0]) {
2046 usedTrack[posFromV0] = kTRUE;
2048 AliESDtrack *esdTrack = esd->GetTrack(posFromV0);
2049 esdTrack->GetPxPyPz(p);
2050 esdTrack->GetXYZ(pos);
2051 esdTrack->GetCovarianceXYZPxPyPz(covTr);
2052 esdTrack->GetESDpid(pid);
2054 vV0->AddDaughter(aodTrack =
2055 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2056 esdTrack->GetLabel(),
2062 (Short_t)esdTrack->GetSign(),
2063 esdTrack->GetITSClusterMap(),
2066 kTRUE, // check if this is right
2067 kFALSE, // check if this is right
2068 AliAODTrack::kSecondary)
2070 aodTrack->ConvertAliPIDtoAODPID();
2073 cerr << "Error: event " << iEvent << " V0 " << nV0
2074 << " track " << posFromV0 << " has already been used!" << endl;
2077 // Add the negative tracks from the V0
2079 if (!usedTrack[negFromV0]) {
2081 usedTrack[negFromV0] = kTRUE;
2083 AliESDtrack *esdTrack = esd->GetTrack(negFromV0);
2084 esdTrack->GetPxPyPz(p);
2085 esdTrack->GetXYZ(pos);
2086 esdTrack->GetCovarianceXYZPxPyPz(covTr);
2087 esdTrack->GetESDpid(pid);
2089 vV0->AddDaughter(aodTrack =
2090 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2091 esdTrack->GetLabel(),
2097 (Short_t)esdTrack->GetSign(),
2098 esdTrack->GetITSClusterMap(),
2101 kTRUE, // check if this is right
2102 kFALSE, // check if this is right
2103 AliAODTrack::kSecondary)
2105 aodTrack->ConvertAliPIDtoAODPID();
2108 cerr << "Error: event " << iEvent << " V0 " << nV0
2109 << " track " << negFromV0 << " has already been used!" << endl;
2112 } // end of the loop on V0s
2114 // Kinks: it is a big mess the access to the information in the kinks
2115 // The loop is on the tracks in order to find the mother and daugther of each kink
2118 for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) {
2121 AliESDtrack * esdTrack = esd->GetTrack(iTrack);
2123 Int_t ikink = esdTrack->GetKinkIndex(0);
2126 // Negative kink index: mother, positive: daughter
2128 // Search for the second track of the kink
2130 for (Int_t jTrack = iTrack+1; jTrack<nTracks; ++jTrack) {
2132 AliESDtrack * esdTrack1 = esd->GetTrack(jTrack);
2134 Int_t jkink = esdTrack1->GetKinkIndex(0);
2136 if ( TMath::Abs(ikink)==TMath::Abs(jkink) ) {
2138 // The two tracks are from the same kink
2140 if (usedKink[TMath::Abs(ikink)-1]) continue; // skip used kinks
2143 Int_t idaughter = -1;
2145 if (ikink<0 && jkink>0) {
2150 else if (ikink>0 && jkink<0) {
2156 cerr << "Error: Wrong combination of kink indexes: "
2157 << ikink << " " << jkink << endl;
2161 // Add the mother track
2163 AliAODTrack * mother = NULL;
2165 if (!usedTrack[imother]) {
2167 usedTrack[imother] = kTRUE;
2169 AliESDtrack *esdTrack = esd->GetTrack(imother);
2170 esdTrack->GetPxPyPz(p);
2171 esdTrack->GetXYZ(pos);
2172 esdTrack->GetCovarianceXYZPxPyPz(covTr);
2173 esdTrack->GetESDpid(pid);
2176 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2177 esdTrack->GetLabel(),
2183 (Short_t)esdTrack->GetSign(),
2184 esdTrack->GetITSClusterMap(),
2187 kTRUE, // check if this is right
2188 kTRUE, // check if this is right
2189 AliAODTrack::kPrimary);
2190 primary->AddDaughter(mother);
2191 mother->ConvertAliPIDtoAODPID();
2194 cerr << "Error: event " << iEvent << " kink " << TMath::Abs(ikink)-1
2195 << " track " << imother << " has already been used!" << endl;
2198 // Add the kink vertex
2199 AliESDkink * kink = esd->GetKink(TMath::Abs(ikink)-1);
2201 AliAODVertex * vkink =
2202 new(vertices[jVertices++]) AliAODVertex(kink->GetPosition(),
2206 AliAODVertex::kKink);
2207 // Add the daughter track
2209 AliAODTrack * daughter = NULL;
2211 if (!usedTrack[idaughter]) {
2213 usedTrack[idaughter] = kTRUE;
2215 AliESDtrack *esdTrack = esd->GetTrack(idaughter);
2216 esdTrack->GetPxPyPz(p);
2217 esdTrack->GetXYZ(pos);
2218 esdTrack->GetCovarianceXYZPxPyPz(covTr);
2219 esdTrack->GetESDpid(pid);
2222 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2223 esdTrack->GetLabel(),
2229 (Short_t)esdTrack->GetSign(),
2230 esdTrack->GetITSClusterMap(),
2233 kTRUE, // check if this is right
2234 kTRUE, // check if this is right
2235 AliAODTrack::kPrimary);
2236 vkink->AddDaughter(daughter);
2237 daughter->ConvertAliPIDtoAODPID();
2240 cerr << "Error: event " << iEvent << " kink " << TMath::Abs(ikink)-1
2241 << " track " << idaughter << " has already been used!" << endl;
2253 // Tracks (primary and orphan)
2255 for (Int_t nTrack = 0; nTrack < nTracks; ++nTrack) {
2258 if (usedTrack[nTrack]) continue;
2260 AliESDtrack *esdTrack = esd->GetTrack(nTrack);
2261 esdTrack->GetPxPyPz(p);
2262 esdTrack->GetXYZ(pos);
2263 esdTrack->GetCovarianceXYZPxPyPz(covTr);
2264 esdTrack->GetESDpid(pid);
2266 Float_t impactXY, impactZ;
2268 esdTrack->GetImpactParameters(impactXY,impactZ);
2271 // track inside the beam pipe
2273 primary->AddDaughter(aodTrack =
2274 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2275 esdTrack->GetLabel(),
2281 (Short_t)esdTrack->GetSign(),
2282 esdTrack->GetITSClusterMap(),
2285 kTRUE, // check if this is right
2286 kTRUE, // check if this is right
2287 AliAODTrack::kPrimary)
2289 aodTrack->ConvertAliPIDtoAODPID();
2292 // outside the beam pipe: orphan track
2294 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2295 esdTrack->GetLabel(),
2301 (Short_t)esdTrack->GetSign(),
2302 esdTrack->GetITSClusterMap(),
2305 kFALSE, // check if this is right
2306 kFALSE, // check if this is right
2307 AliAODTrack::kOrphan);
2308 aodTrack->ConvertAliPIDtoAODPID();
2310 } // end of loop on tracks
2313 Int_t nMuTracks = esd->GetNumberOfMuonTracks();
2314 for (Int_t nMuTrack = 0; nMuTrack < nMuTracks; ++nMuTrack) {
2316 AliESDMuonTrack *esdMuTrack = esd->GetMuonTrack(nMuTrack);
2317 p[0] = esdMuTrack->Px();
2318 p[1] = esdMuTrack->Py();
2319 p[2] = esdMuTrack->Pz();
2320 pos[0] = primary->GetX();
2321 pos[1] = primary->GetY();
2322 pos[2] = primary->GetZ();
2324 // has to be changed once the muon pid is provided by the ESD
2325 for (Int_t i = 0; i < 10; pid[i++] = 0.); pid[AliAODTrack::kMuon]=1.;
2327 primary->AddDaughter(
2328 new(tracks[jTracks++]) AliAODTrack(0, // no ID provided
2329 0, // no label provided
2334 NULL, // no covariance matrix provided
2335 (Short_t)-99, // no charge provided
2336 0, // no ITSClusterMap
2339 kTRUE, // check if this is right
2340 kTRUE, // not used for vertex fit
2341 AliAODTrack::kPrimary)
2345 // Access to the AOD container of clusters
2346 TClonesArray &clusters = *(aod->GetClusters());
2350 Int_t nClusters = esd->GetNumberOfCaloClusters();
2352 for (Int_t iClust=0; iClust<nClusters; ++iClust) {
2354 AliESDCaloCluster * cluster = esd->GetCaloCluster(iClust);
2356 Int_t id = cluster->GetID();
2358 Float_t energy = cluster->E();
2359 cluster->GetPosition(posF);
2360 AliAODVertex *prodVertex = primary;
2361 AliAODTrack *primTrack = NULL;
2362 Char_t ttype=AliAODCluster::kUndef;
2364 if (cluster->IsPHOS()) ttype=AliAODCluster::kPHOSNeutral;
2365 else if (cluster->IsEMCAL()) {
2367 if (cluster->GetClusterType() == AliESDCaloCluster::kPseudoCluster)
2368 ttype = AliAODCluster::kEMCALPseudoCluster;
2370 ttype = AliAODCluster::kEMCALClusterv1;
2374 new(clusters[jClusters++]) AliAODCluster(id,
2378 NULL, // no covariance matrix provided
2379 NULL, // no pid for clusters provided
2384 } // end of loop on calo clusters
2386 delete [] usedTrack;
2390 // fill the tree for this event
2392 } // end of event loop
2394 aodTree->GetUserInfo()->Add(aod);
2399 // write the tree to the specified file
2400 aodFile = aodTree->GetCurrentFile();
2407 void AliReconstruction::WriteAlignmentData(AliESDEvent* esd)
2409 // Write space-points which are then used in the alignment procedures
2410 // For the moment only ITS, TRD and TPC
2412 // Load TOF clusters
2414 fLoader[3]->LoadRecPoints("read");
2415 TTree* tree = fLoader[3]->TreeR();
2417 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[3]));
2420 fTracker[3]->LoadClusters(tree);
2422 Int_t ntracks = esd->GetNumberOfTracks();
2423 for (Int_t itrack = 0; itrack < ntracks; itrack++)
2425 AliESDtrack *track = esd->GetTrack(itrack);
2428 for (Int_t iDet = 3; iDet >= 0; iDet--)
2429 nsp += track->GetNcls(iDet);
2431 AliTrackPointArray *sp = new AliTrackPointArray(nsp);
2432 track->SetTrackPointArray(sp);
2434 for (Int_t iDet = 3; iDet >= 0; iDet--) {
2435 AliTracker *tracker = fTracker[iDet];
2436 if (!tracker) continue;
2437 Int_t nspdet = track->GetNcls(iDet);
2438 if (nspdet <= 0) continue;
2439 track->GetClusters(iDet,idx);
2443 while (isp < nspdet) {
2444 Bool_t isvalid = tracker->GetTrackPoint(idx[isp2],p); isp2++;
2445 const Int_t kNTPCmax = 159;
2446 if (iDet==1 && isp2>kNTPCmax) break; // to be fixed
2447 if (!isvalid) continue;
2448 sp->AddPoint(isptrack,&p); isptrack++; isp++;
2454 fTracker[3]->UnloadClusters();
2455 fLoader[3]->UnloadRecPoints();
2459 //_____________________________________________________________________________
2460 void AliReconstruction::FillRawDataErrorLog(Int_t iEvent, AliESDEvent* esd)
2462 // The method reads the raw-data error log
2463 // accumulated within the rawReader.
2464 // It extracts the raw-data errors related to
2465 // the current event and stores them into
2466 // a TClonesArray inside the esd object.
2468 if (!fRawReader) return;
2470 for(Int_t i = 0; i < fRawReader->GetNumberOfErrorLogs(); i++) {
2472 AliRawDataErrorLog *log = fRawReader->GetErrorLog(i);
2474 if (iEvent != log->GetEventNumber()) continue;
2476 esd->AddRawDataErrorLog(log);
2481 TNamed* AliReconstruction::CopyFileToTNamed(TString fPath,TString fName){
2482 // Dump a file content into a char in TNamed
2484 in.open(fPath.Data(),ios::in | ios::binary|ios::ate);
2485 Int_t kBytes = (Int_t)in.tellg();
2486 printf("Size: %d \n",kBytes);
2489 char* memblock = new char [kBytes];
2490 in.seekg (0, ios::beg);
2491 in.read (memblock, kBytes);
2493 TString fData(memblock,kBytes);
2494 fn = new TNamed(fName,fData);
2495 printf("fData Size: %d \n",fData.Sizeof());
2496 printf("fName Size: %d \n",fName.Sizeof());
2497 printf("fn Size: %d \n",fn->Sizeof());
2501 AliInfo(Form("Could not Open %s\n",fPath.Data()));
2507 void AliReconstruction::TNamedToFile(TTree* fTree, TString fName){
2508 // This is not really needed in AliReconstruction at the moment
2509 // but can serve as a template
2511 TList *fList = fTree->GetUserInfo();
2512 TNamed *fn = (TNamed*)fList->FindObject(fName.Data());
2513 printf("fn Size: %d \n",fn->Sizeof());
2515 TString fTmp(fn->GetName()); // to be 100% sure in principle fName also works
2516 const char* cdata = fn->GetTitle();
2517 printf("fTmp Size %d\n",fTmp.Sizeof());
2519 int size = fn->Sizeof()-fTmp.Sizeof()-sizeof(UChar_t)-sizeof(Int_t); // see dfinition of TString::SizeOf()...
2520 printf("calculated size %d\n",size);
2521 ofstream out(fName.Data(),ios::out | ios::binary);
2522 out.write(cdata,size);