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
 695   // checking the QA of previous steps
 696   CheckQA() ;
 697
 698   for (Int_t iEvent = 0; iEvent < fRunLoader->GetNumberOfEvents(); iEvent++) {
 699     if (fRawReader) fRawReader->NextEvent();
 700     if ((iEvent < fFirstEvent) || ((fLastEvent >= 0) && (iEvent > fLastEvent))) {
 701       // copy old ESD to the new one
 702       if (treeOld) {
 703         esd->ReadFromTree(treeOld);
 704         treeOld->GetEntry(iEvent);
 705       }
 706       tree->Fill();
 707       if (hlttreeOld) {
 708         esd->ReadFromTree(hlttreeOld);
 709         hlttreeOld->GetEntry(iEvent);
 710       }
 711       hlttree->Fill();
 712       continue;
 713     }
 714
 715     AliInfo(Form("processing event %d", iEvent));
 716     fRunLoader->GetEvent(iEvent);
 717
 718     char aFileName[256];
 719     sprintf(aFileName, "ESD_%d.%d_final.root",
 720             fRunLoader->GetHeader()->GetRun(),
 721             fRunLoader->GetHeader()->GetEventNrInRun());
 722     if (!gSystem->AccessPathName(aFileName)) continue;
 723
 724     // local reconstruction
 725     if (!fRunLocalReconstruction.IsNull()) {
 726       if (!RunLocalEventReconstruction(fRunLocalReconstruction)) {
 727         if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
 728       }
 729     }
 730
 731
 732     esd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
 733     hltesd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
 734     esd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
 735     hltesd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
 736
 737     // Set magnetic field from the tracker
 738     esd->SetMagneticField(AliTracker::GetBz());
 739     hltesd->SetMagneticField(AliTracker::GetBz());
 740
 741
 742
 743     // Fill raw-data error log into the ESD
 744     if (fRawReader) FillRawDataErrorLog(iEvent,esd);
 745
 746     // vertex finder
 747     if (fRunVertexFinder) {
 748       if (!ReadESD(esd, "vertex")) {
 749         if (!RunVertexFinder(esd)) {
 750           if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
 751         }
 752         if (fCheckPointLevel > 0) WriteESD(esd, "vertex");
 753       }
 754     }
 755
 756     // HLT tracking
 757     if (!fRunTracking.IsNull()) {
 758       if (fRunHLTTracking) {
 759         hltesd->SetVertex(esd->GetVertex());
 760         if (!RunHLTTracking(hltesd)) {
 761           if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
 762         }
 763       }
 764     }
 765
 766     // Muon tracking
 767     if (!fRunTracking.IsNull()) {
 768       if (fRunMuonTracking) {
 769         if (!RunMuonTracking(esd)) {
 770           if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
 771         }
 772       }
 773     }
 774
 775     // barrel tracking
 776     if (!fRunTracking.IsNull()) {
 777       if (!ReadESD(esd, "tracking")) {
 778         if (!RunTracking(esd)) {
 779           if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
 780         }
 781         if (fCheckPointLevel > 0) WriteESD(esd, "tracking");
 782       }
 783     }
 784
 785     // fill ESD
 786     if (!fFillESD.IsNull()) {
 787       if (!FillESD(esd, fFillESD)) {
 788         if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
 789       }
 790     }
 791
 792     if (!fFillESD.IsNull())
 793     RunQA(fFillESD.Data(), esd);
 794
 795     // fill Event header information from the RawEventHeader
 796     if (fRawReader){FillRawEventHeaderESD(esd);}
 797
 798     // combined PID
 799     AliESDpid::MakePID(esd);
 800     if (fCheckPointLevel > 1) WriteESD(esd, "PID");
 801
 802     if (fFillTriggerESD) {
 803       if (!ReadESD(esd, "trigger")) {
 804         if (!FillTriggerESD(esd)) {
 805           if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
 806         }
 807         if (fCheckPointLevel > 1) WriteESD(esd, "trigger");
 808       }
 809     }
 810
 811     file->cd();
 812
 813     //Try to improve the reconstructed primary vertex position using the tracks
 814     AliESDVertex *pvtx=0;
 815     Bool_t dovertex=kTRUE;
 816     TObject* obj = fOptions.FindObject("ITS");
 817     if (obj) {
 818       TString optITS = obj->GetTitle();
 819       if (optITS.Contains("cosmics") || optITS.Contains("COSMICS"))
 820         dovertex=kFALSE;
 821     }
 822     if(dovertex) pvtx=tVertexer.FindPrimaryVertex(esd);
 823     if(fDiamondProfile) esd->SetDiamond(fDiamondProfile);
 824
 825     if (pvtx)
 826     if (pvtx->GetStatus()) {
 827        // Store the improved primary vertex
 828        esd->SetPrimaryVertex(pvtx);
 829        // Propagate the tracks to the DCA to the improved primary vertex
 830        Double_t somethingbig = 777.;
 831        Double_t bz = esd->GetMagneticField();
 832        Int_t nt=esd->GetNumberOfTracks();
 833        while (nt--) {
 834          AliESDtrack *t = esd->GetTrack(nt);
 835          t->RelateToVertex(pvtx, bz, somethingbig);
 836        }
 837     }
 838
 839     if (fRunV0Finder) {
 840        // V0 finding
 841        AliV0vertexer vtxer;
 842        vtxer.Tracks2V0vertices(esd);
 843
 844        if (fRunCascadeFinder) {
 845           // Cascade finding
 846           AliCascadeVertexer cvtxer;
 847           cvtxer.V0sTracks2CascadeVertices(esd);
 848        }
 849     }
 850
 851     // write ESD
 852     if (fCleanESD) CleanESD(esd);
 853     if (fWriteESDfriend) {
 854       esdf->~AliESDfriend();
 855       new (esdf) AliESDfriend(); // Reset...
 856       esd->GetESDfriend(esdf);
 857     }
 858     tree->Fill();
 859
 860     // write HLT ESD
 861     hlttree->Fill();
 862
 863     if (fCheckPointLevel > 0)  WriteESD(esd, "final");
 864     esd->Reset();
 865     hltesd->Reset();
 866     if (fWriteESDfriend) {
 867       esdf->~AliESDfriend();
 868       new (esdf) AliESDfriend(); // Reset...
 869     }
 870     // esdf->Reset();
 871     // delete esdf; esdf = 0;
 872   // ESD QA
 873
 874     gSystem->GetProcInfo(&ProcInfo);
 875     AliInfo(Form("Event %d -> Current memory usage %d %d",iEvent, ProcInfo.fMemResident, ProcInfo.fMemVirtual));
 876   }
 877
 878   detStr = fFillESD ;
 879   // write quality assurance ESDs data (one entry for all events)
 880   for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
 881         if (!IsSelected(fgkDetectorName[iDet], detStr))
 882                 continue;
 883     AliQADataMaker * qadm = GetQADataMaker(iDet);
 884     if (!qadm) continue;
 885     qadm->EndOfCycle(AliQA::kRECPOINTS);
 886     qadm->EndOfCycle(AliQA::kESDS);
 887     qadm->Finish(AliQA::kRECPOINTS);
 888     qadm->Finish(AliQA::kESDS) ;
 889   }
 890
 891   tree->GetUserInfo()->Add(esd);
 892   hlttree->GetUserInfo()->Add(hltesd);
 893
 894
 895
 896   if(fESDPar.Contains("ESD.par")){
 897     AliInfo("Attaching ESD.par to Tree");
 898     TNamed *fn = CopyFileToTNamed(fESDPar.Data(),"ESD.par");
 899     tree->GetUserInfo()->Add(fn);
 900   }
 901
 902
 903   file->cd();
 904   if (fWriteESDfriend)
 905     tree->SetBranchStatus("ESDfriend*",0);
 906   // we want to have only one tree version number
 907   tree->Write(tree->GetName(),TObject::kOverwrite);
 908   hlttree->Write();
 909
 910   if (fWriteAOD) {
 911     TFile *aodFile = TFile::Open("AliAOD.root", "RECREATE");
 912     ESDFile2AODFile(file, aodFile);
 913     aodFile->Close();
 914   }
 915
 916   gROOT->cd();
 917   CleanUp(file, fileOld);
 918
 919   // Create tags for the events in the ESD tree (the ESD tree is always present)
 920   // In case of empty events the tags will contain dummy values
 921   AliESDTagCreator *esdtagCreator = new AliESDTagCreator();
 922   esdtagCreator->CreateESDTags(fFirstEvent,fLastEvent);
 923   if (fWriteAOD) {
 924     AliAODTagCreator *aodtagCreator = new AliAODTagCreator();
 925     aodtagCreator->CreateAODTags(fFirstEvent,fLastEvent);
 926   }
 927
 928   return kTRUE;
 929 }
 930
 931
 932 //_____________________________________________________________________________
 933 Bool_t AliReconstruction::RunLocalReconstruction(const TString& detectors)
 934 {
 935 // run the local reconstruction
 936   static Int_t eventNr=0;
 937   AliCodeTimerAuto("")
 938
 939  //  AliCDBManager* man = AliCDBManager::Instance();
 940 //   Bool_t origCache = man->GetCacheFlag();
 941
 942 //   TString detStr = detectors;
 943 //   for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
 944 //     if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
 945 //     AliReconstructor* reconstructor = GetReconstructor(iDet);
 946 //     if (!reconstructor) continue;
 947 //     if (reconstructor->HasLocalReconstruction()) continue;
 948
 949 //     AliCodeTimerStart(Form("running reconstruction for %s", fgkDetectorName[iDet]));
 950 //     AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
 951
 952 //     AliCodeTimerStart(Form("Loading calibration data from OCDB for %s", fgkDetectorName[iDet]));
 953 //     AliInfo(Form("Loading calibration data from OCDB for %s", fgkDetectorName[iDet]));
 954
 955 //     man->SetCacheFlag(kTRUE);
 956 //     TString calibPath = Form("%s/Calib/*", fgkDetectorName[iDet]);
 957 //     man->GetAll(calibPath); // entries are cached!
 958
 959 //     AliCodeTimerStop(Form("Loading calibration data from OCDB for %s", fgkDetectorName[iDet]));
 960
 961 //     if (fRawReader) {
 962 //       fRawReader->RewindEvents();
 963 //       reconstructor->Reconstruct(fRunLoader, fRawReader);
 964 //     } else {
 965 //       reconstructor->Reconstruct(fRunLoader);
 966 //     }
 967
 968 //      AliCodeTimerStop(Form("running reconstruction for %s", fgkDetectorName[iDet]));
 969     // AliSysInfo::AddStamp(Form("LRec%s_%d",fgkDetectorName[iDet],eventNr));
 970
 971 //     // unload calibration data
 972 //     man->UnloadFromCache(calibPath);
 973 //     //man->ClearCache();
 974 //   }
 975
 976 //   man->SetCacheFlag(origCache);
 977
 978 //   if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
 979 //     AliError(Form("the following detectors were not found: %s",
 980 //                   detStr.Data()));
 981 //     if (fStopOnError) return kFALSE;
 982 //   }
 983
 984           eventNr++;
 985   return kTRUE;
 986 }
 987
 988 //_____________________________________________________________________________
 989 Bool_t AliReconstruction::RunLocalEventReconstruction(const TString& detectors)
 990 {
 991 // run the local reconstruction
 992   static Int_t eventNr=0;
 993   AliCodeTimerAuto("")
 994
 995   TString detStr = detectors;
 996   for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
 997     if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
 998     AliReconstructor* reconstructor = GetReconstructor(iDet);
 999     if (!reconstructor) continue;
1000     AliLoader* loader = fLoader[iDet];
1001     if (!loader) {
1002       AliWarning(Form("No loader is defined for %s!",fgkDetectorName[iDet]));
1003       continue;
1004     }
1005
1006     // conversion of digits
1007     if (fRawReader && reconstructor->HasDigitConversion()) {
1008       AliInfo(Form("converting raw data digits into root objects for %s",
1009                    fgkDetectorName[iDet]));
1010       AliCodeTimerAuto(Form("converting raw data digits into root objects for %s",
1011                             fgkDetectorName[iDet]));
1012       loader->LoadDigits("update");
1013       loader->CleanDigits();
1014       loader->MakeDigitsContainer();
1015       TTree* digitsTree = loader->TreeD();
1016       reconstructor->ConvertDigits(fRawReader, digitsTree);
1017       loader->WriteDigits("OVERWRITE");
1018       loader->UnloadDigits();
1019     }
1020
1021     // local reconstruction
1022     AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
1023     AliCodeTimerAuto(Form("running reconstruction for %s", fgkDetectorName[iDet]));
1024     loader->LoadRecPoints("update");
1025     loader->CleanRecPoints();
1026     loader->MakeRecPointsContainer();
1027     TTree* clustersTree = loader->TreeR();
1028     if (fRawReader && !reconstructor->HasDigitConversion()) {
1029       reconstructor->Reconstruct(fRawReader, clustersTree);
1030     } else {
1031       loader->LoadDigits("read");
1032       TTree* digitsTree = loader->TreeD();
1033       if (!digitsTree) {
1034         AliError(Form("Can't get the %s digits tree", fgkDetectorName[iDet]));
1035         if (fStopOnError) return kFALSE;
1036       } else {
1037         reconstructor->Reconstruct(digitsTree, clustersTree);
1038       }
1039       loader->UnloadDigits();
1040     }
1041
1042     AliQADataMaker * qadm = GetQADataMaker(iDet);
1043     if (qadm) {
1044       AliCodeTimerStart(Form("running quality assurance data maker for %s", fgkDetectorName[iDet]));
1045       AliInfo(Form("running quality assurance data maker for %s", fgkDetectorName[iDet]));
1046
1047      if (qadm->IsCycleDone() ) {
1048       qadm->EndOfCycle(AliQA::kRECPOINTS) ;
1049           qadm->EndOfCycle(AliQA::kESDS) ;
1050       qadm->StartOfCycle(AliQA::kRECPOINTS) ;
1051           qadm->StartOfCycle(AliQA::kESDS, "same") ;
1052      }
1053       qadm->Exec(AliQA::kRECPOINTS, clustersTree) ;
1054       AliCodeTimerStop(Form("running quality assurance data maker for %s", fgkDetectorName[iDet]));
1055     }
1056
1057     loader->WriteRecPoints("OVERWRITE");
1058     loader->UnloadRecPoints();
1059      AliSysInfo::AddStamp(Form("LRec%s_%d",fgkDetectorName[iDet],eventNr));
1060   }
1061
1062   if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1063     AliError(Form("the following detectors were not found: %s",
1064                   detStr.Data()));
1065     if (fStopOnError) return kFALSE;
1066   }
1067   eventNr++;
1068   return kTRUE;
1069 }
1070
1071 //_____________________________________________________________________________
1072 Bool_t AliReconstruction::RunVertexFinder(AliESDEvent*& esd)
1073 {
1074 // run the barrel tracking
1075
1076   AliCodeTimerAuto("")
1077
1078   AliESDVertex* vertex = NULL;
1079   Double_t vtxPos[3] = {0, 0, 0};
1080   Double_t vtxErr[3] = {0.07, 0.07, 0.1};
1081   TArrayF mcVertex(3);
1082   if (fRunLoader->GetHeader() && fRunLoader->GetHeader()->GenEventHeader()) {
1083     fRunLoader->GetHeader()->GenEventHeader()->PrimaryVertex(mcVertex);
1084     for (Int_t i = 0; i < 3; i++) vtxPos[i] = mcVertex[i];
1085   }
1086
1087   if (fVertexer) {
1088     if(fDiamondProfile) fVertexer->SetVtxStart(fDiamondProfile);
1089     AliInfo("running the ITS vertex finder");
1090     if (fLoader[0]) fLoader[0]->LoadRecPoints();
1091     vertex = fVertexer->FindVertexForCurrentEvent(fRunLoader->GetEventNumber());
1092     if (fLoader[0]) fLoader[0]->UnloadRecPoints();
1093     if(!vertex){
1094       AliWarning("Vertex not found");
1095       vertex = new AliESDVertex();
1096       vertex->SetName("default");
1097     }
1098     else {
1099       vertex->SetName("reconstructed");
1100     }
1101
1102   } else {
1103     AliInfo("getting the primary vertex from MC");
1104     vertex = new AliESDVertex(vtxPos, vtxErr);
1105   }
1106
1107   if (vertex) {
1108     vertex->GetXYZ(vtxPos);
1109     vertex->GetSigmaXYZ(vtxErr);
1110   } else {
1111     AliWarning("no vertex reconstructed");
1112     vertex = new AliESDVertex(vtxPos, vtxErr);
1113   }
1114   esd->SetVertex(vertex);
1115   // if SPD multiplicity has been determined, it is stored in the ESD
1116   AliMultiplicity *mult = fVertexer->GetMultiplicity();
1117   if(mult)esd->SetMultiplicity(mult);
1118
1119   for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1120     if (fTracker[iDet]) fTracker[iDet]->SetVertex(vtxPos, vtxErr);
1121   }
1122   delete vertex;
1123
1124   return kTRUE;
1125 }
1126
1127 //_____________________________________________________________________________
1128 Bool_t AliReconstruction::RunHLTTracking(AliESDEvent*& esd)
1129 {
1130 // run the HLT barrel tracking
1131
1132   AliCodeTimerAuto("")
1133
1134   if (!fRunLoader) {
1135     AliError("Missing runLoader!");
1136     return kFALSE;
1137   }
1138
1139   AliInfo("running HLT tracking");
1140
1141   // Get a pointer to the HLT reconstructor
1142   AliReconstructor *reconstructor = GetReconstructor(fgkNDetectors-1);
1143   if (!reconstructor) return kFALSE;
1144
1145   // TPC + ITS
1146   for (Int_t iDet = 1; iDet >= 0; iDet--) {
1147     TString detName = fgkDetectorName[iDet];
1148     AliDebug(1, Form("%s HLT tracking", detName.Data()));
1149     reconstructor->SetOption(detName.Data());
1150     AliTracker *tracker = reconstructor->CreateTracker();
1151     if (!tracker) {
1152       AliWarning(Form("couldn't create a HLT tracker for %s", detName.Data()));
1153       if (fStopOnError) return kFALSE;
1154       continue;
1155     }
1156     Double_t vtxPos[3];
1157     Double_t vtxErr[3]={0.005,0.005,0.010};
1158     const AliESDVertex *vertex = esd->GetVertex();
1159     vertex->GetXYZ(vtxPos);
1160     tracker->SetVertex(vtxPos,vtxErr);
1161     if(iDet != 1) {
1162       fLoader[iDet]->LoadRecPoints("read");
1163       TTree* tree = fLoader[iDet]->TreeR();
1164       if (!tree) {
1165         AliError(Form("Can't get the %s cluster tree", detName.Data()));
1166         return kFALSE;
1167       }
1168       tracker->LoadClusters(tree);
1169     }
1170     if (tracker->Clusters2Tracks(esd) != 0) {
1171       AliError(Form("HLT %s Clusters2Tracks failed", fgkDetectorName[iDet]));
1172       return kFALSE;
1173     }
1174     if(iDet != 1) {
1175       tracker->UnloadClusters();
1176     }
1177     delete tracker;
1178   }
1179
1180   return kTRUE;
1181 }
1182
1183 //_____________________________________________________________________________
1184 Bool_t AliReconstruction::RunMuonTracking(AliESDEvent*& esd)
1185 {
1186 // run the muon spectrometer tracking
1187
1188   AliCodeTimerAuto("")
1189
1190   if (!fRunLoader) {
1191     AliError("Missing runLoader!");
1192     return kFALSE;
1193   }
1194   Int_t iDet = 7; // for MUON
1195
1196   AliInfo("is running...");
1197
1198   // Get a pointer to the MUON reconstructor
1199   AliReconstructor *reconstructor = GetReconstructor(iDet);
1200   if (!reconstructor) return kFALSE;
1201
1202
1203   TString detName = fgkDetectorName[iDet];
1204   AliDebug(1, Form("%s tracking", detName.Data()));
1205   AliTracker *tracker =  reconstructor->CreateTracker();
1206   if (!tracker) {
1207     AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
1208     return kFALSE;
1209   }
1210
1211   // create Tracks
1212   fLoader[iDet]->LoadTracks("update");
1213   fLoader[iDet]->CleanTracks();
1214   fLoader[iDet]->MakeTracksContainer();
1215
1216   // read RecPoints
1217   fLoader[iDet]->LoadRecPoints("read");
1218   tracker->LoadClusters(fLoader[iDet]->TreeR());
1219
1220   Int_t rv = tracker->Clusters2Tracks(esd);
1221
1222   fLoader[iDet]->UnloadRecPoints();
1223
1224   if ( rv )
1225   {
1226     AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
1227     return kFALSE;
1228   }
1229
1230   tracker->UnloadClusters();
1231
1232   fLoader[iDet]->UnloadRecPoints();
1233
1234   fLoader[iDet]->WriteTracks("OVERWRITE");
1235   fLoader[iDet]->UnloadTracks();
1236
1237   delete tracker;
1238
1239   return kTRUE;
1240 }
1241
1242
1243 //_____________________________________________________________________________
1244 Bool_t AliReconstruction::RunTracking(AliESDEvent*& esd)
1245 {
1246 // run the barrel tracking
1247   static Int_t eventNr=0;
1248   AliCodeTimerAuto("")
1249
1250   AliInfo("running tracking");
1251
1252   //Fill the ESD with the T0 info (will be used by the TOF)
1253   if (fReconstructor[11] && fLoader[11]) {
1254     fLoader[11]->LoadRecPoints("READ");
1255     TTree *treeR = fLoader[11]->TreeR();
1256     GetReconstructor(11)->FillESD((TTree *)NULL,treeR,esd);
1257   }
1258
1259   // pass 1: TPC + ITS inwards
1260   for (Int_t iDet = 1; iDet >= 0; iDet--) {
1261     if (!fTracker[iDet]) continue;
1262     AliDebug(1, Form("%s tracking", fgkDetectorName[iDet]));
1263
1264     // load clusters
1265     fLoader[iDet]->LoadRecPoints("read");
1266     TTree* tree = fLoader[iDet]->TreeR();
1267     if (!tree) {
1268       AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
1269       return kFALSE;
1270     }
1271     fTracker[iDet]->LoadClusters(tree);
1272
1273     // run tracking
1274     if (fTracker[iDet]->Clusters2Tracks(esd) != 0) {
1275       AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
1276       return kFALSE;
1277     }
1278     if (fCheckPointLevel > 1) {
1279       WriteESD(esd, Form("%s.tracking", fgkDetectorName[iDet]));
1280     }
1281     // preliminary PID in TPC needed by the ITS tracker
1282     if (iDet == 1) {
1283       GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
1284       AliESDpid::MakePID(esd);
1285     }
1286      AliSysInfo::AddStamp(Form("Tracking0%s_%d",fgkDetectorName[iDet],eventNr));
1287   }
1288
1289   // pass 2: ALL backwards
1290   for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1291     if (!fTracker[iDet]) continue;
1292     AliDebug(1, Form("%s back propagation", fgkDetectorName[iDet]));
1293
1294     // load clusters
1295     if (iDet > 1) {     // all except ITS, TPC
1296       TTree* tree = NULL;
1297       fLoader[iDet]->LoadRecPoints("read");
1298       tree = fLoader[iDet]->TreeR();
1299       if (!tree) {
1300         AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
1301         return kFALSE;
1302       }
1303       fTracker[iDet]->LoadClusters(tree);
1304        AliSysInfo::AddStamp(Form("LoadCluster0%s_%d",fgkDetectorName[iDet],eventNr));
1305     }
1306
1307     // run tracking
1308     if (fTracker[iDet]->PropagateBack(esd) != 0) {
1309       AliError(Form("%s backward propagation failed", fgkDetectorName[iDet]));
1310       //      return kFALSE;
1311     }
1312     if (fCheckPointLevel > 1) {
1313       WriteESD(esd, Form("%s.back", fgkDetectorName[iDet]));
1314     }
1315
1316     // unload clusters
1317     if (iDet > 2) {     // all except ITS, TPC, TRD
1318       fTracker[iDet]->UnloadClusters();
1319       fLoader[iDet]->UnloadRecPoints();
1320     }
1321     // updated PID in TPC needed by the ITS tracker -MI
1322     if (iDet == 1) {
1323       GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
1324       AliESDpid::MakePID(esd);
1325     }
1326      AliSysInfo::AddStamp(Form("Tracking1%s_%d",fgkDetectorName[iDet],eventNr));
1327   }
1328
1329   // write space-points to the ESD in case alignment data output
1330   // is switched on
1331   if (fWriteAlignmentData)
1332     WriteAlignmentData(esd);
1333
1334   // pass 3: TRD + TPC + ITS refit inwards
1335   for (Int_t iDet = 2; iDet >= 0; iDet--) {
1336     if (!fTracker[iDet]) continue;
1337     AliDebug(1, Form("%s inward refit", fgkDetectorName[iDet]));
1338
1339     // run tracking
1340     if (fTracker[iDet]->RefitInward(esd) != 0) {
1341       AliError(Form("%s inward refit failed", fgkDetectorName[iDet]));
1342       //      return kFALSE;
1343     }
1344     if (fCheckPointLevel > 1) {
1345       WriteESD(esd, Form("%s.refit", fgkDetectorName[iDet]));
1346     }
1347
1348     // unload clusters
1349     fTracker[iDet]->UnloadClusters();
1350     fLoader[iDet]->UnloadRecPoints();
1351      AliSysInfo::AddStamp(Form("Tracking2%s_%d",fgkDetectorName[iDet],eventNr));
1352   }
1353   //
1354   // Propagate track to the vertex - if not done by ITS
1355   //
1356   Int_t ntracks = esd->GetNumberOfTracks();
1357   for (Int_t itrack=0; itrack<ntracks; itrack++){
1358     const Double_t kRadius  = 3;   // beam pipe radius
1359     const Double_t kMaxStep = 5;   // max step
1360     const Double_t kMaxD    = 123456;  // max distance to prim vertex
1361     Double_t       fieldZ   = AliTracker::GetBz();  //
1362     AliESDtrack * track = esd->GetTrack(itrack);
1363     if (!track) continue;
1364     if (track->IsOn(AliESDtrack::kITSrefit)) continue;
1365    AliTracker::PropagateTrackTo(track,kRadius,track->GetMass(),kMaxStep,kTRUE);
1366     track->RelateToVertex(esd->GetVertex(),fieldZ, kMaxD);
1367   }
1368   eventNr++;
1369   return kTRUE;
1370 }
1371
1372 //_____________________________________________________________________________
1373 Bool_t AliReconstruction::CleanESD(AliESDEvent *esd){
1374   //
1375   // Remove the data which are not needed for the physics analysis.
1376   //
1377
1378   AliInfo("Cleaning the ESD...");
1379   Int_t nTracks=esd->GetNumberOfTracks();
1380   AliInfo(Form("Number of ESD tracks before cleaning %d",nTracks));
1381
1382   Float_t cleanPars[]={fDmax,fZmax};
1383   Bool_t rc=esd->Clean(cleanPars);
1384
1385   nTracks=esd->GetNumberOfTracks();
1386   AliInfo(Form("Number of ESD tracks after cleaning %d",nTracks));
1387
1388   return rc;
1389 }
1390
1391 //_____________________________________________________________________________
1392 Bool_t AliReconstruction::FillESD(AliESDEvent*& esd, const TString& detectors)
1393 {
1394 // fill the event summary data
1395
1396   AliCodeTimerAuto("")
1397     static Int_t eventNr=0;
1398   TString detStr = detectors;
1399   for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1400     if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1401     AliReconstructor* reconstructor = GetReconstructor(iDet);
1402     if (!reconstructor) continue;
1403
1404     if (!ReadESD(esd, fgkDetectorName[iDet])) {
1405       AliDebug(1, Form("filling ESD for %s", fgkDetectorName[iDet]));
1406       TTree* clustersTree = NULL;
1407       if (fLoader[iDet]) {
1408         fLoader[iDet]->LoadRecPoints("read");
1409         clustersTree = fLoader[iDet]->TreeR();
1410         if (!clustersTree) {
1411           AliError(Form("Can't get the %s clusters tree",
1412                         fgkDetectorName[iDet]));
1413           if (fStopOnError) return kFALSE;
1414         }
1415       }
1416       if (fRawReader && !reconstructor->HasDigitConversion()) {
1417         reconstructor->FillESD(fRawReader, clustersTree, esd);
1418       } else {
1419         TTree* digitsTree = NULL;
1420         if (fLoader[iDet]) {
1421           fLoader[iDet]->LoadDigits("read");
1422           digitsTree = fLoader[iDet]->TreeD();
1423           if (!digitsTree) {
1424             AliError(Form("Can't get the %s digits tree",
1425                           fgkDetectorName[iDet]));
1426             if (fStopOnError) return kFALSE;
1427           }
1428         }
1429         reconstructor->FillESD(digitsTree, clustersTree, esd);
1430         if (fLoader[iDet]) fLoader[iDet]->UnloadDigits();
1431       }
1432       if (fLoader[iDet]) {
1433         fLoader[iDet]->UnloadRecPoints();
1434       }
1435
1436       if (fCheckPointLevel > 2) WriteESD(esd, fgkDetectorName[iDet]);
1437     }
1438   }
1439
1440   if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1441     AliError(Form("the following detectors were not found: %s",
1442                   detStr.Data()));
1443     if (fStopOnError) return kFALSE;
1444   }
1445    AliSysInfo::AddStamp(Form("FillESD%d",eventNr));
1446   eventNr++;
1447   return kTRUE;
1448 }
1449
1450 //_____________________________________________________________________________
1451 Bool_t AliReconstruction::FillTriggerESD(AliESDEvent*& esd)
1452 {
1453   // Reads the trigger decision which is
1454   // stored in Trigger.root file and fills
1455   // the corresponding esd entries
1456
1457   AliCodeTimerAuto("")
1458
1459   AliInfo("Filling trigger information into the ESD");
1460
1461   if (fRawReader) {
1462     AliCTPRawStream input(fRawReader);
1463     if (!input.Next()) {
1464       AliError("No valid CTP (trigger) DDL raw data is found ! The trigger information is not stored in the ESD !");
1465       return kFALSE;
1466     }
1467     esd->SetTriggerMask(input.GetClassMask());
1468     esd->SetTriggerCluster(input.GetClusterMask());
1469   }
1470   else {
1471     AliRunLoader *runloader = AliRunLoader::GetRunLoader();
1472     if (runloader) {
1473       if (!runloader->LoadTrigger()) {
1474         AliCentralTrigger *aCTP = runloader->GetTrigger();
1475         esd->SetTriggerMask(aCTP->GetClassMask());
1476         esd->SetTriggerCluster(aCTP->GetClusterMask());
1477       }
1478       else {
1479         AliWarning("No trigger can be loaded! The trigger information is not stored in the ESD !");
1480         return kFALSE;
1481       }
1482     }
1483     else {
1484       AliError("No run loader is available! The trigger information is not stored in the ESD !");
1485       return kFALSE;
1486     }
1487   }
1488
1489   return kTRUE;
1490 }
1491
1492
1493
1494
1495
1496 //_____________________________________________________________________________
1497 Bool_t AliReconstruction::FillRawEventHeaderESD(AliESDEvent*& esd)
1498 {
1499   //
1500   // Filling information from RawReader Header
1501   //
1502
1503   AliInfo("Filling information from RawReader Header");
1504   esd->SetBunchCrossNumber(0);
1505   esd->SetOrbitNumber(0);
1506   esd->SetPeriodNumber(0);
1507   esd->SetTimeStamp(0);
1508   esd->SetEventType(0);
1509   const AliRawEventHeaderBase * eventHeader = fRawReader->GetEventHeader();
1510   if (eventHeader){
1511
1512     const UInt_t *id = eventHeader->GetP("Id");
1513     esd->SetBunchCrossNumber((id)[1]&0x00000fff);
1514     esd->SetOrbitNumber((((id)[0]<<20)&0xf00000)|(((id)[1]>>12)&0xfffff));
1515     esd->SetPeriodNumber(((id)[0]>>4)&0x0fffffff);
1516
1517     esd->SetTimeStamp((eventHeader->Get("Timestamp")));
1518     esd->SetEventType((eventHeader->Get("Type")));
1519   }
1520
1521   return kTRUE;
1522 }
1523
1524
1525 //_____________________________________________________________________________
1526 Bool_t AliReconstruction::IsSelected(TString detName, TString& detectors) const
1527 {
1528 // check whether detName is contained in detectors
1529 // if yes, it is removed from detectors
1530
1531   // check if all detectors are selected
1532   if ((detectors.CompareTo("ALL") == 0) ||
1533       detectors.BeginsWith("ALL ") ||
1534       detectors.EndsWith(" ALL") ||
1535       detectors.Contains(" ALL ")) {
1536     detectors = "ALL";
1537     return kTRUE;
1538   }
1539
1540   // search for the given detector
1541   Bool_t result = kFALSE;
1542   if ((detectors.CompareTo(detName) == 0) ||
1543       detectors.BeginsWith(detName+" ") ||
1544       detectors.EndsWith(" "+detName) ||
1545       detectors.Contains(" "+detName+" ")) {
1546     detectors.ReplaceAll(detName, "");
1547     result = kTRUE;
1548   }
1549
1550   // clean up the detectors string
1551   while (detectors.Contains("  ")) detectors.ReplaceAll("  ", " ");
1552   while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
1553   while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
1554
1555   return result;
1556 }
1557
1558 //_____________________________________________________________________________
1559 Bool_t AliReconstruction::InitRunLoader()
1560 {
1561 // get or create the run loader
1562
1563   if (gAlice) delete gAlice;
1564   gAlice = NULL;
1565
1566   if (!gSystem->AccessPathName(fGAliceFileName.Data())) { // galice.root exists
1567     // load all base libraries to get the loader classes
1568     TString libs = gSystem->GetLibraries();
1569     for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1570       TString detName = fgkDetectorName[iDet];
1571       if (detName == "HLT") continue;
1572       if (libs.Contains("lib" + detName + "base.so")) continue;
1573       gSystem->Load("lib" + detName + "base.so");
1574     }
1575     fRunLoader = AliRunLoader::Open(fGAliceFileName.Data());
1576     if (!fRunLoader) {
1577       AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
1578       CleanUp();
1579       return kFALSE;
1580     }
1581     fRunLoader->CdGAFile();
1582     if (gFile->GetKey(AliRunLoader::GetGAliceName())) {
1583       if (fRunLoader->LoadgAlice() == 0) {
1584         gAlice = fRunLoader->GetAliRun();
1585         AliTracker::SetFieldMap(gAlice->Field(),fUniformField);
1586       }
1587     }
1588     if (!gAlice && !fRawReader) {
1589       AliError(Form("no gAlice object found in file %s",
1590                     fGAliceFileName.Data()));
1591       CleanUp();
1592       return kFALSE;
1593     }
1594
1595     //PH This is a temporary fix to give access to the kinematics
1596     //PH that is needed for the labels of ITS clusters
1597     fRunLoader->LoadKinematics();
1598
1599   } else {               // galice.root does not exist
1600     if (!fRawReader) {
1601       AliError(Form("the file %s does not exist", fGAliceFileName.Data()));
1602       CleanUp();
1603       return kFALSE;
1604     }
1605     fRunLoader = AliRunLoader::Open(fGAliceFileName.Data(),
1606                                     AliConfig::GetDefaultEventFolderName(),
1607                                     "recreate");
1608     if (!fRunLoader) {
1609       AliError(Form("could not create run loader in file %s",
1610                     fGAliceFileName.Data()));
1611       CleanUp();
1612       return kFALSE;
1613     }
1614     fRunLoader->MakeTree("E");
1615     Int_t iEvent = 0;
1616     while (fRawReader->NextEvent()) {
1617       fRunLoader->SetEventNumber(iEvent);
1618       fRunLoader->GetHeader()->Reset(fRawReader->GetRunNumber(),
1619                                      iEvent, iEvent);
1620       fRunLoader->MakeTree("H");
1621       fRunLoader->TreeE()->Fill();
1622       iEvent++;
1623     }
1624     fRawReader->RewindEvents();
1625     if (fNumberOfEventsPerFile > 0)
1626       fRunLoader->SetNumberOfEventsPerFile(fNumberOfEventsPerFile);
1627     else
1628       fRunLoader->SetNumberOfEventsPerFile(iEvent);
1629     fRunLoader->WriteHeader("OVERWRITE");
1630     fRunLoader->CdGAFile();
1631     fRunLoader->Write(0, TObject::kOverwrite);
1632 //    AliTracker::SetFieldMap(???);
1633   }
1634
1635   return kTRUE;
1636 }
1637
1638 //_____________________________________________________________________________
1639 AliReconstructor* AliReconstruction::GetReconstructor(Int_t iDet)
1640 {
1641 // get the reconstructor object and the loader for a detector
1642
1643   if (fReconstructor[iDet]) return fReconstructor[iDet];
1644
1645   // load the reconstructor object
1646   TPluginManager* pluginManager = gROOT->GetPluginManager();
1647   TString detName = fgkDetectorName[iDet];
1648   TString recName = "Ali" + detName + "Reconstructor";
1649   if (gAlice && !gAlice->GetDetector(detName) && (detName != "HLT")) return NULL;
1650
1651   AliReconstructor* reconstructor = NULL;
1652   // first check if a plugin is defined for the reconstructor
1653   TPluginHandler* pluginHandler =
1654     pluginManager->FindHandler("AliReconstructor", detName);
1655   // if not, add a plugin for it
1656   if (!pluginHandler) {
1657     AliDebug(1, Form("defining plugin for %s", recName.Data()));
1658     TString libs = gSystem->GetLibraries();
1659     if (libs.Contains("lib" + detName + "base.so") ||
1660         (gSystem->Load("lib" + detName + "base.so") >= 0)) {
1661       pluginManager->AddHandler("AliReconstructor", detName,
1662                                 recName, detName + "rec", recName + "()");
1663     } else {
1664       pluginManager->AddHandler("AliReconstructor", detName,
1665                                 recName, detName, recName + "()");
1666     }
1667     pluginHandler = pluginManager->FindHandler("AliReconstructor", detName);
1668   }
1669   if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
1670     reconstructor = (AliReconstructor*) pluginHandler->ExecPlugin(0);
1671   }
1672   if (reconstructor) {
1673     TObject* obj = fOptions.FindObject(detName.Data());
1674     if (obj) reconstructor->SetOption(obj->GetTitle());
1675     reconstructor->Init();
1676     fReconstructor[iDet] = reconstructor;
1677   }
1678
1679   // get or create the loader
1680   if (detName != "HLT") {
1681     fLoader[iDet] = fRunLoader->GetLoader(detName + "Loader");
1682     if (!fLoader[iDet]) {
1683       AliConfig::Instance()
1684         ->CreateDetectorFolders(fRunLoader->GetEventFolder(),
1685                                 detName, detName);
1686       // first check if a plugin is defined for the loader
1687       pluginHandler =
1688         pluginManager->FindHandler("AliLoader", detName);
1689       // if not, add a plugin for it
1690       if (!pluginHandler) {
1691         TString loaderName = "Ali" + detName + "Loader";
1692         AliDebug(1, Form("defining plugin for %s", loaderName.Data()));
1693         pluginManager->AddHandler("AliLoader", detName,
1694                                   loaderName, detName + "base",
1695                                   loaderName + "(const char*, TFolder*)");
1696         pluginHandler = pluginManager->FindHandler("AliLoader", detName);
1697       }
1698       if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
1699         fLoader[iDet] =
1700           (AliLoader*) pluginHandler->ExecPlugin(2, detName.Data(),
1701                                                  fRunLoader->GetEventFolder());
1702       }
1703       if (!fLoader[iDet]) {   // use default loader
1704         fLoader[iDet] = new AliLoader(detName, fRunLoader->GetEventFolder());
1705       }
1706       if (!fLoader[iDet]) {
1707         AliWarning(Form("couldn't get loader for %s", detName.Data()));
1708         if (fStopOnError) return NULL;
1709       } else {
1710         fRunLoader->AddLoader(fLoader[iDet]);
1711         fRunLoader->CdGAFile();
1712         if (gFile && !gFile->IsWritable()) gFile->ReOpen("UPDATE");
1713         fRunLoader->Write(0, TObject::kOverwrite);
1714       }
1715     }
1716   }
1717
1718   return reconstructor;
1719 }
1720
1721 //_____________________________________________________________________________
1722 Bool_t AliReconstruction::CreateVertexer()
1723 {
1724 // create the vertexer
1725
1726   fVertexer = NULL;
1727   AliReconstructor* itsReconstructor = GetReconstructor(0);
1728   if (itsReconstructor) {
1729     fVertexer = itsReconstructor->CreateVertexer();
1730   }
1731   if (!fVertexer) {
1732     AliWarning("couldn't create a vertexer for ITS");
1733     if (fStopOnError) return kFALSE;
1734   }
1735
1736   return kTRUE;
1737 }
1738
1739 //_____________________________________________________________________________
1740 Bool_t AliReconstruction::CreateTrackers(const TString& detectors)
1741 {
1742 // create the trackers
1743
1744   TString detStr = detectors;
1745   for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1746     if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1747     AliReconstructor* reconstructor = GetReconstructor(iDet);
1748     if (!reconstructor) continue;
1749     TString detName = fgkDetectorName[iDet];
1750     if (detName == "HLT") {
1751       fRunHLTTracking = kTRUE;
1752       continue;
1753     }
1754     if (detName == "MUON") {
1755       fRunMuonTracking = kTRUE;
1756       continue;
1757     }
1758
1759
1760     fTracker[iDet] = reconstructor->CreateTracker();
1761     if (!fTracker[iDet] && (iDet < 7)) {
1762       AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
1763       if (fStopOnError) return kFALSE;
1764     }
1765   }
1766
1767   return kTRUE;
1768 }
1769
1770 //_____________________________________________________________________________
1771 void AliReconstruction::CleanUp(TFile* file, TFile* fileOld)
1772 {
1773 // delete trackers and the run loader and close and delete the file
1774
1775   for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1776     delete fReconstructor[iDet];
1777     fReconstructor[iDet] = NULL;
1778     fLoader[iDet] = NULL;
1779     delete fTracker[iDet];
1780     fTracker[iDet] = NULL;
1781     delete fQADataMaker[iDet];
1782     fQADataMaker[iDet] = NULL;
1783   }
1784   delete fVertexer;
1785   fVertexer = NULL;
1786   delete fDiamondProfile;
1787   fDiamondProfile = NULL;
1788
1789   delete fRunLoader;
1790   fRunLoader = NULL;
1791   delete fRawReader;
1792   fRawReader = NULL;
1793
1794   if (file) {
1795     file->Close();
1796     delete file;
1797   }
1798
1799   if (fileOld) {
1800     fileOld->Close();
1801     delete fileOld;
1802     gSystem->Unlink("AliESDs.old.root");
1803   }
1804 }
1805
1806 //_____________________________________________________________________________
1807
1808 Bool_t AliReconstruction::ReadESD(AliESDEvent*& esd, const char* recStep) const
1809 {
1810 // read the ESD event from a file
1811
1812   if (!esd) return kFALSE;
1813   char fileName[256];
1814   sprintf(fileName, "ESD_%d.%d_%s.root",
1815           esd->GetRunNumber(), esd->GetEventNumberInFile(), recStep);
1816   if (gSystem->AccessPathName(fileName)) return kFALSE;
1817
1818   AliInfo(Form("reading ESD from file %s", fileName));
1819   AliDebug(1, Form("reading ESD from file %s", fileName));
1820   TFile* file = TFile::Open(fileName);
1821   if (!file || !file->IsOpen()) {
1822     AliError(Form("opening %s failed", fileName));
1823     delete file;
1824     return kFALSE;
1825   }
1826
1827   gROOT->cd();
1828   delete esd;
1829   esd = (AliESDEvent*) file->Get("ESD");
1830   file->Close();
1831   delete file;
1832   return kTRUE;
1833
1834 }
1835
1836
1837
1838 //_____________________________________________________________________________
1839 void AliReconstruction::WriteESD(AliESDEvent* esd, const char* recStep) const
1840 {
1841 // write the ESD event to a file
1842
1843   if (!esd) return;
1844   char fileName[256];
1845   sprintf(fileName, "ESD_%d.%d_%s.root",
1846           esd->GetRunNumber(), esd->GetEventNumberInFile(), recStep);
1847
1848   AliDebug(1, Form("writing ESD to file %s", fileName));
1849   TFile* file = TFile::Open(fileName, "recreate");
1850   if (!file || !file->IsOpen()) {
1851     AliError(Form("opening %s failed", fileName));
1852   } else {
1853     esd->Write("ESD");
1854     file->Close();
1855   }
1856   delete file;
1857 }
1858
1859
1860
1861
1862
1863 //_____________________________________________________________________________
1864 void AliReconstruction::ESDFile2AODFile(TFile* esdFile, TFile* aodFile)
1865 {
1866   // write all files from the given esd file to an aod file
1867
1868   // create an AliAOD object
1869   AliAODEvent *aod = new AliAODEvent();
1870   aod->CreateStdContent();
1871
1872   // go to the file
1873   aodFile->cd();
1874
1875   // create the tree
1876   TTree *aodTree = new TTree("aodTree", "AliAOD tree");
1877   aodTree->Branch(aod->GetList());
1878
1879   // connect to ESD
1880   TTree *t = (TTree*) esdFile->Get("esdTree");
1881   AliESDEvent *esd = new AliESDEvent();
1882   esd->ReadFromTree(t);
1883
1884   Int_t nEvents = t->GetEntries();
1885
1886   // set arrays and pointers
1887   Float_t posF[3];
1888   Double_t pos[3];
1889   Double_t p[3];
1890   Double_t covVtx[6];
1891   Double_t covTr[21];
1892   Double_t pid[10];
1893
1894   // loop over events and fill them
1895   for (Int_t iEvent = 0; iEvent < nEvents; ++iEvent) {
1896     t->GetEntry(iEvent);
1897
1898     // Multiplicity information needed by the header (to be revised!)
1899     Int_t nTracks   = esd->GetNumberOfTracks();
1900     Int_t nPosTracks = 0;
1901     for (Int_t iTrack=0; iTrack<nTracks; ++iTrack)
1902       if (esd->GetTrack(iTrack)->Charge()> 0) nPosTracks++;
1903
1904     // Access the header
1905     AliAODHeader *header = aod->GetHeader();
1906
1907     // fill the header
1908     header->SetRunNumber       (esd->GetRunNumber()       );
1909     header->SetBunchCrossNumber(esd->GetBunchCrossNumber());
1910     header->SetOrbitNumber     (esd->GetOrbitNumber()     );
1911     header->SetPeriodNumber    (esd->GetPeriodNumber()    );
1912     header->SetTriggerMask     (esd->GetTriggerMask()     );
1913     header->SetTriggerCluster  (esd->GetTriggerCluster()  );
1914     header->SetEventType       (esd->GetEventType()       );
1915     header->SetMagneticField   (esd->GetMagneticField()   );
1916     header->SetZDCN1Energy     (esd->GetZDCN1Energy()     );
1917     header->SetZDCP1Energy     (esd->GetZDCP1Energy()     );
1918     header->SetZDCN2Energy     (esd->GetZDCN2Energy()     );
1919     header->SetZDCP2Energy     (esd->GetZDCP2Energy()     );
1920     header->SetZDCEMEnergy     (esd->GetZDCEMEnergy()     );
1921     header->SetRefMultiplicity   (nTracks);
1922     header->SetRefMultiplicityPos(nPosTracks);
1923     header->SetRefMultiplicityNeg(nTracks - nPosTracks);
1924     header->SetMuonMagFieldScale(-999.); // FIXME
1925     header->SetCentrality(-999.);        // FIXME
1926
1927     Int_t nV0s      = esd->GetNumberOfV0s();
1928     Int_t nCascades = esd->GetNumberOfCascades();
1929     Int_t nKinks    = esd->GetNumberOfKinks();
1930     Int_t nVertices = nV0s + nCascades + nKinks;
1931
1932     aod->ResetStd(nTracks, nVertices);
1933     AliAODTrack *aodTrack;
1934
1935     // Array to take into account the tracks already added to the AOD
1936     Bool_t * usedTrack = NULL;
1937     if (nTracks>0) {
1938       usedTrack = new Bool_t[nTracks];
1939       for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) usedTrack[iTrack]=kFALSE;
1940     }
1941     // Array to take into account the V0s already added to the AOD
1942     Bool_t * usedV0 = NULL;
1943     if (nV0s>0) {
1944       usedV0 = new Bool_t[nV0s];
1945       for (Int_t iV0=0; iV0<nV0s; ++iV0) usedV0[iV0]=kFALSE;
1946     }
1947     // Array to take into account the kinks already added to the AOD
1948     Bool_t * usedKink = NULL;
1949     if (nKinks>0) {
1950       usedKink = new Bool_t[nKinks];
1951       for (Int_t iKink=0; iKink<nKinks; ++iKink) usedKink[iKink]=kFALSE;
1952     }
1953
1954     // Access to the AOD container of vertices
1955     TClonesArray &vertices = *(aod->GetVertices());
1956     Int_t jVertices=0;
1957
1958     // Access to the AOD container of tracks
1959     TClonesArray &tracks = *(aod->GetTracks());
1960     Int_t jTracks=0;
1961
1962     // Add primary vertex. The primary tracks will be defined
1963     // after the loops on the composite objects (V0, cascades, kinks)
1964     const AliESDVertex *vtx = esd->GetPrimaryVertex();
1965
1966     vtx->GetXYZ(pos); // position
1967     vtx->GetCovMatrix(covVtx); //covariance matrix
1968
1969     AliAODVertex * primary = new(vertices[jVertices++])
1970       AliAODVertex(pos, covVtx, vtx->GetChi2toNDF(), NULL, -1, AliAODVertex::kPrimary);
1971
1972     // Create vertices starting from the most complex objects
1973
1974     // Cascades
1975     for (Int_t nCascade = 0; nCascade < nCascades; ++nCascade) {
1976       AliESDcascade *cascade = esd->GetCascade(nCascade);
1977
1978       cascade->GetXYZcascade(pos[0], pos[1], pos[2]);  // Bo: bug correction
1979       cascade->GetPosCovXi(covVtx);
1980
1981       // Add the cascade vertex
1982       AliAODVertex * vcascade = new(vertices[jVertices++]) AliAODVertex(pos,
1983                                                                         covVtx,
1984                                                                         cascade->GetChi2Xi(), // = chi2/NDF since NDF = 2*2-3
1985                                                                         primary,
1986                                                                         nCascade,
1987                                                                         AliAODVertex::kCascade);
1988
1989       primary->AddDaughter(vcascade);
1990
1991       // Add the V0 from the cascade. The ESD class have to be optimized...
1992       // Now we have to search for the corresponding V0 in the list of V0s
1993       // using the indeces of the positive and negative tracks
1994
1995       Int_t posFromV0 = cascade->GetPindex();
1996       Int_t negFromV0 = cascade->GetNindex();
1997
1998       AliESDv0 * v0 = 0x0;
1999       Int_t indV0 = -1;
2000
2001       for (Int_t iV0=0; iV0<nV0s; ++iV0) {
2002
2003         v0 = esd->GetV0(iV0);
2004         Int_t posV0 = v0->GetPindex();
2005         Int_t negV0 = v0->GetNindex();
2006
2007         if (posV0==posFromV0 && negV0==negFromV0) {
2008           indV0 = iV0;
2009           break;
2010         }
2011       }
2012
2013       AliAODVertex * vV0FromCascade = 0x0;
2014
2015       if (indV0>-1 && !usedV0[indV0] ) {
2016
2017         // the V0 exists in the array of V0s and is not used
2018
2019         usedV0[indV0] = kTRUE;
2020
2021         v0->GetXYZ(pos[0], pos[1], pos[2]);
2022         v0->GetPosCov(covVtx);
2023
2024         vV0FromCascade = new(vertices[jVertices++]) AliAODVertex(pos,
2025                                                                  covVtx,
2026                                                                  v0->GetChi2V0(), // = chi2/NDF since NDF = 2*2-3
2027                                                                  vcascade,
2028                                                                  indV0,
2029                                                                  AliAODVertex::kV0);
2030       } else {
2031
2032         // the V0 doesn't exist in the array of V0s or was used
2033         cerr << "Error: event " << iEvent << " cascade " << nCascade
2034              << " The V0 " << indV0
2035              << " doesn't exist in the array of V0s or was used!" << endl;
2036
2037         cascade->GetXYZ(pos[0], pos[1], pos[2]);
2038         cascade->GetPosCov(covVtx);
2039
2040         vV0FromCascade = new(vertices[jVertices++]) AliAODVertex(pos,
2041                                                                  covVtx,
2042                                                                  v0->GetChi2V0(), // = chi2/NDF since NDF = 2*2-3
2043                                                                  vcascade,
2044                                                                  indV0,
2045                                                                  AliAODVertex::kV0);
2046         vcascade->AddDaughter(vV0FromCascade);
2047       }
2048
2049       // Add the positive tracks from the V0
2050
2051       if (! usedTrack[posFromV0]) {
2052
2053         usedTrack[posFromV0] = kTRUE;
2054
2055         AliESDtrack *esdTrack = esd->GetTrack(posFromV0);
2056         esdTrack->GetPxPyPz(p);
2057         esdTrack->GetXYZ(pos);
2058         esdTrack->GetCovarianceXYZPxPyPz(covTr);
2059         esdTrack->GetESDpid(pid);
2060
2061         vV0FromCascade->AddDaughter(aodTrack =
2062                                     new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2063                                            esdTrack->GetLabel(),
2064                                            p,
2065                                            kTRUE,
2066                                            pos,
2067                                            kFALSE,
2068                                            covTr,
2069                                            (Short_t)esdTrack->Charge(),
2070                                            esdTrack->GetITSClusterMap(),
2071                                            pid,
2072                                            vV0FromCascade,
2073                                            kTRUE,  // check if this is right
2074                                            kFALSE, // check if this is right
2075                                            AliAODTrack::kSecondary)
2076                 );
2077         aodTrack->ConvertAliPIDtoAODPID();
2078       }
2079       else {
2080         cerr << "Error: event " << iEvent << " cascade " << nCascade
2081              << " track " << posFromV0 << " has already been used!" << endl;
2082       }
2083
2084       // Add the negative tracks from the V0
2085
2086       if (!usedTrack[negFromV0]) {
2087
2088         usedTrack[negFromV0] = kTRUE;
2089
2090         AliESDtrack *esdTrack = esd->GetTrack(negFromV0);
2091         esdTrack->GetPxPyPz(p);
2092         esdTrack->GetXYZ(pos);
2093         esdTrack->GetCovarianceXYZPxPyPz(covTr);
2094         esdTrack->GetESDpid(pid);
2095
2096         vV0FromCascade->AddDaughter(aodTrack =
2097                 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2098                                            esdTrack->GetLabel(),
2099                                            p,
2100                                            kTRUE,
2101                                            pos,
2102                                            kFALSE,
2103                                            covTr,
2104                                            (Short_t)esdTrack->Charge(),
2105                                            esdTrack->GetITSClusterMap(),
2106                                            pid,
2107                                            vV0FromCascade,
2108                                            kTRUE,  // check if this is right
2109                                            kFALSE, // check if this is right
2110                                            AliAODTrack::kSecondary)
2111                 );
2112         aodTrack->ConvertAliPIDtoAODPID();
2113       }
2114       else {
2115         cerr << "Error: event " << iEvent << " cascade " << nCascade
2116              << " track " << negFromV0 << " has already been used!" << endl;
2117       }
2118
2119       // Add the bachelor track from the cascade
2120
2121       Int_t bachelor = cascade->GetBindex();
2122
2123       if(!usedTrack[bachelor]) {
2124
2125         usedTrack[bachelor] = kTRUE;
2126
2127         AliESDtrack *esdTrack = esd->GetTrack(bachelor);
2128         esdTrack->GetPxPyPz(p);
2129         esdTrack->GetXYZ(pos);
2130         esdTrack->GetCovarianceXYZPxPyPz(covTr);
2131         esdTrack->GetESDpid(pid);
2132
2133         vcascade->AddDaughter(aodTrack =
2134                 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2135                                            esdTrack->GetLabel(),
2136                                            p,
2137                                            kTRUE,
2138                                            pos,
2139                                            kFALSE,
2140                                            covTr,
2141                                            (Short_t)esdTrack->Charge(),
2142                                            esdTrack->GetITSClusterMap(),
2143                                            pid,
2144                                            vcascade,
2145                                            kTRUE,  // check if this is right
2146                                            kFALSE, // check if this is right
2147                                            AliAODTrack::kSecondary)
2148                 );
2149         aodTrack->ConvertAliPIDtoAODPID();
2150      }
2151       else {
2152         cerr << "Error: event " << iEvent << " cascade " << nCascade
2153              << " track " << bachelor << " has already been used!" << endl;
2154       }
2155
2156       // Add the primary track of the cascade (if any)
2157
2158     } // end of the loop on cascades
2159
2160     // V0s
2161
2162     for (Int_t nV0 = 0; nV0 < nV0s; ++nV0) {
2163
2164       if (usedV0[nV0]) continue; // skip if aready added to the AOD
2165
2166       AliESDv0 *v0 = esd->GetV0(nV0);
2167
2168       v0->GetXYZ(pos[0], pos[1], pos[2]);
2169       v0->GetPosCov(covVtx);
2170
2171       AliAODVertex * vV0 =
2172         new(vertices[jVertices++]) AliAODVertex(pos,
2173                                                 covVtx,
2174                                                 v0->GetChi2V0(), // = chi2/NDF since NDF = 2*2-3
2175                                                 primary,
2176                                                 nV0,
2177                                                 AliAODVertex::kV0);
2178       primary->AddDaughter(vV0);
2179
2180       Int_t posFromV0 = v0->GetPindex();
2181       Int_t negFromV0 = v0->GetNindex();
2182
2183       // Add the positive tracks from the V0
2184
2185       if (!usedTrack[posFromV0]) {
2186
2187         usedTrack[posFromV0] = kTRUE;
2188
2189         AliESDtrack *esdTrack = esd->GetTrack(posFromV0);
2190         esdTrack->GetPxPyPz(p);
2191         esdTrack->GetXYZ(pos);
2192         esdTrack->GetCovarianceXYZPxPyPz(covTr);
2193         esdTrack->GetESDpid(pid);
2194
2195         vV0->AddDaughter(aodTrack =
2196                 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2197                                            esdTrack->GetLabel(),
2198                                            p,
2199                                            kTRUE,
2200                                            pos,
2201                                            kFALSE,
2202                                            covTr,
2203                                            (Short_t)esdTrack->Charge(),
2204                                            esdTrack->GetITSClusterMap(),
2205                                            pid,
2206                                            vV0,
2207                                            kTRUE,  // check if this is right
2208                                            kFALSE, // check if this is right
2209                                            AliAODTrack::kSecondary)
2210                 );
2211         aodTrack->ConvertAliPIDtoAODPID();
2212       }
2213       else {
2214         cerr << "Error: event " << iEvent << " V0 " << nV0
2215              << " track " << posFromV0 << " has already been used!" << endl;
2216       }
2217
2218       // Add the negative tracks from the V0
2219
2220       if (!usedTrack[negFromV0]) {
2221
2222         usedTrack[negFromV0] = kTRUE;
2223
2224         AliESDtrack *esdTrack = esd->GetTrack(negFromV0);
2225         esdTrack->GetPxPyPz(p);
2226         esdTrack->GetXYZ(pos);
2227         esdTrack->GetCovarianceXYZPxPyPz(covTr);
2228         esdTrack->GetESDpid(pid);
2229
2230         vV0->AddDaughter(aodTrack =
2231                 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2232                                            esdTrack->GetLabel(),
2233                                            p,
2234                                            kTRUE,
2235                                            pos,
2236                                            kFALSE,
2237                                            covTr,
2238                                            (Short_t)esdTrack->Charge(),
2239                                            esdTrack->GetITSClusterMap(),
2240                                            pid,
2241                                            vV0,
2242                                            kTRUE,  // check if this is right
2243                                            kFALSE, // check if this is right
2244                                            AliAODTrack::kSecondary)
2245                 );
2246         aodTrack->ConvertAliPIDtoAODPID();
2247       }
2248       else {
2249         cerr << "Error: event " << iEvent << " V0 " << nV0
2250              << " track " << negFromV0 << " has already been used!" << endl;
2251       }
2252
2253     } // end of the loop on V0s
2254
2255     // Kinks: it is a big mess the access to the information in the kinks
2256     // The loop is on the tracks in order to find the mother and daugther of each kink
2257
2258
2259     for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) {
2260
2261
2262       AliESDtrack * esdTrack = esd->GetTrack(iTrack);
2263
2264       Int_t ikink = esdTrack->GetKinkIndex(0);
2265
2266       if (ikink) {
2267         // Negative kink index: mother, positive: daughter
2268
2269         // Search for the second track of the kink
2270
2271         for (Int_t jTrack = iTrack+1; jTrack<nTracks; ++jTrack) {
2272
2273           AliESDtrack * esdTrack1 = esd->GetTrack(jTrack);
2274
2275           Int_t jkink = esdTrack1->GetKinkIndex(0);
2276
2277           if ( TMath::Abs(ikink)==TMath::Abs(jkink) ) {
2278
2279             // The two tracks are from the same kink
2280
2281             if (usedKink[TMath::Abs(ikink)-1]) continue; // skip used kinks
2282
2283             Int_t imother = -1;
2284             Int_t idaughter = -1;
2285
2286             if (ikink<0 && jkink>0) {
2287
2288               imother = iTrack;
2289               idaughter = jTrack;
2290             }
2291             else if (ikink>0 && jkink<0) {
2292
2293               imother = jTrack;
2294               idaughter = iTrack;
2295             }
2296             else {
2297               cerr << "Error: Wrong combination of kink indexes: "
2298               << ikink << " " << jkink << endl;
2299               continue;
2300             }
2301
2302             // Add the mother track
2303
2304             AliAODTrack * mother = NULL;
2305
2306             if (!usedTrack[imother]) {
2307
2308               usedTrack[imother] = kTRUE;
2309
2310               AliESDtrack *esdTrack = esd->GetTrack(imother);
2311               esdTrack->GetPxPyPz(p);
2312               esdTrack->GetXYZ(pos);
2313               esdTrack->GetCovarianceXYZPxPyPz(covTr);
2314               esdTrack->GetESDpid(pid);
2315
2316               mother =
2317                 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2318                                            esdTrack->GetLabel(),
2319                                            p,
2320                                            kTRUE,
2321                                            pos,
2322                                            kFALSE,
2323                                            covTr,
2324                                            (Short_t)esdTrack->Charge(),
2325                                            esdTrack->GetITSClusterMap(),
2326                                            pid,
2327                                            primary,
2328                                            kTRUE, // check if this is right
2329                                            kTRUE, // check if this is right
2330                                            AliAODTrack::kPrimary);
2331               primary->AddDaughter(mother);
2332               mother->ConvertAliPIDtoAODPID();
2333             }
2334             else {
2335               cerr << "Error: event " << iEvent << " kink " << TMath::Abs(ikink)-1
2336               << " track " << imother << " has already been used!" << endl;
2337             }
2338
2339             // Add the kink vertex
2340             AliESDkink * kink = esd->GetKink(TMath::Abs(ikink)-1);
2341
2342             AliAODVertex * vkink =
2343             new(vertices[jVertices++]) AliAODVertex(kink->GetPosition(),
2344                                                     NULL,
2345                                                     0.,
2346                                                     mother,
2347                                                     esdTrack->GetID(), // This is the track ID of the mother's track!
2348                                                     AliAODVertex::kKink);
2349             // Add the daughter track
2350
2351             AliAODTrack * daughter = NULL;
2352
2353             if (!usedTrack[idaughter]) {
2354
2355               usedTrack[idaughter] = kTRUE;
2356
2357               AliESDtrack *esdTrack = esd->GetTrack(idaughter);
2358               esdTrack->GetPxPyPz(p);
2359               esdTrack->GetXYZ(pos);
2360               esdTrack->GetCovarianceXYZPxPyPz(covTr);
2361               esdTrack->GetESDpid(pid);
2362
2363               daughter =
2364                 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2365                                            esdTrack->GetLabel(),
2366                                            p,
2367                                            kTRUE,
2368                                            pos,
2369                                            kFALSE,
2370                                            covTr,
2371                                            (Short_t)esdTrack->Charge(),
2372                                            esdTrack->GetITSClusterMap(),
2373                                            pid,
2374                                            vkink,
2375                                            kTRUE, // check if this is right
2376                                            kTRUE, // check if this is right
2377                                            AliAODTrack::kPrimary);
2378               vkink->AddDaughter(daughter);
2379               daughter->ConvertAliPIDtoAODPID();
2380             }
2381             else {
2382               cerr << "Error: event " << iEvent << " kink " << TMath::Abs(ikink)-1
2383               << " track " << idaughter << " has already been used!" << endl;
2384             }
2385
2386
2387           }
2388         }
2389
2390       }
2391
2392     }
2393
2394
2395     // Tracks (primary and orphan)
2396
2397     for (Int_t nTrack = 0; nTrack < nTracks; ++nTrack) {
2398
2399
2400       if (usedTrack[nTrack]) continue;
2401
2402       AliESDtrack *esdTrack = esd->GetTrack(nTrack);
2403       esdTrack->GetPxPyPz(p);
2404       esdTrack->GetXYZ(pos);
2405       esdTrack->GetCovarianceXYZPxPyPz(covTr);
2406       esdTrack->GetESDpid(pid);
2407
2408       Float_t impactXY, impactZ;
2409
2410       esdTrack->GetImpactParameters(impactXY,impactZ);
2411
2412       if (impactXY<3) {
2413         // track inside the beam pipe
2414
2415         primary->AddDaughter(aodTrack =
2416             new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2417                                          esdTrack->GetLabel(),
2418                                          p,
2419                                          kTRUE,
2420                                          pos,
2421                                          kFALSE,
2422                                          covTr,
2423                                          (Short_t)esdTrack->Charge(),
2424                                          esdTrack->GetITSClusterMap(),
2425                                          pid,
2426                                          primary,
2427                                          kTRUE, // check if this is right
2428                                          kTRUE, // check if this is right
2429                                          AliAODTrack::kPrimary)
2430             );
2431         aodTrack->ConvertAliPIDtoAODPID();
2432       }
2433       else {
2434         // outside the beam pipe: orphan track
2435             aodTrack =
2436             new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2437                                          esdTrack->GetLabel(),
2438                                          p,
2439                                          kTRUE,
2440                                          pos,
2441                                          kFALSE,
2442                                          covTr,
2443                                          (Short_t)esdTrack->Charge(),
2444                                          esdTrack->GetITSClusterMap(),
2445                                          pid,
2446                                          NULL,
2447                                          kFALSE, // check if this is right
2448                                          kFALSE, // check if this is right
2449                                          AliAODTrack::kOrphan);
2450             aodTrack->ConvertAliPIDtoAODPID();
2451       }
2452     } // end of loop on tracks
2453
2454     // muon tracks
2455     Int_t nMuTracks = esd->GetNumberOfMuonTracks();
2456     for (Int_t nMuTrack = 0; nMuTrack < nMuTracks; ++nMuTrack) {
2457
2458       AliESDMuonTrack *esdMuTrack = esd->GetMuonTrack(nMuTrack);
2459       p[0] = esdMuTrack->Px();
2460       p[1] = esdMuTrack->Py();
2461       p[2] = esdMuTrack->Pz();
2462       pos[0] = primary->GetX();
2463       pos[1] = primary->GetY();
2464       pos[2] = primary->GetZ();
2465
2466       // has to be changed once the muon pid is provided by the ESD
2467       for (Int_t i = 0; i < 10; pid[i++] = 0.); pid[AliAODTrack::kMuon]=1.;
2468
2469       primary->AddDaughter(aodTrack =
2470           new(tracks[jTracks++]) AliAODTrack(0, // no ID provided
2471                                              0, // no label provided
2472                                              p,
2473                                              kTRUE,
2474                                              pos,
2475                                              kFALSE,
2476                                              NULL, // no covariance matrix provided
2477                                              esdMuTrack->Charge(),
2478                                              0, // ITSClusterMap is set below
2479                                              pid,
2480                                              primary,
2481                                              kFALSE,  // muon tracks are not used to fit the primary vtx
2482                                              kFALSE,  // not used for vertex fit
2483                                              AliAODTrack::kPrimary)
2484           );
2485
2486       aodTrack->SetHitsPatternInTrigCh(esdMuTrack->GetHitsPatternInTrigCh());
2487       Int_t track2Trigger = esdMuTrack->GetMatchTrigger();
2488       aodTrack->SetMatchTrigger(track2Trigger);
2489       if (track2Trigger)
2490         aodTrack->SetChi2MatchTrigger(esdMuTrack->GetChi2MatchTrigger());
2491       else
2492         aodTrack->SetChi2MatchTrigger(0.);
2493     }
2494
2495     // Access to the AOD container of clusters
2496     TClonesArray &clusters = *(aod->GetClusters());
2497     Int_t jClusters=0;
2498
2499     // Calo Clusters
2500     Int_t nClusters    = esd->GetNumberOfCaloClusters();
2501
2502     for (Int_t iClust=0; iClust<nClusters; ++iClust) {
2503
2504       AliESDCaloCluster * cluster = esd->GetCaloCluster(iClust);
2505
2506       Int_t id = cluster->GetID();
2507       Int_t label = -1;
2508       Float_t energy = cluster->E();
2509       cluster->GetPosition(posF);
2510       AliAODVertex *prodVertex = primary;
2511       AliAODTrack *primTrack = NULL;
2512       Char_t ttype=AliAODCluster::kUndef;
2513
2514       if (cluster->GetClusterType() == AliESDCaloCluster::kPHOSCluster){
2515         ttype = AliAODCluster::kPHOSNeutral;
2516       }
2517       else if (cluster->GetClusterType() == AliESDCaloCluster::kEMCALClusterv1){
2518         ttype = AliAODCluster::kEMCALClusterv1;
2519       }
2520       else if (cluster->GetClusterType() == AliESDCaloCluster::kEMCALPseudoCluster){
2521         ttype = AliAODCluster::kEMCALPseudoCluster;
2522       }
2523
2524       new(clusters[jClusters++]) AliAODCluster(id,
2525                                                label,
2526                                                energy,
2527                                                pos,
2528                                                NULL, // no covariance matrix provided
2529                                                NULL, // no pid for clusters provided
2530                                                prodVertex,
2531                                                primTrack,
2532                                                ttype);
2533
2534     } // end of loop on calo clusters
2535
2536     // tracklets
2537     const AliMultiplicity *mult = esd->GetMultiplicity();
2538     if (mult) {
2539       if (mult->GetNumberOfTracklets()>0) {
2540         aod->GetTracklets()->CreateContainer(mult->GetNumberOfTracklets());
2541
2542         for (Int_t n=0; n<mult->GetNumberOfTracklets(); n++) {
2543           aod->GetTracklets()->SetTracklet(n, mult->GetTheta(n), mult->GetPhi(n), mult->GetDeltaPhi(n), mult->GetLabel(n));
2544         }
2545       }
2546     } else {
2547       Printf("ERROR: AliMultiplicity could not be retrieved from ESD");
2548     }
2549
2550     delete [] usedTrack;
2551     delete [] usedV0;
2552     delete [] usedKink;
2553
2554     // fill the tree for this event
2555     aodTree->Fill();
2556   } // end of event loop
2557
2558   aodTree->GetUserInfo()->Add(aod);
2559
2560   // close ESD file
2561   esdFile->Close();
2562
2563   // write the tree to the specified file
2564   aodFile = aodTree->GetCurrentFile();
2565   aodFile->cd();
2566   aodTree->Write();
2567
2568   return;
2569 }
2570
2571 void AliReconstruction::WriteAlignmentData(AliESDEvent* esd)
2572 {
2573   // Write space-points which are then used in the alignment procedures
2574   // For the moment only ITS, TRD and TPC
2575
2576   // Load TOF clusters
2577   if (fTracker[3]){
2578     fLoader[3]->LoadRecPoints("read");
2579     TTree* tree = fLoader[3]->TreeR();
2580     if (!tree) {
2581       AliError(Form("Can't get the %s cluster tree", fgkDetectorName[3]));
2582       return;
2583     }
2584     fTracker[3]->LoadClusters(tree);
2585   }
2586   Int_t ntracks = esd->GetNumberOfTracks();
2587   for (Int_t itrack = 0; itrack < ntracks; itrack++)
2588     {
2589       AliESDtrack *track = esd->GetTrack(itrack);
2590       Int_t nsp = 0;
2591       Int_t idx[200];
2592       for (Int_t iDet = 3; iDet >= 0; iDet--)
2593         nsp += track->GetNcls(iDet);
2594       if (nsp) {
2595         AliTrackPointArray *sp = new AliTrackPointArray(nsp);
2596         track->SetTrackPointArray(sp);
2597         Int_t isptrack = 0;
2598         for (Int_t iDet = 3; iDet >= 0; iDet--) {
2599           AliTracker *tracker = fTracker[iDet];
2600           if (!tracker) continue;
2601           Int_t nspdet = track->GetNcls(iDet);
2602           if (nspdet <= 0) continue;
2603           track->GetClusters(iDet,idx);
2604           AliTrackPoint p;
2605           Int_t isp = 0;
2606           Int_t isp2 = 0;
2607           while (isp < nspdet) {
2608             Bool_t isvalid;
2609             if(IsSelected(fgkDetectorName[iDet],fUseTrackingErrorsForAlignment)) {
2610               isvalid = tracker->GetTrackPointTrackingError(idx[isp2],p,track);
2611             } else {
2612               isvalid = tracker->GetTrackPoint(idx[isp2],p);
2613             }
2614             isp2++;
2615             const Int_t kNTPCmax = 159;
2616             if (iDet==1 && isp2>kNTPCmax) break;   // to be fixed
2617             if (!isvalid) continue;
2618             sp->AddPoint(isptrack,&p); isptrack++; isp++;
2619           }
2620         }
2621       }
2622     }
2623   if (fTracker[3]){
2624     fTracker[3]->UnloadClusters();
2625     fLoader[3]->UnloadRecPoints();
2626   }
2627 }
2628
2629 //_____________________________________________________________________________
2630 void AliReconstruction::FillRawDataErrorLog(Int_t iEvent, AliESDEvent* esd)
2631 {
2632   // The method reads the raw-data error log
2633   // accumulated within the rawReader.
2634   // It extracts the raw-data errors related to
2635   // the current event and stores them into
2636   // a TClonesArray inside the esd object.
2637
2638   if (!fRawReader) return;
2639
2640   for(Int_t i = 0; i < fRawReader->GetNumberOfErrorLogs(); i++) {
2641
2642     AliRawDataErrorLog *log = fRawReader->GetErrorLog(i);
2643     if (!log) continue;
2644     if (iEvent != log->GetEventNumber()) continue;
2645
2646     esd->AddRawDataErrorLog(log);
2647   }
2648
2649 }
2650
2651 TNamed* AliReconstruction::CopyFileToTNamed(TString fPath,TString fName){
2652   // Dump a file content into a char in TNamed
2653   ifstream in;
2654   in.open(fPath.Data(),ios::in | ios::binary|ios::ate);
2655   Int_t kBytes = (Int_t)in.tellg();
2656   printf("Size: %d \n",kBytes);
2657   TNamed *fn = 0;
2658   if(in.good()){
2659     char* memblock = new char [kBytes];
2660     in.seekg (0, ios::beg);
2661     in.read (memblock, kBytes);
2662     in.close();
2663     TString fData(memblock,kBytes);
2664     fn = new TNamed(fName,fData);
2665     printf("fData Size: %d \n",fData.Sizeof());
2666     printf("fName Size: %d \n",fName.Sizeof());
2667     printf("fn    Size: %d \n",fn->Sizeof());
2668     delete[] memblock;
2669   }
2670   else{
2671     AliInfo(Form("Could not Open %s\n",fPath.Data()));
2672   }
2673
2674   return fn;
2675 }
2676
2677 void AliReconstruction::TNamedToFile(TTree* fTree, TString fName){
2678   // This is not really needed in AliReconstruction at the moment
2679   // but can serve as a template
2680
2681   TList *fList = fTree->GetUserInfo();
2682   TNamed *fn = (TNamed*)fList->FindObject(fName.Data());
2683   printf("fn Size: %d \n",fn->Sizeof());
2684
2685   TString fTmp(fn->GetName()); // to be 100% sure in principle fName also works
2686   const char* cdata = fn->GetTitle();
2687   printf("fTmp Size %d\n",fTmp.Sizeof());
2688
2689   int size = fn->Sizeof()-fTmp.Sizeof()-sizeof(UChar_t)-sizeof(Int_t); // see dfinition of TString::SizeOf()...
2690   printf("calculated size %d\n",size);
2691   ofstream out(fName.Data(),ios::out | ios::binary);
2692   out.write(cdata,size);
2693   out.close();
2694
2695 }
2696
2697 //_____________________________________________________________________________
2698 AliQADataMaker * AliReconstruction::GetQADataMaker(Int_t iDet)
2699 {
2700 // get the quality assurance data maker object and the loader for a detector
2701
2702   if (fQADataMaker[iDet])
2703     return fQADataMaker[iDet];
2704
2705   // load the QA data maker object
2706   TPluginManager* pluginManager = gROOT->GetPluginManager();
2707   TString detName = fgkDetectorName[iDet];
2708   TString qadmName = "Ali" + detName + "QADataMaker";
2709   if (gAlice && !gAlice->GetDetector(detName) && (detName != "HLT"))
2710     return NULL;
2711
2712   AliQADataMaker * qadm = NULL;
2713   // first check if a plugin is defined for the quality assurance data maker
2714   TPluginHandler* pluginHandler = pluginManager->FindHandler("AliQADataMaker", detName);
2715   // if not, add a plugin for it
2716   if (!pluginHandler) {
2717     AliDebug(1, Form("defining plugin for %s", qadmName.Data()));
2718     TString libs = gSystem->GetLibraries();
2719     if (libs.Contains("lib" + detName + "base.so") ||
2720         (gSystem->Load("lib" + detName + "base.so") >= 0)) {
2721       pluginManager->AddHandler("AliQADataMaker", detName,
2722                                 qadmName, detName + "qadm", qadmName + "()");
2723     } else {
2724       pluginManager->AddHandler("AliQADataMaker", detName,
2725                                 qadmName, detName, qadmName + "()");
2726     }
2727     pluginHandler = pluginManager->FindHandler("AliQADataMaker", detName);
2728   }
2729   if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
2730     qadm = (AliQADataMaker *) pluginHandler->ExecPlugin(0);
2731   }
2732   if (qadm) {
2733     AliInfo(Form("Initializing quality assurance data maker for %s", fgkDetectorName[iDet]));
2734     qadm->Init(AliQA::kRECPOINTS, AliCDBManager::Instance()->GetRun(), GetQACycles(fgkDetectorName[iDet]));
2735     qadm->StartOfCycle(AliQA::kRECPOINTS);
2736     qadm->Init(AliQA::kESDS, AliCDBManager::Instance()->GetRun());
2737     qadm->StartOfCycle(AliQA::kESDS, "same") ;
2738     fQADataMaker[iDet] = qadm;
2739   }
2740
2741   return qadm;
2742 }
2743
2744 //_____________________________________________________________________________
2745 Bool_t AliReconstruction::RunQA(const char* detectors, AliESDEvent *& esd)
2746 {
2747   // run the Quality Assurance data producer
2748
2749   AliCodeTimerAuto("")
2750   TString detStr = detectors;
2751   for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2752    if (!IsSelected(fgkDetectorName[iDet], detStr))
2753      continue;
2754    AliQADataMaker * qadm = GetQADataMaker(iDet);
2755    if (!qadm)
2756      continue;
2757    AliCodeTimerStart(Form("running quality assurance data maker for %s", fgkDetectorName[iDet]));
2758    AliInfo(Form("running quality assurance data maker for %s", fgkDetectorName[iDet]));
2759
2760    qadm->Exec(AliQA::kESDS, esd) ;
2761    qadm->Increment() ;
2762
2763    AliCodeTimerStop(Form("running quality assurance data maker for %s", fgkDetectorName[iDet]));
2764  }
2765  if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
2766    AliError(Form("the following detectors were not found: %s",
2767                  detStr.Data()));
2768    if (fStopOnError)
2769      return kFALSE;
2770  }
2771
2772  return kTRUE;
2773
2774 }
2775
2776
2777 //_____________________________________________________________________________
2778 void AliReconstruction::CheckQA()
2779 {
2780 // check the QA of SIM for this run and remove the detectors
2781 // with status Fatal
2782
2783         TString newDetList ;
2784         for (Int_t iDet = 0; iDet < AliQA::kNDET; iDet++) {
2785                 TString detName(AliQA::GetDetName(iDet)) ;
2786                 if ( fRunLocalReconstruction.Contains(AliQA::GetDetName(iDet)) ||
2787                         fRunLocalReconstruction.Contains("ALL") )  {
2788                         AliQA * qa = AliQA::Instance(AliQA::DETECTORINDEX(iDet)) ;
2789                         if ( qa->IsSet(AliQA::DETECTORINDEX(iDet), AliQA::kSIM, AliQA::kFATAL)) {
2790                                 AliInfo(Form("QA status for %s in Hits and/or SDIGITS  and/or Digits was Fatal; No reconstruction performed", detName.Data())) ;
2791                         } else if ( qa->IsSet(AliQA::DETECTORINDEX(iDet), AliQA::kSIM, AliQA::kERROR)) {
2792                                 AliError(Form("QA status for %s in Hits and/or SDIGITS  and/or Digits was ERROR", detName.Data())) ;
2793                                 newDetList += detName ;
2794                                 newDetList += " " ;
2795                         } else if ( qa->IsSet(AliQA::DETECTORINDEX(iDet), AliQA::kSIM, AliQA::kWARNING) ) {
2796                                 AliWarning(Form("QA status for %s in Hits and/or SDIGITS  and/or Digits was WARNING", detName.Data())) ;
2797                                 newDetList += detName ;
2798                                 newDetList += " " ;
2799                         } else if ( qa->IsSet(AliQA::DETECTORINDEX(iDet), AliQA::kSIM, AliQA::kINFO) ) {
2800                                 AliInfo(Form("QA status for %s in Hits and/or SDIGITS  and/or Digits was INFO", detName.Data())) ;
2801                                 newDetList += detName ;
2802                                 newDetList += " " ;
2803                         }
2804                 }
2805         }
2806         fRunLocalReconstruction = newDetList ;
2807 }
2808
2809 //_____________________________________________________________________________
2810 Int_t AliReconstruction::GetDetIndex(const char* detector)
2811 {
2812   // return the detector index corresponding to detector
2813   Int_t index = -1 ;
2814   for (index = 0; index < fgkNDetectors ; index++) {
2815     if ( strcmp(detector, fgkDetectorName[index]) == 0 )
2816         break ;
2817   }
2818   return index ;
2819 }