Coding convention corrections + few minor bug fixes
[u/mrichter/AliRoot.git] / EVGEN / AliGenExtFile.cxx
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
16 /*
17 $Log$
18 Revision 1.9  2000/06/14 15:20:09  morsch
19 Include clean-up (IH)
20
21 Revision 1.8  2000/06/09 20:36:44  morsch
22 All coding rule violations except RS3 corrected
23
24 Revision 1.7  2000/02/16 14:56:27  morsch
25 Convert geant particle code into pdg code before putting particle on the stack.
26
27 Revision 1.6  1999/11/09 07:38:48  fca
28 Changes for compatibility with version 2.23 of ROOT
29
30 Revision 1.5  1999/09/29 09:24:12  fca
31 Introduction of the Copyright and cvs Log
32
33 */
34
35 #include <iostream.h>
36
37 #include "AliGenExtFile.h"
38 #include "AliMC.h"
39 #include "AliRun.h"
40
41 #include <TDirectory.h>
42 #include <TFile.h>
43 #include <stdlib.h>
44
45  ClassImp(AliGenExtFile)
46      AliGenExtFile::AliGenExtFile()
47          :AliGenerator(-1)
48 {
49 //  Constructor
50     fName="ExtFile";
51     fTitle="Primaries from ext. File";
52     fFileName="dtujet93.root";
53     fTreeNtuple=0;
54     fNcurrent=0;
55 //
56 //  Read all particles
57     fNpart=-1;
58 }
59
60 AliGenExtFile::AliGenExtFile(Int_t npart)
61     :AliGenerator(npart)
62 {
63 //  Constructor
64     fName="ExtFile";
65     fTitle="Primaries from ext. File";
66     fFileName="dtujet93.root";
67     fTreeNtuple=0;
68     fNcurrent=0;
69 }
70
71 AliGenExtFile::AliGenExtFile(const AliGenExtFile & ExtFile)
72 {
73 // copy constructor
74 }
75 //____________________________________________________________
76 AliGenExtFile::~AliGenExtFile()
77 {
78 // Destructor
79     delete fTreeNtuple;
80 }
81
82 //____________________________________________________________
83 void AliGenExtFile::NtupleInit() 
84 {
85 //
86 // reset the existing file environment and open a new root file if
87 // the pointer to the Fluka tree is null
88     
89     TFile *pFile=0;
90     if (fTreeNtuple==0) {
91         if (!pFile) {
92             pFile = new TFile(fFileName);
93             pFile->cd();
94             cout<<"I have opened "<<fFileName<<" file "<<endl;
95         }
96 // get the tree address in the Fluka boundary source file
97         fTreeNtuple = (TTree*)gDirectory->Get("h888");
98     } else {
99         pFile = fTreeNtuple->GetCurrentFile();
100         pFile->cd();
101     }
102
103     TTree *h2=fTreeNtuple;
104 //Set branch addresses
105 //Set branch addresses
106     h2->SetBranchAddress("Nihead",&fNihead);
107     h2->SetBranchAddress("Ihead",fIhead);
108     h2->SetBranchAddress("Nrhead",&fNrhead);
109     h2->SetBranchAddress("Rhead",fRhead);
110     h2->SetBranchAddress("Idpart",&fIdpart);
111     h2->SetBranchAddress("Theta",&fTheta);
112     h2->SetBranchAddress("Phi",&fPhi);
113     h2->SetBranchAddress("P",&fP);
114     h2->SetBranchAddress("E",&fE);
115 }
116
117
118 //____________________________________________________________
119 void AliGenExtFile::Generate()
120 {
121 // Generate particles
122
123   Float_t polar[3]= {0,0,0};
124   //
125   Float_t origin[3]={0,0,0};
126   Float_t p[3];
127   Float_t random[6];
128   Float_t prwn;
129   Int_t i, j, nt, nTracks=0;
130   //
131   NtupleInit();
132   TTree *h2=fTreeNtuple;
133   Int_t nentries = (Int_t) h2->GetEntries();
134   // loop over number of particles
135   Int_t nb = (Int_t)h2->GetEvent(fNcurrent);
136   Int_t i5=fIhead[4];
137   Int_t i6=fIhead[5];
138
139   for (j=0;j<3;j++) origin[j]=fOrigin[j];
140   if(fVertexSmear==kPerTrack) {
141     gMC->Rndm(random,6);
142     for (j=0;j<3;j++) {
143         origin[j]+=fOsigma[j]*TMath::Cos(2*random[2*j]*TMath::Pi())*
144             TMath::Sqrt(-2*TMath::Log(random[2*j+1]));
145     }
146   }
147
148   if (fNcurrent >= nentries) {
149       printf("\n No more entries !!! !\n");
150       return;
151   }
152   
153           
154   if (i5==0) {
155       printf("\n This should never happen !\n");
156   } else {
157       printf("\n Next event contains %d tracks! \n", i6);
158       nTracks=i6;
159   }
160   for (i=0; i<nTracks; i++) {
161       fIdpart=gMC->PDGFromId(fIdpart);
162       Double_t amass = TDatabasePDG::Instance()->GetParticle(fIdpart)->Mass();
163       if(fE<=amass) {
164         Warning("Generate","Particle %d no %d E = %f mass = %f\n",
165                 fIdpart,i,fE,amass);
166         prwn=0;
167       } else {
168         prwn=sqrt((fE+amass)*(fE-amass));
169       }
170
171       fTheta *= TMath::Pi()/180.;
172       fPhi    = (fPhi-180)*TMath::Pi()/180.;      
173       if(fTheta<fThetaMin || fTheta>fThetaMax ||
174          fPhi<fPhiMin || fPhi>fPhiMax         ||
175          prwn<fPMin || prwn>fPMax)          
176       {
177           ;
178       } else {
179           p[0]=prwn*TMath::Sin(fTheta)*TMath::Cos(fPhi);
180           p[1]=prwn*TMath::Sin(fTheta)*TMath::Sin(fPhi);      
181           p[2]=prwn*TMath::Cos(fTheta);
182           
183           if(fVertexSmear==kPerTrack) {
184               gMC->Rndm(random,6);
185               for (j=0;j<3;j++) {
186                   origin[j]=fOrigin[j]
187                       +fOsigma[j]*TMath::Cos(2*random[2*j]*TMath::Pi())*
188                       TMath::Sqrt(-2*TMath::Log(random[2*j+1]));
189               }
190           }
191           gAlice->SetTrack(fTrackIt,-1,fIdpart,p,origin,polar,0,"Primary",nt);
192       }
193       fNcurrent++;
194       nb = (Int_t)h2->GetEvent(fNcurrent); 
195   }
196  
197     TFile *pFile=0;
198 // Get AliRun object or create it 
199     if (!gAlice) {
200         gAlice = (AliRun*)pFile->Get("gAlice");
201         if (gAlice) printf("AliRun object found on file\n");
202         if (!gAlice) gAlice = new AliRun("gAlice","Alice test program");
203     }
204     TTree *fAli=gAlice->TreeK();
205     if (fAli) pFile =fAli->GetCurrentFile();
206     pFile->cd();
207 }
208
209
210 AliGenExtFile& AliGenExtFile::operator=(const  AliGenExtFile& rhs)
211 {
212 // Assignment operator
213     return *this;
214 }
215
216
217
218
219
220
221
222