STEER/AliRun.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 //  Control class for Alice C++                                              //
  21 //  Only one single instance of this class exists.                           //
  22 //  The object is created in main program aliroot                            //
  23 //  and is pointed by the global gAlice.                                     //
  24 //                                                                           //
  25 //   -Supports the list of all Alice Detectors (fModules).                 //
  26 //   -Supports the list of particles (fParticles).                           //
  27 //   -Supports the Trees.                                                    //
  28 //   -Supports the geometry.                                                 //
  29 //   -Supports the event display.                                            //
  30 //Begin_Html
  31 /*
  32 <img src="picts/AliRunClass.gif">
  33 */
  34 //End_Html
  35 //Begin_Html
  36 /*
  37 <img src="picts/alirun.gif">
  38 */
  39 //End_Html
  40 //                                                                           //
  41 ///////////////////////////////////////////////////////////////////////////////
  42
  43 #include <Riostream.h>
  44 #include <stdio.h>
  45 #include <stdlib.h>
  46 #include <string.h>
  47
  48 #include <TBRIK.h>
  49 #include <TBrowser.h>
  50 #include <TCint.h>
  51 #include <TFile.h>
  52 #include <TFolder.h>
  53 #include <TGeometry.h>
  54 #include <TKey.h>
  55 #include <TNode.h>
  56 #include <TObjectTable.h>
  57 #include <TParticle.h>
  58 #include <TROOT.h>
  59 #include <TRandom.h>
  60 #include <TRandom3.h>
  61 #include <TSystem.h>
  62 #include <TTree.h>
  63 #include <TVirtualMC.h>
  64
  65 #include "AliConfig.h"
  66 #include "AliDetector.h"
  67 #include "AliDisplay.h"
  68 #include "AliGenEventHeader.h"
  69 #include "AliGenerator.h"
  70 #include "AliHeader.h"
  71 #include "AliHit.h"
  72 #include "AliLego.h"
  73 #include "AliLegoGenerator.h"
  74 #include "AliLoader.h"
  75 #include "AliMCQA.h"
  76 #include "AliMagFC.h"
  77 #include "AliMagFCM.h"
  78 #include "AliMagFDM.h"
  79 #include "AliPDG.h"
  80 #include "AliRun.h"
  81 #include "AliRunLoader.h"
  82 #include "AliStack.h"
  83
  84 AliRun *gAlice;
  85
  86 ClassImp(AliRun)
  87
  88 //_______________________________________________________________________
  89 AliRun::AliRun():
  90   fRun(0),
  91   fEvent(0),
  92   fEventNrInRun(0),
  93   fEventsPerRun(0),
  94   fDebug(0),
  95   fModules(0),
  96   fGeometry(0),
  97   fDisplay(0),
  98   fTimer(),
  99   fField(0),
 100   fMC(0),
 101   fImedia(0),
 102   fNdets(0),
 103   fTrRmax(1.e10),
 104   fTrZmax(1.e10),
 105   fGenerator(0),
 106   fInitDone(kFALSE),
 107   fLego(0),
 108   fPDGDB(0),  //Particle factory object
 109   fHitLists(0),
 110   fEventEnergy(0),
 111   fSummEnergy(0),
 112   fSum2Energy(0),
 113   fConfigFunction("\0"),
 114   fRandom(0),
 115   fMCQA(0),
 116   fTransParName("\0"),
 117   fRunLoader(0x0)
 118 {
 119   //
 120   // Default constructor for AliRun
 121   //
 122   AliConfig::Instance();//skowron 29 Feb 2002
 123                         //ensures that the folder structure is build
 124 }
 125
 126 //_______________________________________________________________________
 127 AliRun::AliRun(const AliRun& arun):
 128   TVirtualMCApplication(arun),
 129   fRun(0),
 130   fEvent(0),
 131   fEventNrInRun(0),
 132   fEventsPerRun(0),
 133   fDebug(0),
 134   fModules(0),
 135   fGeometry(0),
 136   fDisplay(0),
 137   fTimer(),
 138   fField(0),
 139   fMC(0),
 140   fImedia(0),
 141   fNdets(0),
 142   fTrRmax(1.e10),
 143   fTrZmax(1.e10),
 144   fGenerator(0),
 145   fInitDone(kFALSE),
 146   fLego(0),
 147   fPDGDB(0),  //Particle factory object
 148   fHitLists(0),
 149   fEventEnergy(0),
 150   fSummEnergy(0),
 151   fSum2Energy(0),
 152   fConfigFunction("\0"),
 153   fRandom(0),
 154   fMCQA(0),
 155   fTransParName("\0"),
 156   fRunLoader(0x0)
 157 {
 158   //
 159   // Copy constructor for AliRun
 160   //
 161   arun.Copy(*this);
 162 }
 163
 164 //_____________________________________________________________________________
 165 AliRun::AliRun(const char *name, const char *title):
 166   TVirtualMCApplication(name,title),
 167   fRun(0),
 168   fEvent(0),
 169   fEventNrInRun(0),
 170   fEventsPerRun(0),
 171   fDebug(0),
 172   fModules(new TObjArray(77)), // Support list for the Detectors
 173   fGeometry(0),
 174   fDisplay(0),
 175   fTimer(),
 176   fField(0),
 177   fMC(gMC),
 178   fImedia(new TArrayI(1000)),
 179   fNdets(0),
 180   fTrRmax(1.e10),
 181   fTrZmax(1.e10),
 182   fGenerator(0),
 183   fInitDone(kFALSE),
 184   fLego(0),
 185   fPDGDB(TDatabasePDG::Instance()),        //Particle factory object!
 186   fHitLists(new TList()),                  // Create HitLists list
 187   fEventEnergy(0),
 188   fSummEnergy(0),
 189   fSum2Energy(0),
 190   fConfigFunction("Config();"),
 191   fRandom(new TRandom3()),
 192   fMCQA(0),
 193   fTransParName("\0"),
 194   fRunLoader(0x0)
 195 {
 196   //
 197   //  Constructor for the main processor.
 198   //  Creates the geometry
 199   //  Creates the list of Detectors.
 200   //  Creates the list of particles.
 201   //
 202
 203   gAlice     = this;
 204
 205   // Set random number generator
 206   gRandom = fRandom;
 207
 208   if (gSystem->Getenv("CONFIG_SEED")) {
 209      gRandom->SetSeed(static_cast<UInt_t>(atoi(gSystem->Getenv("CONFIG_SEED"))));
 210   }
 211
 212   // Add to list of browsable
 213   gROOT->GetListOfBrowsables()->Add(this,name);
 214   // Create the TNode geometry for the event display
 215   BuildSimpleGeometry();
 216
 217   // Create default mag field
 218   SetField();
 219
 220   // Prepare the tracking medium lists
 221   for(Int_t i=0;i<1000;i++) (*fImedia)[i]=-99;
 222
 223   // Add particle list to configuration
 224   AliConfig::Instance()->Add(fPDGDB);
 225
 226   // Set transport parameters
 227   SetTransPar();
 228 }
 229
 230
 231 //_______________________________________________________________________
 232 AliRun::~AliRun()
 233 {
 234   //
 235   // Default AliRun destructor
 236   //
 237   gROOT->GetListOfBrowsables()->Remove(this);
 238
 239   if (fRunLoader)
 240    {
 241     TFolder* evfold = fRunLoader->GetEventFolder();
 242     TFolder* modfold = dynamic_cast<TFolder*>(evfold->FindObjectAny(AliConfig::GetModulesFolderName()));
 243     TIter next(fModules);
 244     AliModule *mod;
 245     while((mod = (AliModule*)next()))
 246      {
 247        modfold->Remove(mod);
 248      }
 249    }
 250
 251   delete fImedia;
 252   delete fField;
 253   // delete fMC;
 254   delete gMC; gMC=0;
 255   delete fGeometry;
 256   delete fDisplay;
 257   delete fGenerator;
 258   delete fLego;
 259   if (fModules) {
 260     fModules->Delete();
 261     delete fModules;
 262   }
 263
 264   delete fHitLists;
 265   delete fPDGDB;
 266   delete fMCQA;
 267 }
 268
 269 //_______________________________________________________________________
 270 void AliRun::Copy(AliRun &) const
 271 {
 272   Fatal("Copy","Not implemented!\n");
 273 }
 274
 275 //_______________________________________________________________________
 276 void AliRun::AddHit(Int_t id, Int_t track, Int_t *vol, Float_t *hits) const
 277 {
 278   //
 279   //  Add a hit to detector id
 280   //
 281   TObjArray &dets = *fModules;
 282   if(dets[id]) dynamic_cast<AliModule*>(dets[id])->AddHit(track,vol,hits);
 283 }
 284
 285 //_______________________________________________________________________
 286 void AliRun::AddDigit(Int_t id, Int_t *tracks, Int_t *digits) const
 287 {
 288   //
 289   // Add digit to detector id
 290   //
 291   TObjArray &dets = *fModules;
 292   if(dets[id]) dynamic_cast<AliModule*>(dets[id])->AddDigit(tracks,digits);
 293 }
 294
 295 //_______________________________________________________________________
 296 void AliRun::Browse(TBrowser *b)
 297 {
 298   //
 299   // Called when the item "Run" is clicked on the left pane
 300   // of the Root browser.
 301   // It displays the Root Trees and all detectors.
 302   //
 303   //detectors are in folders anyway
 304   b->Add(fMCQA,"AliMCQA");
 305 }
 306
 307 //_______________________________________________________________________
 308 void AliRun::Build()
 309 {
 310   //
 311   // Initialize Alice geometry
 312   // Dummy routine
 313   //
 314 }
 315
 316 //_______________________________________________________________________
 317 void AliRun::BuildSimpleGeometry()
 318 {
 319   //
 320   // Create a simple TNode geometry used by Root display engine
 321   //
 322   // Initialise geometry
 323   //
 324   fGeometry = new TGeometry("AliceGeom","Galice Geometry for Hits");
 325   new TMaterial("void","Vacuum",0,0,0);  //Everything is void
 326   TBRIK *brik = new TBRIK("S_alice","alice volume","void",2000,2000,3000);
 327   brik->SetVisibility(0);
 328   new TNode("alice","alice","S_alice");
 329 }
 330
 331 //_______________________________________________________________________
 332 void AliRun::CleanDetectors()
 333 {
 334   //
 335   // Clean Detectors at the end of event
 336   //
 337    fRunLoader->CleanDetectors();
 338 }
 339
 340 //_______________________________________________________________________
 341 Int_t AliRun::DistancetoPrimitive(Int_t, Int_t) const
 342 {
 343   //
 344   // Return the distance from the mouse to the AliRun object
 345   // Dummy routine
 346   //
 347   return 9999;
 348 }
 349
 350 //_______________________________________________________________________
 351 void AliRun::DumpPart (Int_t i) const
 352 {
 353   //
 354   // Dumps particle i in the stack
 355   //
 356    if (fRunLoader->Stack())
 357     fRunLoader->Stack()->DumpPart(i);
 358 }
 359
 360 //_______________________________________________________________________
 361 void AliRun::DumpPStack () const
 362 {
 363   //
 364   // Dumps the particle stack
 365   //
 366    if (fRunLoader->Stack())
 367     fRunLoader->Stack()->DumpPStack();
 368 }
 369
 370 //_______________________________________________________________________
 371 void  AliRun::SetField(AliMagF* magField)
 372 {
 373     // Set Magnetic Field Map
 374     fField = magField;
 375     fField->ReadField();
 376 }
 377
 378 //_______________________________________________________________________
 379 void AliRun::SetField(Int_t type, Int_t version, Float_t scale,
 380                       Float_t maxField, char* filename)
 381 {
 382   //
 383   //  Set magnetic field parameters
 384   //  type      Magnetic field transport flag 0=no field, 2=helix, 3=Runge Kutta
 385   //  version   Magnetic field map version (only 1 active now)
 386   //  scale     Scale factor for the magnetic field
 387   //  maxField  Maximum value for the magnetic field
 388
 389   //
 390   // --- Sanity check on mag field flags
 391   if(fField) delete fField;
 392   if(version==1) {
 393     fField = new AliMagFC("Map1"," ",type,scale,maxField);
 394   } else if(version<=2) {
 395     fField = new AliMagFCM("Map2-3",filename,type,scale,maxField);
 396     fField->ReadField();
 397   } else if(version==3) {
 398     fField = new AliMagFDM("Map4",filename,type,scale,maxField);
 399     fField->ReadField();
 400   } else {
 401     Warning("SetField","Invalid map %d\n",version);
 402   }
 403 }
 404
 405 //_____________________________________________________________________________
 406
 407 void AliRun::InitLoaders()
 408 {
 409   //creates list of getters
 410   if (GetDebug()) Info("InitLoaders","");
 411   TIter next(fModules);
 412   AliModule *mod;
 413   while((mod = (AliModule*)next()))
 414    {
 415      AliDetector *det = dynamic_cast<AliDetector*>(mod);
 416      if (det)
 417       {
 418         if (GetDebug()) Info("InitLoaders"," Adding %s ",det->GetName());
 419         fRunLoader->AddLoader(det);
 420       }
 421    }
 422   if (GetDebug()) Info("InitLoaders","Done");
 423 }
 424 //_____________________________________________________________________________
 425
 426 void AliRun::FinishRun()
 427 {
 428   //
 429   // Called at the end of the run.
 430   //
 431
 432   if(fLego)
 433    {
 434     if (GetDebug()) Info("FinishRun"," Finish Lego");
 435     fRunLoader->CdGAFile();
 436     fLego->FinishRun();
 437    }
 438
 439   // Clean detector information
 440   TIter next(fModules);
 441   AliModule *detector;
 442   while((detector = dynamic_cast<AliModule*>(next()))) {
 443     if (GetDebug()) Info("FinishRun"," %s->FinishRun()",detector->GetName());
 444     detector->FinishRun();
 445   }
 446
 447   //Output energy summary tables
 448   if (GetDebug()) Info("FinishRun"," EnergySummary()");
 449   EnergySummary();
 450
 451   if (GetDebug()) Info("FinishRun"," fRunLoader->WriteHeader(OVERWRITE)");
 452   fRunLoader->WriteHeader("OVERWRITE");
 453
 454   // Write AliRun info and all detectors parameters
 455   fRunLoader->CdGAFile();
 456   Write(0,TObject::kOverwrite);//write AliRun
 457   fRunLoader->Write(0,TObject::kOverwrite);//write RunLoader itself
 458
 459   // Clean tree information
 460   if (GetDebug()) Info("FinishRun"," fRunLoader->Stack()->FinishRun()");
 461   fRunLoader->Stack()->FinishRun();
 462
 463   // Clean detector information
 464   if (GetDebug()) Info("FinishRun"," fGenerator->FinishRun()");
 465   fGenerator->FinishRun();
 466
 467   fRunLoader->Synchronize();
 468 }
 469
 470 //_______________________________________________________________________
 471 void AliRun::FlagTrack(Int_t track)
 472 {
 473   // Delegate to stack
 474   //
 475     fRunLoader->Stack()->FlagTrack(track);
 476 }
 477
 478 //_______________________________________________________________________
 479 void AliRun::EnergySummary()
 480 {
 481   //
 482   // Print summary of deposited energy
 483   //
 484
 485   Int_t ndep=0;
 486   Float_t edtot=0;
 487   Float_t ed, ed2;
 488   Int_t kn, i, left, j, id;
 489   const Float_t kzero=0;
 490   Int_t ievent=fRunLoader->GetHeader()->GetEvent()+1;
 491   //
 492   // Energy loss information
 493   if(ievent) {
 494     printf("***************** Energy Loss Information per event (GEV) *****************\n");
 495     for(kn=1;kn<fEventEnergy.GetSize();kn++) {
 496       ed=fSummEnergy[kn];
 497       if(ed>0) {
 498         fEventEnergy[ndep]=kn;
 499         if(ievent>1) {
 500           ed=ed/ievent;
 501           ed2=fSum2Energy[kn];
 502           ed2=ed2/ievent;
 503           ed2=100*TMath::Sqrt(TMath::Max(ed2-ed*ed,kzero))/ed;
 504         } else
 505           ed2=99;
 506         fSummEnergy[ndep]=ed;
 507         fSum2Energy[ndep]=TMath::Min(static_cast<Float_t>(99.),TMath::Max(ed2,kzero));
 508         edtot+=ed;
 509         ndep++;
 510       }
 511     }
 512     for(kn=0;kn<(ndep-1)/3+1;kn++) {
 513       left=ndep-kn*3;
 514       for(i=0;i<(3<left?3:left);i++) {
 515         j=kn*3+i;
 516         id=Int_t (fEventEnergy[j]+0.1);
 517         printf(" %s %10.3f +- %10.3f%%;",gMC->VolName(id),fSummEnergy[j],fSum2Energy[j]);
 518       }
 519       printf("\n");
 520     }
 521     //
 522     // Relative energy loss in different detectors
 523     printf("******************** Relative Energy Loss per event ********************\n");
 524     printf("Total energy loss per event %10.3f GeV\n",edtot);
 525     for(kn=0;kn<(ndep-1)/5+1;kn++) {
 526       left=ndep-kn*5;
 527       for(i=0;i<(5<left?5:left);i++) {
 528         j=kn*5+i;
 529         id=Int_t (fEventEnergy[j]+0.1);
 530         printf(" %s %10.3f%%;",gMC->VolName(id),100*fSummEnergy[j]/edtot);
 531       }
 532       printf("\n");
 533     }
 534     for(kn=0;kn<75;kn++) printf("*");
 535     printf("\n");
 536   }
 537   //
 538   // Reset the TArray's
 539   //  fEventEnergy.Set(0);
 540   //  fSummEnergy.Set(0);
 541   //  fSum2Energy.Set(0);
 542 }
 543
 544 //_______________________________________________________________________
 545 void AliRun::Announce() const
 546 {
 547   //
 548   // Announce the current version of AliRoot
 549   //
 550   printf("%70s",
 551          "****************************************************************\n");
 552   printf("%6s","*");printf("%64s","*\n");
 553
 554   printf("%6s","*");
 555   printf("    You are running AliRoot version NewIO\n");
 556
 557   printf("%6s","*");
 558   printf("    The cvs tag for the current program is $Name$\n");
 559
 560   printf("%6s","*");printf("%64s","*\n");
 561   printf("%70s",
 562          "****************************************************************\n");
 563 }
 564
 565 //_______________________________________________________________________
 566 AliModule *AliRun::GetModule(const char *name) const
 567 {
 568   //
 569   // Return pointer to detector from name
 570   //
 571   return dynamic_cast<AliModule*>(fModules->FindObject(name));
 572 }
 573
 574 //_______________________________________________________________________
 575 AliDetector *AliRun::GetDetector(const char *name) const
 576 {
 577   //
 578   // Return pointer to detector from name
 579   //
 580   return dynamic_cast<AliDetector*>(fModules->FindObject(name));
 581 }
 582
 583 //_______________________________________________________________________
 584 Int_t AliRun::GetModuleID(const char *name) const
 585 {
 586   //
 587   // Return galice internal detector identifier from name
 588   //
 589   Int_t i=-1;
 590   TObject *mod=fModules->FindObject(name);
 591   if(mod) i=fModules->IndexOf(mod);
 592   return i;
 593 }
 594
 595 //_______________________________________________________________________
 596 Int_t AliRun::GetEvent(Int_t event)
 597 {
 598 //
 599 // Reloads data containers in folders # event
 600 // Set branch addresses
 601 //
 602   if (fRunLoader == 0x0)
 603    {
 604      Error("GetEvent","RunLoader is not set. Can not load data.");
 605      return -1;
 606    }
 607 /*****************************************/
 608 /****   P R E    R E L O A D I N G    ****/
 609 /*****************************************/
 610 // Reset existing structures
 611   ResetHits();
 612   ResetTrackReferences();
 613   ResetDigits();
 614   ResetSDigits();
 615
 616 /*****************************************/
 617 /****       R  E  L  O  A  D          ****/
 618 /*****************************************/
 619
 620   fRunLoader->GetEvent(event);
 621
 622 /*****************************************/
 623 /****  P O S T    R E L O A D I N G   ****/
 624 /*****************************************/
 625
 626   // Set Trees branch addresses
 627   TIter next(fModules);
 628   AliModule *detector;
 629   while((detector = dynamic_cast<AliModule*>(next())))
 630    {
 631      detector->SetTreeAddress();
 632    }
 633
 634   return fRunLoader->GetHeader()->GetNtrack();
 635 }
 636
 637 //_______________________________________________________________________
 638 TGeometry *AliRun::GetGeometry()
 639 {
 640   //
 641   // Import Alice geometry from current file
 642   // Return pointer to geometry object
 643   //
 644   if (!fGeometry) fGeometry = dynamic_cast<TGeometry*>(gDirectory->Get("AliceGeom"));
 645   //
 646   // Unlink and relink nodes in detectors
 647   // This is bad and there must be a better way...
 648   //
 649
 650   TIter next(fModules);
 651   AliModule *detector;
 652   while((detector = dynamic_cast<AliModule*>(next()))) {
 653     TList *dnodes=detector->Nodes();
 654     Int_t j;
 655     TNode *node, *node1;
 656     for ( j=0; j<dnodes->GetSize(); j++) {
 657       node = dynamic_cast<TNode*>(dnodes->At(j));
 658       node1 = fGeometry->GetNode(node->GetName());
 659       dnodes->Remove(node);
 660       dnodes->AddAt(node1,j);
 661     }
 662   }
 663   return fGeometry;
 664 }
 665
 666 //_______________________________________________________________________
 667 Int_t AliRun::GetPrimary(Int_t track) const
 668 {
 669   //
 670   // return number of primary that has generated track
 671   //
 672     return fRunLoader->Stack()->GetPrimary(track);
 673 }
 674
 675 //_______________________________________________________________________
 676 void AliRun::MediaTable()
 677 {
 678   //
 679   // Built media table to get from the media number to
 680   // the detector id
 681   //
 682
 683   Int_t kz, nz, idt, lz, i, k, ind;
 684   //  Int_t ibeg;
 685   TObjArray &dets = *gAlice->Detectors();
 686   AliModule *det;
 687   //
 688   // For all detectors
 689   for (kz=0;kz<fNdets;kz++) {
 690     // If detector is defined
 691     if((det=dynamic_cast<AliModule*>(dets[kz]))) {
 692         TArrayI &idtmed = *(det->GetIdtmed());
 693         for(nz=0;nz<100;nz++) {
 694         // Find max and min material number
 695         if((idt=idtmed[nz])) {
 696           det->LoMedium() = det->LoMedium() < idt ? det->LoMedium() : idt;
 697           det->HiMedium() = det->HiMedium() > idt ? det->HiMedium() : idt;
 698         }
 699       }
 700       if(det->LoMedium() > det->HiMedium()) {
 701         det->LoMedium() = 0;
 702         det->HiMedium() = 0;
 703       } else {
 704         if(det->HiMedium() > fImedia->GetSize()) {
 705           Error("MediaTable","Increase fImedia from %d to %d",
 706                 fImedia->GetSize(),det->HiMedium());
 707           return;
 708         }
 709         // Tag all materials in rage as belonging to detector kz
 710         for(lz=det->LoMedium(); lz<= det->HiMedium(); lz++) {
 711           (*fImedia)[lz]=kz;
 712         }
 713       }
 714     }
 715   }
 716   //
 717   // Print summary table
 718   printf(" Traking media ranges:\n");
 719   for(i=0;i<(fNdets-1)/6+1;i++) {
 720     for(k=0;k< (6<fNdets-i*6?6:fNdets-i*6);k++) {
 721       ind=i*6+k;
 722       det=dynamic_cast<AliModule*>(dets[ind]);
 723       if(det)
 724         printf(" %6s: %3d -> %3d;",det->GetName(),det->LoMedium(),
 725                det->HiMedium());
 726       else
 727         printf(" %6s: %3d -> %3d;","NULL",0,0);
 728     }
 729     printf("\n");
 730   }
 731 }
 732
 733 //_______________________________________________________________________
 734 void AliRun::SetGenerator(AliGenerator *generator)
 735 {
 736   //
 737   // Load the event generator
 738   //
 739   if(!fGenerator) fGenerator = generator;
 740 }
 741
 742 //_______________________________________________________________________
 743 void AliRun::ResetGenerator(AliGenerator *generator)
 744 {
 745   //
 746   // Load the event generator
 747   //
 748   if(fGenerator)
 749     if(generator)
 750       Warning("ResetGenerator","Replacing generator %s with %s\n",
 751               fGenerator->GetName(),generator->GetName());
 752     else
 753       Warning("ResetGenerator","Replacing generator %s with NULL\n",
 754               fGenerator->GetName());
 755   fGenerator = generator;
 756 }
 757
 758 //_______________________________________________________________________
 759 void AliRun::SetTransPar(const char *filename)
 760 {
 761   //
 762   // Sets the file name for transport parameters
 763   //
 764   fTransParName = filename;
 765 }
 766
 767 //_______________________________________________________________________
 768 void AliRun::SetBaseFile(const char *filename)
 769 {
 770   fBaseFileName = filename;
 771 }
 772
 773 //_______________________________________________________________________
 774 void AliRun::ReadTransPar()
 775 {
 776   //
 777   // Read filename to set the transport parameters
 778   //
 779
 780
 781   const Int_t kncuts=10;
 782   const Int_t knflags=11;
 783   const Int_t knpars=kncuts+knflags;
 784   const char kpars[knpars][7] = {"CUTGAM" ,"CUTELE","CUTNEU","CUTHAD","CUTMUO",
 785                                "BCUTE","BCUTM","DCUTE","DCUTM","PPCUTM","ANNI",
 786                                "BREM","COMP","DCAY","DRAY","HADR","LOSS",
 787                                "MULS","PAIR","PHOT","RAYL"};
 788   char line[256];
 789   char detName[7];
 790   char* filtmp;
 791   Float_t cut[kncuts];
 792   Int_t flag[knflags];
 793   Int_t i, itmed, iret, ktmed, kz;
 794   FILE *lun;
 795   //
 796   // See whether the file is there
 797   filtmp=gSystem->ExpandPathName(fTransParName.Data());
 798   lun=fopen(filtmp,"r");
 799   delete [] filtmp;
 800   if(!lun) {
 801     Warning("ReadTransPar","File %s does not exist!\n",fTransParName.Data());
 802     return;
 803   }
 804   //
 805   if(fDebug) {
 806     printf(" "); for(i=0;i<60;i++) printf("*"); printf("\n");
 807     printf(" *%59s\n","*");
 808     printf(" *       Please check carefully what you are doing!%10s\n","*");
 809     printf(" *%59s\n","*");
 810   }
 811   //
 812   while(1) {
 813     // Initialise cuts and flags
 814     for(i=0;i<kncuts;i++) cut[i]=-99;
 815     for(i=0;i<knflags;i++) flag[i]=-99;
 816     itmed=0;
 817     for(i=0;i<256;i++) line[i]='\0';
 818     // Read up to the end of line excluded
 819     iret=fscanf(lun,"%[^\n]",line);
 820     if(iret<0) {
 821       //End of file
 822       fclose(lun);
 823       if(fDebug){
 824         printf(" *%59s\n","*");
 825         printf(" "); for(i=0;i<60;i++) printf("*"); printf("\n");
 826       }
 827       return;
 828     }
 829     // Read the end of line
 830     fscanf(lun,"%*c");
 831     if(!iret) continue;
 832     if(line[0]=='*') continue;
 833     // Read the numbers
 834     iret=sscanf(line,"%s %d %f %f %f %f %f %f %f %f %f %f %d %d %d %d %d %d %d %d %d %d %d",
 835                 detName,&itmed,&cut[0],&cut[1],&cut[2],&cut[3],&cut[4],&cut[5],&cut[6],&cut[7],&cut[8],
 836                 &cut[9],&flag[0],&flag[1],&flag[2],&flag[3],&flag[4],&flag[5],&flag[6],&flag[7],
 837                 &flag[8],&flag[9],&flag[10]);
 838     if(!iret) continue;
 839     if(iret<0) {
 840       //reading error
 841       Warning("ReadTransPar","Error reading file %s\n",fTransParName.Data());
 842       continue;
 843     }
 844     // Check that the module exist
 845     AliModule *mod = GetModule(detName);
 846     if(mod) {
 847       // Get the array of media numbers
 848       TArrayI &idtmed = *mod->GetIdtmed();
 849       // Check that the tracking medium code is valid
 850       if(0<=itmed && itmed < 100) {
 851         ktmed=idtmed[itmed];
 852         if(!ktmed) {
 853           Warning("ReadTransPar","Invalid tracking medium code %d for %s\n",itmed,mod->GetName());
 854           continue;
 855         }
 856         // Set energy thresholds
 857         for(kz=0;kz<kncuts;kz++) {
 858           if(cut[kz]>=0) {
 859             if(fDebug) printf(" *  %-6s set to %10.3E for tracking medium code %4d for %s\n",
 860                    kpars[kz],cut[kz],itmed,mod->GetName());
 861             gMC->Gstpar(ktmed,kpars[kz],cut[kz]);
 862           }
 863         }
 864         // Set transport mechanisms
 865         for(kz=0;kz<knflags;kz++) {
 866           if(flag[kz]>=0) {
 867             if(fDebug) printf(" *  %-6s set to %10d for tracking medium code %4d for %s\n",
 868                    kpars[kncuts+kz],flag[kz],itmed,mod->GetName());
 869             gMC->Gstpar(ktmed,kpars[kncuts+kz],Float_t(flag[kz]));
 870           }
 871         }
 872       } else {
 873         Warning("ReadTransPar","Invalid medium code %d *\n",itmed);
 874         continue;
 875       }
 876     } else {
 877       if(fDebug) printf("%s::ReadTransParModule: %s not present\n",ClassName(),detName);
 878       continue;
 879     }
 880   }
 881 }
 882 //_____________________________________________________________________________
 883
 884 void AliRun::BeginEvent()
 885 {
 886   //
 887   // Clean-up previous event
 888   // Energy scores
 889   if (GetDebug())
 890    {
 891      Info("BeginEvent",">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>");
 892      Info("BeginEvent",">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>");
 893      Info("BeginEvent",">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>");
 894      Info("BeginEvent","          BEGINNING EVENT               ");
 895      Info("BeginEvent",">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>");
 896      Info("BeginEvent",">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>");
 897      Info("BeginEvent",">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>");
 898    }
 899
 900   /*******************************/
 901   /*   Clean after eventual      */
 902   /*   previous event            */
 903   /*******************************/
 904
 905
 906   //Set the next event in Run Loader -> Cleans trees (TreeK and all trees in detectors),
 907   fRunLoader->SetEventNumber(++fEventNrInRun);// sets new files, cleans the previous event stuff, if necessary, etc.,
 908   if (GetDebug()) Info("BeginEvent","EventNr is %d",fEventNrInRun);
 909
 910   fEventEnergy.Reset();
 911     // Clean detector information
 912
 913   if (fRunLoader->Stack())
 914     fRunLoader->Stack()->Reset();//clean stack -> tree is unloaded
 915   else
 916     fRunLoader->MakeStack();//or make a new one
 917
 918   if (GetDebug()) Info("BeginEvent","  fRunLoader->MakeTree(K)");
 919   fRunLoader->MakeTree("K");
 920   if (GetDebug()) Info("BeginEvent","  gMC->SetStack(fRunLoader->Stack())");
 921   gMC->SetStack(fRunLoader->Stack());//Was in InitMC - but was moved here
 922                                      //because we don't have guarantee that
 923                                      //stack pointer is not going to change from event to event
 924                          //since it bellobgs to header and is obtained via RunLoader
 925   //
 926   //  Reset all Detectors & kinematics & make/reset trees
 927   //
 928
 929   fRunLoader->GetHeader()->Reset(fRun,fEvent,fEventNrInRun);
 930 //  fRunLoader->WriteKinematics("OVERWRITE");  is there any reason to rewrite here since MakeTree does so
 931
 932   if (GetDebug()) Info("BeginEvent","  fRunLoader->MakeTrackRefsContainer()");
 933   fRunLoader->MakeTrackRefsContainer();//for insurance
 934
 935   if (GetDebug()) Info("BeginEvent","  ResetHits()");
 936   ResetHits();
 937   if (GetDebug()) Info("BeginEvent","  fRunLoader->MakeTree(H)");
 938   fRunLoader->MakeTree("H");
 939
 940   //
 941   if(fLego)
 942    {
 943     fLego->BeginEvent();
 944     return;
 945    }
 946
 947   //create new branches and SetAdresses
 948   TIter next(fModules);
 949   AliModule *detector;
 950   while((detector = (AliModule*)next()))
 951    {
 952     if (GetDebug()) Info("BeginEvent","  %s->MakeBranch(H)",detector->GetName());
 953     detector->MakeBranch("H");
 954     if (GetDebug()) Info("BeginEvent","  %s->MakeBranchTR()",detector->GetName());
 955     detector->MakeBranchTR();
 956     if (GetDebug()) Info("BeginEvent","  %s->SetTreeAddress()",detector->GetName());
 957     detector->SetTreeAddress();
 958    }
 959 }
 960
 961 //_______________________________________________________________________
 962 TParticle* AliRun::Particle(Int_t i) const
 963 {
 964   if (fRunLoader)
 965    if (fRunLoader->Stack())
 966     return fRunLoader->Stack()->Particle(i);
 967   return 0x0;
 968 }
 969
 970 //_______________________________________________________________________
 971 void AliRun::ResetDigits()
 972 {
 973   //
 974   //  Reset all Detectors digits
 975   //
 976   TIter next(fModules);
 977   AliModule *detector;
 978   while((detector = dynamic_cast<AliModule*>(next()))) {
 979      detector->ResetDigits();
 980   }
 981 }
 982
 983 //_______________________________________________________________________
 984 void AliRun::ResetSDigits()
 985 {
 986   //
 987   //  Reset all Detectors digits
 988   //
 989   TIter next(fModules);
 990   AliModule *detector;
 991   while((detector = dynamic_cast<AliModule*>(next()))) {
 992      detector->ResetSDigits();
 993   }
 994 }
 995
 996 //_______________________________________________________________________
 997 void AliRun::ResetHits()
 998 {
 999   //
1000   //  Reset all Detectors hits
1001   //
1002   TIter next(fModules);
1003   AliModule *detector;
1004   while((detector = dynamic_cast<AliModule*>(next()))) {
1005      detector->ResetHits();
1006   }
1007 }
1008 //_______________________________________________________________________
1009
1010 void AliRun::ResetTrackReferences()
1011 {
1012   //
1013   //  Reset all Detectors hits
1014   //
1015   TIter next(fModules);
1016   AliModule *detector;
1017   while((detector = dynamic_cast<AliModule*>(next()))) {
1018      detector->ResetTrackReferences();
1019   }
1020 }
1021 //_______________________________________________________________________
1022
1023 void AliRun::ResetPoints()
1024 {
1025   //
1026   // Reset all Detectors points
1027   //
1028   TIter next(fModules);
1029   AliModule *detector;
1030   while((detector = dynamic_cast<AliModule*>(next()))) {
1031      detector->ResetPoints();
1032   }
1033 }
1034 //_______________________________________________________________________
1035
1036 void AliRun::InitMC(const char *setup)
1037 {
1038   //
1039   // Initialize the Alice setup
1040   //
1041   Announce();
1042
1043   if(fInitDone) {
1044     Warning("Init","Cannot initialise AliRun twice!\n");
1045     return;
1046   }
1047
1048   gROOT->LoadMacro(setup);
1049   gInterpreter->ProcessLine(fConfigFunction.Data());
1050
1051   // Register MC in configuration
1052   AliConfig::Instance()->Add(gMC);
1053
1054   InitLoaders();
1055
1056   fRunLoader->MakeTree("E");
1057   fRunLoader->LoadKinematics("RECREATE");
1058   fRunLoader->LoadTrackRefs("RECREATE");
1059   fRunLoader->LoadHits("all","RECREATE");
1060
1061
1062   fRunLoader->CdGAFile();
1063
1064   gMC->DefineParticles();  //Create standard MC particles
1065   AliPDG::AddParticlesToPdgDataBase();
1066
1067   TObject *objfirst, *objlast;
1068
1069   fNdets = fModules->GetLast()+1;
1070
1071   //
1072   //=================Create Materials and geometry
1073   gMC->Init();
1074
1075   // Added also after in case of interactive initialisation of modules
1076   fNdets = fModules->GetLast()+1;
1077
1078   TIter next(fModules);
1079   AliModule *detector;
1080   while((detector = dynamic_cast<AliModule*>(next())))
1081    {
1082      objlast = gDirectory->GetList()->Last();
1083
1084      // Add Detector histograms in Detector list of histograms
1085      if (objlast) objfirst = gDirectory->GetList()->After(objlast);
1086      else         objfirst = gDirectory->GetList()->First();
1087      while (objfirst)
1088       {
1089         detector->Histograms()->Add(objfirst);
1090         objfirst = gDirectory->GetList()->After(objfirst);
1091       }
1092    }
1093    ReadTransPar(); //Read the cuts for all materials
1094
1095    MediaTable(); //Build the special IMEDIA table
1096
1097    //Initialise geometry deposition table
1098    fEventEnergy.Set(gMC->NofVolumes()+1);
1099    fSummEnergy.Set(gMC->NofVolumes()+1);
1100    fSum2Energy.Set(gMC->NofVolumes()+1);
1101
1102    //Compute cross-sections
1103    gMC->BuildPhysics();
1104
1105    //Write Geometry object to current file.
1106    fRunLoader->WriteGeometry();
1107
1108    fInitDone = kTRUE;
1109
1110    fMCQA = new AliMCQA(fNdets);
1111
1112    //
1113    // Save stuff at the beginning of the file to avoid file corruption
1114    Write();
1115    fEventNrInRun = -1; //important - we start Begin event from increasing current number in run
1116 }
1117
1118 //_______________________________________________________________________
1119
1120 void AliRun::RunMC(Int_t nevent, const char *setup)
1121 {
1122   //
1123   // Main function to be called to process a galice run
1124   // example
1125   //    Root > gAlice.Run();
1126   // a positive number of events will cause the finish routine
1127   // to be called
1128   //
1129   fEventsPerRun = nevent;
1130   // check if initialisation has been done
1131   if (!fInitDone) InitMC(setup);
1132
1133   // Create the Root Tree with one branch per detector
1134   //Hits moved to begin event -> now we are crating separate tree for each event
1135
1136   gMC->ProcessRun(nevent);
1137
1138   // End of this run, close files
1139   if(nevent>0) FinishRun();
1140 }
1141
1142 //_______________________________________________________________________
1143 void AliRun::RunReco(const char *selected, Int_t first, Int_t last)
1144 {
1145   //
1146   // Main function to be called to reconstruct Alice event
1147   //
1148    Int_t nev = fRunLoader->GetNumberOfEvents();
1149    if (GetDebug()) Info("RunReco","Found %d events",nev);
1150    Int_t nFirst = first;
1151    Int_t nLast  = (last < 0)? nev : last;
1152
1153    for (Int_t nevent = nFirst; nevent <= nLast; nevent++) {
1154      if (GetDebug()) Info("RunReco","Processing event %d",nevent);
1155      GetEvent(nevent);
1156      Digits2Reco(selected);
1157    }
1158 }
1159
1160 //_______________________________________________________________________
1161
1162 void AliRun::Hits2Digits(const char *selected)
1163 {
1164
1165    // Convert Hits to sumable digits
1166    //
1167    for (Int_t nevent=0; nevent<gAlice->TreeE()->GetEntries(); nevent++) {
1168      GetEvent(nevent);
1169      Hits2SDigits(selected);
1170      SDigits2Digits(selected);
1171    }
1172 }
1173
1174
1175 //_______________________________________________________________________
1176
1177 void AliRun::Tree2Tree(Option_t *option, const char *selected)
1178 {
1179   //
1180   // Function to transform the content of
1181   //
1182   // - TreeH to TreeS (option "S")
1183   // - TreeS to TreeD (option "D")
1184   // - TreeD to TreeR (option "R")
1185   //
1186   // If multiple options are specified ("SDR"), transformation will be done in sequence for
1187   // selected detector and for all detectors if none is selected (detector string
1188   // can contain blank separated list of detector names).
1189
1190
1191    const char *oS = strstr(option,"S");
1192    const char *oD = strstr(option,"D");
1193    const char *oR = strstr(option,"R");
1194
1195    TObjArray *detectors = Detectors();
1196
1197    TIter next(detectors);
1198
1199    AliDetector *detector = 0;
1200
1201    while((detector = dynamic_cast<AliDetector*>(next()))) {
1202      if (selected)
1203        if (strcmp(detector->GetName(),selected)) continue;
1204      if (detector->IsActive())
1205       {
1206
1207        AliLoader* loader = detector->GetLoader();
1208        if (loader == 0x0) continue;
1209
1210        if (oS)
1211         {
1212           if (GetDebug()) Info("Tree2Tree","Processing Hits2SDigits for %s ...",detector->GetName());
1213           loader->LoadHits("read");
1214           if (loader->TreeS() == 0x0) loader->MakeTree("S");
1215           detector->MakeBranch(option);
1216           detector->SetTreeAddress();
1217           detector->Hits2SDigits();
1218           loader->UnloadHits();
1219           loader->UnloadSDigits();
1220         }
1221        if (oD)
1222         {
1223           if (GetDebug()) Info("Tree2Tree","Processing SDigits2Digits for %s ...",detector->GetName());
1224           loader->LoadSDigits("read");
1225           if (loader->TreeD() == 0x0) loader->MakeTree("D");
1226           detector->MakeBranch(option);
1227           detector->SetTreeAddress();
1228           detector->SDigits2Digits();
1229           loader->UnloadSDigits();
1230           loader->UnloadDigits();
1231         }
1232        if (oR)
1233         {
1234           if (GetDebug()) Info("Tree2Tree","Processing Digits2Reco for %s ...",detector->GetName());
1235           loader->LoadDigits("read");
1236           if (loader->TreeR() == 0x0) loader->MakeTree("R");
1237           detector->MakeBranch(option);
1238           detector->SetTreeAddress();
1239           detector->Digits2Reco();
1240           loader->UnloadDigits();
1241           loader->UnloadRecPoints();
1242
1243         }
1244      }
1245    }
1246 }
1247
1248 //_______________________________________________________________________
1249 void AliRun::RunLego(const char *setup, Int_t nc1, Float_t c1min,
1250                      Float_t c1max,Int_t nc2,Float_t c2min,Float_t c2max,
1251                      Float_t rmin,Float_t rmax,Float_t zmax, AliLegoGenerator* gener)
1252 {
1253   //
1254   // Generates lego plots of:
1255   //    - radiation length map phi vs theta
1256   //    - radiation length map phi vs eta
1257   //    - interaction length map
1258   //    - g/cm2 length map
1259   //
1260   //  ntheta    bins in theta, eta
1261   //  themin    minimum angle in theta (degrees)
1262   //  themax    maximum angle in theta (degrees)
1263   //  nphi      bins in phi
1264   //  phimin    minimum angle in phi (degrees)
1265   //  phimax    maximum angle in phi (degrees)
1266   //  rmin      minimum radius
1267   //  rmax      maximum radius
1268   //
1269   //
1270   //  The number of events generated = ntheta*nphi
1271   //  run input parameters in macro setup (default="Config.C")
1272   //
1273   //  Use macro "lego.C" to visualize the 3 lego plots in spherical coordinates
1274   //Begin_Html
1275   /*
1276     <img src="picts/AliRunLego1.gif">
1277   */
1278   //End_Html
1279   //Begin_Html
1280   /*
1281     <img src="picts/AliRunLego2.gif">
1282   */
1283   //End_Html
1284   //Begin_Html
1285   /*
1286     <img src="picts/AliRunLego3.gif">
1287   */
1288   //End_Html
1289   //
1290
1291   // check if initialisation has been done
1292   if (!fInitDone) InitMC(setup);
1293   //Save current generator
1294   AliGenerator *gen=Generator();
1295
1296   // Set new generator
1297   if (!gener) gener  = new AliLegoGenerator();
1298   ResetGenerator(gener);
1299   //
1300   // Configure Generator
1301   gener->SetRadiusRange(rmin, rmax);
1302   gener->SetZMax(zmax);
1303   gener->SetCoor1Range(nc1, c1min, c1max);
1304   gener->SetCoor2Range(nc2, c2min, c2max);
1305
1306
1307   //Create Lego object
1308   fLego = new AliLego("lego",gener);
1309
1310   //Prepare MC for Lego Run
1311   gMC->InitLego();
1312
1313   //Run Lego Object
1314
1315   //gMC->ProcessRun(nc1*nc2+1);
1316   gMC->ProcessRun(nc1*nc2);
1317
1318   // Create only the Root event Tree
1319   fRunLoader->MakeTree("E");
1320
1321   // End of this run, close files
1322   FinishRun();
1323   // Restore current generator
1324   ResetGenerator(gen);
1325   // Delete Lego Object
1326   delete fLego; fLego=0;
1327 }
1328
1329 //_______________________________________________________________________
1330 void AliRun::SetConfigFunction(const char * config)
1331 {
1332   //
1333   // Set the signature of the function contained in Config.C to configure
1334   // the run
1335   //
1336   fConfigFunction=config;
1337 }
1338
1339 //_______________________________________________________________________
1340 void AliRun::SetCurrentTrack(Int_t track)
1341 {
1342   //
1343   // Set current track number
1344   //
1345     fRunLoader->Stack()->SetCurrentTrack(track);
1346 }
1347
1348 //_______________________________________________________________________
1349 void AliRun::PushTrack(Int_t done, Int_t parent, Int_t pdg, Float_t *pmom,
1350                       Float_t *vpos, Float_t *polar, Float_t tof,
1351                       TMCProcess mech, Int_t &ntr, Float_t weight, Int_t is)
1352 {
1353 // Delegate to stack
1354 //
1355     fRunLoader->Stack()->PushTrack(done, parent, pdg, pmom, vpos, polar, tof,
1356                      mech, ntr, weight, is);
1357 }
1358
1359 //_______________________________________________________________________
1360 void AliRun::PushTrack(Int_t done, Int_t parent, Int_t pdg,
1361                       Double_t px, Double_t py, Double_t pz, Double_t e,
1362                       Double_t vx, Double_t vy, Double_t vz, Double_t tof,
1363                       Double_t polx, Double_t poly, Double_t polz,
1364                       TMCProcess mech, Int_t &ntr, Float_t weight, Int_t is)
1365 {
1366   // Delegate to stack
1367   //
1368   fRunLoader->Stack()->PushTrack(done, parent, pdg, px, py, pz, e, vx, vy, vz, tof,
1369                                 polx, poly, polz, mech, ntr, weight, is);
1370 }
1371
1372 //_______________________________________________________________________
1373 void AliRun::SetHighWaterMark(const Int_t nt)
1374 {
1375     //
1376     // Set high water mark for last track in event
1377     fRunLoader->Stack()->SetHighWaterMark(nt);
1378 }
1379
1380 //_______________________________________________________________________
1381 void AliRun::KeepTrack(const Int_t track)
1382 {
1383   //
1384   // Delegate to stack
1385   //
1386     fRunLoader->Stack()->KeepTrack(track);
1387 }
1388
1389 //
1390 // MC Application
1391 //
1392
1393 //_______________________________________________________________________
1394 void  AliRun::ConstructGeometry()
1395 {
1396   //
1397   // Create modules, materials, geometry
1398   //
1399
1400     TStopwatch stw;
1401     TIter next(fModules);
1402     AliModule *detector;
1403     if (GetDebug()) Info("ConstructGeometry","Geometry creation:");
1404     while((detector = dynamic_cast<AliModule*>(next()))) {
1405       stw.Start();
1406       // Initialise detector materials and geometry
1407       detector->CreateMaterials();
1408       detector->CreateGeometry();
1409       printf("%10s R:%.2fs C:%.2fs\n",
1410              detector->GetName(),stw.RealTime(),stw.CpuTime());
1411     }
1412 }
1413
1414 //_______________________________________________________________________
1415 void  AliRun::InitGeometry()
1416 {
1417   //
1418   // Initialize detectors and display geometry
1419   //
1420
1421    printf("Initialisation:\n");
1422     TStopwatch stw;
1423     TIter next(fModules);
1424     AliModule *detector;
1425     while((detector = dynamic_cast<AliModule*>(next()))) {
1426       stw.Start();
1427       // Initialise detector and display geometry
1428       detector->Init();
1429       detector->BuildGeometry();
1430       printf("%10s R:%.2fs C:%.2fs\n",
1431              detector->GetName(),stw.RealTime(),stw.CpuTime());
1432     }
1433
1434 }
1435 //_______________________________________________________________________
1436
1437 void  AliRun::GeneratePrimaries()
1438 {
1439   //
1440   // Generate primary particles and fill them in the stack.
1441   //
1442
1443   Generator()->Generate();
1444 }
1445 //_______________________________________________________________________
1446
1447 void AliRun::BeginPrimary()
1448 {
1449   //
1450   // Called  at the beginning of each primary track
1451   //
1452
1453   // Reset Hits info
1454   gAlice->ResetHits();
1455   gAlice->ResetTrackReferences();
1456
1457 }
1458
1459 //_______________________________________________________________________
1460 void AliRun::PreTrack()
1461 {
1462      TObjArray &dets = *fModules;
1463      AliModule *module;
1464
1465      for(Int_t i=0; i<=fNdets; i++)
1466        if((module = dynamic_cast<AliModule*>(dets[i])))
1467          module->PreTrack();
1468
1469      fMCQA->PreTrack();
1470 }
1471
1472 //_______________________________________________________________________
1473 void AliRun::Stepping()
1474 {
1475   //
1476   // Called at every step during transport
1477   //
1478
1479   Int_t id = DetFromMate(gMC->GetMedium());
1480   if (id < 0) return;
1481
1482   //
1483   // --- If lego option, do it and leave
1484   if (fLego)
1485     fLego->StepManager();
1486   else {
1487     Int_t copy;
1488     //Update energy deposition tables
1489     AddEnergyDeposit(gMC->CurrentVolID(copy),gMC->Edep());
1490
1491     //Call the appropriate stepping routine;
1492     AliModule *det = dynamic_cast<AliModule*>(fModules->At(id));
1493     if(det && det->StepManagerIsEnabled()) {
1494       fMCQA->StepManager(id);
1495       det->StepManager();
1496     }
1497   }
1498 }
1499
1500 //_______________________________________________________________________
1501 void AliRun::PostTrack()
1502 {
1503      TObjArray &dets = *fModules;
1504      AliModule *module;
1505
1506      for(Int_t i=0; i<=fNdets; i++)
1507        if((module = dynamic_cast<AliModule*>(dets[i])))
1508          module->PostTrack();
1509 }
1510
1511 //_______________________________________________________________________
1512 void AliRun::FinishPrimary()
1513 {
1514   //
1515   // Called  at the end of each primary track
1516   //
1517
1518   //  static Int_t count=0;
1519   //  const Int_t times=10;
1520   // This primary is finished, purify stack
1521   fRunLoader->Stack()->PurifyKine();
1522
1523   TIter next(fModules);
1524   AliModule *detector;
1525   while((detector = dynamic_cast<AliModule*>(next()))) {
1526     detector->FinishPrimary();
1527     if(detector->GetLoader())
1528      {
1529        detector->GetLoader()->TreeH()->Fill();
1530      }
1531   }
1532
1533   // Write out track references if any
1534   if (fRunLoader->TreeTR())
1535    {
1536     fRunLoader->TreeTR()->Fill();
1537    }
1538 }
1539
1540 //_______________________________________________________________________
1541 void AliRun::FinishEvent()
1542 {
1543   //
1544   // Called at the end of the event.
1545   //
1546
1547   //
1548   if(fLego) fLego->FinishEvent();
1549
1550   TIter next(fModules);
1551   AliModule *detector;
1552   while((detector = dynamic_cast<AliModule*>(next()))) {
1553     detector->FinishEvent();
1554   }
1555
1556   //Update the energy deposit tables
1557   Int_t i;
1558   for(i=0;i<fEventEnergy.GetSize();i++)
1559    {
1560     fSummEnergy[i]+=fEventEnergy[i];
1561     fSum2Energy[i]+=fEventEnergy[i]*fEventEnergy[i];
1562    }
1563
1564   AliHeader* header = fRunLoader->GetHeader();
1565   AliStack* stack = fRunLoader->Stack();
1566   if ( (header == 0x0) || (stack == 0x0) )
1567    {//check if we got header and stack. If not cry and exit aliroot
1568     Fatal("AliRun","Can not get the stack or header from LOADER");
1569     return;//never reached
1570    }
1571   // Update Header information
1572   header->SetNprimary(stack->GetNprimary());
1573   header->SetNtrack(stack->GetNtrack());
1574
1575
1576   // Write out the kinematics
1577   stack->FinishEvent();
1578
1579   // Write out the event Header information
1580   TTree* treeE = fRunLoader->TreeE();
1581   if (treeE)
1582    {
1583       header->SetStack(stack);
1584       treeE->Fill();
1585    }
1586   else
1587    {
1588     Error("FinishEvent","Can not get TreeE from RL");
1589    }
1590
1591   fRunLoader->WriteKinematics("OVERWRITE");
1592   fRunLoader->WriteTrackRefs("OVERWRITE");
1593   fRunLoader->WriteHits("OVERWRITE");
1594
1595   if (GetDebug())
1596    {
1597      Info("FinishEvent","<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<");
1598      Info("FinishEvent","<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<");
1599      Info("FinishEvent","<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<");
1600      Info("FinishEvent","          FINISHING EVENT               ");
1601      Info("FinishEvent","<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<");
1602      Info("FinishEvent","<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<");
1603      Info("FinishEvent","<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<");
1604    }
1605 }
1606
1607 //_______________________________________________________________________
1608 void AliRun::Field(const Double_t* x, Double_t *b) const
1609 {
1610    Float_t xfloat[3];
1611    for (Int_t i=0; i<3; i++) xfloat[i] = x[i];
1612
1613    if (Field()) {
1614          Float_t bfloat[3];
1615          Field()->Field(xfloat,bfloat);
1616          for (Int_t j=0; j<3; j++) b[j] = bfloat[j];
1617    }
1618    else {
1619          printf("No mag field defined!\n");
1620          b[0]=b[1]=b[2]=0.;
1621    }
1622
1623 }
1624
1625 //
1626 // End of MC Application
1627 //
1628
1629 //_______________________________________________________________________
1630 void AliRun::Streamer(TBuffer &R__b)
1631 {
1632   // Stream an object of class AliRun.
1633
1634   if (R__b.IsReading()) {
1635     if (!gAlice) gAlice = this;
1636     AliRun::Class()->ReadBuffer(R__b, this);
1637     gROOT->GetListOfBrowsables()->Add(this,"Run");
1638
1639     gRandom = fRandom;
1640   } else {
1641     AliRun::Class()->WriteBuffer(R__b, this);
1642   }
1643 }
1644
1645
1646 //_______________________________________________________________________
1647 Int_t AliRun::GetCurrentTrackNumber() const {
1648   //
1649   // Returns current track
1650   //
1651   return fRunLoader->Stack()->GetCurrentTrackNumber();
1652 }
1653
1654 //_______________________________________________________________________
1655 Int_t AliRun::GetNtrack() const {
1656   //
1657   // Returns number of tracks in stack
1658   //
1659   return fRunLoader->Stack()->GetNtrack();
1660 }
1661 //_______________________________________________________________________
1662
1663 //_______________________________________________________________________
1664 TObjArray* AliRun::Particles() const {
1665   //
1666   // Returns pointer to Particles array
1667   //
1668   if (fRunLoader)
1669    if (fRunLoader->Stack())
1670     return fRunLoader->Stack()->Particles();
1671   return 0x0;
1672 }
1673
1674 //___________________________________________________________________________
1675
1676 //_______________________________________________________________________
1677 void AliRun::SetGenEventHeader(AliGenEventHeader* header)
1678 {
1679   fRunLoader->GetHeader()->SetGenEventHeader(header);
1680 }
1681
1682 //___________________________________________________________________________
1683
1684 Int_t AliRun::GetEvNumber() const
1685 {
1686 //Returns number of current event
1687   if (fRunLoader == 0x0)
1688    {
1689      Error("GetEvent","RunLoader is not set. Can not load data.");
1690      return -1;
1691    }
1692
1693   return fRunLoader->GetEventNumber();
1694 }
1695
1696 void AliRun::SetRunLoader(AliRunLoader* rloader)
1697 {
1698   fRunLoader = rloader;
1699   if (fRunLoader == 0x0) return;
1700
1701   TString evfoldname;
1702   TFolder* evfold = fRunLoader->GetEventFolder();
1703   if (evfold) evfoldname = evfold->GetName();
1704   else Warning("SetRunLoader","Did not get Event Folder from Run Loader");
1705
1706   if ( fRunLoader->GetAliRun() )
1707    {//if alrun already exists in folder
1708     if (fRunLoader->GetAliRun() != this )
1709      {//and is different than this - crash
1710        Fatal("AliRun","AliRun is already in Folder and it is not this object");
1711        return;//pro forma
1712      }//else do nothing
1713    }
1714   else
1715    {
1716      evfold->Add(this);//Post this AliRun to Folder
1717    }
1718
1719   TIter next(fModules);
1720   AliModule *module;
1721   while((module = (AliModule*)next()))
1722    {
1723      if (evfold) AliConfig::Instance()->Add(module,evfoldname);
1724      AliDetector* detector = dynamic_cast<AliDetector*>(module);
1725      if (detector)
1726       {
1727         AliLoader* loader = fRunLoader->GetLoader(detector);
1728         if (loader == 0x0)
1729          {
1730            Error("SetRunLoader","Can not get loader for detector %s",detector->GetName());
1731          }
1732         else
1733          {
1734            if (GetDebug()) Info("SetRunLoader","Setting loader for detector %s",detector->GetName());
1735            detector->SetLoader(loader);
1736          }
1737       }
1738    }
1739 }
1740
1741 void AliRun::AddModule(AliModule* mod)
1742 {
1743   if (mod == 0x0) return;
1744   if (strlen(mod->GetName()) == 0) return;
1745   if (GetModuleID(mod->GetName()) >= 0) return;
1746
1747   if (GetDebug()) Info("AddModule","%s",mod->GetName());
1748   if (fRunLoader == 0x0) AliConfig::Instance()->Add(mod);
1749   else AliConfig::Instance()->Add(mod,fRunLoader->GetEventFolder()->GetName());
1750
1751   Modules()->Add(mod);
1752 }