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