X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=STEER%2FAliRun.cxx;h=9a88b014bd3f9d181b8afa9a5f9715eea62a2779;hb=f5bc14853de66388dae5300a43f945c59eece4ba;hp=9893e1eb8c30d96616be77b825cfca98db889c8d;hpb=ee1dd322f107680b3d9181c595a4e08a7521465f;p=u%2Fmrichter%2FAliRoot.git diff --git a/STEER/AliRun.cxx b/STEER/AliRun.cxx index 9893e1eb8c3..9a88b014bd3 100644 --- a/STEER/AliRun.cxx +++ b/STEER/AliRun.cxx @@ -13,28 +13,7 @@ * provided "as is" without express or implied warranty. * **************************************************************************/ -/* -$Log$ -Revision 1.25 2000/02/23 16:25:22 fca -AliVMC and AliGeant3 classes introduced -ReadEuclid moved from AliRun to AliModule - -Revision 1.24 2000/01/19 17:17:20 fca -Introducing a list of lists of hits -- more hits allowed for detector now - -Revision 1.23 1999/12/03 11:14:31 fca -Fixing previous wrong checking - -Revision 1.21 1999/11/25 10:40:08 fca -Fixing daughters information also in primary tracks - -Revision 1.20 1999/10/04 18:08:49 fca -Adding protection against inconsistent Euclid files - -Revision 1.19 1999/09/29 07:50:40 fca -Introduction of the Copyright and cvs Log - -*/ +/* $Id$ */ /////////////////////////////////////////////////////////////////////////////// // // @@ -61,68 +40,190 @@ Introduction of the Copyright and cvs Log // // /////////////////////////////////////////////////////////////////////////////// -#include -#include -#include -#include -#include -#include -#include - -#include "TParticle.h" -#include "AliRun.h" -#include "AliDisplay.h" -#include "AliVMC.h" - -#include #include +#include #include -AliRun *gAlice; +#include "Riostream.h" +#include "TBRIK.h" +#include "TBrowser.h" +#include "TCint.h" +#include "TFile.h" +#include "TFolder.h" +#include "TGeometry.h" +#include "TNode.h" +#include "TParticle.h" +#include "TRandom3.h" +#include "TROOT.h" +#include "TSystem.h" +#include "TTree.h" -static AliHeader *header; +#include "AliConfig.h" +#include "AliDetector.h" +#include "AliDisplay.h" +#include "AliGenerator.h" +#include "AliHeader.h" +#include "AliLego.h" +#include "AliLegoGenerator.h" +#include "AliMCQA.h" +#include "AliMagFC.h" +#include "AliMagFCM.h" +#include "AliMagFDM.h" +#include "AliPDG.h" +#include "AliRun.h" +#include "AliStack.h" -static TArrayF sEventEnergy; -static TArrayF sSummEnergy; -static TArrayF sSum2Energy; +AliRun *gAlice; ClassImp(AliRun) -//_____________________________________________________________________________ -AliRun::AliRun() +//_______________________________________________________________________ +AliRun::AliRun(): + fRun(0), + fEvent(0), + fEventNrInRun(0), + fEventsPerRun(0), + fDebug(0), + fHeader(0), + fTreeD(0), + fTreeS(0), + fTreeH(0), + fTreeTR(0), + fTreeE(0), + fTreeR(0), + fModules(0), + fGeometry(0), + fDisplay(0), + fTimer(), + fField(0), + fMC(0), + fImedia(0), + fNdets(0), + fTrRmax(1.e10), + fTrZmax(1.e10), + fGenerator(0), + fInitDone(kFALSE), + fLego(0), + fPDGDB(0), //Particle factory object + fHitLists(0), + fEventEnergy(0), + fSummEnergy(0), + fSum2Energy(0), + fConfigFunction("\0"), + fRandom(0), + fMCQA(0), + fTransParName("\0"), + fBaseFileName("\0"), + fStack(0), + fTreeDFileName(""), + fTreeDFile(0), + fTreeSFileName(""), + fTreeSFile(0), + fTreeRFileName(""), + fTreeRFile(0) { // // Default constructor for AliRun // - header=&fHeader; - fRun = 0; - fEvent = 0; - fCurrent = -1; - fModules = 0; - fGenerator = 0; - fTreeD = 0; - fTreeK = 0; - fTreeH = 0; - fTreeE = 0; - fTreeR = 0; - fParticles = 0; - fGeometry = 0; - fDisplay = 0; - fField = 0; - fMC = 0; - fNdets = 0; - fImedia = 0; - fTrRmax = 1.e10; - fTrZmax = 1.e10; - fInitDone = kFALSE; - fLego = 0; - fPDGDB = 0; //Particle factory object! - fHitLists = 0; +} + +//_______________________________________________________________________ +AliRun::AliRun(const AliRun& arun): + TVirtualMCApplication(arun), + fRun(0), + fEvent(0), + fEventNrInRun(0), + fEventsPerRun(0), + fDebug(0), + fHeader(0), + fTreeD(0), + fTreeS(0), + fTreeH(0), + fTreeTR(0), + fTreeE(0), + fTreeR(0), + fModules(0), + fGeometry(0), + fDisplay(0), + fTimer(), + fField(0), + fMC(0), + fImedia(0), + fNdets(0), + fTrRmax(1.e10), + fTrZmax(1.e10), + fGenerator(0), + fInitDone(kFALSE), + fLego(0), + fPDGDB(0), //Particle factory object + fHitLists(0), + fEventEnergy(0), + fSummEnergy(0), + fSum2Energy(0), + fConfigFunction("\0"), + fRandom(0), + fMCQA(0), + fTransParName("\0"), + fBaseFileName("\0"), + fStack(0), + fTreeDFileName(""), + fTreeDFile(0), + fTreeSFileName(""), + fTreeSFile(0), + fTreeRFileName(""), + fTreeRFile(0) +{ + // + // Copy constructor for AliRun + // + arun.Copy(*this); } //_____________________________________________________________________________ -AliRun::AliRun(const char *name, const char *title) - : TNamed(name,title) +AliRun::AliRun(const char *name, const char *title): + TVirtualMCApplication(name,title), + fRun(0), + fEvent(0), + fEventNrInRun(0), + fEventsPerRun(0), + fDebug(0), + fHeader(new AliHeader()), + fTreeD(0), + fTreeS(0), + fTreeH(0), + fTreeTR(0), + fTreeE(0), + fTreeR(0), + fModules(new TObjArray(77)), // Support list for the Detectors + fGeometry(0), + fDisplay(0), + fTimer(), + fField(0), + fMC(gMC), + fImedia(new TArrayI(1000)), + fNdets(0), + fTrRmax(1.e10), + fTrZmax(1.e10), + fGenerator(0), + fInitDone(kFALSE), + fLego(0), + fPDGDB(TDatabasePDG::Instance()), //Particle factory object! + fHitLists(new TList()), // Create HitLists list + fEventEnergy(0), + fSummEnergy(0), + fSum2Energy(0), + fConfigFunction("Config();"), + fRandom(new TRandom3()), + fMCQA(0), + fTransParName("\0"), + fBaseFileName("\0"), + fStack(new AliStack(10000)), //Particle stack + fTreeDFileName(""), + fTreeDFile(0), + fTreeSFileName(""), + fTreeSFile(0), + fTreeRFileName(""), + fTreeRFile(0) { // // Constructor for the main processor. @@ -130,109 +231,123 @@ AliRun::AliRun(const char *name, const char *title) // Creates the list of Detectors. // Creates the list of particles. // - Int_t i; - + gAlice = this; - fTreeD = 0; - fTreeK = 0; - fTreeH = 0; - fTreeE = 0; - fTreeR = 0; - fTrRmax = 1.e10; - fTrZmax = 1.e10; - fGenerator = 0; - fInitDone = kFALSE; - fLego = 0; - fField = 0; - + + // Set random number generator + gRandom = fRandom; + + if (gSystem->Getenv("CONFIG_SEED")) { + gRandom->SetSeed(static_cast(atoi(gSystem->Getenv("CONFIG_SEED")))); + } + + // Add to list of browsable gROOT->GetListOfBrowsables()->Add(this,name); - // - // create the support list for the various Detectors - fModules = new TObjArray(77); - // + // Create the TNode geometry for the event display - BuildSimpleGeometry(); - - fNtrack=0; - fHgwmk=0; - fCurrent=-1; - header=&fHeader; - fRun = 0; - fEvent = 0; - // - // Create the particle stack - fParticles = new TClonesArray("TParticle",100); - - fDisplay = 0; - // // Create default mag field SetField(); - // - fMC = gMC; - // + // Prepare the tracking medium lists - fImedia = new TArrayI(1000); - for(i=0;i<1000;i++) (*fImedia)[i]=-99; - // - // Make particles - fPDGDB = TDatabasePDG::Instance(); //Particle factory object! - // - // Create HitLists list - fHitLists = new TList(); + for(Int_t i=0;i<1000;i++) (*fImedia)[i]=-99; + + // Add particle list to configuration + AliConfig::Instance()->Add(fPDGDB); + + // Set transport parameters + SetTransPar(); } -//_____________________________________________________________________________ + +//_______________________________________________________________________ AliRun::~AliRun() { // - // Defaullt AliRun destructor + // Default AliRun destructor // + TFile *curfil =0; + if(fTreeE)curfil=fTreeE->GetCurrentFile(); delete fImedia; delete fField; - delete fMC; + // delete fMC; + delete gMC; gMC=0; delete fGeometry; delete fDisplay; delete fGenerator; delete fLego; delete fTreeD; - delete fTreeK; delete fTreeH; + delete fTreeTR; delete fTreeE; delete fTreeR; + delete fTreeS; if (fModules) { fModules->Delete(); delete fModules; } - if (fParticles) { - fParticles->Delete(); - delete fParticles; - } + delete fStack; delete fHitLists; + delete fPDGDB; + delete fMCQA; + delete fHeader; + // avoid to delete TFile objects not owned by this object + // avoid multiple deletions + if(curfil == fTreeDFile) fTreeDFile=0; + if(curfil == fTreeSFile) fTreeSFile=0; + if(curfil == fTreeRFile) fTreeRFile=0; + if(fTreeSFile == fTreeDFile) fTreeSFile=0; + if(fTreeRFile == fTreeDFile) fTreeRFile=0; + if(fTreeRFile == fTreeSFile) fTreeRFile=0; + if(fTreeDFile){ + if(fTreeDFile->IsOpen())fTreeDFile->Close(); + delete fTreeDFile; + } + if(fTreeSFile){ + if(fTreeSFile->IsOpen())fTreeSFile->Close(); + delete fTreeSFile; + } + if(fTreeRFile){ + if(fTreeRFile->IsOpen())fTreeRFile->Close(); + delete fTreeRFile; + } + if (gROOT->GetListOfBrowsables()) + gROOT->GetListOfBrowsables()->Remove(this); + + gAlice=0; } -//_____________________________________________________________________________ +//_______________________________________________________________________ +void AliRun::Copy(AliRun &) const +{ + // + // Copy method ... not implemented + // + Fatal("Copy","Not implemented!\n"); +} + +//_______________________________________________________________________ void AliRun::AddHit(Int_t id, Int_t track, Int_t *vol, Float_t *hits) const { // // Add a hit to detector id // TObjArray &dets = *fModules; - if(dets[id]) ((AliModule*) dets[id])->AddHit(track,vol,hits); + if(dets[id]) dynamic_cast(dets[id])->AddHit(track,vol,hits); } -//_____________________________________________________________________________ +//_______________________________________________________________________ void AliRun::AddDigit(Int_t id, Int_t *tracks, Int_t *digits) const { // // Add digit to detector id // TObjArray &dets = *fModules; - if(dets[id]) ((AliModule*) dets[id])->AddDigit(tracks,digits); + if(dets[id]) dynamic_cast(dets[id])->AddDigit(tracks,digits); } -//_____________________________________________________________________________ +//_______________________________________________________________________ void AliRun::Browse(TBrowser *b) { // @@ -240,20 +355,26 @@ void AliRun::Browse(TBrowser *b) // of the Root browser. // It displays the Root Trees and all detectors. // - if (fTreeK) b->Add(fTreeK,fTreeK->GetName()); + if(!fStack) fStack=fHeader->Stack(); + TTree* pTreeK = fStack->TreeK(); + + if (pTreeK) b->Add(pTreeK,pTreeK->GetName()); if (fTreeH) b->Add(fTreeH,fTreeH->GetName()); + if (fTreeTR) b->Add(fTreeTR,fTreeH->GetName()); if (fTreeD) b->Add(fTreeD,fTreeD->GetName()); if (fTreeE) b->Add(fTreeE,fTreeE->GetName()); if (fTreeR) b->Add(fTreeR,fTreeR->GetName()); + if (fTreeS) b->Add(fTreeS,fTreeS->GetName()); TIter next(fModules); AliModule *detector; - while((detector = (AliModule*)next())) { + while((detector = dynamic_cast(next()))) { b->Add(detector,detector->GetName()); } + b->Add(fMCQA,"AliMCQA"); } -//_____________________________________________________________________________ +//_______________________________________________________________________ void AliRun::Build() { // @@ -262,7 +383,7 @@ void AliRun::Build() // } -//_____________________________________________________________________________ +//_______________________________________________________________________ void AliRun::BuildSimpleGeometry() { // @@ -277,7 +398,7 @@ void AliRun::BuildSimpleGeometry() new TNode("alice","alice","S_alice"); } -//_____________________________________________________________________________ +//_______________________________________________________________________ void AliRun::CleanDetectors() { // @@ -285,32 +406,13 @@ void AliRun::CleanDetectors() // TIter next(fModules); AliModule *detector; - while((detector = (AliModule*)next())) { + while((detector = dynamic_cast(next()))) { detector->FinishEvent(); } } -//_____________________________________________________________________________ -void AliRun::CleanParents() -{ - // - // Clean Particles stack. - // Set parent/daughter relations - // - TClonesArray &particles = *(gAlice->Particles()); - TParticle *part; - int i; - for(i=0; iTestBit(Daughters_Bit)) { - part->SetFirstDaughter(-1); - part->SetLastDaughter(-1); - } - } -} - -//_____________________________________________________________________________ -Int_t AliRun::DistancetoPrimitive(Int_t, Int_t) +//_______________________________________________________________________ +Int_t AliRun::DistancetoPrimitive(Int_t, Int_t) const { // // Return the distance from the mouse to the AliRun object @@ -319,35 +421,33 @@ Int_t AliRun::DistancetoPrimitive(Int_t, Int_t) return 9999; } -//_____________________________________________________________________________ -void AliRun::DumpPart (Int_t i) +//_______________________________________________________________________ +void AliRun::DumpPart (Int_t i) const { // // Dumps particle i in the stack // - TClonesArray &particles = *fParticles; - ((TParticle*) particles[i])->Print(); + fStack->DumpPart(i); } -//_____________________________________________________________________________ -void AliRun::DumpPStack () +//_______________________________________________________________________ +void AliRun::DumpPStack () const { // // Dumps the particle stack // - TClonesArray &particles = *fParticles; - printf( - "\n\n=======================================================================\n"); - for (Int_t i=0;i %d ",i); ((TParticle*) particles[i])->Print(); - printf("--------------------------------------------------------------\n"); - } - printf( - "\n=======================================================================\n\n"); + fStack->DumpPStack(); } -//_____________________________________________________________________________ +//_______________________________________________________________________ +void AliRun::SetField(AliMagF* magField) +{ + // Set Magnetic Field Map + fField = magField; + fField->ReadField(); +} + +//_______________________________________________________________________ void AliRun::SetField(Int_t type, Int_t version, Float_t scale, Float_t maxField, char* filename) { @@ -360,182 +460,86 @@ void AliRun::SetField(Int_t type, Int_t version, Float_t scale, // // --- Sanity check on mag field flags - if(type<0 || type > 2) { - Warning("SetField", - "Invalid magnetic field flag: %5d; Helix tracking chosen instead\n" - ,type); - type=2; - } if(fField) delete fField; if(version==1) { - fField = new AliMagFC("Map1"," ",type,version,scale,maxField); - } else if(version<=3) { - fField = new AliMagFCM("Map2-3",filename,type,version,scale,maxField); + fField = new AliMagFC("Map1"," ",type,scale,maxField); + } else if(version<=2) { + fField = new AliMagFCM("Map2-3",filename,type,scale,maxField); + fField->ReadField(); + } else if(version==3) { + fField = new AliMagFDM("Map4",filename,type,scale,maxField); fField->ReadField(); } else { Warning("SetField","Invalid map %d\n",version); } } - -//_____________________________________________________________________________ -void AliRun::FillTree() -{ - // - // Fills all AliRun TTrees - // - if (fTreeK) fTreeK->Fill(); - if (fTreeH) fTreeH->Fill(); - if (fTreeD) fTreeD->Fill(); - if (fTreeR) fTreeR->Fill(); -} -//_____________________________________________________________________________ -void AliRun::FinishPrimary() -{ - // - // Called at the end of each primary track - // - - // static Int_t count=0; - // const Int_t times=10; - // This primary is finished, purify stack - gAlice->PurifyKine(); - - // Write out hits if any - if (gAlice->TreeH()) { - gAlice->TreeH()->Fill(); - } - - // Reset Hits info - gAlice->ResetHits(); - - // - // if(++count%times==1) gObjectTable->Print(); -} - -//_____________________________________________________________________________ -void AliRun::FinishEvent() -{ - // - // Called at the end of the event. - // - - //Update the energy deposit tables - Int_t i; - for(i=0;iFill(); - } - - // Write out the digits - if (fTreeD) { - fTreeD->Fill(); - ResetDigits(); - } - - // Write out reconstructed clusters - if (fTreeR) { - fTreeR->Fill(); - } - - // Write out the event Header information - if (fTreeE) fTreeE->Fill(); - - // Reset stack info - ResetStack(); - - // Write Tree headers - // Int_t ievent = fHeader.GetEvent(); - // char hname[30]; - // sprintf(hname,"TreeK%d",ievent); - if (fTreeK) fTreeK->Write(); - // sprintf(hname,"TreeH%d",ievent); - if (fTreeH) fTreeH->Write(); - // sprintf(hname,"TreeD%d",ievent); - if (fTreeD) fTreeD->Write(); - // sprintf(hname,"TreeR%d",ievent); - if (fTreeR) fTreeR->Write(); -} - -//_____________________________________________________________________________ +//_______________________________________________________________________ void AliRun::FinishRun() { // // Called at the end of the run. // + // + if(fLego) fLego->FinishRun(); + // Clean detector information TIter next(fModules); AliModule *detector; - while((detector = (AliModule*)next())) { + while((detector = dynamic_cast(next()))) { detector->FinishRun(); } //Output energy summary tables EnergySummary(); + + TFile *file = fTreeE->GetCurrentFile(); + + file->cd(); - // file is retrieved from whatever tree - TFile *File = 0; - if (fTreeK) File = fTreeK->GetCurrentFile(); - if ((!File) && (fTreeH)) File = fTreeH->GetCurrentFile(); - if ((!File) && (fTreeD)) File = fTreeD->GetCurrentFile(); - if ((!File) && (fTreeE)) File = fTreeE->GetCurrentFile(); - if( NULL==File ) { - Error("FinishRun","There isn't root file!"); - exit(1); - } - File->cd(); - fTreeE->Write(); + fTreeE->Write(0,TObject::kOverwrite); + // Write AliRun info and all detectors parameters + Write(0,TObject::kOverwrite); + // Clean tree information - delete fTreeK; fTreeK = 0; - delete fTreeH; fTreeH = 0; - delete fTreeD; fTreeD = 0; - delete fTreeR; fTreeR = 0; - delete fTreeE; fTreeE = 0; + + fStack->FinishRun(); - // Write AliRun info and all detectors parameters - Write(); + if (fTreeH) { + delete fTreeH; fTreeH = 0; + } + if (fTreeTR) { + delete fTreeTR; fTreeTR = 0; + } + if (fTreeD) { + delete fTreeD; fTreeD = 0; + } + if (fTreeR) { + delete fTreeR; fTreeR = 0; + } +// if (fTreeE) { +// delete fTreeE; fTreeE = 0; +// } + if (fTreeS) { + delete fTreeS; fTreeS = 0; + } + fGenerator->FinishRun(); // Close output file - File->Write(); + file->Write(); } -//_____________________________________________________________________________ +//_______________________________________________________________________ void AliRun::FlagTrack(Int_t track) { + // Delegate to stack // - // Flags a track and all its family tree to be kept - // - int curr; - TParticle *particle; - - curr=track; - while(1) { - particle=(TParticle*)fParticles->UncheckedAt(curr); - - // If the particle is flagged the three from here upward is saved already - if(particle->TestBit(Keep_Bit)) return; - - // Save this particle - particle->SetBit(Keep_Bit); - - // Move to father if any - if((curr=particle->GetFirstMother())==-1) return; - } + fStack->FlagTrack(track); } -//_____________________________________________________________________________ +//_______________________________________________________________________ void AliRun::EnergySummary() { // @@ -546,25 +550,25 @@ void AliRun::EnergySummary() Float_t edtot=0; Float_t ed, ed2; Int_t kn, i, left, j, id; - const Float_t zero=0; - Int_t ievent=fHeader.GetEvent()+1; + const Float_t kzero=0; + Int_t ievent=fHeader->GetEvent()+1; // // Energy loss information if(ievent) { printf("***************** Energy Loss Information per event (GEV) *****************\n"); - for(kn=1;kn0) { - sEventEnergy[ndep]=kn; + fEventEnergy[ndep]=kn; if(ievent>1) { ed=ed/ievent; - ed2=sSum2Energy[kn]; + ed2=fSum2Energy[kn]; ed2=ed2/ievent; - ed2=100*TMath::Sqrt(TMath::Max(ed2-ed*ed,zero))/ed; + ed2=100*TMath::Sqrt(TMath::Max(ed2-ed*ed,kzero))/ed; } else ed2=99; - sSummEnergy[ndep]=ed; - sSum2Energy[ndep]=TMath::Min((Float_t) 99.,TMath::Max(ed2,zero)); + fSummEnergy[ndep]=ed; + fSum2Energy[ndep]=TMath::Min(static_cast(99.),TMath::Max(ed2,kzero)); edtot+=ed; ndep++; } @@ -573,8 +577,8 @@ void AliRun::EnergySummary() left=ndep-kn*3; for(i=0;i<(3VolName(id),sSummEnergy[j],sSum2Energy[j]); + id=Int_t (fEventEnergy[j]+0.1); + printf(" %s %10.3f +- %10.3f%%;",gMC->VolName(id),fSummEnergy[j],fSum2Energy[j]); } printf("\n"); } @@ -586,8 +590,8 @@ void AliRun::EnergySummary() left=ndep-kn*5; for(i=0;i<(5VolName(id),100*sSummEnergy[j]/edtot); + id=Int_t (fEventEnergy[j]+0.1); + printf(" %s %10.3f%%;",gMC->VolName(id),100*fSummEnergy[j]/edtot); } printf("\n"); } @@ -596,31 +600,52 @@ void AliRun::EnergySummary() } // // Reset the TArray's - // sEventEnergy.Set(0); - // sSummEnergy.Set(0); - // sSum2Energy.Set(0); + // fEventEnergy.Set(0); + // fSummEnergy.Set(0); + // fSum2Energy.Set(0); } -//_____________________________________________________________________________ -AliModule *AliRun::GetModule(const char *name) +//_______________________________________________________________________ +void AliRun::Announce() const +{ + // + // Announce the current version of AliRoot + // + printf("%70s", + "****************************************************************\n"); + printf("%6s","*");printf("%64s","*\n"); + + printf("%6s","*"); + printf(" You are running AliRoot version v3-09-07\n"); + + printf("%6s","*"); + printf(" The cvs tag for the current program is $Name$\n"); + + printf("%6s","*");printf("%64s","*\n"); + printf("%70s", + "****************************************************************\n"); +} + +//_______________________________________________________________________ +AliModule *AliRun::GetModule(const char *name) const { // // Return pointer to detector from name // - return (AliModule*)fModules->FindObject(name); + return dynamic_cast(fModules->FindObject(name)); } -//_____________________________________________________________________________ -AliDetector *AliRun::GetDetector(const char *name) +//_______________________________________________________________________ +AliDetector *AliRun::GetDetector(const char *name) const { // // Return pointer to detector from name // - return (AliDetector*)fModules->FindObject(name); + return dynamic_cast(fModules->FindObject(name)); } -//_____________________________________________________________________________ -Int_t AliRun::GetModuleID(const char *name) +//_______________________________________________________________________ +Int_t AliRun::GetModuleID(const char *name) const { // // Return galice internal detector identifier from name @@ -631,7 +656,7 @@ Int_t AliRun::GetModuleID(const char *name) return i; } -//_____________________________________________________________________________ +//_______________________________________________________________________ Int_t AliRun::GetEvent(Int_t event) { // @@ -640,69 +665,156 @@ Int_t AliRun::GetEvent(Int_t event) // // Reset existing structures - ResetStack(); ResetHits(); + ResetTrackReferences(); ResetDigits(); + ResetSDigits(); // Delete Trees already connected - if (fTreeK) delete fTreeK; - if (fTreeH) delete fTreeH; - if (fTreeD) delete fTreeD; - if (fTreeR) delete fTreeR; - - // Get header from file - if(fTreeE) fTreeE->GetEntry(event); - else Error("GetEvent","Cannot file Header Tree\n"); + if (fTreeH) { delete fTreeH; fTreeH = 0;} + if (fTreeTR) { delete fTreeTR; fTreeTR = 0;} + if (fTreeD) { delete fTreeD; fTreeD = 0;} + if (fTreeR) { delete fTreeR; fTreeR = 0;} + if (fTreeS) { delete fTreeS; fTreeS = 0;} + + // Create the particle stack + if (fHeader) delete fHeader; + fHeader = 0; - // Get Kine Tree from file + // Get header from file + if(fTreeE) { + fTreeE->SetBranchAddress("Header", &fHeader); + + if (!fTreeE->GetEntry(event)) { + Error("GetEvent","Cannot find event:%d\n",event); + return -1; + } + } + else { + Error("GetEvent","Cannot find Header Tree (TE)\n"); + return -1; + } + + // Get the stack from the header, set fStack to 0 if it + // fails to get event + // + TFile *file = fTreeE->GetCurrentFile(); char treeName[20]; - sprintf(treeName,"TreeK%d",event); - fTreeK = (TTree*)gDirectory->Get(treeName); - if (fTreeK) fTreeK->SetBranchAddress("Particles", &fParticles); - else Error("GetEvent","cannot find Kine Tree for event:%d\n",event); + file->cd(); + + if (fStack) delete fStack; + fStack = fHeader->Stack(); + if (fStack) { + if (!fStack->GetEvent(event)) fStack = 0; + } + // Get Hits Tree header from file sprintf(treeName,"TreeH%d",event); - fTreeH = (TTree*)gDirectory->Get(treeName); + fTreeH = dynamic_cast(gDirectory->Get(treeName)); if (!fTreeH) { - Error("GetEvent","cannot find Hits Tree for event:%d\n",event); + Warning("GetEvent","cannot find Hits Tree for event:%d\n",event); } - + + // Get TracReferences Tree header from file + sprintf(treeName,"TreeTR%d",event); + fTreeTR = dynamic_cast(gDirectory->Get(treeName)); + if (!fTreeTR) { + Warning("GetEvent","cannot find TrackReferences Tree for event:%d\n",event); + } + + // get current file name and compare with names containing trees S,D,R + TString curfilname=static_cast(fTreeE->GetCurrentFile()->GetName()); + if(fTreeDFileName==curfilname)fTreeDFileName=""; + if(fTreeSFileName==curfilname)fTreeSFileName=""; + if(fTreeRFileName==curfilname)fTreeRFileName=""; + // Get Digits Tree header from file sprintf(treeName,"TreeD%d",event); - fTreeD = (TTree*)gDirectory->Get(treeName); + + if (!fTreeDFile && fTreeDFileName != "") { + InitTreeFile("D",fTreeDFileName); + } + if (fTreeDFile) { + fTreeD = dynamic_cast(fTreeDFile->Get(treeName)); + } else { + fTreeD = dynamic_cast(file->Get(treeName)); + } if (!fTreeD) { - Warning("GetEvent","cannot find Digits Tree for event:%d\n",event); + // Warning("GetEvent","cannot find Digits Tree for event:%d\n",event); } - + if(fTreeDFileName != ""){ + if(fTreeDFileName==fTreeSFileName) { + fTreeSFileName = ""; + fTreeSFile = fTreeDFile; + } + if(fTreeDFileName==fTreeRFileName) { + fTreeRFileName = ""; + fTreeRFile = fTreeDFile; + } + } + + file->cd(); + + // Get SDigits Tree header from file + sprintf(treeName,"TreeS%d",event); + if (!fTreeSFile && fTreeSFileName != "") { + InitTreeFile("S",fTreeSFileName); + } + if (fTreeSFile) { + fTreeS = dynamic_cast(fTreeSFile->Get(treeName)); + } else { + fTreeS = dynamic_cast(gDirectory->Get(treeName)); + } + if (!fTreeS) { + // Warning("GetEvent","cannot find SDigits Tree for event:%d\n",event); + } + + if(fTreeSFileName != ""){ + if(fTreeSFileName==fTreeRFileName){ + fTreeRFileName = ""; + fTreeRFile = fTreeSFile; + } + } + + file->cd(); // Get Reconstruct Tree header from file sprintf(treeName,"TreeR%d",event); - fTreeR = (TTree*)gDirectory->Get(treeName); + if (!fTreeRFile && fTreeRFileName != "") { + InitTreeFile("R",fTreeRFileName); + } + if(fTreeRFile) { + fTreeR = dynamic_cast(fTreeRFile->Get(treeName)); + } else { + fTreeR = dynamic_cast(gDirectory->Get(treeName)); + } if (!fTreeR) { // printf("WARNING: cannot find Reconstructed Tree for event:%d\n",event); } + + file->cd(); // Set Trees branch addresses TIter next(fModules); AliModule *detector; - while((detector = (AliModule*)next())) { + while((detector = dynamic_cast(next()))) { detector->SetTreeAddress(); } - - if (fTreeK) fTreeK->GetEvent(0); - fNtrack = Int_t (fParticles->GetEntries()); - return fNtrack; + + fEvent=event; //MI change + + return fHeader->GetNtrack(); } -//_____________________________________________________________________________ +//_______________________________________________________________________ TGeometry *AliRun::GetGeometry() { // // Import Alice geometry from current file // Return pointer to geometry object // - if (!fGeometry) fGeometry = (TGeometry*)gDirectory->Get("AliceGeom"); + if (!fGeometry) fGeometry = dynamic_cast(gDirectory->Get("AliceGeom")); // // Unlink and relink nodes in detectors // This is bad and there must be a better way... @@ -710,13 +822,12 @@ TGeometry *AliRun::GetGeometry() TIter next(fModules); AliModule *detector; - while((detector = (AliModule*)next())) { - detector->SetTreeAddress(); + while((detector = dynamic_cast(next()))) { TList *dnodes=detector->Nodes(); Int_t j; TNode *node, *node1; for ( j=0; jGetSize(); j++) { - node = (TNode*) dnodes->At(j); + node = dynamic_cast(dnodes->At(j)); node1 = fGeometry->GetNode(node->GetName()); dnodes->Remove(node); dnodes->AddAt(node1,j); @@ -725,147 +836,31 @@ TGeometry *AliRun::GetGeometry() return fGeometry; } -//_____________________________________________________________________________ -void AliRun::GetNextTrack(Int_t &mtrack, Int_t &ipart, Float_t *pmom, - Float_t &e, Float_t *vpos, Float_t *polar, - Float_t &tof) -{ - // - // Return next track from stack of particles - // - TVector3 pol; - fCurrent=-1; - TParticle *track; - for(Int_t i=fNtrack-1; i>=0; i--) { - track=(TParticle*) fParticles->UncheckedAt(i); - if(!track->TestBit(Done_Bit)) { - // - // The track has not yet been processed - fCurrent=i; - ipart=track->GetPdgCode(); - pmom[0]=track->Px(); - pmom[1]=track->Py(); - pmom[2]=track->Pz(); - e =track->Energy(); - vpos[0]=track->Vx(); - vpos[1]=track->Vy(); - vpos[2]=track->Vz(); - track->GetPolarisation(pol); - polar[0]=pol.X(); - polar[1]=pol.Y(); - polar[2]=pol.Z(); - tof=track->T(); - track->SetBit(Done_Bit); - break; - } - } - mtrack=fCurrent; - // - // stop and start timer when we start a primary track - Int_t nprimaries = fHeader.GetNprimary(); - if (fCurrent >= nprimaries) return; - if (fCurrent < nprimaries-1) { - fTimer.Stop(); - track=(TParticle*) fParticles->UncheckedAt(fCurrent+1); - // track->SetProcessTime(fTimer.CpuTime()); - } - fTimer.Start(); -} - -//_____________________________________________________________________________ -Int_t AliRun::GetPrimary(Int_t track) +//_______________________________________________________________________ +Int_t AliRun::GetPrimary(Int_t track) const { // // return number of primary that has generated track // - int current, parent; - TParticle *part; - // - parent=track; - while (1) { - current=parent; - part = (TParticle *)fParticles->UncheckedAt(current); - parent=part->GetFirstMother(); - if(parent<0) return current; - } + return fStack->GetPrimary(track); } -//_____________________________________________________________________________ -void AliRun::Init(const char *setup) +//_______________________________________________________________________ +void AliRun::MediaTable() { // - // Initialize the Alice setup + // Built media table to get from the media number to + // the detector id // - - gROOT->LoadMacro(setup); - gInterpreter->ProcessLine("Config();"); - - gMC->DefineParticles(); //Create standard MC particles - - TObject *objfirst, *objlast; - - fNdets = fModules->GetLast()+1; - - // - //=================Create Materials, geometry, histograms, etc - TIter next(fModules); - AliModule *detector; - while((detector = (AliModule*)next())) { - detector->SetTreeAddress(); - objlast = gDirectory->GetList()->Last(); - - // Initialise detector materials, geometry, histograms,etc - detector->CreateMaterials(); - detector->CreateGeometry(); - detector->BuildGeometry(); - detector->Init(); - - // Add Detector histograms in Detector list of histograms - if (objlast) objfirst = gDirectory->GetList()->After(objlast); - else objfirst = gDirectory->GetList()->First(); - while (objfirst) { - detector->Histograms()->Add(objfirst); - objfirst = gDirectory->GetList()->After(objfirst); - } - } - SetTransPar(); //Read the cuts for all materials - - MediaTable(); //Build the special IMEDIA table - - //Terminate building of geometry - printf("%p\n",gVMC); - gVMC->FinishGeometry(); - - //Initialise geometry deposition table - sEventEnergy.Set(gMC->NofVolumes()+1); - sSummEnergy.Set(gMC->NofVolumes()+1); - sSum2Energy.Set(gMC->NofVolumes()+1); - - //Compute cross-sections - gVMC->BuildPhysics(); - - //Write Geometry object to current file. - fGeometry->Write(); - - fInitDone = kTRUE; -} - -//_____________________________________________________________________________ -void AliRun::MediaTable() -{ - // - // Built media table to get from the media number to - // the detector id - // - Int_t kz, nz, idt, lz, i, k, ind; - // Int_t ibeg; - TObjArray &dets = *gAlice->Detectors(); - AliModule *det; + Int_t kz, nz, idt, lz, i, k, ind; + // Int_t ibeg; + TObjArray &dets = *gAlice->Detectors(); + AliModule *det; // // For all detectors for (kz=0;kz(dets[kz]))) { TArrayI &idtmed = *(det->GetIdtmed()); for(nz=0;nz<100;nz++) { // Find max and min material number @@ -896,7 +891,7 @@ void AliRun::MediaTable() for(i=0;i<(fNdets-1)/6+1;i++) { for(k=0;k< (6(dets[ind]); if(det) printf(" %6s: %3d -> %3d;",det->GetName(),det->LoMedium(), det->HiMedium()); @@ -907,7 +902,7 @@ void AliRun::MediaTable() } } -//____________________________________________________________________________ +//_______________________________________________________________________ void AliRun::SetGenerator(AliGenerator *generator) { // @@ -916,59 +911,80 @@ void AliRun::SetGenerator(AliGenerator *generator) if(!fGenerator) fGenerator = generator; } -//____________________________________________________________________________ +//_______________________________________________________________________ void AliRun::ResetGenerator(AliGenerator *generator) { // // Load the event generator // if(fGenerator) - Warning("ResetGenerator","Replacing generator %s with %s\n", - fGenerator->GetName(),generator->GetName()); + if(generator) + Warning("ResetGenerator","Replacing generator %s with %s\n", + fGenerator->GetName(),generator->GetName()); + else + Warning("ResetGenerator","Replacing generator %s with NULL\n", + fGenerator->GetName()); fGenerator = generator; } -//____________________________________________________________________________ -void AliRun::SetTransPar(char* filename) +//_______________________________________________________________________ +void AliRun::SetTransPar(const char *filename) +{ + // + // Sets the file name for transport parameters + // + fTransParName = filename; +} + +//_______________________________________________________________________ +void AliRun::SetBaseFile(const char *filename) +{ + fBaseFileName = filename; +} + +//_______________________________________________________________________ +void AliRun::ReadTransPar() { // // Read filename to set the transport parameters // - const Int_t ncuts=10; - const Int_t nflags=11; - const Int_t npars=ncuts+nflags; - const char pars[npars][7] = {"CUTGAM" ,"CUTELE","CUTNEU","CUTHAD","CUTMUO", + const Int_t kncuts=10; + const Int_t knflags=11; + const Int_t knpars=kncuts+knflags; + const char kpars[knpars][7] = {"CUTGAM" ,"CUTELE","CUTNEU","CUTHAD","CUTMUO", "BCUTE","BCUTM","DCUTE","DCUTM","PPCUTM","ANNI", "BREM","COMP","DCAY","DRAY","HADR","LOSS", "MULS","PAIR","PHOT","RAYL"}; char line[256]; char detName[7]; char* filtmp; - Float_t cut[ncuts]; - Int_t flag[nflags]; + Float_t cut[kncuts]; + Int_t flag[knflags]; Int_t i, itmed, iret, ktmed, kz; FILE *lun; // // See whether the file is there - filtmp=gSystem->ExpandPathName(filename); + filtmp=gSystem->ExpandPathName(fTransParName.Data()); lun=fopen(filtmp,"r"); delete [] filtmp; if(!lun) { - Warning("SetTransPar","File %s does not exist!\n",filename); + Warning("ReadTransPar","File %s does not exist!\n",fTransParName.Data()); return; } // - printf(" "); for(i=0;i<60;i++) printf("*"); printf("\n"); - printf(" *%59s\n","*"); - printf(" * Please check carefully what you are doing!%10s\n","*"); - printf(" *%59s\n","*"); + if(fDebug) { + printf(" "); for(i=0;i<60;i++) printf("*"); printf("\n"); + printf(" *%59s\n","*"); + printf(" * Please check carefully what you are doing!%10s\n","*"); + printf(" *%59s\n","*"); + } // while(1) { // Initialise cuts and flags - for(i=0;iGetName()); + Warning("ReadTransPar","Invalid tracking medium code %d for %s\n",itmed,mod->GetName()); continue; } // Set energy thresholds - for(kz=0;kz=0) { - printf(" * %-6s set to %10.3E for tracking medium code %4d for %s\n", - pars[kz],cut[kz],itmed,mod->GetName()); - gMC->Gstpar(ktmed,pars[kz],cut[kz]); + if(fDebug) printf(" * %-6s set to %10.3E for tracking medium code %4d for %s\n", + kpars[kz],cut[kz],itmed,mod->GetName()); + gMC->Gstpar(ktmed,kpars[kz],cut[kz]); } } // Set transport mechanisms - for(kz=0;kz=0) { - printf(" * %-6s set to %10d for tracking medium code %4d for %s\n", - pars[ncuts+kz],flag[kz],itmed,mod->GetName()); - gMC->Gstpar(ktmed,pars[ncuts+kz],Float_t(flag[kz])); + if(fDebug) printf(" * %-6s set to %10d for tracking medium code %4d for %s\n", + kpars[kncuts+kz],flag[kz],itmed,mod->GetName()); + gMC->Gstpar(ktmed,kpars[kncuts+kz],Float_t(flag[kz])); } } } else { - Warning("SetTransPar","Invalid medium code %d *\n",itmed); + Warning("ReadTransPar","Invalid medium code %d *\n",itmed); continue; } } else { - Warning("SetTransPar","Module %s not present\n",detName); + if(fDebug) printf("%s::ReadTransParModule: %s not present\n",ClassName(),detName); continue; } } } -//_____________________________________________________________________________ -void AliRun::MakeTree(Option_t *option) + +//_______________________________________________________________________ +void AliRun::MakeTree(Option_t *option, const char *file) { // // Create the ROOT trees @@ -1045,220 +1064,152 @@ void AliRun::MakeTree(Option_t *option) char hname[30]; // // Analyse options - char *K = strstr(option,"K"); - char *H = strstr(option,"H"); - char *E = strstr(option,"E"); - char *D = strstr(option,"D"); - char *R = strstr(option,"R"); - // - if (K && !fTreeK) { - sprintf(hname,"TreeK%d",fEvent); - fTreeK = new TTree(hname,"Kinematics"); - // Create a branch for particles - fTreeK->Branch("Particles",&fParticles,4000); - } - if (H && !fTreeH) { + const char *oK = strstr(option,"K"); + const char *oH = strstr(option,"H"); + const char *oTR = strstr(option,"T"); + const char *oE = strstr(option,"E"); + const char *oD = strstr(option,"D"); + const char *oR = strstr(option,"R"); + const char *oS = strstr(option,"S"); + // + + TDirectory *cwd = gDirectory; + + TBranch *branch = 0; + + if (oK) fStack->MakeTree(fEvent, file); + + if (oE && !fTreeE) { + fTreeE = new TTree("TE","Header"); + // branch = fTreeE->Branch("Header", "AliHeader", &fHeader, 4000, 0); + branch = fTreeE->Branch("Header", "AliHeader", &fHeader, 4000, 0); + branch->SetAutoDelete(kFALSE); + TFolder *folder = dynamic_cast(gROOT->FindObjectAny("/Folders/RunMC/Event/Header")); + if (folder) folder->Add(fHeader); +// branch = fTreeE->Branch("Stack","AliStack", &fStack, 4000, 0); +// branch->SetAutoDelete(kFALSE); +// if (folder) folder->Add(fStack); + fTreeE->Write(0,TObject::kOverwrite); + } + + if (file && branch) { + char * outFile = new char[strlen(gAlice->GetBaseFile())+strlen(file)+2]; + sprintf(outFile,"%s/%s",GetBaseFile(),file); + branch->SetFile(outFile); + TIter next( branch->GetListOfBranches()); + while ((branch=dynamic_cast(next()))) { + branch->SetFile(outFile); + } + if (GetDebug()>1) + printf("* MakeBranch * Diverting Branch %s to file %s\n", branch->GetName(),file); + cwd->cd(); + delete outFile; + } + + if (oH && !fTreeH) { sprintf(hname,"TreeH%d",fEvent); fTreeH = new TTree(hname,"Hits"); fTreeH->SetAutoSave(1000000000); //no autosave + fTreeH->Write(0,TObject::kOverwrite); + } + + if (oTR && !fTreeTR) { + sprintf(hname,"TreeTR%d",fEvent); + fTreeTR = new TTree(hname,"TrackReferences"); + fTreeTR->SetAutoSave(1000000000); //no autosave + fTreeTR->Write(0,TObject::kOverwrite); } - if (D && !fTreeD) { + + if (oD && !fTreeD) { sprintf(hname,"TreeD%d",fEvent); fTreeD = new TTree(hname,"Digits"); + fTreeD->Write(0,TObject::kOverwrite); + } + if (oS && !fTreeS) { + sprintf(hname,"TreeS%d",fEvent); + fTreeS = new TTree(hname,"SDigits"); + fTreeS->Write(0,TObject::kOverwrite); } - if (R && !fTreeR) { + if (oR && !fTreeR) { sprintf(hname,"TreeR%d",fEvent); fTreeR = new TTree(hname,"Reconstruction"); + fTreeR->Write(0,TObject::kOverwrite); } - if (E && !fTreeE) { - fTreeE = new TTree("TE","Header"); - // Create a branch for Header - fTreeE->Branch("Header","AliHeader",&header,4000); - } + // // Create a branch for hits/digits for each detector // Each branch is a TClonesArray. Each data member of the Hits classes // will be in turn a subbranch of the detector master branch TIter next(fModules); AliModule *detector; - while((detector = (AliModule*)next())) { - if (H || D || R) detector->MakeBranch(option); + while((detector = dynamic_cast(next()))) { + if (oH) detector->MakeBranch(option,file); + if (oTR) detector->MakeBranchTR(option,file); } } -//_____________________________________________________________________________ -Int_t AliRun::PurifyKine(Int_t lastSavedTrack, Int_t nofTracks) +//_______________________________________________________________________ +TParticle* AliRun::Particle(Int_t i) const { // - // PurifyKine with external parameters + // Returns particle i on the simulation stack // - fHgwmk = lastSavedTrack; - fNtrack = nofTracks; - PurifyKine(); - return fHgwmk; + return fStack->Particle(i); } -//_____________________________________________________________________________ -void AliRun::PurifyKine() +//_______________________________________________________________________ +void AliRun::ResetDigits() { // - // Compress kinematic tree keeping only flagged particles - // and renaming the particle id's in all the hits + // Reset all Detectors digits // - TClonesArray &particles = *fParticles; - int nkeep=fHgwmk+1, parent, i; - TParticle *part, *partnew, *father; - int *map = new int[particles.GetEntries()]; - - // Save in Header total number of tracks before compression - fHeader.SetNtrack(fHeader.GetNtrack()+fNtrack-fHgwmk); - - // Preset map, to be removed later - for(i=0; iTestBit(Keep_Bit)) { - - // This particle has to be kept - map[i]=nkeep; - if(i!=nkeep) { - - // Old and new are different, have to copy - partnew = (TParticle *)particles.UncheckedAt(nkeep); - *partnew = *part; - } else partnew = part; - - // as the parent is always *before*, it must be already - // in place. This is what we are checking anyway! - if((parent=partnew->GetFirstMother())>fHgwmk) { - if(map[parent]==-99) printf("map[%d] = -99!\n",parent); - partnew->SetFirstMother(map[parent]); - } - nkeep++; - } - } - fNtrack=nkeep; - - // Fix daughters information - for (i=0; iGetFirstMother(); - if(parent>=0) { - father = (TParticle *)particles.UncheckedAt(parent); - if(father->TestBit(Daughters_Bit)) { - - if(iGetFirstDaughter()) father->SetFirstDaughter(i); - if(i>father->GetLastDaughter()) father->SetLastDaughter(i); - } else { - // Iitialise daughters info for first pass - father->SetFirstDaughter(i); - father->SetLastDaughter(i); - father->SetBit(Daughters_Bit); - } - } - } - -#ifdef old - // Now loop on all detectors and reset the hits - AliHit *OneHit; TIter next(fModules); AliModule *detector; - while((detector = (AliModule*)next())) { - if (!detector->Hits()) continue; - TClonesArray &vHits=*(detector->Hits()); - if(vHits.GetEntries() != detector->GetNhits()) - printf("vHits.GetEntries()!=detector->GetNhits(): %d != %d\n", - vHits.GetEntries(),detector->GetNhits()); - for (i=0; iGetNhits(); i++) { - OneHit = (AliHit *)vHits.UncheckedAt(i); - OneHit->SetTrack(map[OneHit->GetTrack()]); - } - } -#else - - // Now loop on all registered hit lists - TIter next(fHitLists); - TCollection *hitList; - while((hitList = (TCollection*)next())) { - TIter nexthit(hitList); - AliHit *hit; - while((hit = (AliHit*)nexthit())) { - hit->SetTrack(map[hit->GetTrack()]); - } + while((detector = dynamic_cast(next()))) { + detector->ResetDigits(); } -#endif - - fHgwmk=nkeep-1; - particles.SetLast(fHgwmk); - delete [] map; } -//_____________________________________________________________________________ -void AliRun::Reset(Int_t run, Int_t idevent) +//_______________________________________________________________________ +void AliRun::ResetSDigits() { // - // Reset all Detectors & kinematics & trees - // - char hname[30]; + // Reset all Detectors digits // - ResetStack(); - ResetHits(); - ResetDigits(); - - // Initialise event header - fHeader.Reset(run,idevent); - - if(fTreeK) { - fTreeK->Reset(); - sprintf(hname,"TreeK%d",idevent); - fTreeK->SetName(hname); - } - if(fTreeH) { - fTreeH->Reset(); - sprintf(hname,"TreeH%d",idevent); - fTreeH->SetName(hname); - } - if(fTreeD) { - fTreeD->Reset(); - sprintf(hname,"TreeD%d",idevent); - fTreeD->SetName(hname); - } - if(fTreeR) { - fTreeR->Reset(); - sprintf(hname,"TreeR%d",idevent); - fTreeR->SetName(hname); + TIter next(fModules); + AliModule *detector; + while((detector = dynamic_cast(next()))) { + detector->ResetSDigits(); } } -//_____________________________________________________________________________ -void AliRun::ResetDigits() +//_______________________________________________________________________ +void AliRun::ResetHits() { // - // Reset all Detectors digits + // Reset all Detectors hits // TIter next(fModules); AliModule *detector; - while((detector = (AliModule*)next())) { - detector->ResetDigits(); + while((detector = dynamic_cast(next()))) { + detector->ResetHits(); } } -//_____________________________________________________________________________ -void AliRun::ResetHits() +//_______________________________________________________________________ +void AliRun::ResetTrackReferences() { // // Reset all Detectors hits // TIter next(fModules); AliModule *detector; - while((detector = (AliModule*)next())) { - detector->ResetHits(); + while((detector = dynamic_cast(next()))) { + detector->ResetTrackReferences(); } } -//_____________________________________________________________________________ +//_______________________________________________________________________ void AliRun::ResetPoints() { // @@ -1266,13 +1217,87 @@ void AliRun::ResetPoints() // TIter next(fModules); AliModule *detector; - while((detector = (AliModule*)next())) { + while((detector = dynamic_cast(next()))) { detector->ResetPoints(); } } -//_____________________________________________________________________________ -void AliRun::Run(Int_t nevent, const char *setup) +//_______________________________________________________________________ +void AliRun::InitMC(const char *setup) +{ + // + // Initialize the Alice setup + // + + Announce(); + + if(fInitDone) { + Warning("Init","Cannot initialise AliRun twice!\n"); + return; + } + + gROOT->LoadMacro(setup); + gInterpreter->ProcessLine(fConfigFunction.Data()); + + // Register MC in configuration + AliConfig::Instance()->Add(gMC); + gMC->SetStack(fStack); + + gMC->DefineParticles(); //Create standard MC particles + AliPDG::AddParticlesToPdgDataBase(); + + TObject *objfirst, *objlast; + + fNdets = fModules->GetLast()+1; + + // + //=================Create Materials and geometry + gMC->Init(); + + // Added also after in case of interactive initialisation of modules + fNdets = fModules->GetLast()+1; + + TIter next(fModules); + AliModule *detector; + while((detector = dynamic_cast(next()))) { + detector->SetTreeAddress(); + objlast = gDirectory->GetList()->Last(); + + // Add Detector histograms in Detector list of histograms + if (objlast) objfirst = gDirectory->GetList()->After(objlast); + else objfirst = gDirectory->GetList()->First(); + while (objfirst) { + detector->Histograms()->Add(objfirst); + objfirst = gDirectory->GetList()->After(objfirst); + } + } + ReadTransPar(); //Read the cuts for all materials + + MediaTable(); //Build the special IMEDIA table + + //Initialise geometry deposition table + fEventEnergy.Set(gMC->NofVolumes()+1); + fSummEnergy.Set(gMC->NofVolumes()+1); + fSum2Energy.Set(gMC->NofVolumes()+1); + + //Compute cross-sections + gMC->BuildPhysics(); + + //Write Geometry object to current file. + fGeometry->Write(); + + fInitDone = kTRUE; + + fMCQA = new AliMCQA(fNdets); + + AliConfig::Instance(); + // + // Save stuff at the beginning of the file to avoid file corruption + Write(); +} + +//_______________________________________________________________________ +void AliRun::RunMC(Int_t nevent, const char *setup) { // // Main function to be called to process a galice run @@ -1281,31 +1306,140 @@ void AliRun::Run(Int_t nevent, const char *setup) // a positive number of events will cause the finish routine // to be called // - - Int_t i, todo; + fEventsPerRun = nevent; // check if initialisation has been done - if (!fInitDone) Init(setup); + if (!fInitDone) InitMC(setup); // Create the Root Tree with one branch per detector - gAlice->MakeTree("KHDER"); - todo = TMath::Abs(nevent); - for (i=0; iReset(fRun, fEvent); - gVMC->ProcessEvent(); - gAlice->FinishEvent(); - fEvent++; + MakeTree("ESDRT"); + + if (gSystem->Getenv("CONFIG_SPLIT_FILE")) { + MakeTree("K","Kine.root"); + MakeTree("H","Hits.root"); + } else { + MakeTree("KH"); } - + + gMC->ProcessRun(nevent); + // End of this run, close files - if(nevent>0) gAlice->FinishRun(); + if(nevent>0) FinishRun(); } -//_____________________________________________________________________________ -void AliRun::RunLego(const char *setup,Int_t ntheta,Float_t themin, - Float_t themax,Int_t nphi,Float_t phimin,Float_t phimax, - Float_t rmin,Float_t rmax,Float_t zmax) +//_______________________________________________________________________ +void AliRun::RunReco(const char *selected, Int_t first, Int_t last) +{ + // + // Main function to be called to reconstruct Alice event + // + cout << "Found "<< gAlice->TreeE()->GetEntries() << "events" << endl; + Int_t nFirst = first; + Int_t nLast = (last < 0)? static_cast(gAlice->TreeE()->GetEntries()) : last; + + for (Int_t nevent = nFirst; nevent <= nLast; nevent++) { + cout << "Processing event "<< nevent << endl; + GetEvent(nevent); + // MakeTree("R"); + Digits2Reco(selected); + } +} + +//_______________________________________________________________________ + +void AliRun::Hits2Digits(const char *selected) +{ + + // Convert Hits to sumable digits + // + for (Int_t nevent=0; neventTreeE()->GetEntries(); nevent++) { + GetEvent(nevent); + // MakeTree("D"); + Hits2SDigits(selected); + SDigits2Digits(selected); + } +} + + +//_______________________________________________________________________ + +void AliRun::Tree2Tree(Option_t *option, const char *selected) +{ + // + // Function to transform the content of + // + // - TreeH to TreeS (option "S") + // - TreeS to TreeD (option "D") + // - TreeD to TreeR (option "R") + // + // If multiple options are specified ("SDR"), transformation will be done in sequence for + // selected detector and for all detectors if none is selected (detector string + // can contain blank separated list of detector names). + + + const char *oS = strstr(option,"S"); + const char *oD = strstr(option,"D"); + const char *oR = strstr(option,"R"); + + TObjArray *detectors = Detectors(); + + TIter next(detectors); + + AliDetector *detector = 0; + + TDirectory *cwd = gDirectory; + + TObject *obj; + + char outFile[32]; + + while((obj = next())) { + if (!dynamic_cast(obj)) + Fatal("Tree2Tree","Wrong type in fModules array\n"); + if (!(detector = dynamic_cast(obj))) continue; + if (selected) + if (strcmp(detector->GetName(),selected)) continue; + if (!detector->IsActive()) continue; + if (gSystem->Getenv("CONFIG_SPLIT_FILE")) { + if (oS) { + sprintf(outFile,"SDigits.%s.root",detector->GetName()); + detector->MakeBranch("S",outFile); + } + if (oD) { + sprintf(outFile,"Digits.%s.root",detector->GetName()); + detector->MakeBranch("D",outFile); + } + if (oR) { + sprintf(outFile,"Reco.%s.root",detector->GetName()); + detector->MakeBranch("R",outFile); + } + } else { + detector->MakeBranch(option); + } + + cwd->cd(); + + if (oS) { + cout << "Hits2SDigits: Processing " << detector->GetName() << "..." << endl; + detector->Hits2SDigits(); + } + if (oD) { + cout << "SDigits2Digits: Processing " << detector->GetName() << "..." << endl; + detector->SDigits2Digits(); + } + if (oR) { + cout << "Digits2Reco: Processing " << detector->GetName() << "..." << endl; + detector->Digits2Reco(); + } + + cwd->cd(); + } +} + +//_______________________________________________________________________ +void AliRun::RunLego(const char *setup, Int_t nc1, Float_t c1min, + Float_t c1max,Int_t nc2,Float_t c2min,Float_t c2max, + Float_t rmin,Float_t rmax,Float_t zmax, AliLegoGenerator* gener) { // // Generates lego plots of: @@ -1346,113 +1480,269 @@ void AliRun::RunLego(const char *setup,Int_t ntheta,Float_t themin, // // check if initialisation has been done - if (!fInitDone) Init(setup); - + if (!fInitDone) InitMC(setup); + //Save current generator + AliGenerator *gen=Generator(); + + // Set new generator + if (!gener) gener = new AliLegoGenerator(); + ResetGenerator(gener); + // + // Configure Generator + gener->SetRadiusRange(rmin, rmax); + gener->SetZMax(zmax); + gener->SetCoor1Range(nc1, c1min, c1max); + gener->SetCoor2Range(nc2, c2min, c2max); + + //Create Lego object - fLego = new AliLego("lego",ntheta,themin,themax,nphi,phimin,phimax,rmin,rmax,zmax); + fLego = new AliLego("lego",gener); + //Prepare MC for Lego Run + gMC->InitLego(); + //Run Lego Object - fLego->Run(); + + //gMC->ProcessRun(nc1*nc2+1); + gMC->ProcessRun(nc1*nc2); // Create only the Root event Tree - gAlice->MakeTree("E"); + MakeTree("E"); // End of this run, close files - gAlice->FinishRun(); + FinishRun(); + // Restore current generator + ResetGenerator(gen); + // Delete Lego Object + delete fLego; fLego=0; } -//_____________________________________________________________________________ +//_______________________________________________________________________ +void AliRun::SetConfigFunction(const char * config) +{ + // + // Set the signature of the function contained in Config.C to configure + // the run + // + fConfigFunction=config; +} + +//_______________________________________________________________________ void AliRun::SetCurrentTrack(Int_t track) { // // Set current track number // - fCurrent = track; + fStack->SetCurrentTrack(track); } -//_____________________________________________________________________________ +//_______________________________________________________________________ void AliRun::SetTrack(Int_t done, Int_t parent, Int_t pdg, Float_t *pmom, Float_t *vpos, Float_t *polar, Float_t tof, - const char *mecha, Int_t &ntr, Float_t weight) + TMCProcess mech, Int_t &ntr, Float_t weight, Int_t is) { +// Delegate to stack +// + + fStack->SetTrack(done, parent, pdg, pmom, vpos, polar, tof, + mech, ntr, weight, is); +} + +//_______________________________________________________________________ +void AliRun::SetTrack(Int_t done, Int_t parent, Int_t pdg, + Double_t px, Double_t py, Double_t pz, Double_t e, + Double_t vx, Double_t vy, Double_t vz, Double_t tof, + Double_t polx, Double_t poly, Double_t polz, + TMCProcess mech, Int_t &ntr, Float_t weight, Int_t is) +{ + // Delegate to stack // - // Load a track on the stack - // - // done 0 if the track has to be transported - // 1 if not - // parent identifier of the parent track. -1 for a primary - // pdg particle code - // pmom momentum GeV/c - // vpos position - // polar polarisation - // tof time of flight in seconds - // mecha production mechanism - // ntr on output the number of the track stored - // - TClonesArray &particles = *fParticles; - TParticle *particle; - Float_t mass; - const Int_t firstdaughter=-1; - const Int_t lastdaughter=-1; - const Int_t KS=0; - // const Float_t tlife=0; - - // - // Here we get the static mass - // For MC is ok, but a more sophisticated method could be necessary - // if the calculated mass is required - // also, this method is potentially dangerous if the mass - // used in the MC is not the same of the PDG database - // - mass = TDatabasePDG::Instance()->GetParticle(pdg)->Mass(); - Float_t e=TMath::Sqrt(mass*mass+pmom[0]*pmom[0]+ - pmom[1]*pmom[1]+pmom[2]*pmom[2]); - - //printf("Loading particle %s mass %f ene %f No %d ip %d pos %f %f %f mom %f %f %f KS %d m %s\n", - //pname,mass,e,fNtrack,pdg,vpos[0],vpos[1],vpos[2],pmom[0],pmom[1],pmom[2],KS,mecha); - - particle=new(particles[fNtrack]) TParticle(pdg,KS,parent,-1,firstdaughter, - lastdaughter,pmom[0],pmom[1],pmom[2], - e,vpos[0],vpos[1],vpos[2],tof); - // polar[0],polar[1],polar[2],tof, - // mecha,weight); - ((TParticle*)particles[fNtrack])->SetPolarisation(TVector3(polar[0],polar[1],polar[2])); - ((TParticle*)particles[fNtrack])->SetWeight(weight); - if(!done) particle->SetBit(Done_Bit); - - if(parent>=0) { - particle=(TParticle*) fParticles->UncheckedAt(parent); - particle->SetLastDaughter(fNtrack); - if(particle->GetFirstDaughter()<0) particle->SetFirstDaughter(fNtrack); - } else { - // - // This is a primary track. Set high water mark for this event - fHgwmk=fNtrack; + fStack->SetTrack(done, parent, pdg, px, py, pz, e, vx, vy, vz, tof, + polx, poly, polz, mech, ntr, weight, is); + +} + +//_______________________________________________________________________ +void AliRun::SetHighWaterMark(const Int_t nt) +{ // - // Set also number if primary tracks - fHeader.SetNprimary(fHgwmk+1); - fHeader.SetNtrack(fHgwmk+1); - } - ntr = fNtrack++; + // Set high water mark for last track in event + fStack->SetHighWaterMark(nt); } -//_____________________________________________________________________________ +//_______________________________________________________________________ void AliRun::KeepTrack(const Int_t track) { // - // flags a track to be kept + // Delegate to stack // - TClonesArray &particles = *fParticles; - ((TParticle*)particles[track])->SetBit(Keep_Bit); + fStack->KeepTrack(track); } -//_____________________________________________________________________________ -void AliRun::StepManager(Int_t id) const +// +// MC Application +// + +//_______________________________________________________________________ +void AliRun::ConstructGeometry() +{ + // + // Create modules, materials, geometry + // + + TStopwatch stw; + TIter next(fModules); + AliModule *detector; + printf("Geometry creation:\n"); + while((detector = dynamic_cast(next()))) { + stw.Start(); + // Initialise detector materials and geometry + detector->CreateMaterials(); + detector->CreateGeometry(); + printf("%10s R:%.2fs C:%.2fs\n", + detector->GetName(),stw.RealTime(),stw.CpuTime()); + } +} + +//_______________________________________________________________________ +void AliRun::InitGeometry() +{ + // + // Initialize detectors and display geometry + // + + printf("Initialisation:\n"); + TStopwatch stw; + TIter next(fModules); + AliModule *detector; + while((detector = dynamic_cast(next()))) { + stw.Start(); + // Initialise detector and display geometry + detector->Init(); + detector->BuildGeometry(); + printf("%10s R:%.2fs C:%.2fs\n", + detector->GetName(),stw.RealTime(),stw.CpuTime()); + } + +} + +//_______________________________________________________________________ +void AliRun::GeneratePrimaries() +{ + // + // Generate primary particles and fill them in the stack. + // + + Generator()->Generate(); +} + +//_______________________________________________________________________ +void AliRun::BeginEvent() +{ + // Clean-up previous event + // Energy scores + fEventEnergy.Reset(); + // Clean detector information + CleanDetectors(); + // Reset stack info + fStack->Reset(); + + + // + // Reset all Detectors & kinematics & trees + // + char hname[30]; + // + // Initialise event header + fHeader->Reset(fRun,fEvent,fEventNrInRun); + // + fStack->BeginEvent(fEvent); + + // + if(fLego) { + fLego->BeginEvent(); + return; + } + + // + + ResetHits(); + ResetTrackReferences(); + ResetDigits(); + ResetSDigits(); + + + if(fTreeH) { + fTreeH->Reset(); + sprintf(hname,"TreeH%d",fEvent); + fTreeH->SetName(hname); + } + + if(fTreeTR) { + fTreeTR->Reset(); + sprintf(hname,"TreeTR%d",fEvent); + fTreeTR->SetName(hname); + } + + if(fTreeD) { + fTreeD->Reset(); + sprintf(hname,"TreeD%d",fEvent); + fTreeD->SetName(hname); + fTreeD->Write(0,TObject::kOverwrite); + } + if(fTreeS) { + fTreeS->Reset(); + sprintf(hname,"TreeS%d",fEvent); + fTreeS->SetName(hname); + fTreeS->Write(0,TObject::kOverwrite); + } + if(fTreeR) { + fTreeR->Reset(); + sprintf(hname,"TreeR%d",fEvent); + fTreeR->SetName(hname); + fTreeR->Write(0,TObject::kOverwrite); + } +} + +//_______________________________________________________________________ +void AliRun::BeginPrimary() +{ + // + // Called at the beginning of each primary track + // + + // Reset Hits info + gAlice->ResetHits(); + gAlice->ResetTrackReferences(); + +} + +//_______________________________________________________________________ +void AliRun::PreTrack() +{ + // + // Method called before each track + // + TObjArray &dets = *fModules; + AliModule *module; + + for(Int_t i=0; i<=fNdets; i++) + if((module = dynamic_cast(dets[i]))) + module->PreTrack(); + + fMCQA->PreTrack(); +} + +//_______________________________________________________________________ +void AliRun::Stepping() { // // Called at every step during transport // + Int_t id = DetFromMate(gMC->GetMedium()); + if (id < 0) return; + // // --- If lego option, do it and leave if (fLego) @@ -1460,62 +1750,450 @@ void AliRun::StepManager(Int_t id) const else { Int_t copy; //Update energy deposition tables - sEventEnergy[gMC->CurrentVolID(copy)]+=gMC->Edep(); + AddEnergyDeposit(gMC->CurrentVolID(copy),gMC->Edep()); //Call the appropriate stepping routine; - AliModule *det = (AliModule*)fModules->At(id); - if(det) det->StepManager(); + AliModule *det = dynamic_cast(fModules->At(id)); + if(det && det->StepManagerIsEnabled()) { + fMCQA->StepManager(id); + det->StepManager(); + } } } -//_____________________________________________________________________________ -void AliRun::Streamer(TBuffer &R__b) +//_______________________________________________________________________ +void AliRun::PostTrack() { // - // Stream an object of class AliRun. + // Called after a track has been trasported + // + TObjArray &dets = *fModules; + AliModule *module; + + for(Int_t i=0; i<=fNdets; i++) + if((module = dynamic_cast(dets[i]))) + module->PostTrack(); +} + +//_______________________________________________________________________ +void AliRun::FinishPrimary() +{ + // + // Called at the end of each primary track + // + + // static Int_t count=0; + // const Int_t times=10; + // This primary is finished, purify stack + fStack->PurifyKine(); + + TIter next(fModules); + AliModule *detector; + while((detector = dynamic_cast(next()))) { + detector->FinishPrimary(); + } + + // Write out hits if any + if (gAlice->TreeH()) { + gAlice->TreeH()->Fill(); + } + + // Write out hits if any + if (gAlice->TreeTR()) { + gAlice->TreeTR()->Fill(); + } + + // + // if(++count%times==1) gObjectTable->Print(); +} + +//_______________________________________________________________________ +void AliRun::FinishEvent() +{ + // + // Called at the end of the event. + // + + // + if(fLego) fLego->FinishEvent(); + + //Update the energy deposit tables + Int_t i; + for(i=0;iSetNprimary(fStack->GetNprimary()); + fHeader->SetNtrack(fStack->GetNtrack()); + + + // Write out the kinematics + fStack->FinishEvent(); + + // Write out the event Header information + if (fTreeE) { + fHeader->SetStack(fStack); + fTreeE->Fill(); + } + + + // Write Tree headers + TTree* pTreeK = fStack->TreeK(); + if (pTreeK) pTreeK->Write(0,TObject::kOverwrite); + if (fTreeH) fTreeH->Write(0,TObject::kOverwrite); + if (fTreeTR) fTreeTR->Write(0,TObject::kOverwrite); + + ++fEvent; + ++fEventNrInRun; +} + +//_______________________________________________________________________ +void AliRun::Field(const Double_t* x, Double_t *b) const +{ + // + // Returns the magnetic field at point x[3] + // Units are kGauss // + Float_t xfloat[3]; + for (Int_t i=0; i<3; i++) xfloat[i] = x[i]; + + if (Field()) { + Float_t bfloat[3]; + Field()->Field(xfloat,bfloat); + for (Int_t j=0; j<3; j++) b[j] = bfloat[j]; + } + else { + printf("No mag field defined!\n"); + b[0]=b[1]=b[2]=0.; + } + +} + +// +// End of MC Application +// + +//_______________________________________________________________________ +void AliRun::Streamer(TBuffer &R__b) +{ + // Stream an object of class AliRun. + if (R__b.IsReading()) { - Version_t R__v = R__b.ReadVersion(); if (R__v) { } - TNamed::Streamer(R__b); if (!gAlice) gAlice = this; + + AliRun::Class()->ReadBuffer(R__b, this); + // gROOT->GetListOfBrowsables()->Add(this,"Run"); - fTreeE = (TTree*)gDirectory->Get("TE"); - if (fTreeE) fTreeE->SetBranchAddress("Header", &header); + + fTreeE = dynamic_cast(gDirectory->Get("TE")); + if (fTreeE) { + fTreeE->SetBranchAddress("Header", &fHeader); + } else Error("Streamer","cannot find Header Tree\n"); - R__b >> fNtrack; - R__b >> fHgwmk; - R__b >> fDebug; - fHeader.Streamer(R__b); - R__b >> fModules; - R__b >> fParticles; - R__b >> fField; - // R__b >> fMC; - R__b >> fNdets; - R__b >> fTrRmax; - R__b >> fTrZmax; - R__b >> fGenerator; - if(R__v>1) { - R__b >> fPDGDB; //Particle factory object! - fTreeE->GetEntry(0); - } else { - fHeader.SetEvent(0); - fPDGDB = TDatabasePDG::Instance(); //Particle factory object! - } + + fTreeE->GetEntry(0); + gRandom = fRandom; } else { - R__b.WriteVersion(AliRun::IsA()); - TNamed::Streamer(R__b); - R__b << fNtrack; - R__b << fHgwmk; - R__b << fDebug; - fHeader.Streamer(R__b); - R__b << fModules; - R__b << fParticles; - R__b << fField; - // R__b << fMC; - R__b << fNdets; - R__b << fTrRmax; - R__b << fTrZmax; - R__b << fGenerator; - R__b << fPDGDB; //Particle factory object! - } -} + AliRun::Class()->WriteBuffer(R__b, this); + } +} + + +//_______________________________________________________________________ +Int_t AliRun::CurrentTrack() const { + // + // Returns current track + // + return fStack->CurrentTrack(); +} + +//_______________________________________________________________________ +Int_t AliRun::GetNtrack() const { + // + // Returns number of tracks in stack + // + return fStack->GetNtrack(); +} + +//_______________________________________________________________________ +TObjArray* AliRun::Particles() const { + // + // Returns pointer to Particles array + // + return fStack->Particles(); +} + +//_______________________________________________________________________ +TTree* AliRun::TreeK() const { + // + // Returns pointer to the TreeK array + // + return fStack->TreeK(); +} + + +//_______________________________________________________________________ +void AliRun::SetGenEventHeader(AliGenEventHeader* header) +{ + fHeader->SetGenEventHeader(header); +} + +//_______________________________________________________________________ +TFile* AliRun::InitFile(TString fileName) +{ +// +// create the file where the whole tree will be saved +// + TDirectory *wd = gDirectory; + TFile* file = TFile::Open(fileName,"update"); + gDirectory = wd; + if (!file->IsOpen()) { + Error("Cannot open file, %s\n",fileName); + return 0; + } + return file; +} + +//_______________________________________________________________________ +TFile* AliRun::InitTreeFile(Option_t *option, TString fileName) +{ + // + // create the file where one of the following trees will be saved + // trees: S,D,R + // WARNING: by default these trees are saved on the file on which + // hits are stored. If you divert one of these trees, you cannot restore + // it to the original file (usually galice.root) in the same aliroot session + Bool_t oS = (strstr(option,"S")!=0); + Bool_t oR = (strstr(option,"R")!=0); + Bool_t oD = (strstr(option,"D")!=0); + Int_t choice[3]; + for (Int_t i=0; i<3; i++) choice[i] = 0; + if(oS)choice[0] = 1; + if(oD)choice[1] = 1; + if(oR)choice[2] = 1; + + TFile *ptr=0; + + if(!(oS || oR || oD))return ptr; + + Int_t active[3]; + for (Int_t i=0; i<3; i++) active[i] = 0; + if(fTreeSFileName != "") active[0] = 1; + if(fTreeDFileName != "") active[1] = 1; + if(fTreeDFileName != "") active[2] = 1; + + Bool_t alreadyopen1 = kFALSE; + Bool_t alreadyopen2 = kFALSE; + + if(oS){ + // if already active and same name with non-null ptr + if(active[0]==1 && fileName == fTreeSFileName && fTreeSFile){ + Warning("InitTreeFile","File %s already opened",fTreeSFileName.Data()); + ptr = fTreeSFile; + } + else { + // if already active with different name with non-null ptr + if(active[0]==1 && fileName != fTreeSFileName && fTreeSFile){ + // close the active files and also the other possible files in option + CloseTreeFile(option); + } + fTreeSFileName = fileName; + alreadyopen1 = + (active[1] == 1 && fTreeDFileName == fTreeSFileName && fTreeDFile); + alreadyopen2 = + (active[2] == 1 && fTreeRFileName == fTreeSFileName && fTreeRFile); + if(!(alreadyopen1 || alreadyopen2)){ + ptr = InitFile(fileName); + fTreeSFile = ptr; + } + else { + if(alreadyopen1){fTreeSFile = fTreeDFile; ptr = fTreeSFile;} + if(alreadyopen2){fTreeSFile = fTreeRFile; ptr = fTreeSFile;} + } + if(choice[1] == 1) { fTreeDFileName = fileName; fTreeDFile = ptr;} + if(choice[2] == 1) { fTreeRFileName = fileName; fTreeRFile = ptr;} + } + return ptr; + } + + if(oD){ + // if already active and same name with non-null ptr + if(active[1]==1 && fileName == fTreeDFileName && fTreeDFile){ + Warning("InitTreeFile","File %s already opened",fTreeDFileName.Data()); + ptr = fTreeDFile; + } + else { + // if already active with different name with non-null ptr + if(active[1]==1 && fileName != fTreeDFileName && fTreeDFile){ + // close the active files and also the other possible files in option + CloseTreeFile(option); + } + fTreeDFileName = fileName; + alreadyopen1 = + (active[0] == 1 && fTreeSFileName == fTreeDFileName && fTreeSFile); + alreadyopen2 = + (active[2] == 1 && fTreeRFileName == fTreeDFileName && fTreeRFile); + if(!(alreadyopen1 || alreadyopen2)){ + ptr = InitFile(fileName); + fTreeDFile = ptr; + } + else { + if(alreadyopen1){fTreeDFile = fTreeSFile; ptr = fTreeDFile;} + if(alreadyopen2){fTreeDFile = fTreeRFile; ptr = fTreeDFile;} + } + if(choice[2] == 1) { fTreeRFileName = fileName; fTreeRFile = ptr;} + } + return ptr; + } + + if(oR){ + // if already active and same name with non-null ptr + if(active[2]==1 && fileName == fTreeRFileName && fTreeRFile){ + Warning("InitTreeFile","File %s already opened",fTreeRFileName.Data()); + ptr = fTreeRFile; + } + else { + // if already active with different name with non-null ptr + if(active[2]==1 && fileName != fTreeRFileName && fTreeRFile){ + // close the active files and also the other possible files in option + CloseTreeFile(option); + } + fTreeRFileName = fileName; + alreadyopen1 = + (active[1] == 1 && fTreeDFileName == fTreeRFileName && fTreeDFile); + alreadyopen2 = + (active[0]== 1 && fTreeSFileName == fTreeRFileName && fTreeSFile); + if(!(alreadyopen1 || alreadyopen2)){ + ptr = InitFile(fileName); + fTreeRFile = ptr; + } + else { + if(alreadyopen1){fTreeRFile = fTreeDFile; ptr = fTreeRFile;} + if(alreadyopen2){fTreeRFile = fTreeSFile; ptr = fTreeRFile;} + } + } + return ptr; + } + return 0; +} + +//_______________________________________________________________________ +void AliRun::PrintTreeFile() +{ + // + // prints the file names and pointer associated to S,D,R trees + // + cout<<"===================================================\n"; + TFile *file = fTreeE->GetCurrentFile(); + TString curfilname=""; + if(file)curfilname=static_cast(file->GetName()); + cout<<" Current tree file name: "<IsOpen()){ + fTreeSFile->Close(); + delete fTreeSFile; + } + } + fTreeSFile = 0; + } + if(oD){ + fTreeDFileName = ""; + if(fTreeDFile){ + if(!((fTreeDFile == fTreeRFile) || (fTreeDFile == fTreeSFile)) && + fTreeDFile->IsOpen()){ + fTreeDFile->Close(); + delete fTreeDFile; + } + } + fTreeDFile = 0; + } + if(oR){ + fTreeRFileName = ""; + if(fTreeRFile){ + if(!((fTreeRFile == fTreeSFile) || (fTreeRFile == fTreeDFile)) && + fTreeRFile->IsOpen()){ + fTreeRFile->Close(); + delete fTreeRFile; + } + } + fTreeRFile = 0; + } +} + +//_______________________________________________________________________ +void AliRun::MakeTree(Option_t *option, TFile *file) +{ + // + // Create some trees in the separate file + // + const char *oD = strstr(option,"D"); + const char *oR = strstr(option,"R"); + const char *oS = strstr(option,"S"); + + TDirectory *cwd = gDirectory; + char hname[30]; + + if (oD) { + delete fTreeD; + sprintf(hname,"TreeD%d",fEvent); + file->cd(); + fTreeD = static_cast(file->Get("hname")); + if (!fTreeD) { + fTreeD = new TTree(hname,"Digits"); + fTreeD->Write(0,TObject::kOverwrite); + } + cwd->cd(); + } + if (oS) { + delete fTreeS; + sprintf(hname,"TreeS%d",fEvent); + file->cd(); + fTreeS = static_cast(file->Get("hname")); + if (!fTreeS) { + fTreeS = new TTree(hname,"SDigits"); + fTreeS->Write(0,TObject::kOverwrite); + } + cwd->cd(); + } + + if (oR) { + delete fTreeR; + sprintf(hname,"TreeR%d",fEvent); + file->cd(); + fTreeR = static_cast(file->Get("hname")); + if (!fTreeR) { + fTreeR = new TTree(hname,"RecPoint"); + fTreeR->Write(0,TObject::kOverwrite); + } + cwd->cd(); + } +}