1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
18 ///////////////////////////////////////////////////////////////////////////////
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. //
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. //
32 <img src="picts/AliRunClass.gif">
37 <img src="picts/alirun.gif">
41 ///////////////////////////////////////////////////////////////////////////////
43 #include <Riostream.h>
53 #include <TGeometry.h>
56 #include <TObjectTable.h>
57 #include <TParticle.h>
63 #include <TVirtualMC.h>
65 #include "AliConfig.h"
66 #include "AliDetector.h"
67 #include "AliDisplay.h"
68 #include "AliGenEventHeader.h"
69 #include "AliGenerator.h"
70 #include "AliHeader.h"
73 #include "AliLegoGenerator.h"
74 #include "AliLoader.h"
77 #include "AliMagFCM.h"
78 #include "AliMagFDM.h"
81 #include "AliRunLoader.h"
88 //_______________________________________________________________________
108 fPDGDB(0), //Particle factory object
113 fConfigFunction("\0"),
120 // Default constructor for AliRun
122 AliConfig::Instance();//skowron 29 Feb 2002
123 //ensures that the folder structure is build
126 //_______________________________________________________________________
127 AliRun::AliRun(const AliRun& arun):
128 TVirtualMCApplication(arun),
147 fPDGDB(0), //Particle factory object
152 fConfigFunction("\0"),
159 // Copy constructor for AliRun
164 //_____________________________________________________________________________
165 AliRun::AliRun(const char *name, const char *title):
166 TVirtualMCApplication(name,title),
172 fModules(new TObjArray(77)), // Support list for the Detectors
178 fImedia(new TArrayI(1000)),
185 fPDGDB(TDatabasePDG::Instance()), //Particle factory object!
186 fHitLists(new TList()), // Create HitLists list
190 fConfigFunction("Config();"),
191 fRandom(new TRandom3()),
197 // Constructor for the main processor.
198 // Creates the geometry
199 // Creates the list of Detectors.
200 // Creates the list of particles.
205 // Set random number generator
208 if (gSystem->Getenv("CONFIG_SEED")) {
209 gRandom->SetSeed(static_cast<UInt_t>(atoi(gSystem->Getenv("CONFIG_SEED"))));
212 // Add to list of browsable
213 gROOT->GetListOfBrowsables()->Add(this,name);
214 // Create the TNode geometry for the event display
215 BuildSimpleGeometry();
217 // Create default mag field
220 // Prepare the tracking medium lists
221 for(Int_t i=0;i<1000;i++) (*fImedia)[i]=-99;
223 // Add particle list to configuration
224 AliConfig::Instance()->Add(fPDGDB);
226 // Set transport parameters
231 //_______________________________________________________________________
235 // Default AliRun destructor
237 gROOT->GetListOfBrowsables()->Remove(this);
241 TFolder* evfold = fRunLoader->GetEventFolder();
242 TFolder* modfold = dynamic_cast<TFolder*>(evfold->FindObjectAny(AliConfig::GetModulesFolderName()));
243 TIter next(fModules);
245 while((mod = (AliModule*)next()))
247 modfold->Remove(mod);
269 //_______________________________________________________________________
270 void AliRun::Copy(AliRun &) const
272 Fatal("Copy","Not implemented!\n");
275 //_______________________________________________________________________
276 void AliRun::AddHit(Int_t id, Int_t track, Int_t *vol, Float_t *hits) const
279 // Add a hit to detector id
281 TObjArray &dets = *fModules;
282 if(dets[id]) dynamic_cast<AliModule*>(dets[id])->AddHit(track,vol,hits);
285 //_______________________________________________________________________
286 void AliRun::AddDigit(Int_t id, Int_t *tracks, Int_t *digits) const
289 // Add digit to detector id
291 TObjArray &dets = *fModules;
292 if(dets[id]) dynamic_cast<AliModule*>(dets[id])->AddDigit(tracks,digits);
295 //_______________________________________________________________________
296 void AliRun::Browse(TBrowser *b)
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.
303 //detectors are in folders anyway
304 b->Add(fMCQA,"AliMCQA");
307 //_______________________________________________________________________
311 // Initialize Alice geometry
316 //_______________________________________________________________________
317 void AliRun::BuildSimpleGeometry()
320 // Create a simple TNode geometry used by Root display engine
322 // Initialise geometry
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");
331 //_______________________________________________________________________
332 void AliRun::CleanDetectors()
335 // Clean Detectors at the end of event
337 fRunLoader->CleanDetectors();
340 //_______________________________________________________________________
341 Int_t AliRun::DistancetoPrimitive(Int_t, Int_t) const
344 // Return the distance from the mouse to the AliRun object
350 //_______________________________________________________________________
351 void AliRun::DumpPart (Int_t i) const
354 // Dumps particle i in the stack
356 if (fRunLoader->Stack())
357 fRunLoader->Stack()->DumpPart(i);
360 //_______________________________________________________________________
361 void AliRun::DumpPStack () const
364 // Dumps the particle stack
366 if (fRunLoader->Stack())
367 fRunLoader->Stack()->DumpPStack();
370 //_______________________________________________________________________
371 void AliRun::SetField(AliMagF* magField)
373 // Set Magnetic Field Map
378 //_______________________________________________________________________
379 void AliRun::SetField(Int_t type, Int_t version, Float_t scale,
380 Float_t maxField, char* filename)
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
390 // --- Sanity check on mag field flags
391 if(fField) delete fField;
393 fField = new AliMagFC("Map1"," ",type,scale,maxField);
394 } else if(version<=2) {
395 fField = new AliMagFCM("Map2-3",filename,type,scale,maxField);
397 } else if(version==3) {
398 fField = new AliMagFDM("Map4",filename,type,scale,maxField);
401 Warning("SetField","Invalid map %d\n",version);
405 //_____________________________________________________________________________
407 void AliRun::InitLoaders()
409 //creates list of getters
410 if (GetDebug()) Info("InitLoaders","");
411 TIter next(fModules);
413 while((mod = (AliModule*)next()))
415 AliDetector *det = dynamic_cast<AliDetector*>(mod);
418 if (GetDebug()) Info("InitLoaders"," Adding %s ",det->GetName());
419 fRunLoader->AddLoader(det);
422 if (GetDebug()) Info("InitLoaders","Done");
424 //_____________________________________________________________________________
426 void AliRun::FinishRun()
429 // Called at the end of the run.
434 if (GetDebug()) Info("FinishRun"," Finish Lego");
435 fRunLoader->CdGAFile();
439 // Clean detector information
440 TIter next(fModules);
442 while((detector = dynamic_cast<AliModule*>(next()))) {
443 if (GetDebug()) Info("FinishRun"," %s->FinishRun()",detector->GetName());
444 detector->FinishRun();
447 //Output energy summary tables
448 if (GetDebug()) Info("FinishRun"," EnergySummary()");
451 if (GetDebug()) Info("FinishRun"," fRunLoader->WriteHeader(OVERWRITE)");
452 fRunLoader->WriteHeader("OVERWRITE");
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
459 // Clean tree information
460 if (GetDebug()) Info("FinishRun"," fRunLoader->Stack()->FinishRun()");
461 fRunLoader->Stack()->FinishRun();
463 // Clean detector information
464 if (GetDebug()) Info("FinishRun"," fGenerator->FinishRun()");
465 fGenerator->FinishRun();
467 fRunLoader->Synchronize();
470 //_______________________________________________________________________
471 void AliRun::FlagTrack(Int_t track)
475 fRunLoader->Stack()->FlagTrack(track);
478 //_______________________________________________________________________
479 void AliRun::EnergySummary()
482 // Print summary of deposited energy
488 Int_t kn, i, left, j, id;
489 const Float_t kzero=0;
490 Int_t ievent=fRunLoader->GetHeader()->GetEvent()+1;
492 // Energy loss information
494 printf("***************** Energy Loss Information per event (GEV) *****************\n");
495 for(kn=1;kn<fEventEnergy.GetSize();kn++) {
498 fEventEnergy[ndep]=kn;
503 ed2=100*TMath::Sqrt(TMath::Max(ed2-ed*ed,kzero))/ed;
506 fSummEnergy[ndep]=ed;
507 fSum2Energy[ndep]=TMath::Min(static_cast<Float_t>(99.),TMath::Max(ed2,kzero));
512 for(kn=0;kn<(ndep-1)/3+1;kn++) {
514 for(i=0;i<(3<left?3:left);i++) {
516 id=Int_t (fEventEnergy[j]+0.1);
517 printf(" %s %10.3f +- %10.3f%%;",gMC->VolName(id),fSummEnergy[j],fSum2Energy[j]);
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++) {
527 for(i=0;i<(5<left?5:left);i++) {
529 id=Int_t (fEventEnergy[j]+0.1);
530 printf(" %s %10.3f%%;",gMC->VolName(id),100*fSummEnergy[j]/edtot);
534 for(kn=0;kn<75;kn++) printf("*");
538 // Reset the TArray's
539 // fEventEnergy.Set(0);
540 // fSummEnergy.Set(0);
541 // fSum2Energy.Set(0);
544 //_______________________________________________________________________
545 void AliRun::Announce() const
548 // Announce the current version of AliRoot
551 "****************************************************************\n");
552 printf("%6s","*");printf("%64s","*\n");
555 printf(" You are running AliRoot version NewIO\n");
558 printf(" The cvs tag for the current program is $Name$\n");
560 printf("%6s","*");printf("%64s","*\n");
562 "****************************************************************\n");
565 //_______________________________________________________________________
566 AliModule *AliRun::GetModule(const char *name) const
569 // Return pointer to detector from name
571 return dynamic_cast<AliModule*>(fModules->FindObject(name));
574 //_______________________________________________________________________
575 AliDetector *AliRun::GetDetector(const char *name) const
578 // Return pointer to detector from name
580 return dynamic_cast<AliDetector*>(fModules->FindObject(name));
583 //_______________________________________________________________________
584 Int_t AliRun::GetModuleID(const char *name) const
587 // Return galice internal detector identifier from name
590 TObject *mod=fModules->FindObject(name);
591 if(mod) i=fModules->IndexOf(mod);
595 //_______________________________________________________________________
596 Int_t AliRun::GetEvent(Int_t event)
599 // Reloads data containers in folders # event
600 // Set branch addresses
602 if (fRunLoader == 0x0)
604 Error("GetEvent","RunLoader is not set. Can not load data.");
607 /*****************************************/
608 /**** P R E R E L O A D I N G ****/
609 /*****************************************/
610 // Reset existing structures
612 ResetTrackReferences();
616 /*****************************************/
617 /**** R E L O A D ****/
618 /*****************************************/
620 fRunLoader->GetEvent(event);
622 /*****************************************/
623 /**** P O S T R E L O A D I N G ****/
624 /*****************************************/
626 // Set Trees branch addresses
627 TIter next(fModules);
629 while((detector = dynamic_cast<AliModule*>(next())))
631 detector->SetTreeAddress();
634 return fRunLoader->GetHeader()->GetNtrack();
637 //_______________________________________________________________________
638 TGeometry *AliRun::GetGeometry()
641 // Import Alice geometry from current file
642 // Return pointer to geometry object
644 if (!fGeometry) fGeometry = dynamic_cast<TGeometry*>(gDirectory->Get("AliceGeom"));
646 // Unlink and relink nodes in detectors
647 // This is bad and there must be a better way...
650 TIter next(fModules);
652 while((detector = dynamic_cast<AliModule*>(next()))) {
653 TList *dnodes=detector->Nodes();
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);
666 //_______________________________________________________________________
667 Int_t AliRun::GetPrimary(Int_t track) const
670 // return number of primary that has generated track
672 return fRunLoader->Stack()->GetPrimary(track);
675 //_______________________________________________________________________
676 void AliRun::MediaTable()
679 // Built media table to get from the media number to
683 Int_t kz, nz, idt, lz, i, k, ind;
685 TObjArray &dets = *gAlice->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;
700 if(det->LoMedium() > det->HiMedium()) {
704 if(det->HiMedium() > fImedia->GetSize()) {
705 Error("MediaTable","Increase fImedia from %d to %d",
706 fImedia->GetSize(),det->HiMedium());
709 // Tag all materials in rage as belonging to detector kz
710 for(lz=det->LoMedium(); lz<= det->HiMedium(); lz++) {
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++) {
722 det=dynamic_cast<AliModule*>(dets[ind]);
724 printf(" %6s: %3d -> %3d;",det->GetName(),det->LoMedium(),
727 printf(" %6s: %3d -> %3d;","NULL",0,0);
733 //_______________________________________________________________________
734 void AliRun::SetGenerator(AliGenerator *generator)
737 // Load the event generator
739 if(!fGenerator) fGenerator = generator;
742 //_______________________________________________________________________
743 void AliRun::ResetGenerator(AliGenerator *generator)
746 // Load the event generator
750 Warning("ResetGenerator","Replacing generator %s with %s\n",
751 fGenerator->GetName(),generator->GetName());
753 Warning("ResetGenerator","Replacing generator %s with NULL\n",
754 fGenerator->GetName());
755 fGenerator = generator;
758 //_______________________________________________________________________
759 void AliRun::SetTransPar(const char *filename)
762 // Sets the file name for transport parameters
764 fTransParName = filename;
767 //_______________________________________________________________________
768 void AliRun::SetBaseFile(const char *filename)
770 fBaseFileName = filename;
773 //_______________________________________________________________________
774 void AliRun::ReadTransPar()
777 // Read filename to set the transport parameters
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"};
793 Int_t i, itmed, iret, ktmed, kz;
796 // See whether the file is there
797 filtmp=gSystem->ExpandPathName(fTransParName.Data());
798 lun=fopen(filtmp,"r");
801 Warning("ReadTransPar","File %s does not exist!\n",fTransParName.Data());
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","*");
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;
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);
824 printf(" *%59s\n","*");
825 printf(" "); for(i=0;i<60;i++) printf("*"); printf("\n");
829 // Read the end of line
832 if(line[0]=='*') continue;
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]);
841 Warning("ReadTransPar","Error reading file %s\n",fTransParName.Data());
844 // Check that the module exist
845 AliModule *mod = GetModule(detName);
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) {
853 Warning("ReadTransPar","Invalid tracking medium code %d for %s\n",itmed,mod->GetName());
856 // Set energy thresholds
857 for(kz=0;kz<kncuts;kz++) {
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]);
864 // Set transport mechanisms
865 for(kz=0;kz<knflags;kz++) {
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]));
873 Warning("ReadTransPar","Invalid medium code %d *\n",itmed);
877 if(fDebug) printf("%s::ReadTransParModule: %s not present\n",ClassName(),detName);
882 //_____________________________________________________________________________
884 void AliRun::BeginEvent()
887 // Clean-up previous event
891 Info("BeginEvent",">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>");
892 Info("BeginEvent",">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>");
893 Info("BeginEvent",">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>");
894 Info("BeginEvent"," BEGINNING EVENT ");
895 Info("BeginEvent",">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>");
896 Info("BeginEvent",">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>");
897 Info("BeginEvent",">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>");
900 /*******************************/
901 /* Clean after eventual */
903 /*******************************/
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);
910 fEventEnergy.Reset();
911 // Clean detector information
913 if (fRunLoader->Stack())
914 fRunLoader->Stack()->Reset();//clean stack -> tree is unloaded
916 fRunLoader->MakeStack();//or make a new one
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
926 // Reset all Detectors & kinematics & make/reset trees
929 fRunLoader->GetHeader()->Reset(fRun,fEvent,fEventNrInRun);
930 // fRunLoader->WriteKinematics("OVERWRITE"); is there any reason to rewrite here since MakeTree does so
932 if (GetDebug()) Info("BeginEvent"," fRunLoader->MakeTrackRefsContainer()");
933 fRunLoader->MakeTrackRefsContainer();//for insurance
935 if (GetDebug()) Info("BeginEvent"," ResetHits()");
937 if (GetDebug()) Info("BeginEvent"," fRunLoader->MakeTree(H)");
938 fRunLoader->MakeTree("H");
947 //create new branches and SetAdresses
948 TIter next(fModules);
950 while((detector = (AliModule*)next()))
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();
961 //_______________________________________________________________________
962 TParticle* AliRun::Particle(Int_t i) const
965 if (fRunLoader->Stack())
966 return fRunLoader->Stack()->Particle(i);
970 //_______________________________________________________________________
971 void AliRun::ResetDigits()
974 // Reset all Detectors digits
976 TIter next(fModules);
978 while((detector = dynamic_cast<AliModule*>(next()))) {
979 detector->ResetDigits();
983 //_______________________________________________________________________
984 void AliRun::ResetSDigits()
987 // Reset all Detectors digits
989 TIter next(fModules);
991 while((detector = dynamic_cast<AliModule*>(next()))) {
992 detector->ResetSDigits();
996 //_______________________________________________________________________
997 void AliRun::ResetHits()
1000 // Reset all Detectors hits
1002 TIter next(fModules);
1003 AliModule *detector;
1004 while((detector = dynamic_cast<AliModule*>(next()))) {
1005 detector->ResetHits();
1008 //_______________________________________________________________________
1010 void AliRun::ResetTrackReferences()
1013 // Reset all Detectors hits
1015 TIter next(fModules);
1016 AliModule *detector;
1017 while((detector = dynamic_cast<AliModule*>(next()))) {
1018 detector->ResetTrackReferences();
1021 //_______________________________________________________________________
1023 void AliRun::ResetPoints()
1026 // Reset all Detectors points
1028 TIter next(fModules);
1029 AliModule *detector;
1030 while((detector = dynamic_cast<AliModule*>(next()))) {
1031 detector->ResetPoints();
1034 //_______________________________________________________________________
1036 void AliRun::InitMC(const char *setup)
1039 // Initialize the Alice setup
1044 Warning("Init","Cannot initialise AliRun twice!\n");
1048 gROOT->LoadMacro(setup);
1049 gInterpreter->ProcessLine(fConfigFunction.Data());
1051 // Register MC in configuration
1052 AliConfig::Instance()->Add(gMC);
1056 fRunLoader->MakeTree("E");
1057 fRunLoader->LoadKinematics("RECREATE");
1058 fRunLoader->LoadTrackRefs("RECREATE");
1059 fRunLoader->LoadHits("all","RECREATE");
1062 fRunLoader->CdGAFile();
1064 gMC->DefineParticles(); //Create standard MC particles
1065 AliPDG::AddParticlesToPdgDataBase();
1067 TObject *objfirst, *objlast;
1069 fNdets = fModules->GetLast()+1;
1072 //=================Create Materials and geometry
1075 // Added also after in case of interactive initialisation of modules
1076 fNdets = fModules->GetLast()+1;
1078 TIter next(fModules);
1079 AliModule *detector;
1080 while((detector = dynamic_cast<AliModule*>(next())))
1082 objlast = gDirectory->GetList()->Last();
1084 // Add Detector histograms in Detector list of histograms
1085 if (objlast) objfirst = gDirectory->GetList()->After(objlast);
1086 else objfirst = gDirectory->GetList()->First();
1089 detector->Histograms()->Add(objfirst);
1090 objfirst = gDirectory->GetList()->After(objfirst);
1093 ReadTransPar(); //Read the cuts for all materials
1095 MediaTable(); //Build the special IMEDIA table
1097 //Initialise geometry deposition table
1098 fEventEnergy.Set(gMC->NofVolumes()+1);
1099 fSummEnergy.Set(gMC->NofVolumes()+1);
1100 fSum2Energy.Set(gMC->NofVolumes()+1);
1102 //Compute cross-sections
1103 gMC->BuildPhysics();
1105 //Write Geometry object to current file.
1106 fRunLoader->WriteGeometry();
1110 fMCQA = new AliMCQA(fNdets);
1113 // Save stuff at the beginning of the file to avoid file corruption
1115 fEventNrInRun = -1; //important - we start Begin event from increasing current number in run
1118 //_______________________________________________________________________
1120 void AliRun::RunMC(Int_t nevent, const char *setup)
1123 // Main function to be called to process a galice run
1125 // Root > gAlice.Run();
1126 // a positive number of events will cause the finish routine
1129 fEventsPerRun = nevent;
1130 // check if initialisation has been done
1131 if (!fInitDone) InitMC(setup);
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
1136 gMC->ProcessRun(nevent);
1138 // End of this run, close files
1139 if(nevent>0) FinishRun();
1142 //_______________________________________________________________________
1143 void AliRun::RunReco(const char *selected, Int_t first, Int_t last)
1146 // Main function to be called to reconstruct Alice event
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;
1153 for (Int_t nevent = nFirst; nevent <= nLast; nevent++) {
1154 if (GetDebug()) Info("RunReco","Processing event %d",nevent);
1156 Digits2Reco(selected);
1160 //_______________________________________________________________________
1162 void AliRun::Hits2Digits(const char *selected)
1165 // Convert Hits to sumable digits
1167 for (Int_t nevent=0; nevent<gAlice->TreeE()->GetEntries(); nevent++) {
1169 Hits2SDigits(selected);
1170 SDigits2Digits(selected);
1175 //_______________________________________________________________________
1177 void AliRun::Tree2Tree(Option_t *option, const char *selected)
1180 // Function to transform the content of
1182 // - TreeH to TreeS (option "S")
1183 // - TreeS to TreeD (option "D")
1184 // - TreeD to TreeR (option "R")
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).
1191 const char *oS = strstr(option,"S");
1192 const char *oD = strstr(option,"D");
1193 const char *oR = strstr(option,"R");
1195 TObjArray *detectors = Detectors();
1197 TIter next(detectors);
1199 AliDetector *detector = 0;
1201 while((detector = dynamic_cast<AliDetector*>(next()))) {
1203 if (strcmp(detector->GetName(),selected)) continue;
1204 if (detector->IsActive())
1207 AliLoader* loader = detector->GetLoader();
1208 if (loader == 0x0) continue;
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();
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();
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();
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)
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
1260 // ntheta bins in theta, eta
1261 // themin minimum angle in theta (degrees)
1262 // themax maximum angle in theta (degrees)
1264 // phimin minimum angle in phi (degrees)
1265 // phimax maximum angle in phi (degrees)
1266 // rmin minimum radius
1267 // rmax maximum radius
1270 // The number of events generated = ntheta*nphi
1271 // run input parameters in macro setup (default="Config.C")
1273 // Use macro "lego.C" to visualize the 3 lego plots in spherical coordinates
1276 <img src="picts/AliRunLego1.gif">
1281 <img src="picts/AliRunLego2.gif">
1286 <img src="picts/AliRunLego3.gif">
1291 // check if initialisation has been done
1292 if (!fInitDone) InitMC(setup);
1293 //Save current generator
1294 AliGenerator *gen=Generator();
1296 // Set new generator
1297 if (!gener) gener = new AliLegoGenerator();
1298 ResetGenerator(gener);
1300 // Configure Generator
1301 gener->SetRadiusRange(rmin, rmax);
1302 gener->SetZMax(zmax);
1303 gener->SetCoor1Range(nc1, c1min, c1max);
1304 gener->SetCoor2Range(nc2, c2min, c2max);
1307 //Create Lego object
1308 fLego = new AliLego("lego",gener);
1310 //Prepare MC for Lego Run
1315 //gMC->ProcessRun(nc1*nc2+1);
1316 gMC->ProcessRun(nc1*nc2);
1318 // Create only the Root event Tree
1319 fRunLoader->MakeTree("E");
1321 // End of this run, close files
1323 // Restore current generator
1324 ResetGenerator(gen);
1325 // Delete Lego Object
1326 delete fLego; fLego=0;
1329 //_______________________________________________________________________
1330 void AliRun::SetConfigFunction(const char * config)
1333 // Set the signature of the function contained in Config.C to configure
1336 fConfigFunction=config;
1339 //_______________________________________________________________________
1340 void AliRun::SetCurrentTrack(Int_t track)
1343 // Set current track number
1345 fRunLoader->Stack()->SetCurrentTrack(track);
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)
1353 // Delegate to stack
1355 fRunLoader->Stack()->PushTrack(done, parent, pdg, pmom, vpos, polar, tof,
1356 mech, ntr, weight, is);
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)
1366 // Delegate to stack
1368 fRunLoader->Stack()->PushTrack(done, parent, pdg, px, py, pz, e, vx, vy, vz, tof,
1369 polx, poly, polz, mech, ntr, weight, is);
1372 //_______________________________________________________________________
1373 void AliRun::SetHighWaterMark(const Int_t nt)
1376 // Set high water mark for last track in event
1377 fRunLoader->Stack()->SetHighWaterMark(nt);
1380 //_______________________________________________________________________
1381 void AliRun::KeepTrack(const Int_t track)
1384 // Delegate to stack
1386 fRunLoader->Stack()->KeepTrack(track);
1393 //_______________________________________________________________________
1394 void AliRun::ConstructGeometry()
1397 // Create modules, materials, geometry
1401 TIter next(fModules);
1402 AliModule *detector;
1403 if (GetDebug()) Info("ConstructGeometry","Geometry creation:");
1404 while((detector = dynamic_cast<AliModule*>(next()))) {
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());
1414 //_______________________________________________________________________
1415 void AliRun::InitGeometry()
1418 // Initialize detectors and display geometry
1421 printf("Initialisation:\n");
1423 TIter next(fModules);
1424 AliModule *detector;
1425 while((detector = dynamic_cast<AliModule*>(next()))) {
1427 // Initialise detector and display geometry
1429 detector->BuildGeometry();
1430 printf("%10s R:%.2fs C:%.2fs\n",
1431 detector->GetName(),stw.RealTime(),stw.CpuTime());
1435 //_______________________________________________________________________
1437 void AliRun::GeneratePrimaries()
1440 // Generate primary particles and fill them in the stack.
1443 Generator()->Generate();
1445 //_______________________________________________________________________
1447 void AliRun::BeginPrimary()
1450 // Called at the beginning of each primary track
1454 gAlice->ResetHits();
1455 gAlice->ResetTrackReferences();
1459 //_______________________________________________________________________
1460 void AliRun::PreTrack()
1462 TObjArray &dets = *fModules;
1465 for(Int_t i=0; i<=fNdets; i++)
1466 if((module = dynamic_cast<AliModule*>(dets[i])))
1472 //_______________________________________________________________________
1473 void AliRun::Stepping()
1476 // Called at every step during transport
1479 Int_t id = DetFromMate(gMC->GetMedium());
1483 // --- If lego option, do it and leave
1485 fLego->StepManager();
1488 //Update energy deposition tables
1489 AddEnergyDeposit(gMC->CurrentVolID(copy),gMC->Edep());
1491 //Call the appropriate stepping routine;
1492 AliModule *det = dynamic_cast<AliModule*>(fModules->At(id));
1493 if(det && det->StepManagerIsEnabled()) {
1494 fMCQA->StepManager(id);
1500 //_______________________________________________________________________
1501 void AliRun::PostTrack()
1503 TObjArray &dets = *fModules;
1506 for(Int_t i=0; i<=fNdets; i++)
1507 if((module = dynamic_cast<AliModule*>(dets[i])))
1508 module->PostTrack();
1511 //_______________________________________________________________________
1512 void AliRun::FinishPrimary()
1515 // Called at the end of each primary track
1518 // static Int_t count=0;
1519 // const Int_t times=10;
1520 // This primary is finished, purify stack
1521 fRunLoader->Stack()->PurifyKine();
1523 TIter next(fModules);
1524 AliModule *detector;
1525 while((detector = dynamic_cast<AliModule*>(next()))) {
1526 detector->FinishPrimary();
1527 if(detector->GetLoader())
1529 detector->GetLoader()->TreeH()->Fill();
1533 // Write out track references if any
1534 if (fRunLoader->TreeTR())
1536 fRunLoader->TreeTR()->Fill();
1540 //_______________________________________________________________________
1541 void AliRun::FinishEvent()
1544 // Called at the end of the event.
1548 if(fLego) fLego->FinishEvent();
1550 TIter next(fModules);
1551 AliModule *detector;
1552 while((detector = dynamic_cast<AliModule*>(next()))) {
1553 detector->FinishEvent();
1556 //Update the energy deposit tables
1558 for(i=0;i<fEventEnergy.GetSize();i++)
1560 fSummEnergy[i]+=fEventEnergy[i];
1561 fSum2Energy[i]+=fEventEnergy[i]*fEventEnergy[i];
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
1571 // Update Header information
1572 header->SetNprimary(stack->GetNprimary());
1573 header->SetNtrack(stack->GetNtrack());
1576 // Write out the kinematics
1577 stack->FinishEvent();
1579 // Write out the event Header information
1580 TTree* treeE = fRunLoader->TreeE();
1583 header->SetStack(stack);
1588 Error("FinishEvent","Can not get TreeE from RL");
1591 fRunLoader->WriteKinematics("OVERWRITE");
1592 fRunLoader->WriteTrackRefs("OVERWRITE");
1593 fRunLoader->WriteHits("OVERWRITE");
1597 Info("FinishEvent","<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<");
1598 Info("FinishEvent","<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<");
1599 Info("FinishEvent","<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<");
1600 Info("FinishEvent"," FINISHING EVENT ");
1601 Info("FinishEvent","<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<");
1602 Info("FinishEvent","<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<");
1603 Info("FinishEvent","<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<");
1607 //_______________________________________________________________________
1608 void AliRun::Field(const Double_t* x, Double_t *b) const
1611 for (Int_t i=0; i<3; i++) xfloat[i] = x[i];
1615 Field()->Field(xfloat,bfloat);
1616 for (Int_t j=0; j<3; j++) b[j] = bfloat[j];
1619 printf("No mag field defined!\n");
1626 // End of MC Application
1629 //_______________________________________________________________________
1630 void AliRun::Streamer(TBuffer &R__b)
1632 // Stream an object of class AliRun.
1634 if (R__b.IsReading()) {
1635 if (!gAlice) gAlice = this;
1636 AliRun::Class()->ReadBuffer(R__b, this);
1637 gROOT->GetListOfBrowsables()->Add(this,"Run");
1641 AliRun::Class()->WriteBuffer(R__b, this);
1646 //_______________________________________________________________________
1647 Int_t AliRun::GetCurrentTrackNumber() const {
1649 // Returns current track
1651 return fRunLoader->Stack()->GetCurrentTrackNumber();
1654 //_______________________________________________________________________
1655 Int_t AliRun::GetNtrack() const {
1657 // Returns number of tracks in stack
1659 return fRunLoader->Stack()->GetNtrack();
1661 //_______________________________________________________________________
1663 //_______________________________________________________________________
1664 TObjArray* AliRun::Particles() const {
1666 // Returns pointer to Particles array
1669 if (fRunLoader->Stack())
1670 return fRunLoader->Stack()->Particles();
1674 //___________________________________________________________________________
1676 //_______________________________________________________________________
1677 void AliRun::SetGenEventHeader(AliGenEventHeader* header)
1679 fRunLoader->GetHeader()->SetGenEventHeader(header);
1682 //___________________________________________________________________________
1684 Int_t AliRun::GetEvNumber() const
1686 //Returns number of current event
1687 if (fRunLoader == 0x0)
1689 Error("GetEvent","RunLoader is not set. Can not load data.");
1693 return fRunLoader->GetEventNumber();
1696 void AliRun::SetRunLoader(AliRunLoader* rloader)
1698 fRunLoader = rloader;
1699 if (fRunLoader == 0x0) return;
1702 TFolder* evfold = fRunLoader->GetEventFolder();
1703 if (evfold) evfoldname = evfold->GetName();
1704 else Warning("SetRunLoader","Did not get Event Folder from Run Loader");
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");
1716 evfold->Add(this);//Post this AliRun to Folder
1719 TIter next(fModules);
1721 while((module = (AliModule*)next()))
1723 if (evfold) AliConfig::Instance()->Add(module,evfoldname);
1724 AliDetector* detector = dynamic_cast<AliDetector*>(module);
1727 AliLoader* loader = fRunLoader->GetLoader(detector);
1730 Error("SetRunLoader","Can not get loader for detector %s",detector->GetName());
1734 if (GetDebug()) Info("SetRunLoader","Setting loader for detector %s",detector->GetName());
1735 detector->SetLoader(loader);
1741 void AliRun::AddModule(AliModule* mod)
1743 if (mod == 0x0) return;
1744 if (strlen(mod->GetName()) == 0) return;
1745 if (GetModuleID(mod->GetName()) >= 0) return;
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());
1751 Modules()->Add(mod);