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