14ab5373 |
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 | |
88cb7938 |
16 | /* $Id$ */ |
14ab5373 |
17 | |
18 | // Read background particles from a boundary source |
19 | // Very specialized generator to simulate background from beam halo. |
20 | // The input file is a text file specially prepared |
21 | // for this purpose. |
22 | // Author: andreas.morsch@cern.ch |
23 | |
116cbefd |
24 | #include <stdlib.h> |
14ab5373 |
25 | |
26 | #include <TDatabasePDG.h> |
116cbefd |
27 | #include <TPDGCode.h> |
eac2c8a3 |
28 | #include <TSystem.h> |
116cbefd |
29 | |
30 | #include "AliGenHaloProtvino.h" |
31 | #include "AliRun.h" |
14ab5373 |
32 | |
33 | ClassImp(AliGenHaloProtvino) |
34 | AliGenHaloProtvino::AliGenHaloProtvino() |
35 | :AliGenerator(-1) |
36 | { |
37 | // Constructor |
154d524b |
38 | |
39 | fName = "HaloProtvino"; |
40 | fTitle = "Halo from LHC Tunnel"; |
14ab5373 |
41 | // |
42 | // Read all particles |
154d524b |
43 | fNpart = -1; |
44 | fFile = 0; |
45 | fSide = 1; |
eac2c8a3 |
46 | // |
47 | SetRunPeriod(); |
7f98ef21 |
48 | SetTimePerEvent(); |
49 | SetAnalog(0); |
14ab5373 |
50 | } |
51 | |
52 | AliGenHaloProtvino::AliGenHaloProtvino(Int_t npart) |
53 | :AliGenerator(npart) |
54 | { |
55 | // Constructor |
154d524b |
56 | fName = "Halo"; |
57 | fTitle= "Halo from LHC Tunnel"; |
14ab5373 |
58 | // |
154d524b |
59 | fNpart = npart; |
60 | fFile = 0; |
61 | fSide = 1; |
eac2c8a3 |
62 | // |
63 | SetRunPeriod(); |
7f98ef21 |
64 | SetTimePerEvent(); |
65 | SetAnalog(0); |
14ab5373 |
66 | } |
67 | |
68 | AliGenHaloProtvino::AliGenHaloProtvino(const AliGenHaloProtvino & HaloProtvino) |
69 | { |
70 | // copy constructor |
71 | } |
72 | |
73 | |
74 | //____________________________________________________________ |
75 | AliGenHaloProtvino::~AliGenHaloProtvino() |
76 | { |
77 | // Destructor |
78 | } |
79 | |
80 | //____________________________________________________________ |
81 | void AliGenHaloProtvino::Init() |
82 | { |
83 | // Initialisation |
154d524b |
84 | fFile = fopen(fFileName,"r"); |
85 | if (fFile) { |
86 | printf("\n File %s opened for reading, %p ! \n ", fFileName.Data(), fFile); |
87 | } else { |
88 | printf("\n Opening of file %s failed, %p ! \n ", fFileName.Data(), fFile); |
89 | } |
eac2c8a3 |
90 | // |
91 | // |
92 | // |
93 | // Read file with gas pressure values |
94 | char *name; |
95 | |
96 | name = gSystem->ExpandPathName("$(ALICE_ROOT)/LHC/gasPressure.dat" ); |
97 | FILE* file = fopen(name, "r"); |
98 | Float_t z; |
99 | Int_t i, j; |
100 | |
101 | // |
102 | // Ring 1 |
103 | // |
7f98ef21 |
104 | for (i = 0; i < 21; i++) |
eac2c8a3 |
105 | { |
106 | fscanf(file, "%f %f %f %f %f %f", &z, |
107 | &fG1[i][0], &fG1[i][1], &fG1[i][2], &fG1[i][3], &fG1[i][4]); |
108 | if (i > 0) { |
109 | fZ1[i] = fZ1[i-1] + z; |
110 | } else { |
111 | fZ1[i] = 20.; |
112 | } |
113 | } |
114 | // |
115 | // Ring 2 |
116 | // |
117 | for (i = 0; i < 21; i++) |
118 | { |
119 | fscanf(file, "%f %f %f %f %f %f", &z, |
120 | &fG2[i][0], &fG2[i][1], &fG2[i][2], &fG2[i][3], &fG2[i][4]); |
121 | if (i > 0) { |
122 | fZ2[i] = fZ2[i-1] + z; |
123 | } else { |
124 | fZ2[i] = 20.; |
125 | } |
126 | } |
127 | // |
128 | // Transform into interaction rates |
129 | // |
130 | const Float_t crossSection = 0.094e-28; // m^2 |
718e5512 |
131 | Float_t pFlux[5] = {0.2, 0.2, 0.3, 0.3, 1.0}; |
eac2c8a3 |
132 | |
133 | for (j = 0; j < 5; j++) { |
718e5512 |
134 | pFlux[j] *= 1.e11/25.e-9; |
eac2c8a3 |
135 | for (i = 0; i < 21; i++) |
136 | { |
718e5512 |
137 | fG1[i][j] = fG1[i][j] * crossSection * pFlux[j]; // 1/m/s |
138 | fG2[i][j] = fG2[i][j] * crossSection * pFlux[j]; // 1/m/s |
eac2c8a3 |
139 | } |
140 | } |
7f98ef21 |
141 | |
eac2c8a3 |
142 | |
7f98ef21 |
143 | Float_t sum1 = 0.; |
144 | Float_t sum2 = 0.; |
eac2c8a3 |
145 | |
7f98ef21 |
146 | for (Int_t i = 0; i < 300; i++) { |
eac2c8a3 |
147 | Float_t z = 20.+i*1.; |
148 | z*=100; |
7f98ef21 |
149 | Float_t wgt1 = GassPressureWeight(z); |
150 | Float_t wgt2 = GassPressureWeight(-z); |
151 | // printf("weight: %f %f %f %f %f \n", z, wgt1, wgt2, fZ1[20], fZ2[20]); |
152 | sum1 += wgt1; |
153 | sum2 += wgt2; |
eac2c8a3 |
154 | } |
7f98ef21 |
155 | sum1/=250.; |
156 | sum2/=250.; |
157 | printf("\n %f %f \n \n", sum1, sum2); |
14ab5373 |
158 | } |
159 | |
160 | //____________________________________________________________ |
161 | void AliGenHaloProtvino::Generate() |
162 | { |
163 | // Generate from input file |
14ab5373 |
164 | |
165 | Float_t polar[3]= {0,0,0}; |
166 | Float_t origin[3]; |
167 | Float_t p[3], p0; |
7f98ef21 |
168 | Float_t tz, txy; |
14ab5373 |
169 | Float_t amass; |
14ab5373 |
170 | // |
7f98ef21 |
171 | Int_t ncols, nt; |
05df024b |
172 | Int_t nskip = 0; |
154d524b |
173 | Int_t nread = 0; |
7f98ef21 |
174 | |
175 | Float_t* zPrimary = new Float_t [fNpart]; |
176 | Int_t * inuc = new Int_t [fNpart]; |
177 | Int_t * ipart = new Int_t [fNpart]; |
178 | Float_t* wgt = new Float_t [fNpart]; |
179 | Float_t* ekin = new Float_t [fNpart]; |
180 | Float_t* vx = new Float_t [fNpart]; |
181 | Float_t* vy = new Float_t [fNpart]; |
182 | Float_t* tx = new Float_t [fNpart]; |
183 | Float_t* ty = new Float_t [fNpart]; |
184 | |
185 | Float_t zVertexOld = -1.e10; |
186 | Int_t nInt = 0; // Counts number of interactions |
187 | |
14ab5373 |
188 | while(1) { |
7f98ef21 |
189 | // |
190 | // Load event into array |
191 | // |
154d524b |
192 | ncols = fscanf(fFile,"%f %d %d %f %f %f %f %f %f", |
7f98ef21 |
193 | &zPrimary[nread], &inuc[nread], &ipart[nread], &wgt[nread], |
194 | &ekin[nread], &vx[nread], &vy[nread], |
195 | &tx[nread], &ty[nread]); |
05df024b |
196 | |
14ab5373 |
197 | if (ncols < 0) break; |
7f98ef21 |
198 | // Skip fNskip events |
05df024b |
199 | nskip++; |
200 | if (fNpart !=-1 && nskip <= fNskip) continue; |
7f98ef21 |
201 | // Count interactions |
202 | if (zPrimary[nread] != zVertexOld) { |
203 | nInt++; |
204 | zVertexOld = zPrimary[nread]; |
205 | } |
206 | // Count tracks |
05df024b |
207 | nread++; |
14ab5373 |
208 | if (fNpart !=-1 && nread > fNpart) break; |
7f98ef21 |
209 | } |
210 | // |
211 | // Mean time between interactions |
212 | // |
213 | Float_t dT = fTimePerEvent/nInt; // sec |
214 | Float_t t = 0; // sec |
215 | |
216 | // |
217 | // Loop over primaries |
218 | // |
219 | zVertexOld = -1.e10; |
220 | Double_t arg = 0.; |
221 | |
222 | for (Int_t nprim = 0; nprim < fNpart; nprim++) |
223 | { |
224 | amass = TDatabasePDG::Instance()->GetParticle(ipart[nprim])->Mass(); |
14ab5373 |
225 | |
226 | // |
227 | // Momentum vector |
228 | // |
7f98ef21 |
229 | p0=sqrt(ekin[nprim]*ekin[nprim] + 2.*amass*ekin[nprim]); |
14ab5373 |
230 | |
7f98ef21 |
231 | txy=TMath::Sqrt(tx[nprim]*tx[nprim]+ty[nprim]*ty[nprim]); |
14ab5373 |
232 | if (txy == 1.) { |
233 | tz=0; |
234 | } else { |
235 | tz=-TMath::Sqrt(1.-txy); |
236 | } |
237 | |
7f98ef21 |
238 | p[0] = p0*tx[nprim]; |
239 | p[1] = p0*ty[nprim]; |
240 | p[2] =-p0*tz; |
05df024b |
241 | |
7f98ef21 |
242 | origin[0] = vx[nprim]; |
243 | origin[1] = vy[nprim]; |
05df024b |
244 | origin[2] = -2196.5; |
14ab5373 |
245 | |
246 | // |
247 | // |
248 | // Particle weight |
249 | |
154d524b |
250 | Float_t originP[3] = {0., 0., 0.}; |
7f98ef21 |
251 | originP[2] = zPrimary[nprim]; |
154d524b |
252 | |
253 | Float_t pP[3] = {0., 0., 0.}; |
254 | Int_t ntP; |
255 | |
256 | if (fSide == -1) { |
7f98ef21 |
257 | originP[2] = -zPrimary[nprim]; |
154d524b |
258 | origin[2] = -origin[2]; |
259 | p[2] = -p[2]; |
260 | } |
05df024b |
261 | |
14ab5373 |
262 | // |
7f98ef21 |
263 | // Time |
264 | // |
265 | if (zPrimary[nprim] != zVertexOld) { |
266 | while(arg==0.) arg = gRandom->Rndm(); |
267 | t -= dT*TMath::Log(arg); // (sec) |
268 | zVertexOld = zPrimary[nprim]; |
269 | } |
270 | |
271 | // Get statistical weight according to local gas-pressure |
272 | // |
273 | fParentWeight=wgt[nprim]*GassPressureWeight(zPrimary[nprim]); |
274 | |
275 | if (!fAnalog || gRandom->Rndm() < fParentWeight) { |
276 | // Pass parent particle |
277 | // |
642f15cf |
278 | PushTrack(0,-1,kProton,pP,originP,polar,t,kPNoProcess,ntP, fParentWeight); |
7f98ef21 |
279 | KeepTrack(ntP); |
642f15cf |
280 | PushTrack(fTrackIt,ntP,ipart[nprim],p,origin,polar,t,kPNoProcess,nt,fParentWeight); |
7f98ef21 |
281 | } |
282 | |
283 | // |
284 | // Both sides are considered |
285 | // |
154d524b |
286 | |
05df024b |
287 | if (fSide > 1) { |
7f98ef21 |
288 | fParentWeight=wgt[nprim]*GassPressureWeight(-zPrimary[nprim]); |
289 | if (!fAnalog || gRandom->Rndm() < fParentWeight) { |
290 | origin[2] = -origin[2]; |
291 | originP[2] = -originP[2]; |
292 | p[2]=-p[2]; |
642f15cf |
293 | PushTrack(0,-1,kProton,pP,originP,polar,t,kPNoProcess,ntP, fParentWeight); |
7f98ef21 |
294 | KeepTrack(ntP); |
642f15cf |
295 | PushTrack(fTrackIt,ntP,ipart[nprim],p,origin,polar,t,kPNoProcess,nt,fParentWeight); |
7f98ef21 |
296 | } |
05df024b |
297 | } |
5bfdba97 |
298 | SetHighWaterMark(nt); |
14ab5373 |
299 | } |
7f98ef21 |
300 | delete [] zPrimary; |
301 | delete [] inuc; |
302 | delete [] ipart; |
303 | delete [] wgt; |
304 | delete [] ekin; |
305 | delete [] vx; |
306 | delete [] vy; |
307 | delete [] tx; |
308 | delete [] ty; |
14ab5373 |
309 | } |
310 | |
311 | |
312 | AliGenHaloProtvino& AliGenHaloProtvino::operator=(const AliGenHaloProtvino& rhs) |
313 | { |
314 | // Assignment operator |
315 | return *this; |
316 | } |
317 | |
318 | |
eac2c8a3 |
319 | |
14ab5373 |
320 | Float_t AliGenHaloProtvino::GassPressureWeight(Float_t zPrimary) |
321 | { |
eac2c8a3 |
322 | // |
323 | // Return z-dependent gasspressure weight = interaction rate [1/m/s]. |
324 | // |
325 | Float_t weight = 0.; |
326 | zPrimary/=100.; // m |
327 | Float_t zAbs = TMath::Abs(zPrimary); |
7f98ef21 |
328 | if (zPrimary > 0.) |
eac2c8a3 |
329 | { |
7f98ef21 |
330 | if (zAbs > fZ1[20]) { |
eac2c8a3 |
331 | weight = 2.e4; |
332 | } else { |
7f98ef21 |
333 | for (Int_t i = 1; i < 21; i++) { |
eac2c8a3 |
334 | if (zAbs < fZ1[i]) { |
335 | weight = fG1[i][fRunPeriod]; |
336 | break; |
337 | } |
338 | } |
339 | } |
340 | } else { |
341 | if (zAbs > fZ2[20]) { |
342 | weight = 2.e4; |
343 | } else { |
344 | for (Int_t i = 1; i < 21; i++) { |
345 | if (zAbs < fZ2[i]) { |
346 | weight = fG2[i][fRunPeriod]; |
347 | break; |
348 | } |
349 | } |
350 | } |
351 | } |
eac2c8a3 |
352 | return weight; |
14ab5373 |
353 | } |
354 | |
355 | /* |
356 | # Title: README file for the sources of IR8 machine induced background |
357 | # Author: Vadim Talanov <Vadim.Talanov@cern.ch> |
358 | # Modified: 12-12-2000 |
359 | |
360 | 0. Overview |
361 | |
362 | There are three files, named ring.one.beta.[01,10,50].m, which |
363 | contain the lists of background particles, induced by proton losses |
364 | upstream of IP8 in the LHC ring one, for the beta* values of 1, 10 |
365 | and 50 m, respectively. |
366 | |
367 | 1. File contents |
368 | |
369 | Each line in the files contains the coordinates of particle track |
370 | crossing with the infinite plane, positioned at z=-1m, together with |
371 | the physical properties of corresponding particle, namely: |
372 | |
373 | S - S coordinate of the primary interaction vertex, cm; |
374 | N - type of the gas nuclei at interaction, 1 is H, 2 - C and 3 - O; |
375 | I - particle ID in PDG particle numbering scheme; |
376 | W - particle weight; |
377 | E - particle kinetic energy, GeV; |
378 | X - x coordinate of the crossing point, cm; |
379 | Y - y coordinate of the crossing point, cm; |
380 | Dx - x direction cosine; |
381 | Dy - y direction cosine. |
382 | |
383 | 2. Normalisation |
384 | |
385 | Each file is given per unity of linear density of proton inelastic |
386 | interactions with the gas nuclei, [1 inelastic interaction/m]. |
387 | |
388 | # ~/vtalanov/public/README.mib: the end. |
389 | |
390 | */ |
391 | |
392 | |
393 | |
394 | |