STEER/AliReconstruction.cxx

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