Use gMC and not pMC everywhere
[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   Float_t polar[3]= {0,0,0};
82   //
83   Float_t origin[3]={0,0,0};
84   Float_t p[3];
85   Float_t random[6];
86   Float_t prwn;
87   char name[100];
88   Float_t amass, charge, tlife;
89   Int_t itrtyp;
90   Int_t i, j, nt, Ntracks=0;
91   //
92   NtupleInit();
93   TTree *h2=fTreeNtuple;
94   Int_t nentries = (Int_t) h2->GetEntries();
95   // loop over number of particles
96   Int_t nb = (Int_t)h2->GetEvent(fNcurrent);
97   Int_t i5=Ihead[4];
98   Int_t i6=Ihead[5];
99
100   for (j=0;j<3;j++) origin[j]=fOrigin[j];
101   if(fVertexSmear==perEvent) {
102     gMC->Rndm(random,6);
103     for (j=0;j<3;j++) {
104         origin[j]+=fOsigma[j]*TMath::Cos(2*random[2*j]*TMath::Pi())*
105             TMath::Sqrt(-2*TMath::Log(random[2*j+1]));
106     }
107   }
108
109   if (fNcurrent >= nentries) {
110       printf("\n No more entries !!! !\n");
111       return;
112   }
113   
114           
115   if (i5==0) {
116       printf("\n This should never happen !\n");
117   } else {
118       printf("\n Next event contains %d tracks! \n", i6);
119       Ntracks=i6;
120   }
121   for (i=0; i<Ntracks; i++) {
122
123       gMC->Gfpart(Idpart, name, itrtyp,amass, charge, tlife); 
124       prwn=sqrt((E+amass)*(E-amass));
125
126       Theta *= TMath::Pi()/180.;
127       Phi    = (Phi-180)*TMath::Pi()/180.;      
128       if(Theta<fThetaMin || Theta>fThetaMax ||
129          Phi<fPhiMin || Phi>fPhiMax         ||
130          prwn<fPMin || prwn>fPMax)          
131       {
132           ;
133       } else {
134           p[0]=prwn*TMath::Sin(Theta)*TMath::Cos(Phi);
135           p[1]=prwn*TMath::Sin(Theta)*TMath::Sin(Phi);      
136           p[2]=prwn*TMath::Cos(Theta);
137           
138           if(fVertexSmear==perTrack) {
139               gMC->Rndm(random,6);
140               for (j=0;j<3;j++) {
141                   origin[j]=fOrigin[j]
142                       +fOsigma[j]*TMath::Cos(2*random[2*j]*TMath::Pi())*
143                       TMath::Sqrt(-2*TMath::Log(random[2*j+1]));
144               }
145           }
146           gAlice->SetTrack(fTrackIt,-1,Idpart,p,origin,polar,0,"Primary",nt);
147       }
148       fNcurrent++;
149       nb = (Int_t)h2->GetEvent(fNcurrent); 
150   }
151  
152     TFile *File=0;
153 // Get AliRun object or create it 
154     if (!gAlice) {
155         gAlice = (AliRun*)File->Get("gAlice");
156         if (gAlice) printf("AliRun object found on file\n");
157         if (!gAlice) gAlice = new AliRun("gAlice","Alice test program");
158     }
159     TTree *fAli=gAlice->TreeK();
160     if (fAli) File =fAli->GetCurrentFile();
161     File->cd();
162 }
163
164
165
166
167
168
169
170
171