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 // The name of the galice file can be changed from the default               //
  30 // "galice.root" by passing it as argument to the AliReconstruction          //
  31 // constructor or by                                                         //
  32 //                                                                           //
  33 //   rec.SetGAliceFile("...");                                               //
  34 //                                                                           //
  35 // The reconstruction can be switched on or off for individual detectors by  //
  36 //                                                                           //
  37 //   rec.SetRunReconstruction("...");                                        //
  38 //                                                                           //
  39 // The argument is a (case sensitive) string with the names of the           //
  40 // detectors separated by a space. The special string "ALL" selects all      //
  41 // available detectors. This is the default.                                 //
  42 //                                                                           //
  43 // The tracking in ITS, TPC and TRD and the creation of ESD tracks can be    //
  44 // switched off by                                                           //
  45 //                                                                           //
  46 //   rec.SetRunTracking(kFALSE);                                             //
  47 //                                                                           //
  48 // The filling of additional ESD information can be steered by               //
  49 //                                                                           //
  50 //   rec.SetFillESD("...");                                                  //
  51 //                                                                           //
  52 // Again, the string specifies the list of detectors. The default is "ALL".  //
  53 //                                                                           //
  54 // The reconstruction requires digits as input. For the creation of digits   //
  55 // have a look at the class AliSimulation.                                   //
  56 //                                                                           //
  57 ///////////////////////////////////////////////////////////////////////////////
  58
  59
  60 #include "AliReconstruction.h"
  61 #include "AliRunLoader.h"
  62 #include "AliRun.h"
  63 #include "AliModule.h"
  64 #include "AliDetector.h"
  65 #include "AliTracker.h"
  66 #include "AliESD.h"
  67 #include "AliHeader.h"
  68 #include "AliGenEventHeader.h"
  69 #include "AliESDpid.h"
  70 #include "AliMagF.h"
  71 #include <TArrayF.h>
  72
  73
  74 ClassImp(AliReconstruction)
  75
  76
  77 //_____________________________________________________________________________
  78 AliReconstruction::AliReconstruction(const char* gAliceFilename,
  79                                      const char* name, const char* title) :
  80   TNamed(name, title),
  81
  82   fRunReconstruction("ALL"),
  83   fRunTracking(kTRUE),
  84   fFillESD("ALL"),
  85   fGAliceFileName(gAliceFilename),
  86   fStopOnError(kFALSE),
  87
  88   fRunLoader(NULL),
  89   fITSLoader(NULL),
  90   fITSTracker(NULL),
  91   fTPCLoader(NULL),
  92   fTPCTracker(NULL),
  93   fTRDLoader(NULL),
  94   fTRDTracker(NULL),
  95   fTOFLoader(NULL),
  96   fTOFTracker(NULL)
  97 {
  98 // create reconstruction object with default parameters
  99
 100 }
 101
 102 //_____________________________________________________________________________
 103 AliReconstruction::AliReconstruction(const AliReconstruction& rec) :
 104   TNamed(rec),
 105
 106   fRunReconstruction(rec.fRunReconstruction),
 107   fRunTracking(rec.fRunTracking),
 108   fFillESD(rec.fFillESD),
 109   fGAliceFileName(rec.fGAliceFileName),
 110   fStopOnError(rec.fStopOnError),
 111
 112   fRunLoader(NULL),
 113   fITSLoader(NULL),
 114   fITSTracker(NULL),
 115   fTPCLoader(NULL),
 116   fTPCTracker(NULL),
 117   fTRDLoader(NULL),
 118   fTRDTracker(NULL),
 119   fTOFLoader(NULL),
 120   fTOFTracker(NULL)
 121 {
 122 // copy constructor
 123
 124 }
 125
 126 //_____________________________________________________________________________
 127 AliReconstruction& AliReconstruction::operator = (const AliReconstruction& rec)
 128 {
 129 // assignment operator
 130
 131   this->~AliReconstruction();
 132   new(this) AliReconstruction(rec);
 133   return *this;
 134 }
 135
 136 //_____________________________________________________________________________
 137 AliReconstruction::~AliReconstruction()
 138 {
 139 // clean up
 140
 141   CleanUp();
 142 }
 143
 144
 145 //_____________________________________________________________________________
 146 void AliReconstruction::SetGAliceFile(const char* fileName)
 147 {
 148 // set the name of the galice file
 149
 150   fGAliceFileName = fileName;
 151 }
 152
 153
 154 //_____________________________________________________________________________
 155 Bool_t AliReconstruction::Run()
 156 {
 157 // run the reconstruction
 158
 159   // open the run loader
 160   fRunLoader = AliRunLoader::Open(fGAliceFileName.Data());
 161   if (!fRunLoader) {
 162     Error("Run", "no run loader found in file %s",
 163           fGAliceFileName.Data());
 164     CleanUp();
 165     return kFALSE;
 166   }
 167   fRunLoader->LoadgAlice();
 168   AliRun* aliRun = fRunLoader->GetAliRun();
 169   if (!aliRun) {
 170     Error("Run", "no gAlice object found in file %s",
 171           fGAliceFileName.Data());
 172     CleanUp();
 173     return kFALSE;
 174   }
 175   gAlice = aliRun;
 176
 177   // local reconstruction
 178   if (!fRunReconstruction.IsNull()) {
 179     if (!RunReconstruction(fRunReconstruction)) {
 180       if (fStopOnError) {CleanUp(); return kFALSE;}
 181     }
 182   }
 183   if (!fRunTracking && fFillESD.IsNull()) return kTRUE;
 184
 185   if (fRunTracking) {
 186     // get loaders and trackers
 187     fITSLoader = fRunLoader->GetLoader("ITSLoader");
 188     if (!fITSLoader) {
 189       Error("Run", "no ITS loader found");
 190       if (fStopOnError) {CleanUp(); return kFALSE;}
 191     }
 192     fITSTracker = NULL;
 193     if (aliRun->GetDetector("ITS")) {
 194       fITSTracker = aliRun->GetDetector("ITS")->CreateTracker();
 195     }
 196     if (!fITSTracker) {
 197       Error("Run", "couldn't create a tracker for ITS");
 198       if (fStopOnError) {CleanUp(); return kFALSE;}
 199     }
 200
 201     fTPCLoader = fRunLoader->GetLoader("TPCLoader");
 202     if (!fTPCLoader) {
 203       Error("Run", "no TPC loader found");
 204       if (fStopOnError) {CleanUp(); return kFALSE;}
 205     }
 206     fTPCTracker = NULL;
 207     if (aliRun->GetDetector("TPC")) {
 208       fTPCTracker = aliRun->GetDetector("TPC")->CreateTracker();
 209     }
 210     if (!fTPCTracker) {
 211       Error("Run", "couldn't create a tracker for TPC");
 212       if (fStopOnError) {CleanUp(); return kFALSE;}
 213     }
 214
 215     fTRDLoader = fRunLoader->GetLoader("TRDLoader");
 216     if (!fTRDLoader) {
 217       Error("Run", "no TRD loader found");
 218       if (fStopOnError) {CleanUp(); return kFALSE;}
 219     }
 220     fTRDTracker = NULL;
 221     if (aliRun->GetDetector("TRD")) {
 222       fTRDTracker = aliRun->GetDetector("TRD")->CreateTracker();
 223     }
 224     if (!fTRDTracker) {
 225       Error("Run", "couldn't create a tracker for TRD");
 226       if (fStopOnError) {CleanUp(); return kFALSE;}
 227     }
 228
 229     fTOFLoader = fRunLoader->GetLoader("TOFLoader");
 230     if (!fTOFLoader) {
 231       Error("Run", "no TOF loader found");
 232       if (fStopOnError) {CleanUp(); return kFALSE;}
 233     }
 234     fTOFTracker = NULL;
 235     if (aliRun->GetDetector("TOF")) {
 236       fTOFTracker = aliRun->GetDetector("TOF")->CreateTracker();
 237     }
 238     if (!fTOFTracker) {
 239       Error("Run", "couldn't create a tracker for TOF");
 240       if (fStopOnError) {CleanUp(); return kFALSE;}
 241     }
 242   }
 243   // create the ESD output file
 244   TFile* file = TFile::Open("AliESDs.root", "RECREATE");
 245   if (!file->IsOpen()) {
 246     Error("Run", "opening AliESDs.root failed");
 247     if (fStopOnError) {CleanUp(file); return kFALSE;}
 248   }
 249
 250   // loop over events
 251   for (Int_t iEvent = 0; iEvent < fRunLoader->GetNumberOfEvents(); iEvent++) {
 252     Info("Run", "processing event %d", iEvent);
 253     AliESD* esd = new AliESD;
 254     fRunLoader->GetEvent(iEvent);
 255     esd->SetRunNumber(aliRun->GetRunNumber());
 256     esd->SetEventNumber(aliRun->GetEvNumber());
 257     esd->SetMagneticField(aliRun->Field()->SolenoidField());
 258
 259     // barrel tracking
 260     if (fRunTracking) {
 261       if (!RunTracking(esd)) {
 262         if (fStopOnError) {CleanUp(file); return kFALSE;}
 263       }
 264     }
 265
 266     // fill ESD
 267     if (!fFillESD.IsNull()) {
 268       if (!FillESD(esd, fFillESD)) {
 269         if (fStopOnError) {CleanUp(file); return kFALSE;}
 270       }
 271     }
 272
 273     // combined PID
 274     AliESDpid::MakePID(esd);
 275
 276     // write ESD
 277     char name[100];
 278     sprintf(name, "ESD%d", iEvent);
 279     file->cd();
 280     if (!esd->Write(name)) {
 281       Error("Run", "writing ESD failed");
 282       if (fStopOnError) {CleanUp(file); return kFALSE;}
 283     }
 284     delete esd;
 285   }
 286
 287   CleanUp(file);
 288
 289   return kTRUE;
 290 }
 291
 292
 293 //_____________________________________________________________________________
 294 Bool_t AliReconstruction::RunReconstruction(const TString& detectors)
 295 {
 296 // run the reconstruction
 297
 298   TStopwatch stopwatch;
 299   stopwatch.Start();
 300
 301   TString detStr = detectors;
 302   TObjArray* detArray = fRunLoader->GetAliRun()->Detectors();
 303   for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
 304     AliModule* det = (AliModule*) detArray->At(iDet);
 305     if (!det || !det->IsActive()) continue;
 306     if (IsSelected(det->GetName(), detStr)) {
 307       Info("RunReconstruction", "running reconstruction for %s",
 308            det->GetName());
 309       TStopwatch stopwatchDet;
 310       stopwatchDet.Start();
 311       det->Reconstruct();
 312       Info("RunReconstruction", "execution time for %s:", det->GetName());
 313       stopwatchDet.Print();
 314     }
 315   }
 316
 317   if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
 318     Error("RunReconstruction", "the following detectors were not found: %s",
 319           detStr.Data());
 320     if (fStopOnError) return kFALSE;
 321   }
 322
 323   Info("RunReconstruction", "execution time:");
 324   stopwatch.Print();
 325
 326   return kTRUE;
 327 }
 328
 329 //_____________________________________________________________________________
 330 Bool_t AliReconstruction::RunTracking(AliESD* esd)
 331 {
 332 // run the barrel tracking
 333
 334   TStopwatch stopwatch;
 335   stopwatch.Start();
 336
 337   // get the primary vertex (from MC for the moment)
 338   TArrayF vertex(3);
 339   fRunLoader->GetHeader()->GenEventHeader()->PrimaryVertex(vertex);
 340   Double_t vtxPos[3] = {vertex[0], vertex[1], vertex[2]};
 341   Double_t vtxCov[6] = {
 342     0.005,
 343     0.000, 0.005,
 344     0.000, 0.000, 0.010
 345   };
 346   Double_t vtxErr[3] = {vtxCov[0], vtxCov[2], vtxCov[5]}; // diag. elements
 347   esd->SetVertex(vtxPos, vtxCov);
 348   fITSTracker->SetVertex(vtxPos, vtxErr);
 349   fTPCTracker->SetVertex(vtxPos, vtxErr);
 350   fTRDTracker->SetVertex(vtxPos, vtxErr);
 351
 352   // TPC tracking
 353   Info("RunTracking", "TPC tracking");
 354   fTPCLoader->LoadRecPoints("read");
 355   TTree* tpcTree = fTPCLoader->TreeR();
 356   if (!tpcTree) {
 357     Error("RunTracking", "Can't get the TPC cluster tree");
 358     return kFALSE;
 359   }
 360   fTPCTracker->LoadClusters(tpcTree);
 361   if (fTPCTracker->Clusters2Tracks(esd) != 0) {
 362     Error("RunTracking", "TPC Clusters2Tracks failed");
 363     return kFALSE;
 364   }
 365
 366   fRunLoader->GetAliRun()->GetDetector("TPC")->FillESD(esd); // preliminary PID
 367   AliESDpid::MakePID(esd);                  // for the ITS tracker
 368
 369   // ITS tracking
 370   Info("RunTracking", "ITS tracking");
 371   fITSLoader->LoadRecPoints("read");
 372   TTree* itsTree = fITSLoader->TreeR();
 373   if (!itsTree) {
 374     Error("RunTracking", "Can't get the ITS cluster tree");
 375     return kFALSE;
 376   }
 377   fITSTracker->LoadClusters(itsTree);
 378   if (fITSTracker->Clusters2Tracks(esd) != 0) {
 379     Error("RunTracking", "ITS Clusters2Tracks failed");
 380     return kFALSE;
 381   }
 382
 383   // ITS back propagation
 384   Info("RunTracking", "ITS back propagation");
 385   if (fITSTracker->PropagateBack(esd) != 0) {
 386     Error("RunTracking", "ITS backward propagation failed");
 387     return kFALSE;
 388   }
 389
 390   // TPC back propagation
 391   Info("RunTracking", "TPC back propagation");
 392   if (fTPCTracker->PropagateBack(esd) != 0) {
 393     Error("RunTracking", "TPC backward propagation failed");
 394     return kFALSE;
 395   }
 396
 397   // TRD back propagation
 398   Info("RunTracking", "TRD back propagation");
 399   fTRDLoader->LoadRecPoints("read");
 400   TTree* trdTree = fTRDLoader->TreeR();
 401   if (!trdTree) {
 402     Error("RunTracking", "Can't get the TRD cluster tree");
 403     return kFALSE;
 404   }
 405   fTRDTracker->LoadClusters(trdTree);
 406   if (fTRDTracker->PropagateBack(esd) != 0) {
 407     Error("RunTracking", "TRD backward propagation failed");
 408     return kFALSE;
 409   }
 410
 411   // TOF back propagation
 412   Info("RunTracking", "TOF back propagation");
 413   fTOFLoader->LoadDigits("read");
 414   TTree* tofTree = fTOFLoader->TreeD();
 415   if (!tofTree) {
 416     Error("RunTracking", "Can't get the TOF digits tree");
 417     return kFALSE;
 418   }
 419   fTOFTracker->LoadClusters(tofTree);
 420   if (fTOFTracker->PropagateBack(esd) != 0) {
 421     Error("RunTracking", "TOF backward propagation failed");
 422     return kFALSE;
 423   }
 424   fTOFTracker->UnloadClusters();
 425   fTOFLoader->UnloadDigits();
 426
 427   // TRD inward refit
 428   Info("RunTracking", "TRD inward refit");
 429   if (fTRDTracker->RefitInward(esd) != 0) {
 430     Error("RunTracking", "TRD inward refit failed");
 431     return kFALSE;
 432   }
 433   fTRDTracker->UnloadClusters();
 434   fTRDLoader->UnloadRecPoints();
 435
 436   // TPC inward refit
 437   Info("RunTracking", "TPC inward refit");
 438   if (fTPCTracker->RefitInward(esd) != 0) {
 439     Error("RunTracking", "TPC inward refit failed");
 440     return kFALSE;
 441   }
 442   fTPCTracker->UnloadClusters();
 443   fTPCLoader->UnloadRecPoints();
 444
 445   // ITS inward refit
 446   Info("RunTracking", "ITS inward refit");
 447   if (fITSTracker->RefitInward(esd) != 0) {
 448     Error("RunTracking", "ITS inward refit failed");
 449     return kFALSE;
 450   }
 451   fITSTracker->UnloadClusters();
 452   fITSLoader->UnloadRecPoints();
 453
 454   Info("RunTracking", "execution time:");
 455   stopwatch.Print();
 456
 457   return kTRUE;
 458 }
 459
 460 //_____________________________________________________________________________
 461 Bool_t AliReconstruction::FillESD(AliESD* esd, const TString& detectors)
 462 {
 463 // fill the event summary data
 464
 465   TStopwatch stopwatch;
 466   stopwatch.Start();
 467
 468   TString detStr = detectors;
 469   TObjArray* detArray = fRunLoader->GetAliRun()->Detectors();
 470   for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
 471     AliModule* det = (AliModule*) detArray->At(iDet);
 472     if (!det || !det->IsActive()) continue;
 473     if (IsSelected(det->GetName(), detStr)) {
 474       Info("FillESD", "filling ESD for %s",
 475            det->GetName());
 476       det->FillESD(esd);
 477     }
 478   }
 479
 480   if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
 481     Error("FillESD", "the following detectors were not found: %s",
 482           detStr.Data());
 483     if (fStopOnError) return kFALSE;
 484   }
 485
 486   Info("FillESD", "execution time:");
 487   stopwatch.Print();
 488
 489   return kTRUE;
 490 }
 491
 492
 493 //_____________________________________________________________________________
 494 Bool_t AliReconstruction::IsSelected(TString detName, TString& detectors) const
 495 {
 496 // check whether detName is contained in detectors
 497 // if yes, it is removed from detectors
 498
 499   // check if all detectors are selected
 500   if ((detectors.CompareTo("ALL") == 0) ||
 501       detectors.BeginsWith("ALL ") ||
 502       detectors.EndsWith(" ALL") ||
 503       detectors.Contains(" ALL ")) {
 504     detectors = "ALL";
 505     return kTRUE;
 506   }
 507
 508   // search for the given detector
 509   Bool_t result = kFALSE;
 510   if ((detectors.CompareTo(detName) == 0) ||
 511       detectors.BeginsWith(detName+" ") ||
 512       detectors.EndsWith(" "+detName) ||
 513       detectors.Contains(" "+detName+" ")) {
 514     detectors.ReplaceAll(detName, "");
 515     result = kTRUE;
 516   }
 517
 518   // clean up the detectors string
 519   while (detectors.Contains("  ")) detectors.ReplaceAll("  ", " ");
 520   while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
 521   while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
 522
 523   return result;
 524 }
 525
 526 //_____________________________________________________________________________
 527 void AliReconstruction::CleanUp(TFile* file)
 528 {
 529 // delete trackers and the run loader and close and delete the file
 530
 531   delete fITSTracker;
 532   fITSTracker = NULL;
 533   delete fTPCTracker;
 534   fTPCTracker = NULL;
 535   delete fTRDTracker;
 536   fTRDTracker = NULL;
 537   delete fTOFTracker;
 538   fTOFTracker = NULL;
 539
 540   delete fRunLoader;
 541   fRunLoader = NULL;
 542
 543   if (file) {
 544     file->Close();
 545     delete file;
 546   }
 547 }