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