f6c206816af15f5dfedcdaa7f930b0d984dea259
[u/mrichter/AliRoot.git] / STEER / AliReconstruction.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15
16 /* $Id$ */
17
18 ///////////////////////////////////////////////////////////////////////////////
19 //                                                                           //
20 // class for running the reconstruction                                      //
21 //                                                                           //
22 // Clusters and tracks are created for all detectors and all events by       //
23 // typing:                                                                   //
24 //                                                                           //
25 //   AliReconstruction rec;                                                  //
26 //   rec.Run();                                                              //
27 //                                                                           //
28 // The Run method returns kTRUE in case of successful execution.             //
29 //                                                                           //
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   //
32 //                                                                           //
33 //   rec.SetInput("...");                                                    //
34 //                                                                           //
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                               //
40 //                                                                           //
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                 //
44 //                                                                           //
45 //   rec.SetEventRange(..., ...);                                            //
46 //                                                                           //
47 // The index -1 (default) can be used for the last event to indicate no      //
48 // upper limit of the event range.                                           //
49 //                                                                           //
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):                     //
55 //                                                                           //
56 //   rec.SetNumberOfEventsPerFile(...);                                      //
57 //                                                                           //
58 //                                                                           //
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                                                         //
62 //                                                                           //
63 //   rec.SetGAliceFile("...");                                               //
64 //                                                                           //
65 // The local reconstruction can be switched on or off for individual         //
66 // detectors by                                                              //
67 //                                                                           //
68 //   rec.SetRunLocalReconstruction("...");                                   //
69 //                                                                           //
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.                                 //
73 //                                                                           //
74 // The reconstruction of the primary vertex position can be switched off by  //
75 //                                                                           //
76 //   rec.SetRunVertexFinder(kFALSE);                                         //
77 //                                                                           //
78 // The tracking and the creation of ESD tracks can be switched on for        //
79 // selected detectors by                                                     //
80 //                                                                           //
81 //   rec.SetRunTracking("...");                                              //
82 //                                                                           //
83 // Uniform/nonuniform field tracking switches (default: uniform field)       //
84 //                                                                           //
85 //   rec.SetUniformFieldTracking(); ( rec.SetUniformFieldTracking(kFALSE); ) //
86 //                                                                           //
87 // The filling of additional ESD information can be steered by               //
88 //                                                                           //
89 //   rec.SetFillESD("...");                                                  //
90 //                                                                           //
91 // Again, for both methods the string specifies the list of detectors.       //
92 // The default is "ALL".                                                     //
93 //                                                                           //
94 // The call of the shortcut method                                           //
95 //                                                                           //
96 //   rec.SetRunReconstruction("...");                                        //
97 //                                                                           //
98 // is equivalent to calling SetRunLocalReconstruction, SetRunTracking and    //
99 // SetFillESD with the same detector selecting string as argument.           //
100 //                                                                           //
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.            //
103 //                                                                           //
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.     //
114 //                                                                           //
115 ///////////////////////////////////////////////////////////////////////////////
116
117 #include <TArrayF.h>
118 #include <TFile.h>
119 #include <TList.h>
120 #include <TSystem.h>
121 #include <TROOT.h>
122 #include <TPluginManager.h>
123 #include <TGeoManager.h>
124 #include <TLorentzVector.h>
125 #include <TArrayS.h>
126 #include <TArrayD.h>
127
128 #include "AliReconstruction.h"
129 #include "AliCodeTimer.h"
130 #include "AliReconstructor.h"
131 #include "AliLog.h"
132 #include "AliRunLoader.h"
133 #include "AliRun.h"
134 #include "AliRawReaderFile.h"
135 #include "AliRawReaderDate.h"
136 #include "AliRawReaderRoot.h"
137 #include "AliRawEventHeaderBase.h"
138 #include "AliESDEvent.h"
139 #include "AliESDMuonTrack.h"
140 #include "AliESDfriend.h"
141 #include "AliESDVertex.h"
142 #include "AliESDcascade.h"
143 #include "AliESDkink.h"
144 #include "AliESDtrack.h"
145 #include "AliESDCaloCluster.h"
146 #include "AliMultiplicity.h"
147 #include "AliTracker.h"
148 #include "AliVertexer.h"
149 #include "AliVertexerTracks.h"
150 #include "AliV0vertexer.h"
151 #include "AliCascadeVertexer.h"
152 #include "AliHeader.h"
153 #include "AliGenEventHeader.h"
154 #include "AliPID.h"
155 #include "AliESDpid.h"
156 #include "AliESDtrack.h"
157 #include "AliESDPmdTrack.h"
158
159 #include "AliESDTagCreator.h"
160 #include "AliAODTagCreator.h"
161
162 #include "AliGeomManager.h"
163 #include "AliTrackPointArray.h"
164 #include "AliCDBManager.h"
165 #include "AliCDBEntry.h"
166 #include "AliAlignObj.h"
167
168 #include "AliCentralTrigger.h"
169 #include "AliCTPRawStream.h"
170
171 #include "AliAODEvent.h"
172 #include "AliAODHeader.h"
173 #include "AliAODTrack.h"
174 #include "AliAODVertex.h"
175 #include "AliAODv0.h"
176 #include "AliAODJet.h"
177 #include "AliAODCaloCells.h"
178 #include "AliAODCaloCluster.h"
179 #include "AliAODPmdCluster.h"
180 #include "AliAODFmdCluster.h"
181 #include "AliAODTracklets.h"
182
183 //#include "AliQADataMaker.h" 
184 #include "AliQA.h"
185 #include "AliQADataMakerSteer.h"
186
187 #include "AliSysInfo.h" // memory snapshots
188
189
190 ClassImp(AliReconstruction)
191
192
193 //_____________________________________________________________________________
194 const char* AliReconstruction::fgkDetectorName[AliReconstruction::fgkNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "HMPID", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "T0", "VZERO", "ACORDE", "HLT"};
195
196 //_____________________________________________________________________________
197 AliReconstruction::AliReconstruction(const char* gAliceFilename, const char* cdbUri,
198                                      const char* name, const char* title) :
199   TNamed(name, title),
200
201   fUniformField(kTRUE),
202   fRunVertexFinder(kTRUE),
203   fRunHLTTracking(kFALSE),
204   fRunMuonTracking(kFALSE),
205   fRunV0Finder(kTRUE),
206   fRunCascadeFinder(kTRUE),
207   fStopOnError(kFALSE),
208   fWriteAlignmentData(kFALSE),
209   fWriteESDfriend(kFALSE),
210   fWriteAOD(kFALSE),
211   fFillTriggerESD(kTRUE),
212
213   fCleanESD(kTRUE),
214   fV0DCAmax(3.),
215   fV0CsPmin(0.),
216   fDmax(50.),
217   fZmax(50.),
218
219   fRunLocalReconstruction("ALL"),
220   fRunTracking("ALL"),
221   fFillESD("ALL"),
222   fUseTrackingErrorsForAlignment(""),
223   fGAliceFileName(gAliceFilename),
224   fInput(""),
225   fEquipIdMap(""),
226   fFirstEvent(0),
227   fLastEvent(-1),
228   fNumberOfEventsPerFile(1),
229   fCheckPointLevel(0),
230   fOptions(),
231   fLoadAlignFromCDB(kTRUE),
232   fLoadAlignData("ALL"),
233   fESDPar(""),
234
235   fRunLoader(NULL),
236   fRawReader(NULL),
237
238   fVertexer(NULL),
239   fDiamondProfile(NULL),
240
241   fGRPList(NULL),
242
243   fAlignObjArray(NULL),
244   fCDBUri(cdbUri),
245   fRemoteCDBUri(""),
246   fSpecCDBUri(), 
247   fRunQA(kTRUE) 
248
249 {
250 // create reconstruction object with default parameters
251   
252   for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
253     fReconstructor[iDet] = NULL;
254     fLoader[iDet] = NULL;
255     fTracker[iDet] = NULL;
256 //    fQADataMaker[iDet] = NULL;
257 //      fQACycles[iDet] = 999999;       
258   }
259   AliPID pid;
260 }
261
262 //_____________________________________________________________________________
263 AliReconstruction::AliReconstruction(const AliReconstruction& rec) :
264   TNamed(rec),
265
266   fUniformField(rec.fUniformField),
267   fRunVertexFinder(rec.fRunVertexFinder),
268   fRunHLTTracking(rec.fRunHLTTracking),
269   fRunMuonTracking(rec.fRunMuonTracking),
270   fRunV0Finder(rec.fRunV0Finder),
271   fRunCascadeFinder(rec.fRunCascadeFinder),
272   fStopOnError(rec.fStopOnError),
273   fWriteAlignmentData(rec.fWriteAlignmentData),
274   fWriteESDfriend(rec.fWriteESDfriend),
275   fWriteAOD(rec.fWriteAOD),
276   fFillTriggerESD(rec.fFillTriggerESD),
277
278   fCleanESD(rec.fCleanESD),
279   fV0DCAmax(rec.fV0DCAmax),
280   fV0CsPmin(rec.fV0CsPmin),
281   fDmax(rec.fDmax),
282   fZmax(rec.fZmax),
283
284   fRunLocalReconstruction(rec.fRunLocalReconstruction),
285   fRunTracking(rec.fRunTracking),
286   fFillESD(rec.fFillESD),
287   fUseTrackingErrorsForAlignment(rec.fUseTrackingErrorsForAlignment),
288   fGAliceFileName(rec.fGAliceFileName),
289   fInput(rec.fInput),
290   fEquipIdMap(rec.fEquipIdMap),
291   fFirstEvent(rec.fFirstEvent),
292   fLastEvent(rec.fLastEvent),
293   fNumberOfEventsPerFile(rec.fNumberOfEventsPerFile),
294   fCheckPointLevel(0),
295   fOptions(),
296   fLoadAlignFromCDB(rec.fLoadAlignFromCDB),
297   fLoadAlignData(rec.fLoadAlignData),
298   fESDPar(rec.fESDPar),
299
300   fRunLoader(NULL),
301   fRawReader(NULL),
302
303   fVertexer(NULL),
304   fDiamondProfile(NULL),
305
306   fGRPList(NULL),
307
308   fAlignObjArray(rec.fAlignObjArray),
309   fCDBUri(rec.fCDBUri),
310   fRemoteCDBUri(rec.fRemoteCDBUri),
311   fSpecCDBUri(), 
312   fRunQA(kTRUE)
313 {
314 // copy constructor
315
316   for (Int_t i = 0; i < rec.fOptions.GetEntriesFast(); i++) {
317     if (rec.fOptions[i]) fOptions.Add(rec.fOptions[i]->Clone());
318   }
319   for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
320     fReconstructor[iDet] = NULL;
321     fLoader[iDet] = NULL;
322     fTracker[iDet] = NULL;
323 //    fQADataMaker[iDet] = NULL;
324 //      fQACycles[iDet] = rec.fQACycles[iDet];  
325   }
326   for (Int_t i = 0; i < rec.fSpecCDBUri.GetEntriesFast(); i++) {
327     if (rec.fSpecCDBUri[i]) fSpecCDBUri.Add(rec.fSpecCDBUri[i]->Clone());
328   }
329 }
330
331 //_____________________________________________________________________________
332 AliReconstruction& AliReconstruction::operator = (const AliReconstruction& rec)
333 {
334 // assignment operator
335
336   this->~AliReconstruction();
337   new(this) AliReconstruction(rec);
338   return *this;
339 }
340
341 //_____________________________________________________________________________
342 AliReconstruction::~AliReconstruction()
343 {
344 // clean up
345
346   CleanUp();
347   fOptions.Delete();
348   fSpecCDBUri.Delete();
349
350   AliCodeTimer::Instance()->Print();
351 }
352
353 //_____________________________________________________________________________
354 void AliReconstruction::InitCDBStorage()
355 {
356 // activate a default CDB storage
357 // First check if we have any CDB storage set, because it is used 
358 // to retrieve the calibration and alignment constants
359
360   AliCDBManager* man = AliCDBManager::Instance();
361   if (man->IsDefaultStorageSet())
362   {
363     AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
364     AliWarning("Default CDB storage has been already set !");
365     AliWarning(Form("Ignoring the default storage declared in AliReconstruction: %s",fCDBUri.Data()));
366     AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
367     fCDBUri = "";
368   }
369   else {
370     AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
371     AliDebug(2, Form("Default CDB storage is set to: %s",fCDBUri.Data()));
372     AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
373     man->SetDefaultStorage(fCDBUri);
374   }
375
376   // Remote storage (the Grid storage) is used if it is activated
377   // and if the object is not found in the default storage
378   // OBSOLETE: Removed
379   //  if (man->IsRemoteStorageSet())
380   //  {
381   //    AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
382   //    AliWarning("Remote CDB storage has been already set !");
383   //    AliWarning(Form("Ignoring the remote storage declared in AliReconstruction: %s",fRemoteCDBUri.Data()));
384   //    AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
385   //    fRemoteCDBUri = "";
386   //  }
387   //  else {
388   //    AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
389   //    AliDebug(2, Form("Remote CDB storage is set to: %s",fRemoteCDBUri.Data()));
390   //    AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
391   //    man->SetRemoteStorage(fRemoteCDBUri);
392   //  }
393
394   // Now activate the detector specific CDB storage locations
395   for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
396     TObject* obj = fSpecCDBUri[i];
397     if (!obj) continue;
398     AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
399     AliDebug(2, Form("Specific CDB storage for %s is set to: %s",obj->GetName(),obj->GetTitle()));
400     AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
401     man->SetSpecificStorage(obj->GetName(), obj->GetTitle());
402   }
403   man->Print();
404 }
405
406 //_____________________________________________________________________________
407 void AliReconstruction::SetDefaultStorage(const char* uri) {
408 // Store the desired default CDB storage location
409 // Activate it later within the Run() method
410
411   fCDBUri = uri;
412
413 }
414
415 //_____________________________________________________________________________
416 void AliReconstruction::SetRemoteStorage(const char* uri) {
417 // Store the desired remote CDB storage location
418 // Activate it later within the Run() method
419 // Remote storage (the Grid storage) is used if it is activated
420 // and if the object is not found in the default storage
421
422   fRemoteCDBUri = uri;
423
424 }
425
426 //_____________________________________________________________________________
427 void AliReconstruction::SetSpecificStorage(const char* calibType, const char* uri) {
428 // Store a detector-specific CDB storage location
429 // Activate it later within the Run() method
430
431   AliCDBPath aPath(calibType);
432   if(!aPath.IsValid()){
433         // if calibType is not wildcard but it is a valid detector, add "/*" to make it a valid path
434         for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
435                 if(!strcmp(calibType, fgkDetectorName[iDet])) {
436                         aPath.SetPath(Form("%s/*", calibType));
437                         AliInfo(Form("Path for specific storage set to %s", aPath.GetPath().Data()));
438                         break;
439                 }
440         }
441         if(!aPath.IsValid()){
442                 AliError(Form("Not a valid path or detector: %s", calibType));
443                 return;
444         }
445   }
446
447 //  // check that calibType refers to a "valid" detector name
448 //  Bool_t isDetector = kFALSE;
449 //  for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
450 //    TString detName = fgkDetectorName[iDet];
451 //    if(aPath.GetLevel0() == detName) {
452 //      isDetector = kTRUE;
453 //      break;
454 //    }
455 //  }
456 //
457 //  if(!isDetector) {
458 //      AliError(Form("Not a valid detector: %s", aPath.GetLevel0().Data()));
459 //      return;
460 //  }
461
462   TObject* obj = fSpecCDBUri.FindObject(aPath.GetPath().Data());
463   if (obj) fSpecCDBUri.Remove(obj);
464   fSpecCDBUri.Add(new TNamed(aPath.GetPath().Data(), uri));
465
466 }
467
468
469
470
471 //_____________________________________________________________________________
472 Bool_t AliReconstruction::SetRunNumber()
473 {
474   // The method is called in Run() in order
475   // to set a correct run number.
476   // In case of raw data reconstruction the
477   // run number is taken from the raw data header
478
479   if(AliCDBManager::Instance()->GetRun() < 0) {
480     if (!fRunLoader) {
481       AliError("No run loader is found !"); 
482       return kFALSE;
483     }
484     // read run number from gAlice
485     if(fRunLoader->GetAliRun())
486       AliCDBManager::Instance()->SetRun(fRunLoader->GetHeader()->GetRun());
487     else {
488       if(fRawReader) {
489         if(fRawReader->NextEvent()) {
490           AliCDBManager::Instance()->SetRun(fRawReader->GetRunNumber());
491           fRawReader->RewindEvents();
492         }
493         else {
494           AliError("No raw-data events found !");
495           return kFALSE;
496         }
497       }
498       else {
499         AliError("Neither gAlice nor RawReader objects are found !");
500         return kFALSE;
501       }
502     }
503     AliInfo(Form("CDB Run number: %d",AliCDBManager::Instance()->GetRun()));
504   }
505   return kTRUE;
506 }
507
508 //_____________________________________________________________________________
509 Bool_t AliReconstruction::MisalignGeometry(const TString& detectors)
510 {
511   // Read the alignment objects from CDB.
512   // Each detector is supposed to have the
513   // alignment objects in DET/Align/Data CDB path.
514   // All the detector objects are then collected,
515   // sorted by geometry level (starting from ALIC) and
516   // then applied to the TGeo geometry.
517   // Finally an overlaps check is performed.
518
519   // Load alignment data from CDB and fill fAlignObjArray 
520   if(fLoadAlignFromCDB){
521         
522     TString detStr = detectors;
523     TString loadAlObjsListOfDets = "";
524     
525     for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
526       if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
527       loadAlObjsListOfDets += fgkDetectorName[iDet];
528       loadAlObjsListOfDets += " ";
529     } // end loop over detectors
530     loadAlObjsListOfDets.Prepend("GRP "); //add alignment objects for non-sensitive modules
531     AliGeomManager::ApplyAlignObjsFromCDB(loadAlObjsListOfDets.Data());
532   }else{
533     // Check if the array with alignment objects was
534     // provided by the user. If yes, apply the objects
535     // to the present TGeo geometry
536     if (fAlignObjArray) {
537       if (gGeoManager && gGeoManager->IsClosed()) {
538         if (AliGeomManager::ApplyAlignObjsToGeom(*fAlignObjArray) == kFALSE) {
539           AliError("The misalignment of one or more volumes failed!"
540                    "Compare the list of simulated detectors and the list of detector alignment data!");
541           return kFALSE;
542         }
543       }
544       else {
545         AliError("Can't apply the misalignment! gGeoManager doesn't exist or it is still opened!");
546         return kFALSE;
547       }
548     }
549   }
550   
551   delete fAlignObjArray; fAlignObjArray=0;
552
553   return kTRUE;
554 }
555
556 //_____________________________________________________________________________
557 void AliReconstruction::SetGAliceFile(const char* fileName)
558 {
559 // set the name of the galice file
560
561   fGAliceFileName = fileName;
562 }
563
564 //_____________________________________________________________________________
565 void AliReconstruction::SetOption(const char* detector, const char* option)
566 {
567 // set options for the reconstruction of a detector
568
569   TObject* obj = fOptions.FindObject(detector);
570   if (obj) fOptions.Remove(obj);
571   fOptions.Add(new TNamed(detector, option));
572 }
573
574
575 //_____________________________________________________________________________
576 Bool_t AliReconstruction::Run(const char* input)
577 {
578 // run the reconstruction
579
580   AliCodeTimerAuto("")
581   
582   // set the input
583   if (!input) input = fInput.Data();
584      
585   TString fileName(input);
586   if (fileName.EndsWith("/")) {
587     fRawReader = new AliRawReaderFile(fileName);
588   } else if (fileName.EndsWith(".root")) {
589     fRawReader = new AliRawReaderRoot(fileName);
590   } else if (!fileName.IsNull()) {
591     fRawReader = new AliRawReaderDate(fileName);
592     fRawReader->SelectEvents(7);
593   }
594   if (!fEquipIdMap.IsNull() && fRawReader)
595     fRawReader->LoadEquipmentIdsMap(fEquipIdMap);
596
597    AliSysInfo::AddStamp("Start");
598   // get the run loader
599   if (!InitRunLoader()) return kFALSE;
600    AliSysInfo::AddStamp("LoadLoader");
601
602   // Initialize the CDB storage
603   InitCDBStorage();
604    AliSysInfo::AddStamp("LoadCDB");
605
606   // Set run number in CDBManager (if it is not already set by the user)
607   if (!SetRunNumber()) if (fStopOnError) return kFALSE;
608
609   // Import ideal TGeo geometry and apply misalignment
610   if (!gGeoManager) {
611     TString geom(gSystem->DirName(fGAliceFileName));
612     geom += "/geometry.root";
613     AliGeomManager::LoadGeometry(geom.Data());
614     if (!gGeoManager) if (fStopOnError) return kFALSE;
615   }
616
617   if (!MisalignGeometry(fLoadAlignData)) if (fStopOnError) return kFALSE;
618    AliSysInfo::AddStamp("LoadGeom");
619
620   //QA 
621   AliQADataMakerSteer qas ; 
622   if ( fRunQA && fRawReader) 
623                 qas.Run("ALL", fRawReader) ; 
624  
625  // checking the QA of previous steps
626   CheckQA() ; 
627  
628   // local reconstruction
629   if (!fRunLocalReconstruction.IsNull()) {
630     if (!RunLocalReconstruction(fRunLocalReconstruction)) {
631       if (fStopOnError) {CleanUp(); return kFALSE;}
632     }
633   }
634 //  if (!fRunVertexFinder && fRunTracking.IsNull() && 
635 //      fFillESD.IsNull()) return kTRUE;
636
637   // get vertexer
638   if (fRunVertexFinder && !CreateVertexer()) {
639     if (fStopOnError) {
640       CleanUp(); 
641       return kFALSE;
642     }
643   }
644    AliSysInfo::AddStamp("Vertexer");
645
646   // get trackers
647   if (!fRunTracking.IsNull() && !CreateTrackers(fRunTracking)) {
648     if (fStopOnError) {
649       CleanUp(); 
650       return kFALSE;
651     }      
652   }
653    AliSysInfo::AddStamp("LoadTrackers");
654
655   // get the possibly already existing ESD file and tree
656   AliESDEvent* esd = new AliESDEvent(); AliESDEvent* hltesd = new AliESDEvent();
657   TFile* fileOld = NULL;
658   TTree* treeOld = NULL; TTree *hlttreeOld = NULL;
659   if (!gSystem->AccessPathName("AliESDs.root")){
660     gSystem->CopyFile("AliESDs.root", "AliESDs.old.root", kTRUE);
661     fileOld = TFile::Open("AliESDs.old.root");
662     if (fileOld && fileOld->IsOpen()) {
663       treeOld = (TTree*) fileOld->Get("esdTree");
664       if (treeOld)esd->ReadFromTree(treeOld);
665       hlttreeOld = (TTree*) fileOld->Get("HLTesdTree");
666       if (hlttreeOld)   hltesd->ReadFromTree(hlttreeOld);
667     }
668   }
669
670   // create the ESD output file and tree
671   TFile* file = TFile::Open("AliESDs.root", "RECREATE");
672   file->SetCompressionLevel(2);
673   if (!file->IsOpen()) {
674     AliError("opening AliESDs.root failed");
675     if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}    
676   }
677
678   TTree* tree = new TTree("esdTree", "Tree with ESD objects");
679   esd = new AliESDEvent();
680   esd->CreateStdContent();
681   esd->WriteToTree(tree);
682
683   TTree* hlttree = new TTree("HLTesdTree", "Tree with HLT ESD objects");
684   hltesd = new AliESDEvent();
685   hltesd->CreateStdContent();
686   hltesd->WriteToTree(hlttree);
687
688   /* CKB Why?
689   delete esd; delete hltesd;
690   esd = NULL; hltesd = NULL;
691   */
692   // create the branch with ESD additions
693
694
695
696   AliESDfriend *esdf = 0; 
697   if (fWriteESDfriend) {
698     esdf = new AliESDfriend();
699     TBranch *br=tree->Branch("ESDfriend.","AliESDfriend", &esdf);
700     br->SetFile("AliESDfriends.root");
701     esd->AddObject(esdf);
702   }
703
704   
705   // Get the GRP CDB entry
706   AliCDBEntry* entryGRP = AliCDBManager::Instance()->Get("GRP/GRP/Data");
707         
708   if(entryGRP) {
709         fGRPList = dynamic_cast<TList*> (entryGRP->GetObject());  
710   } else {
711         AliError("No GRP entry found in OCDB!");
712   }
713
714   // Get the diamond profile from OCDB
715   AliCDBEntry* entry = AliCDBManager::Instance()
716         ->Get("GRP/Calib/MeanVertex");
717         
718   if(entry) {
719         fDiamondProfile = dynamic_cast<AliESDVertex*> (entry->GetObject());  
720   } else {
721         AliError("No diamond profile found in OCDB!");
722   }
723
724   AliVertexerTracks tVertexer(AliTracker::GetBz());
725   if(fDiamondProfile) tVertexer.SetVtxStart(fDiamondProfile);
726
727   // loop over events
728  
729   if (fRawReader) fRawReader->RewindEvents();
730    TString detStr(fFillESD) ; 
731
732   ProcInfo_t ProcInfo;
733   gSystem->GetProcInfo(&ProcInfo);
734   AliInfo(Form("Current memory usage %d %d", ProcInfo.fMemResident, ProcInfo.fMemVirtual));
735   
736   for (Int_t iEvent = 0; iEvent < fRunLoader->GetNumberOfEvents(); iEvent++) {
737     if (fRawReader) fRawReader->NextEvent();
738     if ((iEvent < fFirstEvent) || ((fLastEvent >= 0) && (iEvent > fLastEvent))) {
739       // copy old ESD to the new one
740       if (treeOld) {
741         esd->ReadFromTree(treeOld);
742         treeOld->GetEntry(iEvent);
743       }
744       tree->Fill();
745       if (hlttreeOld) {
746         esd->ReadFromTree(hlttreeOld);
747         hlttreeOld->GetEntry(iEvent);
748       }
749       hlttree->Fill();
750       continue;
751     }
752     
753     AliInfo(Form("processing event %d", iEvent));
754     fRunLoader->GetEvent(iEvent);
755
756     char aFileName[256];
757     sprintf(aFileName, "ESD_%d.%d_final.root", 
758             fRunLoader->GetHeader()->GetRun(), 
759             fRunLoader->GetHeader()->GetEventNrInRun());
760     if (!gSystem->AccessPathName(aFileName)) continue;
761
762     // local reconstruction
763     if (!fRunLocalReconstruction.IsNull()) {
764       if (!RunLocalEventReconstruction(fRunLocalReconstruction)) {
765         if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
766       }
767     }
768
769   
770     esd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
771     hltesd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
772     esd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
773     hltesd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
774     
775     // Set magnetic field from the tracker
776     esd->SetMagneticField(AliTracker::GetBz());
777     hltesd->SetMagneticField(AliTracker::GetBz());
778
779     
780     
781     // Fill raw-data error log into the ESD
782     if (fRawReader) FillRawDataErrorLog(iEvent,esd);
783
784     // vertex finder
785     if (fRunVertexFinder) {
786       if (!ReadESD(esd, "vertex")) {
787         if (!RunVertexFinder(esd)) {
788           if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
789         }
790         if (fCheckPointLevel > 0) WriteESD(esd, "vertex");
791       }
792     }
793
794     // HLT tracking
795     if (!fRunTracking.IsNull()) {
796       if (fRunHLTTracking) {
797         hltesd->SetVertex(esd->GetVertex());
798         if (!RunHLTTracking(hltesd)) {
799           if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
800         }
801       }
802     }
803
804     // Muon tracking
805     if (!fRunTracking.IsNull()) {
806       if (fRunMuonTracking) {
807         if (!RunMuonTracking(esd)) {
808           if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
809         }
810       }
811     }
812
813     // barrel tracking
814     if (!fRunTracking.IsNull()) {
815       if (!ReadESD(esd, "tracking")) {
816         if (!RunTracking(esd)) {
817           if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
818         }
819         if (fCheckPointLevel > 0) WriteESD(esd, "tracking");
820       }
821     }
822
823     // fill ESD
824     if (!fFillESD.IsNull()) {
825       if (!FillESD(esd, fFillESD)) {
826         if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
827       }
828     }
829   
830 //    if (!fFillESD.IsNull()) 
831 //    RunQA(fFillESD.Data(), esd);
832
833     // fill Event header information from the RawEventHeader
834     if (fRawReader){FillRawEventHeaderESD(esd);}
835
836     // combined PID
837     AliESDpid::MakePID(esd);
838     if (fCheckPointLevel > 1) WriteESD(esd, "PID");
839
840     if (fFillTriggerESD) {
841       if (!ReadESD(esd, "trigger")) {
842         if (!FillTriggerESD(esd)) {
843           if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
844         }
845         if (fCheckPointLevel > 1) WriteESD(esd, "trigger");
846       }
847     }
848
849     file->cd();
850
851     //Try to improve the reconstructed primary vertex position using the tracks
852     AliESDVertex *pvtx=0;
853     Bool_t dovertex=kTRUE;
854     TObject* obj = fOptions.FindObject("ITS");
855     if (obj) {
856       TString optITS = obj->GetTitle();
857       if (optITS.Contains("cosmics") || optITS.Contains("COSMICS")) 
858         dovertex=kFALSE;
859     }
860     if(dovertex) pvtx=tVertexer.FindPrimaryVertex(esd);
861     if(fDiamondProfile) esd->SetDiamond(fDiamondProfile);
862     
863     if (pvtx)
864     if (pvtx->GetStatus()) {
865        // Store the improved primary vertex
866        esd->SetPrimaryVertex(pvtx);
867        // Propagate the tracks to the DCA to the improved primary vertex
868        Double_t somethingbig = 777.;
869        Double_t bz = esd->GetMagneticField();
870        Int_t nt=esd->GetNumberOfTracks();
871        while (nt--) {
872          AliESDtrack *t = esd->GetTrack(nt);
873          t->RelateToVertex(pvtx, bz, somethingbig);
874        } 
875     }
876
877     if (fRunV0Finder) {
878        // V0 finding
879        AliV0vertexer vtxer;
880        vtxer.Tracks2V0vertices(esd);
881
882        if (fRunCascadeFinder) {
883           // Cascade finding
884           AliCascadeVertexer cvtxer;
885           cvtxer.V0sTracks2CascadeVertices(esd);
886        }
887     }
888  
889     // write ESD
890     if (fCleanESD) CleanESD(esd);
891     if (fWriteESDfriend) {
892       esdf->~AliESDfriend();
893       new (esdf) AliESDfriend(); // Reset...
894       esd->GetESDfriend(esdf);
895     }
896     tree->Fill();
897
898     // write HLT ESD
899     hlttree->Fill();
900
901     if (fCheckPointLevel > 0)  WriteESD(esd, "final"); 
902     esd->Reset();
903     hltesd->Reset();
904     if (fWriteESDfriend) {
905       esdf->~AliESDfriend();
906       new (esdf) AliESDfriend(); // Reset...
907     }
908     // esdf->Reset();
909     // delete esdf; esdf = 0;
910   // ESD QA 
911  
912     gSystem->GetProcInfo(&ProcInfo);
913     AliInfo(Form("Event %d -> Current memory usage %d %d",iEvent, ProcInfo.fMemResident, ProcInfo.fMemVirtual));
914   }
915   
916   detStr = fFillESD ; 
917 //  // write quality assurance ESDs data (one entry for all events)
918 //  for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
919 //      if (!IsSelected(fgkDetectorName[iDet], detStr)) 
920 //              continue;
921 //    AliQADataMaker * qadm = GetQADataMaker(iDet);
922 //    if (!qadm) continue;
923 //    qadm->EndOfCycle(AliQA::kRECPOINTS);
924 //    qadm->EndOfCycle(AliQA::kESDS);
925 //    qadm->Finish();
926 //  }
927
928   tree->GetUserInfo()->Add(esd);
929   hlttree->GetUserInfo()->Add(hltesd);
930   
931   const TMap *cdbMap = AliCDBManager::Instance()->GetStorageMap();       
932   const TList *cdbList = AliCDBManager::Instance()->GetRetrievedIds();   
933                  
934    TMap *cdbMapCopy = new TMap(cdbMap->GetEntries());    
935    cdbMapCopy->SetOwner(1);      
936    cdbMapCopy->SetName("cdbMap");        
937    TIter iter(cdbMap->GetTable());       
938          
939    TPair* pair = 0;      
940    while((pair = dynamic_cast<TPair*> (iter.Next()))){   
941          TObjString* keyStr = dynamic_cast<TObjString*> (pair->Key());   
942          TObjString* valStr = dynamic_cast<TObjString*> (pair->Value());         
943          cdbMapCopy->Add(new TObjString(keyStr->GetName()), new TObjString(valStr->GetName()));  
944    }     
945          
946    TList *cdbListCopy = new TList();     
947    cdbListCopy->SetOwner(1);     
948    cdbListCopy->SetName("cdbList");      
949          
950    TIter iter2(cdbList);         
951          
952    AliCDBId* id=0;       
953    while((id = dynamic_cast<AliCDBId*> (iter2.Next()))){         
954          cdbListCopy->Add(id->Clone());  
955    }     
956          
957    tree->GetUserInfo()->Add(cdbMapCopy);         
958    tree->GetUserInfo()->Add(cdbListCopy);
959
960
961   if(fESDPar.Contains("ESD.par")){
962     AliInfo("Attaching ESD.par to Tree");
963     TNamed *fn = CopyFileToTNamed(fESDPar.Data(),"ESD.par");
964     tree->GetUserInfo()->Add(fn);
965   }
966
967
968   file->cd();
969   if (fWriteESDfriend)
970     tree->SetBranchStatus("ESDfriend*",0);
971   // we want to have only one tree version number
972   tree->Write(tree->GetName(),TObject::kOverwrite);
973   hlttree->Write();
974
975   if (fWriteAOD) {
976     TFile *aodFile = TFile::Open("AliAOD.root", "RECREATE");
977     ESDFile2AODFile(file, aodFile);
978     aodFile->Close();
979   }
980
981   gROOT->cd();
982   CleanUp(file, fileOld);
983   
984   // Create tags for the events in the ESD tree (the ESD tree is always present)
985   // In case of empty events the tags will contain dummy values
986   AliESDTagCreator *esdtagCreator = new AliESDTagCreator();
987   esdtagCreator->CreateESDTags(fFirstEvent,fLastEvent,fGRPList);
988   if (fWriteAOD) {
989     AliAODTagCreator *aodtagCreator = new AliAODTagCreator();
990     aodtagCreator->CreateAODTags(fFirstEvent,fLastEvent,fGRPList);
991   }
992
993   //QA 
994   if ( fRunQA ) {
995         qas.Run(fRunLocalReconstruction.Data(), AliQA::kRECPOINTS) ;
996 //      qas.Reset() ;
997         qas.Run(fRunTracking.Data(), AliQA::kESDS) ;
998   }
999   return kTRUE;
1000 }
1001
1002
1003 //_____________________________________________________________________________
1004 Bool_t AliReconstruction::RunLocalReconstruction(const TString& detectors)
1005 {
1006 // run the local reconstruction
1007   static Int_t eventNr=0;
1008   AliCodeTimerAuto("")
1009
1010  //  AliCDBManager* man = AliCDBManager::Instance();
1011 //   Bool_t origCache = man->GetCacheFlag();
1012
1013 //   TString detStr = detectors;
1014 //   for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1015 //     if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1016 //     AliReconstructor* reconstructor = GetReconstructor(iDet);
1017 //     if (!reconstructor) continue;
1018 //     if (reconstructor->HasLocalReconstruction()) continue;
1019
1020 //     AliCodeTimerStart(Form("running reconstruction for %s", fgkDetectorName[iDet]));
1021 //     AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
1022     
1023 //     AliCodeTimerStart(Form("Loading calibration data from OCDB for %s", fgkDetectorName[iDet]));                          
1024 //     AliInfo(Form("Loading calibration data from OCDB for %s", fgkDetectorName[iDet]));
1025
1026 //     man->SetCacheFlag(kTRUE);
1027 //     TString calibPath = Form("%s/Calib/*", fgkDetectorName[iDet]);
1028 //     man->GetAll(calibPath); // entries are cached!
1029
1030 //     AliCodeTimerStop(Form("Loading calibration data from OCDB for %s", fgkDetectorName[iDet]));
1031      
1032 //     if (fRawReader) {
1033 //       fRawReader->RewindEvents();
1034 //       reconstructor->Reconstruct(fRunLoader, fRawReader);
1035 //     } else {
1036 //       reconstructor->Reconstruct(fRunLoader);
1037 //     }
1038      
1039 //      AliCodeTimerStop(Form("running reconstruction for %s", fgkDetectorName[iDet]));
1040     // AliSysInfo::AddStamp(Form("LRec%s_%d",fgkDetectorName[iDet],eventNr));
1041
1042 //     // unload calibration data
1043 //     man->UnloadFromCache(calibPath);
1044 //     //man->ClearCache();
1045 //   }
1046
1047 //   man->SetCacheFlag(origCache);
1048
1049 //   if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1050 //     AliError(Form("the following detectors were not found: %s",
1051 //                   detStr.Data()));
1052 //     if (fStopOnError) return kFALSE;
1053 //   }
1054
1055           eventNr++;
1056   return kTRUE;
1057 }
1058
1059 //_____________________________________________________________________________
1060 Bool_t AliReconstruction::RunLocalEventReconstruction(const TString& detectors)
1061 {
1062 // run the local reconstruction
1063   static Int_t eventNr=0;
1064   AliCodeTimerAuto("")
1065
1066   TString detStr = detectors;
1067   for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1068     if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1069     AliReconstructor* reconstructor = GetReconstructor(iDet);
1070     if (!reconstructor) continue;
1071     AliLoader* loader = fLoader[iDet];
1072     if (!loader) {
1073       AliWarning(Form("No loader is defined for %s!",fgkDetectorName[iDet]));
1074       continue;
1075     }
1076
1077     // conversion of digits
1078     if (fRawReader && reconstructor->HasDigitConversion()) {
1079       AliInfo(Form("converting raw data digits into root objects for %s", 
1080                    fgkDetectorName[iDet]));
1081       AliCodeTimerAuto(Form("converting raw data digits into root objects for %s", 
1082                             fgkDetectorName[iDet]));
1083       loader->LoadDigits("update");
1084       loader->CleanDigits();
1085       loader->MakeDigitsContainer();
1086       TTree* digitsTree = loader->TreeD();
1087       reconstructor->ConvertDigits(fRawReader, digitsTree);
1088       loader->WriteDigits("OVERWRITE");
1089       loader->UnloadDigits();
1090     }
1091
1092     // local reconstruction
1093     AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
1094     AliCodeTimerAuto(Form("running reconstruction for %s", fgkDetectorName[iDet]));
1095     loader->LoadRecPoints("update");
1096     loader->CleanRecPoints();
1097     loader->MakeRecPointsContainer();
1098     TTree* clustersTree = loader->TreeR();
1099     if (fRawReader && !reconstructor->HasDigitConversion()) {
1100       reconstructor->Reconstruct(fRawReader, clustersTree);
1101     } else {
1102       loader->LoadDigits("read");
1103       TTree* digitsTree = loader->TreeD();
1104       if (!digitsTree) {
1105         AliError(Form("Can't get the %s digits tree", fgkDetectorName[iDet]));
1106         if (fStopOnError) return kFALSE;
1107       } else {
1108         reconstructor->Reconstruct(digitsTree, clustersTree);
1109       }
1110       loader->UnloadDigits();
1111     }
1112
1113 //    AliQADataMaker * qadm = GetQADataMaker(iDet);
1114 //    if (qadm) {
1115 //      AliCodeTimerStart(Form("running quality assurance data maker for %s", fgkDetectorName[iDet]));
1116 //      AliInfo(Form("running quality assurance data maker for %s", fgkDetectorName[iDet]));
1117 //        
1118 //     if (qadm->IsCycleDone() ) {
1119 //      qadm->EndOfCycle(AliQA::kRECPOINTS) ; 
1120 //        qadm->EndOfCycle(AliQA::kESDS) ; 
1121 //      qadm->StartOfCycle(AliQA::kRECPOINTS) ; 
1122 //        qadm->StartOfCycle(AliQA::kESDS, "same") ; 
1123 //     }
1124 //      qadm->Exec(AliQA::kRECPOINTS, clustersTree) ; 
1125 //      AliCodeTimerStop(Form("running quality assurance data maker for %s", fgkDetectorName[iDet]));
1126 //    }
1127
1128     loader->WriteRecPoints("OVERWRITE");
1129     loader->UnloadRecPoints();
1130     AliSysInfo::AddStamp(Form("LRec%s_%d",fgkDetectorName[iDet],eventNr), iDet,1,eventNr);
1131   }
1132
1133   if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1134     AliError(Form("the following detectors were not found: %s",
1135                   detStr.Data()));
1136     if (fStopOnError) return kFALSE;
1137   }
1138   eventNr++;
1139   return kTRUE;
1140 }
1141
1142 //_____________________________________________________________________________
1143 Bool_t AliReconstruction::RunVertexFinder(AliESDEvent*& esd)
1144 {
1145 // run the barrel tracking
1146
1147   AliCodeTimerAuto("")
1148
1149   AliESDVertex* vertex = NULL;
1150   Double_t vtxPos[3] = {0, 0, 0};
1151   Double_t vtxErr[3] = {0.07, 0.07, 0.1};
1152   TArrayF mcVertex(3); 
1153   if (fRunLoader->GetHeader() && fRunLoader->GetHeader()->GenEventHeader()) {
1154     fRunLoader->GetHeader()->GenEventHeader()->PrimaryVertex(mcVertex);
1155     for (Int_t i = 0; i < 3; i++) vtxPos[i] = mcVertex[i];
1156   }
1157
1158   if (fVertexer) {
1159     if(fDiamondProfile) fVertexer->SetVtxStart(fDiamondProfile);
1160     AliInfo("running the ITS vertex finder");
1161     if (fLoader[0]) fLoader[0]->LoadRecPoints();
1162     vertex = fVertexer->FindVertexForCurrentEvent(fRunLoader->GetEventNumber());
1163     if (fLoader[0]) fLoader[0]->UnloadRecPoints();
1164     if(!vertex){
1165       AliWarning("Vertex not found");
1166       vertex = new AliESDVertex();
1167       vertex->SetName("default");
1168     }
1169     else {
1170       vertex->SetName("reconstructed");
1171     }
1172
1173   } else {
1174     AliInfo("getting the primary vertex from MC");
1175     vertex = new AliESDVertex(vtxPos, vtxErr);
1176   }
1177
1178   if (vertex) {
1179     vertex->GetXYZ(vtxPos);
1180     vertex->GetSigmaXYZ(vtxErr);
1181   } else {
1182     AliWarning("no vertex reconstructed");
1183     vertex = new AliESDVertex(vtxPos, vtxErr);
1184   }
1185   esd->SetVertex(vertex);
1186   // if SPD multiplicity has been determined, it is stored in the ESD
1187   AliMultiplicity *mult = fVertexer->GetMultiplicity();
1188   if(mult)esd->SetMultiplicity(mult);
1189
1190   for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1191     if (fTracker[iDet]) fTracker[iDet]->SetVertex(vtxPos, vtxErr);
1192   }  
1193   delete vertex;
1194
1195   return kTRUE;
1196 }
1197
1198 //_____________________________________________________________________________
1199 Bool_t AliReconstruction::RunHLTTracking(AliESDEvent*& esd)
1200 {
1201 // run the HLT barrel tracking
1202
1203   AliCodeTimerAuto("")
1204
1205   if (!fRunLoader) {
1206     AliError("Missing runLoader!");
1207     return kFALSE;
1208   }
1209
1210   AliInfo("running HLT tracking");
1211
1212   // Get a pointer to the HLT reconstructor
1213   AliReconstructor *reconstructor = GetReconstructor(fgkNDetectors-1);
1214   if (!reconstructor) return kFALSE;
1215
1216   // TPC + ITS
1217   for (Int_t iDet = 1; iDet >= 0; iDet--) {
1218     TString detName = fgkDetectorName[iDet];
1219     AliDebug(1, Form("%s HLT tracking", detName.Data()));
1220     reconstructor->SetOption(detName.Data());
1221     AliTracker *tracker = reconstructor->CreateTracker();
1222     if (!tracker) {
1223       AliWarning(Form("couldn't create a HLT tracker for %s", detName.Data()));
1224       if (fStopOnError) return kFALSE;
1225       continue;
1226     }
1227     Double_t vtxPos[3];
1228     Double_t vtxErr[3]={0.005,0.005,0.010};
1229     const AliESDVertex *vertex = esd->GetVertex();
1230     vertex->GetXYZ(vtxPos);
1231     tracker->SetVertex(vtxPos,vtxErr);
1232     if(iDet != 1) {
1233       fLoader[iDet]->LoadRecPoints("read");
1234       TTree* tree = fLoader[iDet]->TreeR();
1235       if (!tree) {
1236         AliError(Form("Can't get the %s cluster tree", detName.Data()));
1237         return kFALSE;
1238       }
1239       tracker->LoadClusters(tree);
1240     }
1241     if (tracker->Clusters2Tracks(esd) != 0) {
1242       AliError(Form("HLT %s Clusters2Tracks failed", fgkDetectorName[iDet]));
1243       return kFALSE;
1244     }
1245     if(iDet != 1) {
1246       tracker->UnloadClusters();
1247     }
1248     delete tracker;
1249   }
1250
1251   return kTRUE;
1252 }
1253
1254 //_____________________________________________________________________________
1255 Bool_t AliReconstruction::RunMuonTracking(AliESDEvent*& esd)
1256 {
1257 // run the muon spectrometer tracking
1258
1259   AliCodeTimerAuto("")
1260
1261   if (!fRunLoader) {
1262     AliError("Missing runLoader!");
1263     return kFALSE;
1264   }
1265   Int_t iDet = 7; // for MUON
1266
1267   AliInfo("is running...");
1268
1269   // Get a pointer to the MUON reconstructor
1270   AliReconstructor *reconstructor = GetReconstructor(iDet);
1271   if (!reconstructor) return kFALSE;
1272
1273   
1274   TString detName = fgkDetectorName[iDet];
1275   AliDebug(1, Form("%s tracking", detName.Data()));
1276   AliTracker *tracker =  reconstructor->CreateTracker();
1277   if (!tracker) {
1278     AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
1279     return kFALSE;
1280   }
1281      
1282   // create Tracks
1283   fLoader[iDet]->LoadTracks("update");
1284   fLoader[iDet]->CleanTracks();
1285   fLoader[iDet]->MakeTracksContainer();
1286
1287   // read RecPoints
1288   fLoader[iDet]->LoadRecPoints("read");  
1289   tracker->LoadClusters(fLoader[iDet]->TreeR());
1290   
1291   Int_t rv = tracker->Clusters2Tracks(esd);
1292   
1293   fLoader[iDet]->UnloadRecPoints();
1294   
1295   if ( rv )
1296   {
1297     AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
1298     return kFALSE;
1299   }
1300   
1301   tracker->UnloadClusters();
1302   
1303   fLoader[iDet]->UnloadRecPoints();
1304
1305   fLoader[iDet]->WriteTracks("OVERWRITE");
1306   fLoader[iDet]->UnloadTracks();
1307
1308   delete tracker;
1309   
1310   return kTRUE;
1311 }
1312
1313
1314 //_____________________________________________________________________________
1315 Bool_t AliReconstruction::RunTracking(AliESDEvent*& esd)
1316 {
1317 // run the barrel tracking
1318   static Int_t eventNr=0;
1319   AliCodeTimerAuto("")
1320
1321   AliInfo("running tracking");
1322
1323   //Fill the ESD with the T0 info (will be used by the TOF) 
1324   if (fReconstructor[11] && fLoader[11]) {
1325     fLoader[11]->LoadRecPoints("READ");
1326     TTree *treeR = fLoader[11]->TreeR();
1327     GetReconstructor(11)->FillESD((TTree *)NULL,treeR,esd);
1328   }
1329
1330   // pass 1: TPC + ITS inwards
1331   for (Int_t iDet = 1; iDet >= 0; iDet--) {
1332     if (!fTracker[iDet]) continue;
1333     AliDebug(1, Form("%s tracking", fgkDetectorName[iDet]));
1334
1335     // load clusters
1336     fLoader[iDet]->LoadRecPoints("read");
1337     AliSysInfo::AddStamp(Form("RLoadCluster%s_%d",fgkDetectorName[iDet],eventNr),iDet,1, eventNr);
1338     TTree* tree = fLoader[iDet]->TreeR();
1339     if (!tree) {
1340       AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
1341       return kFALSE;
1342     }
1343     fTracker[iDet]->LoadClusters(tree);
1344     AliSysInfo::AddStamp(Form("TLoadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,2, eventNr);
1345     // run tracking
1346     if (fTracker[iDet]->Clusters2Tracks(esd) != 0) {
1347       AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
1348       return kFALSE;
1349     }
1350     if (fCheckPointLevel > 1) {
1351       WriteESD(esd, Form("%s.tracking", fgkDetectorName[iDet]));
1352     }
1353     // preliminary PID in TPC needed by the ITS tracker
1354     if (iDet == 1) {
1355       GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
1356       AliESDpid::MakePID(esd);
1357     } 
1358     AliSysInfo::AddStamp(Form("Tracking0%s_%d",fgkDetectorName[iDet],eventNr), iDet,3,eventNr);
1359   }
1360
1361   // pass 2: ALL backwards
1362   for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1363     if (!fTracker[iDet]) continue;
1364     AliDebug(1, Form("%s back propagation", fgkDetectorName[iDet]));
1365
1366     // load clusters
1367     if (iDet > 1) {     // all except ITS, TPC
1368       TTree* tree = NULL;
1369       fLoader[iDet]->LoadRecPoints("read");
1370       AliSysInfo::AddStamp(Form("RLoadCluster0%s_%d",fgkDetectorName[iDet],eventNr), iDet,1, eventNr);
1371       tree = fLoader[iDet]->TreeR();
1372       if (!tree) {
1373         AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
1374         return kFALSE;
1375       }
1376       fTracker[iDet]->LoadClusters(tree); 
1377       AliSysInfo::AddStamp(Form("TLoadCluster0%s_%d",fgkDetectorName[iDet],eventNr), iDet,2, eventNr);
1378     }
1379
1380     // run tracking
1381     if (fTracker[iDet]->PropagateBack(esd) != 0) {
1382       AliError(Form("%s backward propagation failed", fgkDetectorName[iDet]));
1383       //      return kFALSE;
1384     }
1385     if (fCheckPointLevel > 1) {
1386       WriteESD(esd, Form("%s.back", fgkDetectorName[iDet]));
1387     }
1388
1389     // unload clusters
1390     if (iDet > 2) {     // all except ITS, TPC, TRD
1391       fTracker[iDet]->UnloadClusters();
1392       fLoader[iDet]->UnloadRecPoints();
1393     }
1394     // updated PID in TPC needed by the ITS tracker -MI
1395     if (iDet == 1) {
1396       GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
1397       AliESDpid::MakePID(esd);
1398     }
1399     AliSysInfo::AddStamp(Form("Tracking1%s_%d",fgkDetectorName[iDet],eventNr), iDet,3, eventNr);
1400   }
1401
1402   // write space-points to the ESD in case alignment data output
1403   // is switched on
1404   if (fWriteAlignmentData)
1405     WriteAlignmentData(esd);
1406
1407   // pass 3: TRD + TPC + ITS refit inwards
1408   for (Int_t iDet = 2; iDet >= 0; iDet--) {
1409     if (!fTracker[iDet]) continue;
1410     AliDebug(1, Form("%s inward refit", fgkDetectorName[iDet]));
1411
1412     // run tracking
1413     if (fTracker[iDet]->RefitInward(esd) != 0) {
1414       AliError(Form("%s inward refit failed", fgkDetectorName[iDet]));
1415       //      return kFALSE;
1416     }
1417     if (fCheckPointLevel > 1) {
1418       WriteESD(esd, Form("%s.refit", fgkDetectorName[iDet]));
1419     }
1420     AliSysInfo::AddStamp(Form("Tracking2%s_%d",fgkDetectorName[iDet],eventNr), iDet,3, eventNr);
1421     // unload clusters
1422     fTracker[iDet]->UnloadClusters();
1423     AliSysInfo::AddStamp(Form("TUnloadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,4, eventNr);
1424     fLoader[iDet]->UnloadRecPoints();
1425     AliSysInfo::AddStamp(Form("RUnloadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,5, eventNr);
1426   }
1427   //
1428   // Propagate track to the vertex - if not done by ITS
1429   //
1430   Int_t ntracks = esd->GetNumberOfTracks();
1431   for (Int_t itrack=0; itrack<ntracks; itrack++){
1432     const Double_t kRadius  = 3;   // beam pipe radius
1433     const Double_t kMaxStep = 5;   // max step
1434     const Double_t kMaxD    = 123456;  // max distance to prim vertex
1435     Double_t       fieldZ   = AliTracker::GetBz();  //
1436     AliESDtrack * track = esd->GetTrack(itrack);
1437     if (!track) continue;
1438     if (track->IsOn(AliESDtrack::kITSrefit)) continue;
1439    AliTracker::PropagateTrackTo(track,kRadius,track->GetMass(),kMaxStep,kTRUE);
1440     track->RelateToVertex(esd->GetVertex(),fieldZ, kMaxD);
1441   }
1442   eventNr++;
1443   return kTRUE;
1444 }
1445
1446 //_____________________________________________________________________________
1447 Bool_t AliReconstruction::CleanESD(AliESDEvent *esd){
1448   //
1449   // Remove the data which are not needed for the physics analysis.
1450   //
1451
1452   AliInfo("Cleaning the ESD...");
1453   Int_t nTracks=esd->GetNumberOfTracks();
1454   AliInfo(Form("Number of ESD tracks before cleaning %d",nTracks));
1455   Int_t nV0s=esd->GetNumberOfV0s();
1456   AliInfo(Form("Number of ESD V0s before cleaning %d",nV0s));
1457
1458   Float_t cleanPars[]={fV0DCAmax,fV0CsPmin,fDmax,fZmax};
1459   Bool_t rc=esd->Clean(cleanPars);
1460
1461   nTracks=esd->GetNumberOfTracks();
1462   AliInfo(Form("Number of ESD tracks after cleaning %d",nTracks));
1463   nV0s=esd->GetNumberOfV0s();
1464   AliInfo(Form("Number of ESD V0s after cleaning %d",nV0s));
1465
1466   return rc;
1467 }
1468
1469 //_____________________________________________________________________________
1470 Bool_t AliReconstruction::FillESD(AliESDEvent*& esd, const TString& detectors)
1471 {
1472 // fill the event summary data
1473
1474   AliCodeTimerAuto("")
1475     static Int_t eventNr=0; 
1476   TString detStr = detectors;
1477   
1478   for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1479   if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1480     AliReconstructor* reconstructor = GetReconstructor(iDet);
1481     if (!reconstructor) continue;
1482     if (!ReadESD(esd, fgkDetectorName[iDet])) {
1483       AliDebug(1, Form("filling ESD for %s", fgkDetectorName[iDet]));
1484       TTree* clustersTree = NULL;
1485       if (fLoader[iDet]) {
1486         fLoader[iDet]->LoadRecPoints("read");
1487         clustersTree = fLoader[iDet]->TreeR();
1488         if (!clustersTree) {
1489           AliError(Form("Can't get the %s clusters tree", 
1490                         fgkDetectorName[iDet]));
1491           if (fStopOnError) return kFALSE;
1492         }
1493       }
1494       if (fRawReader && !reconstructor->HasDigitConversion()) {
1495         reconstructor->FillESD(fRawReader, clustersTree, esd);
1496       } else {
1497         TTree* digitsTree = NULL;
1498         if (fLoader[iDet]) {
1499           fLoader[iDet]->LoadDigits("read");
1500           digitsTree = fLoader[iDet]->TreeD();
1501           if (!digitsTree) {
1502             AliError(Form("Can't get the %s digits tree", 
1503                           fgkDetectorName[iDet]));
1504             if (fStopOnError) return kFALSE;
1505           }
1506         }
1507         reconstructor->FillESD(digitsTree, clustersTree, esd);
1508         if (fLoader[iDet]) fLoader[iDet]->UnloadDigits();
1509       }
1510       if (fLoader[iDet]) {
1511         fLoader[iDet]->UnloadRecPoints();
1512       }
1513
1514       if (fCheckPointLevel > 2) WriteESD(esd, fgkDetectorName[iDet]);
1515     }
1516   }
1517
1518   if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1519     AliError(Form("the following detectors were not found: %s", 
1520                   detStr.Data()));
1521     if (fStopOnError) return kFALSE;
1522   }
1523   AliSysInfo::AddStamp(Form("FillESD%d",eventNr), 0,1, eventNr);
1524   eventNr++;
1525   return kTRUE;
1526 }
1527
1528 //_____________________________________________________________________________
1529 Bool_t AliReconstruction::FillTriggerESD(AliESDEvent*& esd)
1530 {
1531   // Reads the trigger decision which is
1532   // stored in Trigger.root file and fills
1533   // the corresponding esd entries
1534
1535   AliCodeTimerAuto("")
1536   
1537   AliInfo("Filling trigger information into the ESD");
1538
1539   if (fRawReader) {
1540     AliCTPRawStream input(fRawReader);
1541     if (!input.Next()) {
1542       AliError("No valid CTP (trigger) DDL raw data is found ! The trigger information is not stored in the ESD !");
1543       return kFALSE;
1544     }
1545     esd->SetTriggerMask(input.GetClassMask());
1546     esd->SetTriggerCluster(input.GetClusterMask());
1547   }
1548   else {
1549     AliRunLoader *runloader = AliRunLoader::GetRunLoader();
1550     if (runloader) {
1551       if (!runloader->LoadTrigger()) {
1552         AliCentralTrigger *aCTP = runloader->GetTrigger();
1553         esd->SetTriggerMask(aCTP->GetClassMask());
1554         esd->SetTriggerCluster(aCTP->GetClusterMask());
1555       }
1556       else {
1557         AliWarning("No trigger can be loaded! The trigger information is not stored in the ESD !");
1558         return kFALSE;
1559       }
1560     }
1561     else {
1562       AliError("No run loader is available! The trigger information is not stored in the ESD !");
1563       return kFALSE;
1564     }
1565   }
1566
1567   return kTRUE;
1568 }
1569
1570
1571
1572
1573
1574 //_____________________________________________________________________________
1575 Bool_t AliReconstruction::FillRawEventHeaderESD(AliESDEvent*& esd)
1576 {
1577   // 
1578   // Filling information from RawReader Header
1579   // 
1580
1581   AliInfo("Filling information from RawReader Header");
1582   esd->SetBunchCrossNumber(0);
1583   esd->SetOrbitNumber(0);
1584   esd->SetPeriodNumber(0);
1585   esd->SetTimeStamp(0);
1586   esd->SetEventType(0);
1587   const AliRawEventHeaderBase * eventHeader = fRawReader->GetEventHeader();
1588   if (eventHeader){
1589
1590     const UInt_t *id = eventHeader->GetP("Id");
1591     esd->SetBunchCrossNumber((id)[1]&0x00000fff);
1592     esd->SetOrbitNumber((((id)[0]<<20)&0xf00000)|(((id)[1]>>12)&0xfffff));
1593     esd->SetPeriodNumber(((id)[0]>>4)&0x0fffffff);
1594
1595     esd->SetTimeStamp((eventHeader->Get("Timestamp")));  
1596     esd->SetEventType((eventHeader->Get("Type")));
1597   }
1598
1599   return kTRUE;
1600 }
1601
1602
1603 //_____________________________________________________________________________
1604 Bool_t AliReconstruction::IsSelected(TString detName, TString& detectors) const
1605 {
1606 // check whether detName is contained in detectors
1607 // if yes, it is removed from detectors
1608
1609   // check if all detectors are selected
1610   if ((detectors.CompareTo("ALL") == 0) ||
1611       detectors.BeginsWith("ALL ") ||
1612       detectors.EndsWith(" ALL") ||
1613       detectors.Contains(" ALL ")) {
1614     detectors = "ALL";
1615     return kTRUE;
1616   }
1617
1618   // search for the given detector
1619   Bool_t result = kFALSE;
1620   if ((detectors.CompareTo(detName) == 0) ||
1621       detectors.BeginsWith(detName+" ") ||
1622       detectors.EndsWith(" "+detName) ||
1623       detectors.Contains(" "+detName+" ")) {
1624     detectors.ReplaceAll(detName, "");
1625     result = kTRUE;
1626   }
1627
1628   // clean up the detectors string
1629   while (detectors.Contains("  ")) detectors.ReplaceAll("  ", " ");
1630   while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
1631   while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
1632
1633   return result;
1634 }
1635
1636 //_____________________________________________________________________________
1637 Bool_t AliReconstruction::InitRunLoader()
1638 {
1639 // get or create the run loader
1640
1641   if (gAlice) delete gAlice;
1642   gAlice = NULL;
1643
1644   if (!gSystem->AccessPathName(fGAliceFileName.Data())) { // galice.root exists
1645     // load all base libraries to get the loader classes
1646     TString libs = gSystem->GetLibraries();
1647     for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1648       TString detName = fgkDetectorName[iDet];
1649       if (detName == "HLT") continue;
1650       if (libs.Contains("lib" + detName + "base.so")) continue;
1651       gSystem->Load("lib" + detName + "base.so");
1652     }
1653     fRunLoader = AliRunLoader::Open(fGAliceFileName.Data());
1654     if (!fRunLoader) {
1655       AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
1656       CleanUp();
1657       return kFALSE;
1658     }
1659     fRunLoader->CdGAFile();
1660     if (gFile->GetKey(AliRunLoader::GetGAliceName())) {
1661       if (fRunLoader->LoadgAlice() == 0) {
1662         gAlice = fRunLoader->GetAliRun();
1663         AliTracker::SetFieldMap(gAlice->Field(),fUniformField);
1664       }
1665     }
1666     if (!gAlice && !fRawReader) {
1667       AliError(Form("no gAlice object found in file %s",
1668                     fGAliceFileName.Data()));
1669       CleanUp();
1670       return kFALSE;
1671     }
1672
1673     //PH This is a temporary fix to give access to the kinematics
1674     //PH that is needed for the labels of ITS clusters
1675     fRunLoader->LoadHeader();
1676     fRunLoader->LoadKinematics();
1677
1678   } else {               // galice.root does not exist
1679     if (!fRawReader) {
1680       AliError(Form("the file %s does not exist", fGAliceFileName.Data()));
1681       CleanUp();
1682       return kFALSE;
1683     }
1684     fRunLoader = AliRunLoader::Open(fGAliceFileName.Data(),
1685                                     AliConfig::GetDefaultEventFolderName(),
1686                                     "recreate");
1687     if (!fRunLoader) {
1688       AliError(Form("could not create run loader in file %s", 
1689                     fGAliceFileName.Data()));
1690       CleanUp();
1691       return kFALSE;
1692     }
1693     fRunLoader->MakeTree("E");
1694     Int_t iEvent = 0;
1695     while (fRawReader->NextEvent()) {
1696       fRunLoader->SetEventNumber(iEvent);
1697       fRunLoader->GetHeader()->Reset(fRawReader->GetRunNumber(), 
1698                                      iEvent, iEvent);
1699       fRunLoader->MakeTree("H");
1700       fRunLoader->TreeE()->Fill();
1701       iEvent++;
1702     }
1703     fRawReader->RewindEvents();
1704     if (fNumberOfEventsPerFile > 0)
1705       fRunLoader->SetNumberOfEventsPerFile(fNumberOfEventsPerFile);
1706     else
1707       fRunLoader->SetNumberOfEventsPerFile(iEvent);
1708     fRunLoader->WriteHeader("OVERWRITE");
1709     fRunLoader->CdGAFile();
1710     fRunLoader->Write(0, TObject::kOverwrite);
1711 //    AliTracker::SetFieldMap(???);
1712   }
1713
1714   return kTRUE;
1715 }
1716
1717 //_____________________________________________________________________________
1718 AliReconstructor* AliReconstruction::GetReconstructor(Int_t iDet)
1719 {
1720 // get the reconstructor object and the loader for a detector
1721
1722   if (fReconstructor[iDet]) return fReconstructor[iDet];
1723
1724   // load the reconstructor object
1725   TPluginManager* pluginManager = gROOT->GetPluginManager();
1726   TString detName = fgkDetectorName[iDet];
1727   TString recName = "Ali" + detName + "Reconstructor";
1728   if (gAlice && !gAlice->GetDetector(detName) && (detName != "HLT")) return NULL;
1729
1730   AliReconstructor* reconstructor = NULL;
1731   // first check if a plugin is defined for the reconstructor
1732   TPluginHandler* pluginHandler = 
1733     pluginManager->FindHandler("AliReconstructor", detName);
1734   // if not, add a plugin for it
1735   if (!pluginHandler) {
1736     AliDebug(1, Form("defining plugin for %s", recName.Data()));
1737     TString libs = gSystem->GetLibraries();
1738     if (libs.Contains("lib" + detName + "base.so") ||
1739         (gSystem->Load("lib" + detName + "base.so") >= 0)) {
1740       pluginManager->AddHandler("AliReconstructor", detName, 
1741                                 recName, detName + "rec", recName + "()");
1742     } else {
1743       pluginManager->AddHandler("AliReconstructor", detName, 
1744                                 recName, detName, recName + "()");
1745     }
1746     pluginHandler = pluginManager->FindHandler("AliReconstructor", detName);
1747   }
1748   if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
1749     reconstructor = (AliReconstructor*) pluginHandler->ExecPlugin(0);
1750   }
1751   if (reconstructor) {
1752     TObject* obj = fOptions.FindObject(detName.Data());
1753     if (obj) reconstructor->SetOption(obj->GetTitle());
1754     reconstructor->Init();
1755     fReconstructor[iDet] = reconstructor;
1756   }
1757
1758   // get or create the loader
1759   if (detName != "HLT") {
1760     fLoader[iDet] = fRunLoader->GetLoader(detName + "Loader");
1761     if (!fLoader[iDet]) {
1762       AliConfig::Instance()
1763         ->CreateDetectorFolders(fRunLoader->GetEventFolder(), 
1764                                 detName, detName);
1765       // first check if a plugin is defined for the loader
1766       pluginHandler = 
1767         pluginManager->FindHandler("AliLoader", detName);
1768       // if not, add a plugin for it
1769       if (!pluginHandler) {
1770         TString loaderName = "Ali" + detName + "Loader";
1771         AliDebug(1, Form("defining plugin for %s", loaderName.Data()));
1772         pluginManager->AddHandler("AliLoader", detName, 
1773                                   loaderName, detName + "base", 
1774                                   loaderName + "(const char*, TFolder*)");
1775         pluginHandler = pluginManager->FindHandler("AliLoader", detName);
1776       }
1777       if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
1778         fLoader[iDet] = 
1779           (AliLoader*) pluginHandler->ExecPlugin(2, detName.Data(), 
1780                                                  fRunLoader->GetEventFolder());
1781       }
1782       if (!fLoader[iDet]) {   // use default loader
1783         fLoader[iDet] = new AliLoader(detName, fRunLoader->GetEventFolder());
1784       }
1785       if (!fLoader[iDet]) {
1786         AliWarning(Form("couldn't get loader for %s", detName.Data()));
1787         if (fStopOnError) return NULL;
1788       } else {
1789         fRunLoader->AddLoader(fLoader[iDet]);
1790         fRunLoader->CdGAFile();
1791         if (gFile && !gFile->IsWritable()) gFile->ReOpen("UPDATE");
1792         fRunLoader->Write(0, TObject::kOverwrite);
1793       }
1794     }
1795   }
1796       
1797   return reconstructor;
1798 }
1799
1800 //_____________________________________________________________________________
1801 Bool_t AliReconstruction::CreateVertexer()
1802 {
1803 // create the vertexer
1804
1805   fVertexer = NULL;
1806   AliReconstructor* itsReconstructor = GetReconstructor(0);
1807   if (itsReconstructor) {
1808     fVertexer = itsReconstructor->CreateVertexer();
1809   }
1810   if (!fVertexer) {
1811     AliWarning("couldn't create a vertexer for ITS");
1812     if (fStopOnError) return kFALSE;
1813   }
1814
1815   return kTRUE;
1816 }
1817
1818 //_____________________________________________________________________________
1819 Bool_t AliReconstruction::CreateTrackers(const TString& detectors)
1820 {
1821 // create the trackers
1822
1823   TString detStr = detectors;
1824   for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1825     if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1826     AliReconstructor* reconstructor = GetReconstructor(iDet);
1827     if (!reconstructor) continue;
1828     TString detName = fgkDetectorName[iDet];
1829     if (detName == "HLT") {
1830       fRunHLTTracking = kTRUE;
1831       continue;
1832     }
1833     if (detName == "MUON") {
1834       fRunMuonTracking = kTRUE;
1835       continue;
1836     }
1837
1838
1839     fTracker[iDet] = reconstructor->CreateTracker();
1840     if (!fTracker[iDet] && (iDet < 7)) {
1841       AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
1842       if (fStopOnError) return kFALSE;
1843     }
1844     AliSysInfo::AddStamp(Form("LTracker%s",fgkDetectorName[iDet]), iDet,0);
1845   }
1846
1847   return kTRUE;
1848 }
1849
1850 //_____________________________________________________________________________
1851 void AliReconstruction::CleanUp(TFile* file, TFile* fileOld)
1852 {
1853 // delete trackers and the run loader and close and delete the file
1854
1855   for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1856     delete fReconstructor[iDet];
1857     fReconstructor[iDet] = NULL;
1858     fLoader[iDet] = NULL;
1859     delete fTracker[iDet];
1860     fTracker[iDet] = NULL;
1861 //    delete fQADataMaker[iDet];
1862 //    fQADataMaker[iDet] = NULL;
1863   }
1864   delete fVertexer;
1865   fVertexer = NULL;
1866   delete fDiamondProfile;
1867   fDiamondProfile = NULL;
1868
1869   delete fGRPList;
1870   fGRPList = NULL;
1871
1872   delete fRunLoader;
1873   fRunLoader = NULL;
1874   delete fRawReader;
1875   fRawReader = NULL;
1876
1877   if (file) {
1878     file->Close();
1879     delete file;
1880   }
1881
1882   if (fileOld) {
1883     fileOld->Close();
1884     delete fileOld;
1885     gSystem->Unlink("AliESDs.old.root");
1886   }
1887 }
1888
1889 //_____________________________________________________________________________
1890
1891 Bool_t AliReconstruction::ReadESD(AliESDEvent*& esd, const char* recStep) const
1892 {
1893 // read the ESD event from a file
1894
1895   if (!esd) return kFALSE;
1896   char fileName[256];
1897   sprintf(fileName, "ESD_%d.%d_%s.root", 
1898           esd->GetRunNumber(), esd->GetEventNumberInFile(), recStep);
1899   if (gSystem->AccessPathName(fileName)) return kFALSE;
1900
1901   AliInfo(Form("reading ESD from file %s", fileName));
1902   AliDebug(1, Form("reading ESD from file %s", fileName));
1903   TFile* file = TFile::Open(fileName);
1904   if (!file || !file->IsOpen()) {
1905     AliError(Form("opening %s failed", fileName));
1906     delete file;
1907     return kFALSE;
1908   }
1909
1910   gROOT->cd();
1911   delete esd;
1912   esd = (AliESDEvent*) file->Get("ESD");
1913   file->Close();
1914   delete file;
1915   return kTRUE;
1916
1917 }
1918
1919
1920
1921 //_____________________________________________________________________________
1922 void AliReconstruction::WriteESD(AliESDEvent* esd, const char* recStep) const
1923 {
1924 // write the ESD event to a file
1925
1926   if (!esd) return;
1927   char fileName[256];
1928   sprintf(fileName, "ESD_%d.%d_%s.root", 
1929           esd->GetRunNumber(), esd->GetEventNumberInFile(), recStep);
1930
1931   AliDebug(1, Form("writing ESD to file %s", fileName));
1932   TFile* file = TFile::Open(fileName, "recreate");
1933   if (!file || !file->IsOpen()) {
1934     AliError(Form("opening %s failed", fileName));
1935   } else {
1936     esd->Write("ESD");
1937     file->Close();
1938   }
1939   delete file;
1940 }
1941
1942
1943
1944
1945
1946 //_____________________________________________________________________________
1947 void AliReconstruction::ESDFile2AODFile(TFile* esdFile, TFile* aodFile)
1948 {
1949   // write all files from the given esd file to an aod file
1950
1951   // create an AliAOD object 
1952   AliAODEvent *aod = new AliAODEvent();
1953   aod->CreateStdContent();
1954   
1955   // go to the file
1956   aodFile->cd();
1957   
1958   // create the tree
1959   TTree *aodTree = new TTree("aodTree", "AliAOD tree");
1960   aodTree->Branch(aod->GetList());
1961
1962   // connect to ESD
1963   TTree *t = (TTree*) esdFile->Get("esdTree");
1964   AliESDEvent *esd = new AliESDEvent();
1965   esd->ReadFromTree(t);
1966
1967   Int_t nEvents = t->GetEntries();
1968
1969   // set arrays and pointers
1970   Float_t posF[3];
1971   Double_t pos[3];
1972   Double_t p[3];
1973   Double_t p_pos[3];
1974   Double_t p_neg[3];
1975   Double_t covVtx[6];
1976   Double_t covTr[21];
1977   Double_t pid[10];
1978
1979   // loop over events and fill them
1980   for (Int_t iEvent = 0; iEvent < nEvents; ++iEvent) {
1981     //cout << "event: " << iEvent << endl;
1982     t->GetEntry(iEvent);
1983
1984     // Multiplicity information needed by the header (to be revised!)
1985     Int_t nTracks   = esd->GetNumberOfTracks();
1986     Int_t nPosTracks = 0;
1987     for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) 
1988       if (esd->GetTrack(iTrack)->Charge()> 0) nPosTracks++;
1989
1990     // Access the header
1991     AliAODHeader *header = aod->GetHeader();
1992
1993     // fill the header
1994     header->SetRunNumber       (esd->GetRunNumber()       );
1995     header->SetBunchCrossNumber(esd->GetBunchCrossNumber());
1996     header->SetOrbitNumber     (esd->GetOrbitNumber()     );
1997     header->SetPeriodNumber    (esd->GetPeriodNumber()    );
1998     header->SetTriggerMask     (esd->GetTriggerMask()     ); 
1999     header->SetTriggerCluster  (esd->GetTriggerCluster()  );
2000     header->SetEventType       (esd->GetEventType()       );
2001     header->SetMagneticField   (esd->GetMagneticField()   );
2002     header->SetZDCN1Energy     (esd->GetZDCN1Energy()     );
2003     header->SetZDCP1Energy     (esd->GetZDCP1Energy()     );
2004     header->SetZDCN2Energy     (esd->GetZDCN2Energy()     );
2005     header->SetZDCP2Energy     (esd->GetZDCP2Energy()     );
2006     header->SetZDCEMEnergy     (esd->GetZDCEMEnergy()     );
2007     header->SetRefMultiplicity   (nTracks);
2008     header->SetRefMultiplicityPos(nPosTracks);
2009     header->SetRefMultiplicityNeg(nTracks - nPosTracks);
2010     header->SetMuonMagFieldScale(-999.); // FIXME
2011     header->SetCentrality(-999.);        // FIXME
2012
2013     Int_t nV0s      = esd->GetNumberOfV0s();
2014     Int_t nCascades = esd->GetNumberOfCascades();
2015     Int_t nKinks    = esd->GetNumberOfKinks();
2016         Int_t nVertices = nV0s + 2*nCascades /*could lead to two vertices, one V0 and the Xi */+ nKinks + 1 /* = prim. vtx*/;    
2017         Int_t nJets     = 0;
2018     Int_t nCaloClus = esd->GetNumberOfCaloClusters();
2019     Int_t nFmdClus  = 0;
2020     Int_t nPmdClus  = esd->GetNumberOfPmdTracks();
2021    
2022     aod->ResetStd(nTracks, nVertices, nV0s+nCascades, nJets, nCaloClus, nFmdClus, nPmdClus);
2023     
2024     // Array to take into account the tracks already added to the AOD
2025     Bool_t * usedTrack = NULL;
2026     if (nTracks>0) {
2027       usedTrack = new Bool_t[nTracks];
2028       for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) usedTrack[iTrack]=kFALSE;
2029     }
2030     // Array to take into account the V0s already added to the AOD
2031     Bool_t * usedV0 = NULL;
2032     if (nV0s>0) {
2033       usedV0 = new Bool_t[nV0s];
2034       for (Int_t iV0=0; iV0<nV0s; ++iV0) usedV0[iV0]=kFALSE;
2035     }
2036     // Array to take into account the kinks already added to the AOD
2037     Bool_t * usedKink = NULL;
2038     if (nKinks>0) {
2039       usedKink = new Bool_t[nKinks];
2040       for (Int_t iKink=0; iKink<nKinks; ++iKink) usedKink[iKink]=kFALSE;
2041     }
2042     
2043     // Access to the AOD container of vertices
2044     TClonesArray &vertices = *(aod->GetVertices());
2045     Int_t jVertices=0;
2046
2047     // Access to the AOD container of tracks
2048     TClonesArray &tracks = *(aod->GetTracks());
2049     Int_t jTracks=0; 
2050    
2051     // Access to the AOD container of V0s
2052     TClonesArray &V0s = *(aod->GetV0s());
2053     Int_t jV0s=0;
2054     
2055     // Add primary vertex. The primary tracks will be defined
2056     // after the loops on the composite objects (V0, cascades, kinks)
2057     const AliESDVertex *vtx = esd->GetPrimaryVertex();
2058       
2059     vtx->GetXYZ(pos); // position
2060     vtx->GetCovMatrix(covVtx); //covariance matrix
2061
2062     AliAODVertex * primary = new(vertices[jVertices++])
2063       AliAODVertex(pos, covVtx, vtx->GetChi2toNDF(), NULL, -1, AliAODVertex::kPrimary);
2064          
2065
2066     AliAODTrack *aodTrack = 0x0;
2067     
2068     // Create vertices starting from the most complex objects
2069
2070     // Cascades
2071     for (Int_t nCascade = 0; nCascade < nCascades; ++nCascade) {
2072       AliESDcascade *cascade = esd->GetCascade(nCascade);
2073       
2074       cascade->GetXYZ(pos[0], pos[1], pos[2]);
2075       cascade->GetPosCovXi(covVtx);
2076      
2077       // Add the cascade vertex
2078       AliAODVertex * vcascade = new(vertices[jVertices++]) AliAODVertex(pos,
2079                                                                         covVtx,
2080                                                                         cascade->GetChi2Xi(), // = chi2/NDF since NDF = 2*2-3
2081                                                                         primary,
2082                                                                         nCascade,
2083                                                                         AliAODVertex::kCascade);
2084
2085       primary->AddDaughter(vcascade); // the cascade 'particle' (represented by a vertex) is added as a daughter to the primary vertex
2086
2087       // Add the V0 from the cascade. The ESD class have to be optimized...
2088       // Now we have to search for the corresponding V0 in the list of V0s
2089       // using the indeces of the positive and negative tracks
2090
2091       Int_t posFromV0 = cascade->GetPindex();
2092       Int_t negFromV0 = cascade->GetNindex();
2093
2094
2095       AliESDv0 * v0 = 0x0;
2096       Int_t indV0 = -1;
2097
2098       for (Int_t iV0=0; iV0<nV0s; ++iV0) {
2099
2100         v0 = esd->GetV0(iV0);
2101         Int_t posV0 = v0->GetPindex();
2102         Int_t negV0 = v0->GetNindex();
2103
2104         if (posV0==posFromV0 && negV0==negFromV0) {
2105           indV0 = iV0;
2106           break;
2107         }
2108       }
2109
2110       AliAODVertex * vV0FromCascade = 0x0;
2111
2112       if (indV0>-1 && !usedV0[indV0]) {
2113         
2114         // the V0 exists in the array of V0s and is not used
2115
2116         usedV0[indV0] = kTRUE;
2117         
2118         v0->GetXYZ(pos[0], pos[1], pos[2]);
2119         v0->GetPosCov(covVtx);
2120         
2121         vV0FromCascade = new(vertices[jVertices++]) AliAODVertex(pos,
2122                                                                  covVtx,
2123                                                                  v0->GetChi2V0(), // = chi2/NDF since NDF = 2*2-3
2124                                                                  vcascade,
2125                                                                  indV0,
2126                                                                  AliAODVertex::kV0);
2127       } else {
2128
2129         // the V0 doesn't exist in the array of V0s or was used
2130         cerr << "Error: event " << iEvent << " cascade " << nCascade
2131              << " The V0 " << indV0 
2132              << " doesn't exist in the array of V0s or was used!" << endl;
2133
2134         cascade->GetXYZ(pos[0], pos[1], pos[2]);
2135         cascade->GetPosCov(covVtx);
2136       
2137         vV0FromCascade = new(vertices[jVertices++]) AliAODVertex(pos,
2138                                                                  covVtx,
2139                                                                  v0->GetChi2V0(), // = chi2/NDF since NDF = 2*2-3
2140                                                                  vcascade,
2141                                                                  indV0,
2142                                                                  AliAODVertex::kV0);
2143         vcascade->AddDaughter(vV0FromCascade);
2144
2145       }
2146
2147       // Add the positive tracks from the V0
2148
2149       if (! usedTrack[posFromV0]) {
2150
2151         usedTrack[posFromV0] = kTRUE;
2152
2153         AliESDtrack *esdTrack = esd->GetTrack(posFromV0);
2154         esdTrack->GetPxPyPz(p_pos);
2155         esdTrack->GetXYZ(pos);
2156         esdTrack->GetCovarianceXYZPxPyPz(covTr);
2157         esdTrack->GetESDpid(pid);
2158         
2159         vV0FromCascade->AddDaughter(aodTrack =
2160                                     new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2161                                            esdTrack->GetLabel(), 
2162                                            p_pos, 
2163                                            kTRUE,
2164                                            pos,
2165                                            kFALSE,
2166                                            covTr, 
2167                                            (Short_t)esdTrack->Charge(),
2168                                            esdTrack->GetITSClusterMap(), 
2169                                            pid,
2170                                            vV0FromCascade,
2171                                            kTRUE,  // check if this is right
2172                                            kFALSE, // check if this is right
2173                                            AliAODTrack::kSecondary)
2174                 );
2175         aodTrack->ConvertAliPIDtoAODPID();
2176       }
2177       else {
2178         cerr << "Error: event " << iEvent << " cascade " << nCascade
2179              << " track " << posFromV0 << " has already been used!" << endl;
2180       }
2181
2182       // Add the negative tracks from the V0
2183
2184       if (!usedTrack[negFromV0]) {
2185         
2186         usedTrack[negFromV0] = kTRUE;
2187         
2188         AliESDtrack *esdTrack = esd->GetTrack(negFromV0);
2189         esdTrack->GetPxPyPz(p_neg);
2190         esdTrack->GetXYZ(pos);
2191         esdTrack->GetCovarianceXYZPxPyPz(covTr);
2192         esdTrack->GetESDpid(pid);
2193         
2194         vV0FromCascade->AddDaughter(aodTrack =
2195                 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2196                                            esdTrack->GetLabel(),
2197                                            p_neg,
2198                                            kTRUE,
2199                                            pos,
2200                                            kFALSE,
2201                                            covTr, 
2202                                            (Short_t)esdTrack->Charge(),
2203                                            esdTrack->GetITSClusterMap(), 
2204                                            pid,
2205                                            vV0FromCascade,
2206                                            kTRUE,  // check if this is right
2207                                            kFALSE, // check if this is right
2208                                            AliAODTrack::kSecondary)
2209                 );
2210         aodTrack->ConvertAliPIDtoAODPID();
2211       }
2212       else {
2213         cerr << "Error: event " << iEvent << " cascade " << nCascade
2214              << " track " << negFromV0 << " has already been used!" << endl;
2215       }
2216
2217       // add it to the V0 array as well
2218       Double_t d0[2] = { -999., -99.};
2219       // counting is probably wrong
2220       new(V0s[jV0s++]) AliAODv0(vV0FromCascade, -999., -99., p_pos, p_neg, d0); // to be refined
2221
2222       // Add the bachelor track from the cascade
2223
2224       Int_t bachelor = cascade->GetBindex();
2225       
2226       if(!usedTrack[bachelor]) {
2227       
2228         usedTrack[bachelor] = kTRUE;
2229         
2230         AliESDtrack *esdTrack = esd->GetTrack(bachelor);
2231         esdTrack->GetPxPyPz(p);
2232         esdTrack->GetXYZ(pos);
2233         esdTrack->GetCovarianceXYZPxPyPz(covTr);
2234         esdTrack->GetESDpid(pid);
2235
2236         vcascade->AddDaughter(aodTrack =
2237                 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2238                                            esdTrack->GetLabel(),
2239                                            p,
2240                                            kTRUE,
2241                                            pos,
2242                                            kFALSE,
2243                                            covTr, 
2244                                            (Short_t)esdTrack->Charge(),
2245                                            esdTrack->GetITSClusterMap(), 
2246                                            pid,
2247                                            vcascade,
2248                                            kTRUE,  // check if this is right
2249                                            kFALSE, // check if this is right
2250                                            AliAODTrack::kSecondary)
2251                 );
2252         aodTrack->ConvertAliPIDtoAODPID();
2253      }
2254       else {
2255         cerr << "Error: event " << iEvent << " cascade " << nCascade
2256              << " track " << bachelor << " has already been used!" << endl;
2257       }
2258       
2259       // Add the primary track of the cascade (if any)
2260       
2261     } // end of the loop on cascades
2262  
2263     // V0s
2264         
2265     for (Int_t nV0 = 0; nV0 < nV0s; ++nV0) {
2266
2267       if (usedV0[nV0]) continue; // skip if aready added to the AOD
2268
2269       AliESDv0 *v0 = esd->GetV0(nV0); 
2270      
2271       v0->GetXYZ(pos[0], pos[1], pos[2]);
2272       v0->GetPosCov(covVtx);
2273
2274       AliAODVertex * vV0 = 
2275         new(vertices[jVertices++]) AliAODVertex(pos,
2276                                                 covVtx,
2277                                                 v0->GetChi2V0(), // = chi2/NDF since NDF = 2*2-3
2278                                                 primary,
2279                                                 nV0,
2280                                                 AliAODVertex::kV0);
2281       primary->AddDaughter(vV0);
2282
2283       Int_t posFromV0 = v0->GetPindex();
2284       Int_t negFromV0 = v0->GetNindex();
2285       
2286       // Add the positive tracks from the V0
2287
2288       if (!usedTrack[posFromV0]) {
2289         
2290         usedTrack[posFromV0] = kTRUE;
2291
2292         AliESDtrack *esdTrack = esd->GetTrack(posFromV0);
2293         esdTrack->GetPxPyPz(p_pos);
2294         esdTrack->GetXYZ(pos);
2295         esdTrack->GetCovarianceXYZPxPyPz(covTr);
2296         esdTrack->GetESDpid(pid);
2297         
2298         vV0->AddDaughter(aodTrack =
2299                 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2300                                            esdTrack->GetLabel(), 
2301                                            p_pos, 
2302                                            kTRUE,
2303                                            pos,
2304                                            kFALSE,
2305                                            covTr, 
2306                                            (Short_t)esdTrack->Charge(),
2307                                            esdTrack->GetITSClusterMap(), 
2308                                            pid,
2309                                            vV0,
2310                                            kTRUE,  // check if this is right
2311                                            kFALSE, // check if this is right
2312                                            AliAODTrack::kSecondary)
2313                 );
2314         aodTrack->ConvertAliPIDtoAODPID();
2315       }
2316       else {
2317         cerr << "Error: event " << iEvent << " V0 " << nV0
2318              << " track " << posFromV0 << " has already been used!" << endl;
2319       }
2320
2321       // Add the negative tracks from the V0
2322
2323       if (!usedTrack[negFromV0]) {
2324
2325         usedTrack[negFromV0] = kTRUE;
2326
2327         AliESDtrack *esdTrack = esd->GetTrack(negFromV0);
2328         esdTrack->GetPxPyPz(p_neg);
2329         esdTrack->GetXYZ(pos);
2330         esdTrack->GetCovarianceXYZPxPyPz(covTr);
2331         esdTrack->GetESDpid(pid);
2332
2333         vV0->AddDaughter(aodTrack =
2334                 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2335                                            esdTrack->GetLabel(),
2336                                            p_neg,
2337                                            kTRUE,
2338                                            pos,
2339                                            kFALSE,
2340                                            covTr, 
2341                                            (Short_t)esdTrack->Charge(),
2342                                            esdTrack->GetITSClusterMap(), 
2343                                            pid,
2344                                            vV0,
2345                                            kTRUE,  // check if this is right
2346                                            kFALSE, // check if this is right
2347                                            AliAODTrack::kSecondary)
2348                 );
2349         aodTrack->ConvertAliPIDtoAODPID();
2350       }
2351       else {
2352         cerr << "Error: event " << iEvent << " V0 " << nV0
2353              << " track " << negFromV0 << " has already been used!" << endl;
2354       }
2355
2356       // add it to the V0 array as well
2357       Double_t d0[2] = { 999., 99.};
2358       new(V0s[jV0s++]) AliAODv0(vV0, 999., 99., p_pos, p_neg, d0); // to be refined
2359     }
2360         V0s.Expand(jV0s);        
2361     // end of the loop on V0s
2362     
2363     // Kinks: it is a big mess the access to the information in the kinks
2364     // The loop is on the tracks in order to find the mother and daugther of each kink
2365
2366
2367     for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) {
2368
2369       AliESDtrack * esdTrack = esd->GetTrack(iTrack);
2370
2371       Int_t ikink = esdTrack->GetKinkIndex(0);
2372
2373       if (ikink) {
2374         // Negative kink index: mother, positive: daughter
2375
2376         // Search for the second track of the kink
2377
2378         for (Int_t jTrack = iTrack+1; jTrack<nTracks; ++jTrack) {
2379
2380           AliESDtrack * esdTrack1 = esd->GetTrack(jTrack);
2381
2382           Int_t jkink = esdTrack1->GetKinkIndex(0);
2383
2384           if ( TMath::Abs(ikink)==TMath::Abs(jkink) ) {
2385
2386             // The two tracks are from the same kink
2387           
2388             if (usedKink[TMath::Abs(ikink)-1]) continue; // skip used kinks
2389
2390             Int_t imother = -1;
2391             Int_t idaughter = -1;
2392
2393             if (ikink<0 && jkink>0) {
2394
2395               imother = iTrack;
2396               idaughter = jTrack;
2397             }
2398             else if (ikink>0 && jkink<0) {
2399
2400               imother = jTrack;
2401               idaughter = iTrack;
2402             }
2403             else {
2404               cerr << "Error: Wrong combination of kink indexes: "
2405               << ikink << " " << jkink << endl;
2406               continue;
2407             }
2408
2409             // Add the mother track
2410
2411             AliAODTrack * mother = NULL;
2412
2413             if (!usedTrack[imother]) {
2414         
2415               usedTrack[imother] = kTRUE;
2416         
2417               AliESDtrack *esdTrack = esd->GetTrack(imother);
2418               esdTrack->GetPxPyPz(p);
2419               esdTrack->GetXYZ(pos);
2420               esdTrack->GetCovarianceXYZPxPyPz(covTr);
2421               esdTrack->GetESDpid(pid);
2422
2423               mother = 
2424                 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2425                                            esdTrack->GetLabel(),
2426                                            p,
2427                                            kTRUE,
2428                                            pos,
2429                                            kFALSE,
2430                                            covTr, 
2431                                            (Short_t)esdTrack->Charge(),
2432                                            esdTrack->GetITSClusterMap(), 
2433                                            pid,
2434                                            primary,
2435                                            kTRUE, // check if this is right
2436                                            kTRUE, // check if this is right
2437                                            AliAODTrack::kPrimary);
2438               primary->AddDaughter(mother);
2439               mother->ConvertAliPIDtoAODPID();
2440             }
2441             else {
2442               cerr << "Error: event " << iEvent << " kink " << TMath::Abs(ikink)-1
2443               << " track " << imother << " has already been used!" << endl;
2444             }
2445
2446             // Add the kink vertex
2447             AliESDkink * kink = esd->GetKink(TMath::Abs(ikink)-1);
2448
2449             AliAODVertex * vkink = 
2450             new(vertices[jVertices++]) AliAODVertex(kink->GetPosition(),
2451                                                     NULL,
2452                                                     0.,
2453                                                     mother,
2454                                                     esdTrack->GetID(), // This is the track ID of the mother's track!
2455                                                     AliAODVertex::kKink);
2456             // Add the daughter track
2457
2458             AliAODTrack * daughter = NULL;
2459
2460             if (!usedTrack[idaughter]) {
2461         
2462               usedTrack[idaughter] = kTRUE;
2463         
2464               AliESDtrack *esdTrack = esd->GetTrack(idaughter);
2465               esdTrack->GetPxPyPz(p);
2466               esdTrack->GetXYZ(pos);
2467               esdTrack->GetCovarianceXYZPxPyPz(covTr);
2468               esdTrack->GetESDpid(pid);
2469
2470               daughter = 
2471                 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2472                                            esdTrack->GetLabel(),
2473                                            p,
2474                                            kTRUE,
2475                                            pos,
2476                                            kFALSE,
2477                                            covTr, 
2478                                            (Short_t)esdTrack->Charge(),
2479                                            esdTrack->GetITSClusterMap(), 
2480                                            pid,
2481                                            vkink,
2482                                            kTRUE, // check if this is right
2483                                            kTRUE, // check if this is right
2484                                            AliAODTrack::kPrimary);
2485               vkink->AddDaughter(daughter);
2486               daughter->ConvertAliPIDtoAODPID();
2487             }
2488             else {
2489               cerr << "Error: event " << iEvent << " kink " << TMath::Abs(ikink)-1
2490               << " track " << idaughter << " has already been used!" << endl;
2491             }
2492           }
2493         }
2494       }
2495     }
2496     vertices.Expand(jVertices);
2497
2498     // Tracks (primary and orphan)
2499     for (Int_t nTrack = 0; nTrack < nTracks; ++nTrack) {
2500
2501       if (usedTrack[nTrack]) continue;
2502
2503       AliESDtrack *esdTrack = esd->GetTrack(nTrack);
2504       esdTrack->GetPxPyPz(p);
2505       esdTrack->GetXYZ(pos);
2506       esdTrack->GetCovarianceXYZPxPyPz(covTr);
2507       esdTrack->GetESDpid(pid);
2508
2509       Float_t impactXY, impactZ;
2510
2511       esdTrack->GetImpactParameters(impactXY,impactZ);
2512
2513       if (impactXY<3.) {
2514         // track inside the beam pipe
2515       
2516         primary->AddDaughter(aodTrack =
2517             new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2518                                          esdTrack->GetLabel(),
2519                                          p,
2520                                          kTRUE,
2521                                          pos,
2522                                          kFALSE,
2523                                          covTr, 
2524                                          (Short_t)esdTrack->Charge(),
2525                                          esdTrack->GetITSClusterMap(), 
2526                                          pid,
2527                                          primary,
2528                                          kTRUE, // check if this is right
2529                                          kTRUE, // check if this is right
2530                                          AliAODTrack::kPrimary)
2531             );
2532         aodTrack->ConvertAliPIDtoAODPID();
2533       }
2534       else {
2535         // outside the beam pipe: orphan track
2536         // Don't write them anymore!
2537         continue;
2538       } 
2539     } // end of loop on tracks
2540     
2541     // muon tracks
2542     Int_t nMuTracks = esd->GetNumberOfMuonTracks();
2543     for (Int_t nMuTrack = 0; nMuTrack < nMuTracks; ++nMuTrack) {
2544       
2545       AliESDMuonTrack *esdMuTrack = esd->GetMuonTrack(nMuTrack);     
2546       p[0] = esdMuTrack->Px(); 
2547       p[1] = esdMuTrack->Py(); 
2548       p[2] = esdMuTrack->Pz();
2549       pos[0] = primary->GetX(); 
2550       pos[1] = primary->GetY(); 
2551       pos[2] = primary->GetZ();
2552       
2553       // has to be changed once the muon pid is provided by the ESD
2554       for (Int_t i = 0; i < 10; pid[i++] = 0.); pid[AliAODTrack::kMuon]=1.;
2555       
2556       primary->AddDaughter(aodTrack =
2557           new(tracks[jTracks++]) AliAODTrack(0, // no ID provided
2558                                              0, // no label provided
2559                                              p,
2560                                              kTRUE,
2561                                              pos,
2562                                              kFALSE,
2563                                              NULL, // no covariance matrix provided
2564                                              esdMuTrack->Charge(),
2565                                              0, // ITSClusterMap is set below
2566                                              pid,
2567                                              primary,
2568                                              kFALSE,  // muon tracks are not used to fit the primary vtx
2569                                              kFALSE,  // not used for vertex fit
2570                                              AliAODTrack::kPrimary)
2571           );
2572
2573       aodTrack->SetHitsPatternInTrigCh(esdMuTrack->GetHitsPatternInTrigCh());
2574       Int_t track2Trigger = esdMuTrack->GetMatchTrigger();
2575       aodTrack->SetMatchTrigger(track2Trigger);
2576       if (track2Trigger) 
2577         aodTrack->SetChi2MatchTrigger(esdMuTrack->GetChi2MatchTrigger());
2578       else 
2579         aodTrack->SetChi2MatchTrigger(0.);
2580     }
2581     tracks.Expand(jTracks); // remove 'empty slots' due to unwritten tracks
2582         
2583     // Access to the AOD container of PMD clusters
2584     TClonesArray &pmdClusters = *(aod->GetPmdClusters());
2585     Int_t jPmdClusters=0;
2586   
2587     for (Int_t iPmd = 0; iPmd < nPmdClus; ++iPmd) {
2588       // file pmd clusters, to be revised!
2589       AliESDPmdTrack *pmdTrack = esd->GetPmdTrack(iPmd);
2590       Int_t nLabel = 0;
2591       Int_t *label = 0x0;
2592       Double_t pos[3] = { pmdTrack->GetClusterX(), pmdTrack->GetClusterY(), pmdTrack->GetClusterZ() };
2593       Double_t pid[9] = { 0., 0., 0., 0., 0., 0., 0., 0., 0. }; // to be revised!
2594       // type not set!
2595       // assoc cluster not set
2596       new(pmdClusters[jPmdClusters++]) AliAODPmdCluster(iPmd, nLabel, label, pmdTrack->GetClusterADC(), pos, pid);
2597     }
2598
2599     // Access to the AOD container of clusters
2600     TClonesArray &caloClusters = *(aod->GetCaloClusters());
2601     Int_t jClusters=0;
2602
2603     // Calo Clusters
2604     TArrayS EMCCellNumber(15000);
2605     TArrayD EMCCellAmplitude(15000);
2606     Int_t nEMCCells = 0;
2607     const Float_t fEMCAmpScale = 1./500;
2608  
2609     for (Int_t iClust=0; iClust<nCaloClus; ++iClust) {
2610
2611       AliESDCaloCluster * cluster = esd->GetCaloCluster(iClust);
2612
2613       Int_t id = cluster->GetID();
2614       Int_t nLabel = 0;
2615       Int_t *label = 0x0;
2616       Float_t energy = cluster->E();
2617       cluster->GetPosition(posF);
2618       Char_t ttype=AliAODCluster::kUndef;
2619
2620       if (cluster->GetClusterType() == AliESDCaloCluster::kPHOSCluster) {
2621         ttype=AliAODCluster::kPHOSNeutral;
2622       } 
2623       else if (cluster->GetClusterType() == AliESDCaloCluster::kEMCALClusterv1) {
2624         ttype = AliAODCluster::kEMCALClusterv1;
2625       }
2626       else if (cluster->GetClusterType() == AliESDCaloCluster::kEMCALPseudoCluster) {
2627         // Collect raw tower info
2628         for (Int_t iDig = 0; iDig < cluster->GetNumberOfDigits(); iDig++) {
2629           EMCCellNumber[nEMCCells] = cluster->GetDigitIndex()->At(iDig);
2630           EMCCellAmplitude[nEMCCells] = fEMCAmpScale*cluster->GetDigitAmplitude()->At(iDig);
2631           nEMCCells++;
2632         }
2633         // don't write cluster data (it's just a pseudo cluster, holding the tower information)
2634         continue; 
2635       }
2636       
2637       AliAODCaloCluster *caloCluster = new(caloClusters[jClusters++]) AliAODCaloCluster(id,
2638                                                                                         nLabel,
2639                                                                                         label,
2640                                                                                         energy,
2641                                                                                         pos,
2642                                                                                         NULL,
2643                                                                                         ttype);
2644       
2645       caloCluster->SetCaloCluster(); // to be refined!
2646
2647     } 
2648         caloClusters.Expand(jClusters); // resize TObjArray to 'remove' slots for pseudo clusters        
2649         // end of loop on calo clusters
2650
2651     // fill EMC cell info
2652     AliAODCaloCells &EMCCells = *(aod->GetCaloCells());
2653     EMCCells.CreateContainer(nEMCCells);
2654     EMCCells.SetType(AliAODCaloCells::kEMCAL);
2655     for (Int_t iCell = 0; iCell < nEMCCells; iCell++) {      
2656       EMCCells.SetCell(iCell,EMCCellNumber[iCell],EMCCellAmplitude[iCell]);
2657     }
2658     EMCCells.Sort();
2659
2660     // tracklets    
2661     AliAODTracklets &SPDTracklets = *(aod->GetTracklets());
2662     const AliMultiplicity *mult = esd->GetMultiplicity();
2663     if (mult) {
2664       if (mult->GetNumberOfTracklets()>0) {
2665         SPDTracklets.CreateContainer(mult->GetNumberOfTracklets());
2666
2667         for (Int_t n=0; n<mult->GetNumberOfTracklets(); n++) {
2668           SPDTracklets.SetTracklet(n, mult->GetTheta(n), mult->GetPhi(n), mult->GetDeltaPhi(n), mult->GetLabel(n));
2669         }
2670       }
2671     } else {
2672       Printf("ERROR: AliMultiplicity could not be retrieved from ESD");
2673     }
2674
2675     delete [] usedTrack;
2676     delete [] usedV0;
2677     delete [] usedKink;
2678
2679     // fill the tree for this event
2680     aodTree->Fill();
2681   } // end of event loop
2682
2683   aodTree->GetUserInfo()->Add(aod);
2684
2685   // write the tree to the specified file
2686   aodFile = aodTree->GetCurrentFile();
2687   aodFile->cd();
2688   aodTree->Write();
2689
2690   return;
2691 }
2692
2693 void AliReconstruction::WriteAlignmentData(AliESDEvent* esd)
2694 {
2695   // Write space-points which are then used in the alignment procedures
2696   // For the moment only ITS, TRD and TPC
2697
2698   // Load TOF clusters
2699   if (fTracker[3]){
2700     fLoader[3]->LoadRecPoints("read");
2701     TTree* tree = fLoader[3]->TreeR();
2702     if (!tree) {
2703       AliError(Form("Can't get the %s cluster tree", fgkDetectorName[3]));
2704       return;
2705     }
2706     fTracker[3]->LoadClusters(tree);
2707   }
2708   Int_t ntracks = esd->GetNumberOfTracks();
2709   for (Int_t itrack = 0; itrack < ntracks; itrack++)
2710     {
2711       AliESDtrack *track = esd->GetTrack(itrack);
2712       Int_t nsp = 0;
2713       Int_t idx[200];
2714       for (Int_t iDet = 3; iDet >= 0; iDet--)
2715         nsp += track->GetNcls(iDet);
2716       if (nsp) {
2717         AliTrackPointArray *sp = new AliTrackPointArray(nsp);
2718         track->SetTrackPointArray(sp);
2719         Int_t isptrack = 0;
2720         for (Int_t iDet = 3; iDet >= 0; iDet--) {
2721           AliTracker *tracker = fTracker[iDet];
2722           if (!tracker) continue;
2723           Int_t nspdet = track->GetNcls(iDet);
2724           if (nspdet <= 0) continue;
2725           track->GetClusters(iDet,idx);
2726           AliTrackPoint p;
2727           Int_t isp = 0;
2728           Int_t isp2 = 0;
2729           while (isp < nspdet) {
2730             Bool_t isvalid;
2731             if(IsSelected(fgkDetectorName[iDet],fUseTrackingErrorsForAlignment)) {
2732               isvalid = tracker->GetTrackPointTrackingError(idx[isp2],p,track);
2733             } else {
2734               isvalid = tracker->GetTrackPoint(idx[isp2],p); 
2735             } 
2736             isp2++;
2737             const Int_t kNTPCmax = 159;
2738             if (iDet==1 && isp2>kNTPCmax) break;   // to be fixed
2739             if (!isvalid) continue;
2740             sp->AddPoint(isptrack,&p); isptrack++; isp++;
2741           }
2742         }       
2743       }
2744     }
2745   if (fTracker[3]){
2746     fTracker[3]->UnloadClusters();
2747     fLoader[3]->UnloadRecPoints();
2748   }
2749 }
2750
2751 //_____________________________________________________________________________
2752 void AliReconstruction::FillRawDataErrorLog(Int_t iEvent, AliESDEvent* esd)
2753 {
2754   // The method reads the raw-data error log
2755   // accumulated within the rawReader.
2756   // It extracts the raw-data errors related to
2757   // the current event and stores them into
2758   // a TClonesArray inside the esd object.
2759
2760   if (!fRawReader) return;
2761
2762   for(Int_t i = 0; i < fRawReader->GetNumberOfErrorLogs(); i++) {
2763
2764     AliRawDataErrorLog *log = fRawReader->GetErrorLog(i);
2765     if (!log) continue;
2766     if (iEvent != log->GetEventNumber()) continue;
2767
2768     esd->AddRawDataErrorLog(log);
2769   }
2770
2771 }
2772
2773 TNamed* AliReconstruction::CopyFileToTNamed(TString fPath,TString fName){
2774   // Dump a file content into a char in TNamed
2775   ifstream in;
2776   in.open(fPath.Data(),ios::in | ios::binary|ios::ate);
2777   Int_t kBytes = (Int_t)in.tellg();
2778   printf("Size: %d \n",kBytes);
2779   TNamed *fn = 0;
2780   if(in.good()){
2781     char* memblock = new char [kBytes];
2782     in.seekg (0, ios::beg);
2783     in.read (memblock, kBytes);
2784     in.close();
2785     TString fData(memblock,kBytes);
2786     fn = new TNamed(fName,fData);
2787     printf("fData Size: %d \n",fData.Sizeof());
2788     printf("fName Size: %d \n",fName.Sizeof());
2789     printf("fn    Size: %d \n",fn->Sizeof());
2790     delete[] memblock;
2791   }
2792   else{
2793     AliInfo(Form("Could not Open %s\n",fPath.Data()));
2794   }
2795
2796   return fn;
2797 }
2798
2799 void AliReconstruction::TNamedToFile(TTree* fTree, TString fName){
2800   // This is not really needed in AliReconstruction at the moment
2801   // but can serve as a template
2802
2803   TList *fList = fTree->GetUserInfo();
2804   TNamed *fn = (TNamed*)fList->FindObject(fName.Data());
2805   printf("fn Size: %d \n",fn->Sizeof());
2806
2807   TString fTmp(fn->GetName()); // to be 100% sure in principle fName also works
2808   const char* cdata = fn->GetTitle();
2809   printf("fTmp Size %d\n",fTmp.Sizeof());
2810
2811   int size = fn->Sizeof()-fTmp.Sizeof()-sizeof(UChar_t)-sizeof(Int_t); // see dfinition of TString::SizeOf()...
2812   printf("calculated size %d\n",size);
2813   ofstream out(fName.Data(),ios::out | ios::binary);
2814   out.write(cdata,size);
2815   out.close();
2816
2817 }
2818   
2819
2820
2821 //_____________________________________________________________________________
2822 //AliQADataMaker * AliReconstruction::GetQADataMaker(Int_t iDet)
2823 //{
2824 // get the quality assurance data maker object and the loader for a detector
2825 //
2826 //  if (fQADataMaker[iDet]) 
2827 //    return fQADataMaker[iDet];
2828 //
2829 //  // load the QA data maker object
2830 //  TPluginManager* pluginManager = gROOT->GetPluginManager();
2831 //  TString detName = fgkDetectorName[iDet];
2832 //  TString qadmName = "Ali" + detName + "QADataMaker";
2833 //  if (gAlice && !gAlice->GetDetector(detName) && (detName != "HLT")) 
2834 //    return NULL;
2835 //
2836 //  AliQADataMaker * qadm = NULL;
2837 //  // first check if a plugin is defined for the quality assurance data maker
2838 //  TPluginHandler* pluginHandler = pluginManager->FindHandler("AliQADataMaker", detName);
2839 //  // if not, add a plugin for it
2840 //  if (!pluginHandler) {
2841 //    AliDebug(1, Form("defining plugin for %s", qadmName.Data()));
2842 //    TString libs = gSystem->GetLibraries();
2843 //    if (libs.Contains("lib" + detName + "base.so") ||
2844 //      (gSystem->Load("lib" + detName + "base.so") >= 0)) {
2845 //      pluginManager->AddHandler("AliQADataMaker", detName, 
2846 //                              qadmName, detName + "qadm", qadmName + "()");
2847 //    } else {
2848 //      pluginManager->AddHandler("AliQADataMaker", detName, 
2849 //                              qadmName, detName, qadmName + "()");
2850 //    }
2851 //    pluginHandler = pluginManager->FindHandler("AliQADataMaker", detName);
2852 //  }
2853 //  if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
2854 //    qadm = (AliQADataMaker *) pluginHandler->ExecPlugin(0);
2855 //  }
2856 //  if (qadm) {
2857 //    AliInfo(Form("Initializing quality assurance data maker for %s", fgkDetectorName[iDet]));
2858 //    qadm->Init(AliQA::kRECPOINTS, AliCDBManager::Instance()->GetRun(), GetQACycles(fgkDetectorName[iDet]));
2859 //    qadm->StartOfCycle(AliQA::kRECPOINTS);
2860 //    qadm->Init(AliQA::kESDS, AliCDBManager::Instance()->GetRun());
2861 //    qadm->StartOfCycle(AliQA::kESDS, "same") ;        
2862 //    fQADataMaker[iDet] = qadm;
2863 //  }
2864 //
2865 //  return qadm;
2866 //}
2867 //
2868 //_____________________________________________________________________________
2869 //Bool_t AliReconstruction::RunQA(const char* detectors, AliESDEvent *& esd)
2870 //{
2871 //  // run the Quality Assurance data producer
2872 //
2873 //  AliCodeTimerAuto("")
2874 //  TString detStr = detectors;
2875 //  for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2876 //   if (!IsSelected(fgkDetectorName[iDet], detStr)) 
2877 //     continue;
2878 //   AliQADataMaker * qadm = GetQADataMaker(iDet);
2879 //   if (!qadm) 
2880 //     continue;
2881 //   AliCodeTimerStart(Form("running quality assurance data maker for %s", fgkDetectorName[iDet]));
2882 //   AliInfo(Form("running quality assurance data maker for %s", fgkDetectorName[iDet]));
2883 //    
2884 //   qadm->Exec(AliQA::kESDS, esd) ; 
2885 //   qadm->Increment() ; 
2886 //
2887 //   AliCodeTimerStop(Form("running quality assurance data maker for %s", fgkDetectorName[iDet]));
2888 // }
2889 // if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
2890 //   AliError(Form("the following detectors were not found: %s",
2891 //               detStr.Data()));
2892 //   if (fStopOnError) 
2893 //     return kFALSE;
2894 // }
2895 // 
2896 // return kTRUE;
2897 //  
2898 //}
2899 //
2900
2901 //_____________________________________________________________________________
2902 void AliReconstruction::CheckQA()
2903 {
2904 // check the QA of SIM for this run and remove the detectors 
2905 // with status Fatal
2906   
2907         TString newRunLocalReconstruction ; 
2908         TString newRunTracking ;
2909         TString newFillESD ;
2910          
2911         for (Int_t iDet = 0; iDet < AliQA::kNDET; iDet++) {
2912                 TString detName(AliQA::GetDetName(iDet)) ;
2913                 AliQA * qa = AliQA::Instance(AliQA::DETECTORINDEX(iDet)) ; 
2914                 if ( qa->IsSet(AliQA::DETECTORINDEX(iDet), AliQA::kSIM, AliQA::kFATAL)) {
2915                                 AliInfo(Form("QA status for %s in Hits and/or SDIGITS  and/or Digits was Fatal; No reconstruction performed", detName.Data())) ;
2916                 } else {
2917                         if ( fRunLocalReconstruction.Contains(AliQA::GetDetName(iDet)) || 
2918                                         fRunLocalReconstruction.Contains("ALL") )  {
2919                                 newRunLocalReconstruction += detName ; 
2920                                 newRunLocalReconstruction += " " ;                      
2921                         }
2922                         if ( fRunTracking.Contains(AliQA::GetDetName(iDet)) || 
2923                                         fRunTracking.Contains("ALL") )  {
2924                                 newRunTracking += detName ; 
2925                                 newRunTracking += " " ;                         
2926                         }
2927                         if ( fFillESD.Contains(AliQA::GetDetName(iDet)) || 
2928                                         fFillESD.Contains("ALL") )  {
2929                                 newFillESD += detName ; 
2930                                 newFillESD += " " ;                     
2931                         }
2932                 }
2933         }
2934         fRunLocalReconstruction = newRunLocalReconstruction ; 
2935         fRunTracking            = newRunTracking ; 
2936         fFillESD                = newFillESD ; 
2937 }
2938
2939 //_____________________________________________________________________________
2940 Int_t AliReconstruction::GetDetIndex(const char* detector)
2941 {
2942   // return the detector index corresponding to detector
2943   Int_t index = -1 ; 
2944   for (index = 0; index < fgkNDetectors ; index++) {
2945     if ( strcmp(detector, fgkDetectorName[index]) == 0 )
2946         break ; 
2947   }     
2948   return index ; 
2949 }