Problem in memory management fixed
[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.24  2002/10/21 09:20:51  alibrary
19 Introduce Riostream.h and remove unused variables
20
21 Revision 1.23  2002/10/14 14:55:35  hristov
22 Merging the VirtualMC branch to the main development branch (HEAD)
23
24 Revision 1.20.4.1  2002/06/10 14:57:41  hristov
25 Merged with v3-08-02
26
27 Revision 1.22  2002/05/15 11:59:49  morsch
28 CdEventFile() method added.
29
30 Revision 1.21  2002/04/26 10:39:31  morsch
31 AliGenExtFile derives from AliGenMC. Generate() uses methods from AliGenMC (N. Carrer)
32
33 Revision 1.20  2002/03/22 09:43:28  morsch
34 Don't delete the TParticle.
35
36 Revision 1.19  2002/02/08 16:50:50  morsch
37 Add name and title in constructor.
38
39 Revision 1.18  2001/11/12 14:31:00  morsch
40 Memory leaks fixed. (M. Bondila)
41
42 Revision 1.17  2001/11/09 09:12:58  morsch
43 Generalization by using AliGenReader object to read particles from file.
44
45 Revision 1.16  2001/07/27 17:09:35  morsch
46 Use local SetTrack, KeepTrack and SetHighWaterMark methods
47 to delegate either to local stack or to stack owned by AliRun.
48 (Piotr Skowronski, A.M.)
49
50 Revision 1.15  2001/01/23 13:29:37  morsch
51 Add method SetParticleCode and enum type Code_t to handle both PDG (new ntuples)
52 and GEANT3 codes (old ntuples) in input file.
53
54 Revision 1.14  2000/12/21 16:24:06  morsch
55 Coding convention clean-up
56
57 Revision 1.13  2000/11/30 07:12:50  alibrary
58 Introducing new Rndm and QA classes
59
60 Revision 1.12  2000/10/27 13:54:45  morsch
61 Remove explicite reference to input file from constuctor.
62
63 Revision 1.11  2000/10/02 21:28:06  fca
64 Removal of useless dependecies via forward declarations
65
66 Revision 1.10  2000/07/11 18:24:55  fca
67 Coding convention corrections + few minor bug fixes
68
69 Revision 1.9  2000/06/14 15:20:09  morsch
70 Include clean-up (IH)
71
72 Revision 1.8  2000/06/09 20:36:44  morsch
73 All coding rule violations except RS3 corrected
74
75 Revision 1.7  2000/02/16 14:56:27  morsch
76 Convert geant particle code into pdg code before putting particle on the stack.
77
78 Revision 1.6  1999/11/09 07:38:48  fca
79 Changes for compatibility with version 2.23 of ROOT
80
81 Revision 1.5  1999/09/29 09:24:12  fca
82 Introduction of the Copyright and cvs Log
83
84 */
85
86
87 // Event generator that using an instance of type AliGenReader
88 // reads particles from a file and applies cuts. 
89
90 #include <Riostream.h>
91
92 #include "AliGenExtFile.h"
93 #include "AliRun.h"
94
95 #include <TParticle.h>
96 #include <TFile.h>
97 #include <TTree.h>
98
99
100  ClassImp(AliGenExtFile)
101
102 AliGenExtFile::AliGenExtFile()
103   :AliGenMC()
104 {
105 //  Constructor
106 //
107 //  Read all particles
108     fNpart  =- 1;
109     fReader =  0;
110 }
111
112 AliGenExtFile::AliGenExtFile(Int_t npart)
113     :AliGenMC(npart)
114 {
115 //  Constructor
116     fName   = "ExtFile";
117     fTitle  = "Primaries from ext. File";
118     fReader = 0;
119 }
120
121 AliGenExtFile::AliGenExtFile(const AliGenExtFile & ExtFile)
122 {
123 // copy constructor
124 }
125 //____________________________________________________________
126 AliGenExtFile::~AliGenExtFile()
127 {
128 // Destructor
129     delete fReader;
130 }
131
132 //___________________________________________________________
133 void AliGenExtFile::Init()
134 {
135 // Initialize
136     if (fReader) fReader->Init();
137 }
138
139     
140 void AliGenExtFile::Generate()
141 {
142 // Generate particles
143
144   Float_t polar[3]  = {0,0,0};
145   //
146   Float_t origin[3] = {0,0,0};
147   Float_t p[3];
148   Float_t random[6];
149   Int_t i, j, nt;
150   //
151   for (j=0;j<3;j++) origin[j]=fOrigin[j];
152   if(fVertexSmear == kPerTrack) {
153     Rndm(random,6);
154     for (j = 0; j < 3; j++) {
155         origin[j] += fOsigma[j]*TMath::Cos(2*random[2*j]*TMath::Pi())*
156             TMath::Sqrt(-2*TMath::Log(random[2*j+1]));
157     }
158   }
159
160   while(1) {
161     Int_t nTracks = fReader->NextEvent();       
162     if (nTracks == 0) {
163       // printf("\n No more events !!! !\n");
164       Warning("AliGenExtFile::Generate","\nNo more events in external file!!!\nLast event may be empty or incomplete.\n");
165       return;
166     }
167
168     //
169     // Particle selection loop
170     //
171     // The selction criterium for the external file generator is as follows:
172     //
173     // 1) All tracs are subjects to the cuts defined by AliGenerator, i.e.
174     //    fThetaMin, fThetaMax, fPhiMin, fPhiMax, fPMin, fPMax, fPtMin, fPtMax,
175     //    fYMin, fYMax.
176     //    If the particle does not satisfy these cuts, it is not put on the
177     //    stack.
178     // 2) If fCutOnChild and some specific child is selected (e.g. if
179     //    fForceDecay==kSemiElectronic) the event is rejected if NOT EVEN ONE
180     //    child falls into the child-cuts.
181     if(fCutOnChild) {
182       // Count the selected children
183       Int_t nSelected = 0;
184
185       for (i = 0; i < nTracks; i++) {
186         TParticle* iparticle = fReader->NextParticle();
187         Int_t kf = CheckPDGCode(iparticle->GetPdgCode());
188         kf = TMath::Abs(kf);
189         if (ChildSelected(kf) && KinematicSelection(iparticle, 1)) {
190           nSelected++;
191         }
192       }
193       if (!nSelected) continue;    // No particle selected:  Go to next event
194       fReader->RewindEvent();
195     }
196
197     //
198     // Stack filling loop
199     //
200     for (i = 0; i < nTracks; i++) {
201
202       TParticle* iparticle = fReader->NextParticle();
203       if (!KinematicSelection(iparticle,0)) {
204         Double_t  pz   = iparticle->Pz();
205         Double_t  e    = iparticle->Energy();
206         Double_t  y;
207         if ((e-pz) == 0) {
208           y = 20.;
209         } else if ((e+pz) == 0.) {
210           y = -20.;
211         } else {
212           y = 0.5*TMath::Log((e+pz)/(e-pz));      
213         }
214         printf("\n Not selected %d %f %f %f %f %f", i,
215                iparticle->Theta(),
216                iparticle->Phi(),
217                iparticle->P(),
218                iparticle->Pt(),
219                y);
220         //PH    delete iparticle;
221         continue;
222       }
223       p[0] = iparticle->Px();
224       p[1] = iparticle->Py();
225       p[2] = iparticle->Pz();
226       Int_t idpart = iparticle->GetPdgCode();
227       if(fVertexSmear==kPerTrack) {
228           Rndm(random,6);
229           for (j = 0; j < 3; j++) {
230               origin[j]=fOrigin[j]
231                   +fOsigma[j]*TMath::Cos(2*random[2*j]*TMath::Pi())*
232                   TMath::Sqrt(-2*TMath::Log(random[2*j+1]));
233           }
234       }
235       Int_t decayed = iparticle->GetFirstDaughter();
236       Int_t doTracking = fTrackIt && (decayed < 0) &&
237                          (TMath::Abs(idpart) > 10);
238       // printf("*** pdg, first daughter, trk = %d, %d, %d\n",
239       //   idpart,decayed, doTracking);
240       SetTrack(doTracking,-1,idpart,p,origin,polar,0,kPPrimary,nt);
241       KeepTrack(nt);
242     } // track loop
243
244     break;
245
246   } // event loop
247
248   SetHighWaterMark(nt);
249   CdEventFile();
250 }
251
252 void AliGenExtFile::CdEventFile()
253 {
254 // CD back to the event file
255   TFile *pFile=0;
256   if (!gAlice) {
257       gAlice = (AliRun*)pFile->Get("gAlice");
258       if (gAlice) printf("AliRun object found on file\n");
259       if (!gAlice) gAlice = new AliRun("gAlice","Alice test program");
260   }
261   TTree *fAli=gAlice->TreeK();
262   if (fAli) pFile =fAli->GetCurrentFile();
263   pFile->cd();    
264 }
265
266
267 AliGenExtFile& AliGenExtFile::operator=(const  AliGenExtFile& rhs)
268 {
269 // Assignment operator
270     return *this;
271 }
272
273
274
275
276
277
278
279