X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;ds=sidebyside;f=EVGEN%2FAliGenExtFile.cxx;h=7382b1950481632481b4812d75263ffe04cb0581;hb=88cb7938ca21d4a80991d4e7aa564008c29340f7;hp=364c02b804b64b966f185071d9ee60b8ba24b66a;hpb=cfce887023b4d349cb7b743880a8f9086a600d9d;p=u%2Fmrichter%2FAliRoot.git diff --git a/EVGEN/AliGenExtFile.cxx b/EVGEN/AliGenExtFile.cxx index 364c02b804b..7382b195048 100644 --- a/EVGEN/AliGenExtFile.cxx +++ b/EVGEN/AliGenExtFile.cxx @@ -1,166 +1,221 @@ +/************************************************************************** + * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * + * * + * Author: The ALICE Off-line Project. * + * Contributors are mentioned in the code where appropriate. * + * * + * Permission to use, copy, modify and distribute this software and its * + * documentation strictly for non-commercial purposes is hereby granted * + * without fee, provided that the above copyright notice appears in all * + * copies and that both the copyright notice and this permission notice * + * appear in the supporting documentation. The authors make no claims * + * about the suitability of this software for any purpose. It is * + * provided "as is" without express or implied warranty. * + **************************************************************************/ + +/* $Id$ */ + +// Event generator that using an instance of type AliGenReader +// reads particles from a file and applies cuts. +// Example: In your Config.C you can include the following lines +// AliGenExtFile *gener = new AliGenExtFile(-1); +// gener->SetMomentumRange(0,999); +// gener->SetPhiRange(-180.,180.); +// gener->SetThetaRange(0,180); +// gener->SetYRange(-999,999); +// AliGenReaderTreeK * reader = new AliGenReaderTreeK(); +// reader->SetFileName("myFileWithTreeK.root"); +// gener->SetReader(reader); +// gener->Init(); + + +#include + #include "AliGenExtFile.h" -#include "AliGenMUONlib.h" -#include "AliMC.h" #include "AliRun.h" -#include + +#include #include #include -#include + + ClassImp(AliGenExtFile) - AliGenExtFile::AliGenExtFile() - :AliGenerator(-1) + +AliGenExtFile::AliGenExtFile() + :AliGenMC() { - // - fName="ExtFile"; - fTitle="Primaries from ext. File"; - fFileName="dtujet93.root"; - fTreeNtuple=0; - fNcurrent=0; +// Constructor // // Read all particles - fNpart=-1; + fNpart =- 1; + fReader = 0; } AliGenExtFile::AliGenExtFile(Int_t npart) - :AliGenerator(npart) + :AliGenMC(npart) { - // - fName="ExtFile"; - fTitle="Primaries from ext. File"; - fFileName="dtujet93.root"; - fTreeNtuple=0; - fNcurrent=0; +// Constructor + fName = "ExtFile"; + fTitle = "Primaries from ext. File"; + fReader = 0; } +AliGenExtFile::AliGenExtFile(const AliGenExtFile & ExtFile) +{ +// copy constructor +} //____________________________________________________________ AliGenExtFile::~AliGenExtFile() { - delete fTreeNtuple; +// Destructor + delete fReader; } -//____________________________________________________________ -void AliGenExtFile::NtupleInit() +//___________________________________________________________ +void AliGenExtFile::Init() { -// -// reset the existing file environment and open a new root file if -// the pointer to the Fluka tree is null - - TFile *File=0; - if (fTreeNtuple==0) { - if (!File) { - File = new TFile(fFileName); - File->cd(); - cout<<"I have opened "<Get("h888"); - } else { - File = fTreeNtuple->GetCurrentFile(); - File->cd(); - } - - TTree *h2=fTreeNtuple; -//Set branch addresses -//Set branch addresses - h2->SetBranchAddress("Nihead",&Nihead); - h2->SetBranchAddress("Ihead",Ihead); - h2->SetBranchAddress("Nrhead",&Nrhead); - h2->SetBranchAddress("Rhead",Rhead); - h2->SetBranchAddress("Idpart",&Idpart); - h2->SetBranchAddress("Theta",&Theta); - h2->SetBranchAddress("Phi",&Phi); - h2->SetBranchAddress("P",&P); - h2->SetBranchAddress("E",&E); +// Initialize + if (fReader) fReader->Init(); } - -//____________________________________________________________ + void AliGenExtFile::Generate() { +// Generate particles - Float_t polar[3]= {0,0,0}; + Float_t polar[3] = {0,0,0}; // - Float_t origin[3]={0,0,0}; + Float_t origin[3] = {0,0,0}; Float_t p[3]; Float_t random[6]; - Float_t prwn; - char name[100]; - Float_t amass, charge, tlife; - Int_t itrtyp; - Int_t i, j, nt, Ntracks=0; + Int_t i = 0, j, nt; // - NtupleInit(); - TTree *h2=fTreeNtuple; - Int_t nentries = (Int_t) h2->GetEntries(); - // loop over number of particles - Int_t nb = (Int_t)h2->GetEvent(fNcurrent); - Int_t i5=Ihead[4]; - Int_t i6=Ihead[5]; - for (j=0;j<3;j++) origin[j]=fOrigin[j]; - if(fVertexSmear==perEvent) { - gMC->Rndm(random,6); - for (j=0;j<3;j++) { - origin[j]+=fOsigma[j]*TMath::Cos(2*random[2*j]*TMath::Pi())* + if(fVertexSmear == kPerTrack) { + Rndm(random,6); + for (j = 0; j < 3; j++) { + origin[j] += fOsigma[j]*TMath::Cos(2*random[2*j]*TMath::Pi())* TMath::Sqrt(-2*TMath::Log(random[2*j+1])); } } - if (fNcurrent >= nentries) { - printf("\n No more entries !!! !\n"); + while(1) { + Int_t nTracks = fReader->NextEvent(); + if (nTracks == 0) { + // printf("\n No more events !!! !\n"); + Warning("AliGenExtFile::Generate","\nNo more events in external file!!!\nLast event may be empty or incomplete.\n"); return; - } - - - if (i5==0) { - printf("\n This should never happen !\n"); - } else { - printf("\n Next event contains %d tracks! \n", i6); - Ntracks=i6; - } - for (i=0; iGfpart(Idpart, name, itrtyp,amass, charge, tlife); - prwn=sqrt((E+amass)*(E-amass)); - - Theta *= TMath::Pi()/180.; - Phi = (Phi-180)*TMath::Pi()/180.; - if(ThetafThetaMax || - PhifPhiMax || - prwnfPMax) - { - ; - } else { - p[0]=prwn*TMath::Sin(Theta)*TMath::Cos(Phi); - p[1]=prwn*TMath::Sin(Theta)*TMath::Sin(Phi); - p[2]=prwn*TMath::Cos(Theta); - - if(fVertexSmear==perTrack) { - gMC->Rndm(random,6); - for (j=0;j<3;j++) { - origin[j]=fOrigin[j] - +fOsigma[j]*TMath::Cos(2*random[2*j]*TMath::Pi())* - TMath::Sqrt(-2*TMath::Log(random[2*j+1])); - } - } - gAlice->SetTrack(fTrackIt,-1,Idpart,p,origin,polar,0,"Primary",nt); + } + + // + // Particle selection loop + // + // The selction criterium for the external file generator is as follows: + // + // 1) All tracks are subjects to the cuts defined by AliGenerator, i.e. + // fThetaMin, fThetaMax, fPhiMin, fPhiMax, fPMin, fPMax, fPtMin, fPtMax, + // fYMin, fYMax. + // If the particle does not satisfy these cuts, it is not put on the + // stack. + // 2) If fCutOnChild and some specific child is selected (e.g. if + // fForceDecay==kSemiElectronic) the event is rejected if NOT EVEN ONE + // child falls into the child-cuts. + TParticle* iparticle = 0x0; + + if(fCutOnChild) { + // Count the selected children + Int_t nSelected = 0; + while ( (iparticle=fReader->NextParticle()) ) { + ; + Int_t kf = CheckPDGCode(iparticle->GetPdgCode()); + kf = TMath::Abs(kf); + if (ChildSelected(kf) && KinematicSelection(iparticle, 1)) { + nSelected++; + } } - fNcurrent++; - nb = (Int_t)h2->GetEvent(fNcurrent); - } - - TFile *File=0; -// Get AliRun object or create it - if (!gAlice) { - gAlice = (AliRun*)File->Get("gAlice"); - if (gAlice) printf("AliRun object found on file\n"); - if (!gAlice) gAlice = new AliRun("gAlice","Alice test program"); + if (!nSelected) continue; // No particle selected: Go to next event + fReader->RewindEvent(); } - TTree *fAli=gAlice->TreeK(); - if (fAli) File =fAli->GetCurrentFile(); - File->cd(); + + // + // Stack filling loop + // + for (i = 0; i < nTracks; i++) { + + TParticle* iparticle = fReader->NextParticle(); + Bool_t selected = KinematicSelection(iparticle,0); + if (!selected) { + Double_t pz = iparticle->Pz(); + Double_t e = iparticle->Energy(); + Double_t y; + if ((e-pz) == 0) { + y = 20.; + } else if ((e+pz) == 0.) { + y = -20.; + } else { + y = 0.5*TMath::Log((e+pz)/(e-pz)); + } + printf("\n Not selected %d %f %f %f %f %f", i, + iparticle->Theta(), + iparticle->Phi(), + iparticle->P(), + iparticle->Pt(), + y); + //PH delete iparticle; + // continue; + } + p[0] = iparticle->Px(); + p[1] = iparticle->Py(); + p[2] = iparticle->Pz(); + Int_t idpart = iparticle->GetPdgCode(); + if(fVertexSmear==kPerTrack) + { + Rndm(random,6); + for (j = 0; j < 3; j++) { + origin[j]=fOrigin[j]+ + fOsigma[j]*TMath::Cos(2*random[2*j]*TMath::Pi())* + TMath::Sqrt(-2*TMath::Log(random[2*j+1])); + } + } + Int_t decayed = iparticle->GetFirstDaughter(); + Int_t doTracking = fTrackIt && (decayed < 0) && + (TMath::Abs(idpart) > 10) && selected; + // printf("*** pdg, first daughter, trk = %d, %d, %d\n", + // idpart,decayed, doTracking); + //PH SetTrack(doTracking,-1,idpart,p,origin,polar,0,kPPrimary,nt); + Int_t parent = iparticle->GetFirstMother(); + SetTrack(doTracking,parent,idpart,p,origin,polar,0,kPPrimary,nt); + KeepTrack(nt); + } // track loop + + break; + + } // event loop + + SetHighWaterMark(nt); + CdEventFile(); } +void AliGenExtFile::CdEventFile() +{ +// CD back to the event file + TFile *pFile=0; + if (!gAlice) { + gAlice = (AliRun*)pFile->Get("gAlice"); + if (gAlice) printf("AliRun object found on file\n"); + if (!gAlice) gAlice = new AliRun("gAlice","Alice test program"); + } + TTree *fAli=gAlice->TreeK(); + if (fAli) pFile =fAli->GetCurrentFile(); + pFile->cd(); +} + + +AliGenExtFile& AliGenExtFile::operator=(const AliGenExtFile& rhs) +{ +// Assignment operator + return *this; +}