]> git.uio.no Git - u/mrichter/AliRoot.git/blame - EVGEN/AliGenFLUKAsource.cxx
Default for decays is kAll.
[u/mrichter/AliRoot.git] / EVGEN / AliGenFLUKAsource.cxx
CommitLineData
4c039060 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$
a99cf51f 18Revision 1.14 2001/03/21 11:28:20 morsch
19Use enum constants for particle selection.
20
09d89034 21Revision 1.13 2000/12/21 16:24:06 morsch
22Coding convention clean-up
23
675e9664 24Revision 1.12 2000/11/30 07:12:50 alibrary
25Introducing new Rndm and QA classes
26
65fb704d 27Revision 1.11 2000/06/14 15:20:40 morsch
28Include clean-up (IH)
29
5c3fd7ea 30Revision 1.10 2000/06/09 20:31:34 morsch
31All coding rule violations except RS3 corrected
32
f87cfe57 33Revision 1.9 2000/03/07 13:52:54 morsch
34static Int_t irwn=0;
35
e7bafa6d 36Revision 1.8 2000/02/14 14:49:38 morsch
37Correct particle type for gamma and neutrons
38More consistent calculation of momentum from kin. energy and mass
39
02769dba 40Revision 1.7 1999/11/03 17:43:20 fca
41New version from G.Martinez & A.Morsch
42
886b6f73 43Revision 1.6 1999/09/29 09:24:12 fca
44Introduction of the Copyright and cvs Log
45
4c039060 46*/
47
675e9664 48
49
50// Read background particles from a FLUKA boundary source file
51// This is a very special generator that works for background studies for the muon-spectrometer.
52// The input files come from FLUKA simulations.
53// Files can be chained.
54// Author: andreas.morsch@cern.ch
55
fe4da5cc 56#include "AliGenFLUKAsource.h"
fe4da5cc 57#include "AliMC.h"
58#include "AliRun.h"
c5970876 59#include "AliPDG.h"
5c3fd7ea 60
61
fe4da5cc 62#include <TFile.h>
63#include <TTree.h>
f87cfe57 64#include <TChain.h>
fe4da5cc 65#include <stdlib.h>
66 ClassImp(AliGenFLUKAsource)
67 AliGenFLUKAsource::AliGenFLUKAsource()
68 :AliGenerator(-1)
69{
f87cfe57 70 // Constructor
fe4da5cc 71 fName="FLUKA";
72 fTitle="FLUKA Boundary Source";
73 // Read in all particle types by default
74 fIkine=6;
75 // Set maximum admitted age of particles to 1.e-05 by default
76 fAgeMax=1.e-05;
77 // Do not add weight
78 fAddWeight=1.;
79 // Shift the z-coordinate of the impact point by 4.5 cm only if it reads
80 // from specific boundary source to the chamber (fZshift=4.5;),else there
81 // is no need to shift as it reads boundary source for the whole volume of
82 // the Muon Arm; the default value corresponds to boundary source for the
83 // whole volume of the MUON Arm
84 fZshift=0;
85 // Set the default file
09d89034 86 fFileName="";
fe4da5cc 87
88 fTreeFluka=0;
886b6f73 89 fTreeChain = new TChain("h1");
fe4da5cc 90//
91// Read all particles
92 fNpart=-1;
fe4da5cc 93}
94
95AliGenFLUKAsource::AliGenFLUKAsource(Int_t npart)
96 :AliGenerator(npart)
97{
f87cfe57 98 // Constructor
675e9664 99 fName = "FLUKA";
100 fTitle = "FLUKA Boundary Source";
fe4da5cc 101 // Read in all particle types by default
102 fIkine=6;
103 // Set maximum admitted age of particles to 1.e-05 by default
104 fAgeMax=1.e-05;
105 // Do not add weight
106 fAddWeight=1.;
107 // Shift the z-coordinate of the impact point by 4.5 cm only if it reads
108 // from specific boundary source to the chamber (fZshift=4.5;),else there
109 // is no need to shift as it reads boundary source for the whole volume of
110 // the Muon Arm; the default value corresponds to boundary source for the
111 // whole volume of the MUON Arm
112 fZshift=0;
113 // Set the default file
09d89034 114 fFileName="";
fe4da5cc 115
116 fTreeFluka=0;
886b6f73 117 fTreeChain = new TChain("h1");
f87cfe57 118 fSourceId=-1;
fe4da5cc 119}
120
f87cfe57 121AliGenFLUKAsource::AliGenFLUKAsource(const AliGenFLUKAsource & FLUKAsource)
122{
123// copy constructor
124}
125
126
fe4da5cc 127//____________________________________________________________
128AliGenFLUKAsource::~AliGenFLUKAsource()
129{
f87cfe57 130// Destructor
886b6f73 131 if (fTreeFluka) delete fTreeFluka;
132 if (fTreeChain) delete fTreeChain;
fe4da5cc 133}
134
886b6f73 135void AliGenFLUKAsource::AddFile(const Text_t *filename)
fe4da5cc 136{
f87cfe57 137// Add a file to the chain
886b6f73 138 fTreeChain->Add(filename);
fe4da5cc 139
886b6f73 140}
141
fe4da5cc 142
886b6f73 143//____________________________________________________________
144void AliGenFLUKAsource::FlukaInit()
145{
f87cfe57 146// Set branch addresses of data entries
886b6f73 147 TChain *h2=fTreeChain;
fe4da5cc 148//Set branch addresses
f87cfe57 149 h2->SetBranchAddress("Ip",&fIp);
150 h2->SetBranchAddress("Ipp",&fIpp);
151 h2->SetBranchAddress("Xi",&fXi);
152 h2->SetBranchAddress("Yi",&fYi);
153 h2->SetBranchAddress("Zi",&fZi);
154 h2->SetBranchAddress("Px",&fPx);
155 h2->SetBranchAddress("Py",&fPy);
156 h2->SetBranchAddress("Pz",&fPz);
157 h2->SetBranchAddress("Ekin",&fEkin);
158 h2->SetBranchAddress("Zv",&fZv);
159 h2->SetBranchAddress("Rv",&fRv);
160 h2->SetBranchAddress("Itra",&fItra);
161 h2->SetBranchAddress("Igas",&fIgas);
162 h2->SetBranchAddress("Wgt",&fWgt);
163 h2->SetBranchAddress("Etag",&fEtag);
164 h2->SetBranchAddress("Ptg",&fPtg);
165 h2->SetBranchAddress("Age",&fAge);
fe4da5cc 166}
167
168//____________________________________________________________
169void AliGenFLUKAsource::Generate()
f87cfe57 170{
171// Generate one event
886b6f73 172 AliMC* gMC = AliMC::GetMC();
fe4da5cc 173
f87cfe57 174 const Int_t kIfluge[28]={kProton, kProtonBar, kElectron, kPositron,
1578254f 175 kNuE, kNuEBar, kGamma, kNeutron, kNeutronBar,
176 kMuonPlus, kMuonMinus, kK0Long , kPiPlus, kPiMinus,
177 kKPlus, kKMinus, kLambda0, kLambda0Bar, kK0Short,
178 kSigmaMinus, kSigmaPlus, kSigma0, kPi0, kK0, kK0Bar,
179 0,kNuMu,kNuMuBar};
886b6f73 180 Float_t polar[3]= {0,0,0};
fe4da5cc 181 //
886b6f73 182 Float_t origin[3];
183 Float_t p[3];
184 Float_t prwn;
185 Float_t wgt, fwgt;
186 Float_t phi;
187 char name[100];
188 Float_t amass, charge, tlife;
189 Int_t itrtyp;
190 Int_t iwgt;
191 Int_t i, j, part, nt;
e7bafa6d 192 static Int_t irwn=0;
886b6f73 193 //
194 Float_t random[2];
fe4da5cc 195 //
886b6f73 196 FlukaInit();
197 TChain *h2=fTreeChain;
198 Int_t nentries = (Int_t) h2->GetEntries();
199 if (fNpart == -1) fNpart=Int_t(nentries*fFrac);
09d89034 200
201
fe4da5cc 202 // loop over number of particles
886b6f73 203 Int_t nb=0;
cfce8870 204 Int_t ev=gMC->CurrentEvent();
886b6f73 205 for (i=0; i<fNpart;i++) {
206 Int_t entry=fNpart*(ev-1)+i;
207 nb = (Int_t)h2->GetEvent(entry);
208 if (irwn > nentries) {
209 printf("No more entries in the FLUKA boundary source file\n");
f87cfe57 210 TFile *pFile=0;
886b6f73 211 // Get AliRun object or create it
212 if (!gAlice) {
f87cfe57 213 gAlice = (AliRun*)pFile->Get("gAlice");
886b6f73 214 if (gAlice) printf("AliRun object found on file\n");
215 if (!gAlice) gAlice = new AliRun("gAlice","Alice test program");
216 }
217 TTree *fAli=gAlice->TreeK();
f87cfe57 218 if (fAli) pFile =fAli->GetCurrentFile();
219 pFile->cd();
886b6f73 220 printf("Generate - I'm out \n");
221 return;
222 }
09d89034 223
224 Int_t ifip = Int_t(fIp);
225
f87cfe57 226
227 if (fSourceId != -1 && fIgas !=fSourceId) {
228 irwn++;
229 continue;
230 }
231
09d89034 232 if (ifip > 28 || ifip < 0) {
886b6f73 233 irwn++;
234 continue;
235 }
236
09d89034 237 if ((ifip != fIkine && fIkine != kAll && fIkine != kCharged
238 && fIkine != 10) || fAge > fAgeMax){
886b6f73 239 irwn++;
240 continue;
09d89034 241 } else if (fIkine == kCharged) {
242 if (ifip == 7 || ifip == 8 || fAge > fAgeMax) {
886b6f73 243 irwn++;
244 continue;
245 }
09d89034 246 } else if (fIkine == kNoNeutron) {
247 if (ifip == 8 || fAge > fAgeMax) {
886b6f73 248 irwn++;
249 continue;
250 }
251 }
fe4da5cc 252
fe4da5cc 253
886b6f73 254 irwn++;
02769dba 255//
09d89034 256// PDG code from FLUKA particle type (ifip)
257 part=kIfluge[int(ifip)-1];
02769dba 258//
259// Calculate momentum from kinetic energy and mass of the particle
260 gMC->Gfpart(part, name, itrtyp,
261 amass, charge, tlife);
f87cfe57 262 prwn=fEkin*sqrt(1. + 2.*amass/fEkin);
02769dba 263
264
f87cfe57 265 origin[0]=fYi;
266 origin[1]=fXi;
267 origin[2]=fZi;
886b6f73 268
f87cfe57 269 p[0]=fPy*prwn;
270 p[1]=fPx*prwn;
271 p[2]=fPz*prwn;
fe4da5cc 272
886b6f73 273 //handle particle weight correctly
f87cfe57 274 wgt = (part == 13) ? fWgt*fAddWeight : fWgt;
886b6f73 275 iwgt=Int_t(wgt);
276 fwgt=wgt-Float_t(iwgt);
65fb704d 277 Rndm(random,2);
886b6f73 278 if (random[0] < fwgt) iwgt++;
279 if (part==1 && iwgt>100) iwgt=100;
280 Int_t nstack=0;
281 for (j=0; j<iwgt; j++) {
a99cf51f 282 SetTrack(fTrackIt,-1,part,p,origin,polar,fAge,kPPrimary,nt);
65fb704d 283 Rndm(random,2);
886b6f73 284 phi=2*random[1]*TMath::Pi();
285 Float_t pn1=p[0]*TMath::Sin(phi) - p[1]*TMath::Cos(phi);
286 Float_t pn2=p[0]*TMath::Cos(phi) + p[1]*TMath::Sin(phi);
287 p[0]=pn1;
288 p[1]=pn2;
289 Float_t on1=origin[0]*TMath::Sin(phi)-origin[1]*TMath::Cos(phi);
290 Float_t on2=origin[0]*TMath::Cos(phi)+origin[1]*TMath::Sin(phi);
291 origin[0]=on1;
292 origin[1]=on2;
293 nstack++;
294 }
295 if (nstack == 0) continue;
fe4da5cc 296 }
886b6f73 297
f87cfe57 298 TFile *pFile=0;
fe4da5cc 299// Get AliRun object or create it
300 if (!gAlice) {
f87cfe57 301 gAlice = (AliRun*)pFile->Get("gAlice");
fe4da5cc 302 if (gAlice) printf("AliRun object found on file\n");
303 if (!gAlice) gAlice = new AliRun("gAlice","Alice test program");
304 }
305 TTree *fAli=gAlice->TreeK();
f87cfe57 306 if (fAli) pFile =fAli->GetCurrentFile();
307 pFile->cd();
fe4da5cc 308}
309
310
f87cfe57 311AliGenFLUKAsource& AliGenFLUKAsource::operator=(const AliGenFLUKAsource& rhs)
312{
313// Assignment operator
314 return *this;
315}
fe4da5cc 316
317
318
319
320
321