STEER/AliReconstruction.cxx

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