Adding Generation via an external file.
[u/mrichter/AliRoot.git] / EVGEN / AliGenExtFile.cxx
1 #include "AliGenExtFile.h"
2 #include "AliGenMUONlib.h"
3 #include "AliMC.h"
4 #include "AliRun.h"
5 #include <TDirectory.h>
6 #include <TFile.h>
7 #include <TTree.h>
8 #include <stdlib.h>
9  ClassImp(AliGenExtFile)
10      AliGenExtFile::AliGenExtFile()
11          :AliGenerator(-1)
12 {
13     //
14     fName="ExtFile";
15     fTitle="Primaries from ext. File";
16     fFileName="dtujet93.root";
17     fTreeNtuple=0;
18     fNcurrent=0;
19 //
20 //  Read all particles
21     fNpart=-1;
22 }
23
24 AliGenExtFile::AliGenExtFile(Int_t npart)
25     :AliGenerator(npart)
26 {
27     //
28     fName="ExtFile";
29     fTitle="Primaries from ext. File";
30     fFileName="dtujet93.root";
31     fTreeNtuple=0;
32     fNcurrent=0;
33 }
34
35 //____________________________________________________________
36 AliGenExtFile::~AliGenExtFile()
37 {
38     delete fTreeNtuple;
39 }
40
41 //____________________________________________________________
42 void AliGenExtFile::NtupleInit() 
43 {
44 //
45 // reset the existing file environment and open a new root file if
46 // the pointer to the Fluka tree is null
47     
48     TFile *File=0;
49     if (fTreeNtuple==0) {
50         if (!File) {
51             File = new TFile(fFileName);
52             File->cd();
53             cout<<"I have opened "<<fFileName<<" file "<<endl;
54         }
55 // get the tree address in the Fluka boundary source file
56         fTreeNtuple = (TTree*)gDirectory->Get("h888");
57     } else {
58         File = fTreeNtuple->GetCurrentFile();
59         File->cd();
60     }
61
62     TTree *h2=fTreeNtuple;
63 //Set branch addresses
64 //Set branch addresses
65     h2->SetBranchAddress("Nihead",&Nihead);
66     h2->SetBranchAddress("Ihead",Ihead);
67     h2->SetBranchAddress("Nrhead",&Nrhead);
68     h2->SetBranchAddress("Rhead",Rhead);
69     h2->SetBranchAddress("Idpart",&Idpart);
70     h2->SetBranchAddress("Theta",&Theta);
71     h2->SetBranchAddress("Phi",&Phi);
72     h2->SetBranchAddress("P",&P);
73     h2->SetBranchAddress("E",&E);
74 }
75
76
77 //____________________________________________________________
78 void AliGenExtFile::Generate()
79 {
80
81   AliMC* pMC = AliMC::GetMC();
82
83   Float_t polar[3]= {0,0,0};
84   //
85   Float_t origin[3]={0,0,0};
86   Float_t p[3];
87   Float_t random[6];
88   Float_t prwn;
89   char name[100];
90   Float_t amass, charge, tlife;
91   Int_t itrtyp;
92   Int_t i, j, nt, Ntracks;
93   //
94   NtupleInit();
95   TTree *h2=fTreeNtuple;
96   Int_t nentries = (Int_t) h2->GetEntries();
97   // loop over number of particles
98   Int_t nb = (Int_t)h2->GetEvent(fNcurrent);
99   Int_t i5=Ihead[4];
100   Int_t i6=Ihead[5];
101
102   for (j=0;j<3;j++) origin[j]=fOrigin[j];
103   if(fVertexSmear==perEvent) {
104     pMC->Rndm(random,6);
105     for (j=0;j<3;j++) {
106         origin[j]+=fOsigma[j]*TMath::Cos(2*random[2*j]*TMath::Pi())*
107             TMath::Sqrt(-2*TMath::Log(random[2*j+1]));
108     }
109   }
110
111   if (fNcurrent >= nentries) {
112       printf("\n No more entries !!! !\n");
113       return;
114   }
115   
116           
117   if (i5==0) {
118       printf("\n This should never happen !\n");
119   } else {
120       printf("\n Next event contains %d tracks! \n", i6);
121       Ntracks=i6;
122   }
123   for (i=0; i<Ntracks; i++) {
124
125       pMC->Gfpart(Idpart, name, itrtyp,amass, charge, tlife); 
126       prwn=sqrt((E+amass)*(E-amass));
127
128       Theta *= TMath::Pi()/180.;
129       Phi    = (Phi-180)*TMath::Pi()/180.;      
130       if(Theta<fThetaMin || Theta>fThetaMax ||
131          Phi<fPhiMin || Phi>fPhiMax         ||
132          prwn<fPMin || prwn>fPMax)          
133       {
134           ;
135       } else {
136           p[0]=prwn*TMath::Sin(Theta)*TMath::Cos(Phi);
137           p[1]=prwn*TMath::Sin(Theta)*TMath::Sin(Phi);      
138           p[2]=prwn*TMath::Cos(Theta);
139           
140           if(fVertexSmear==perTrack) {
141               pMC->Rndm(random,6);
142               for (j=0;j<3;j++) {
143                   origin[j]=fOrigin[j]
144                       +fOsigma[j]*TMath::Cos(2*random[2*j]*TMath::Pi())*
145                       TMath::Sqrt(-2*TMath::Log(random[2*j+1]));
146               }
147           }
148           gAlice->SetTrack(fTrackIt,-1,Idpart,p,origin,polar,0,"Primary",nt);
149       }
150       fNcurrent++;
151       nb = (Int_t)h2->GetEvent(fNcurrent); 
152   }
153  
154     TFile *File=0;
155 // Get AliRun object or create it 
156     if (!gAlice) {
157         gAlice = (AliRun*)File->Get("gAlice");
158         if (gAlice) printf("AliRun object found on file\n");
159         if (!gAlice) gAlice = new AliRun("gAlice","Alice test program");
160     }
161     TTree *fAli=gAlice->TreeK();
162     if (fAli) File =fAli->GetCurrentFile();
163     File->cd();
164 }
165
166
167
168
169
170
171
172
173