COSMOS v7.655  COSMOSv7655
(AirShowerMC)
interface.f
Go to the documentation of this file.
1 ! FNODATDEF if undef, ascii main output; else binary output:
2 ! execSSHtemplate.sh
3 ! or execSGEtemplate.sh must be modified.
4 !
5 #define FNODATDEF 33
6  subroutine xbgrun
7  implicit none
8 #include "Zmaxdef.h"
9 #include "Zmanagerp.h"
10 #include "Ztrack.h"
11 #include "Ztrackv.h"
12 #include "Ztrackp.h"
13 #include "Zcode.h"
14 #include "Zheavyp.h"
15 #include "Zobs.h"
16 #include "Zobsp.h"
17 #include "Zobsv.h"
18 #include "Zstackv.h"
19 #include "Zincidentv.h"
20 #include "Zprivate.h"
21 #include "Zprivate1.h"
22 #include "Zprivate2.h"
23 #include "Zprivate3.h"
24 
25  type(track):: incident
26  type(coord):: AngleAtObs
27 
28  logical HEobs ! if T, currently observing
29  common /zheobs/ heobs ! particle is the one obsrved at skeelton making time
30 
31 
32 
33 
34  save rspec, lossrspec, arspec, respec
35  save rzspec, zfspec, rtspec1, rtspec2, retspec1, retspec2
36  save rezspec, rzfspec, rfspec, efspec, refspec
37 
38 
39  integer id ! input. 1 ==> aTrack is going out from
40 ! outer boundery.
41 ! 2 ==> reached at an observation level
42 ! 3 ==> reached at inner boundery.
43  type(track):: aTrack
44 
45  type(track):: inci
46  type(coord):: angle
47  type(coord):: tetafai
48 
49  character*128 input
50  character*64 dirstr
51  real sr, dr
52  integer i, j, k, icon, ir, l
53  integer NN
54  integer iij, code, codex
55  integer i1, i2, ic
56 
57  real*8 r, Eloss, rinmu, cosang
58  real*8 dedt, rho, dist, disto, dedxmu
59  real*8 fai, rminm
60 
61  real*8 u
62  logical accept
63  real*8 wx, wy, wz, temp
64  real za
65  real de, Ek, f, molu
66  real*8 cvh2den
67  integer ldep, ridx, faiidx
68  integer lengenv
69  integer ncpu, mcpu ! no. of smashed skeletons, and actully used skeletonsn
70  integer margin
71  real*4 enhance ! since we use only mcpu, the result must be enhanced
72  ! by a factor of ncpu/mcpu
73  integer binw
74  character*9 ptcl(4)
75  data ptcl/"Photons", "Electrons","Muons", "hadrons"/
76  character*9 ptcl2(3)
77  data ptcl2/"Electrons", "Muons","All"/
78  real power(3)
79  integer nstr
80  real E0, cosz, limit, limite
81  data power/2.,2.,1./
82  real power2(3)
83  data power2/2.,1.,2./
84  character*128 title
85  character*96 evid(nsites)
86  real*8 cog, cog2, sumne, obstimes, Savederg(5)
87  real*8 firstz, dd
88  logical dosort
89  character*2 kd(3)
90  save
91 
92 
93 !
94 ! ***********************
95 #include "interface2.f"
96 ! *********************
97 ! example
98 ! histdep: 2 5 6 7 10 /
99 ! depth 1000 2000 3000 4000 5000 6000 7000 8000 9000 10000
100 ! ansites = 5
101 ! w2hl: 0 1 0 0 2 3 4 0 0 5 0 0 ...
102 
103 
104  limit = 20000
105 
106  lengenv = kgetenv2("NCPU", input)
107  read(input(1:lengenv),*) ncpu
108  lengenv = kgetenv2("MCPU", input)
109  read(input(1:lengenv),*) mcpu
110  lengenv = kgetenv2("MARGIN", input)
111  read(input(1:lengenv),*) margin
112  enhance = ncpu
113  enhance = enhance/mcpu
114  limite = limit * enhance ! this must be enhanced from the first.
115  write(0,*) ' *** ncpu=',ncpu,
116  * ' mcpu=',mcpu, ' enhance factor=',enhance
117  write(0,*) ' margin=',margin
118  do i = 1, nsites
119  w2hl(i) = 0
120  w2il(i) = 0
121  enddo
122  do i = 1, nsites
123  if(histdep(i) .eq. 0) exit
124  hnsites = i
125  w2hl( histdep(i)) = i
126  enddo
127 
128  do i = 1, nsites
129  if( indivdep(i) .eq. 0) exit
130  ansites = i
131  w2il( indivdep(i) ) = i
132  enddo
133 
134 ! now hnsites is the actual number of sites where statistics
135 ! by histogram is needed.
136 ! mapping; depth index to array index; where2loc
137 !
138 
139 
140  call kwhistso( binw ) ! specify binary or ascii write of histogaram
141  ! 1--> ascii 2--> binary
142 
143 
144  rminm = rmin/10.**(bin/2.0d0)
145  r=rmin
146  dr = 10.**bin
147  rbin(1) = 0.
148  do i = 2, nrbin
149  rbin(i) = r
150  r = r* dr
151  enddo
152 
153  return
154 ! ******************
155  entry ihist
156 !
157 !
158 
159 ! histogram: instanciate
160 ! rspec (lateral):
161  if(tklat) then
162  do i = 1, hnsites
163  do j = 1, 4 ! g,e,mu,had
164  call kwhisti(rspec(j, i),
165  * rmin, bin, nrbin, b'00011' )
166  call kwhistai(rspec(j,i),
167  * "Lateral Dist. of "//ptcl(j),
168  * "lat", "ptcls", .true., power(j), "r", "m.u")
169  enddo
170  enddo
171  endif
172 
173 
174 ! elosrspec (energy loss lateral) 10m-10km log bin 0.1. e+mu, e,mu
175 !
176 !
177  if(tkelosslat) then
178  do i = 1, hnsites
179  do j = 1, 3 ! e, mu, e+mu
180  call kwhisti( lossrspec(j, i),
181  * rmin, bin, nrbin, b'00011')
182  call kwhistai(lossrspec(j, i),
183  * "dE/dx lateral dist. of "//ptcl2(j),
184  * "dEdxLat", "GeV/(g/cm^2)", .true., power2(j),
185  * "r", "m.u")
186  enddo
187  enddo
188  endif
189 
190 
191  if(tkarspec) then
192  do i = 1, hnsites
193  do j=1, 4
194  call kwhisti2(arspec(j, i),
195  * 0., 30.0, 12, b'00010',
196  * rmin, bin, nrbin, b'00011' )
197  call kwhistai2(arspec(j,i),
198  * "Lateral dits. of "//ptcl(j)//" with given Fai bin",
199  * "ar", "ptcls", .true., power(j),
200  * "azimuth", "deg", "r", "m.u")
201  enddo
202  enddo
203  endif
204 
205 
206 !
207  if(tkrespec) then
208  do i = 1, hnsites
209  do j= 1, 2
210  call kwhisti2( respec(j, i),
211  * 0.01, 0.2, 20, b'00011',
212  * 500.e-6, 0.1, 50, b'00001')
213  call kwhiststep2(respec(j, i), 2)
214  call kwhistai2( respec(j, i),
215  * "Energy Spec. of "//ptcl(j)//" at diff. r",
216  * "re", "ptcls", .true., 1.,
217  * "r", "m.u", "E", "GeV")
218  enddo
219 ! mu
220  call kwhisti2( respec(3, i),
221  * 0.01, 0.2, 20 , b'00011',
222  * 0.031627, 0.1, 38, b'00011' )
223  call kwhiststep2(respec(3, i), 2)
224  call kwhistai2( respec(3, i),
225  * "Energy Spec. of mu at diff. r",
226  * "re", "ptcls", .true., 0.,
227  * "r", "m.u", "E", "GeV")
228  enddo
229  endif
230 
231 
232  if(tkrzspec) then
233  do i = 1, hnsites
234  do j= 1, 2
235  call kwhisti2( rzspec(j, i),
236  * 0.01, 0.2, 20, b'00011',
237  * 0., 1.0, 20, b'10000')
238  call kwhiststep2(rzspec(j, i), 2)
239  call kwhistai2( rzspec(j, i),
240  * "Zenith angle dist. of "//ptcl(j)//" at diff. r",
241  * "rz", "ptcls", .true., 0.,
242  * "r", "m.u", "cosz", " ")
243  enddo
244  call kwhisti2(rzspec(3, i),
245  * 0.01, 0.2, 20, b'00011',
246  * 0., 1.0, 20, b'10000' )
247  call kwhiststep2(rzspec(3, i), 2)
248  call kwhistai2( rzspec(3, i),
249  * "Zenith angle dist. of m with diff. r",
250  * "rz", "ptcls", .true., 0.,
251  * "r", "m.u", "cosz", " ")
252  enddo
253  endif
254 
255 
256  if(tkzfspec) then
257  do i = 1, hnsites
258  do j= 1, 3
259  call kwhisti2( zfspec(j, i),
260  * 0.0, 1.0, 10, b'10000',
261  * -1.0, 1.0, 50, b'10000')
262  call kwhiststep2(zfspec(j, i), 2)
263  call kwhistai2( zfspec(j, i),
264  * "f=(wx,wy)*(x,y) spectrum of "//ptcl(j)//
265  * " with diff. cosz",
266  * "zf", "ptcls", .true., 0.,
267  * "cosz", " ", "f", " ")
268  enddo
269  enddo
270  endif
271 
272 
273  if(tkrfspec) then
274  do i = 1, hnsites
275  do j= 1, 3
276  call kwhisti2( rfspec(j, i),
277  * 0.01, 0.2, 20, b'00011',
278  * -1.0, 1.0, 50, b'10000')
279  call kwhiststep2(rfspec(j, i), 2)
280  call kwhistai2( rfspec(j, i),
281  * "f spectrum of "//ptcl(j)//" with diff. r",
282  * "rf", "ptcls", .true., 0.,
283  * "r", "m.u", "f", " ")
284  enddo
285  enddo
286  endif
287 
288 
289  if(tkefspec) then
290  do i = 1, hnsites
291  do j = 1, 3
292  call kwhisti2( efspec(j, i),
293  * 500.e-6, 0.2, 20, b'00001',
294  * -1.0, 1.0, 50, b'10000')
295  call kwhiststep2(efspec(j, i), 2)
296  call kwhistai2( efspec(j, i),
297  * "f spectrum of "//ptcl(j)//" with diff. E",
298  * "ef", "ptcls",.true., 0.,
299  * "E", "GeV", "f", " ")
300  enddo
301  enddo
302  endif
303 
304  if(tkrtspec) then
305  do i = 1, hnsites
306  do j = 1, 3
307  call kwhisti2( rtspec1(j, i),
308  * 1.0, 0.2, 10, b'00011',
309  * 10., 0.1, 55, b'00011' )
310  call kwhiststep2(rtspec1(j, i), 2)
311  call kwhistai2( rtspec1(j, i),
312  * "Arrival time dist. of "//ptcl(j)//" at diff. r",
313  * "rt", "ptcls", .true., 0.,
314  * "r", "m.u", "time", "ns")
315 
316  call kwhisti2( rtspec2(j, i),
317  * 0., 1.0, 10, b'10000',
318  * 0., 0.25, 300, b'00000' )
319  call kwhiststep2(rtspec2(j, i), 2)
320  call kwhistai2( rtspec2(j, i),
321  * "Arrival time dist. of "//ptcl(j)//" at diff. r",
322  * "rt", "ptcls", .true., 0.,
323  * "r", "m.u", "time", "ns")
324  enddo
325  enddo
326  endif
327 
328  if(tkretspec) then
329  do i = 1, hnsites
330  do j = 1, 2
331 ! g,e
332  call kwhisti3(retspec1(j, i),
333  * 1.0, 0.2, 10, b'00011',
334  * 500.e-6, 0.25, 15, b'01001',
335  * 10.0, 0.1, 55, b'00011' )
336  call kwhiststep3(retspec1(j, i), 2, 2)
337  call kwhistai3(retspec1(j, i),
338  * "Arrival time dist. of "//ptcl(j)//" with diff. r&E",
339  * "ret", "ptcls", .true., 0.,
340  * "r", "m.u", "E", "GeV", "Time", "ns")
341 
342  call kwhisti3(retspec2(j, i),
343  * 0.0, 1.0, 10, b'10000',
344  * 500.e-6, 0.25, 12, b'01001',
345  * 0.0, 0.25, 200, b'00000' )
346  call kwhiststep3(retspec2(j, i), 2,2)
347  call kwhistai3(retspec2(j, i),
348  * "Arrival time dist. of "//ptcl(j)//" with diff. r&E",
349  * "ret2", "ptcls", .true., 0.,
350  * "r", "m.u", "E", "GeV", "Time", "ns")
351  enddo
352 ! mu
353  call kwhisti3(retspec1(3, i),
354  * 1.0, 0.2, 10, b'00011',
355  * 100.e-3, 0.25, 10, b'01011',
356  * 10.0, 0.1, 72, b'00011' )
357  call kwhiststep3(retspec1(3, i), 2, 2)
358  call kwhistai3(retspec1(3, i),
359  * "Arrival time dist. of "//ptcl(3)//" with diff. r&E",
360  * "ret", "ptcls", .true., 0.,
361  * "r", "m.u", "E", "GeV", "Time", "ns")
362 
363 
364  call kwhisti3(retspec2(3, i),
365  * 0.0, 1.0, 10, b'10000',
366  * 100.e-3, 0.25, 10, b'01011',
367  * 0.0, 0.25, 200, b'00000' )
368  call kwhiststep3(retspec2(3, i), 2,2)
369  call kwhistai3(retspec2(3, i),
370  * "Arrival time dist. of "//ptcl(3)//" with diff. r&E",
371  * "ret2", "ptcls", .true., 0.,
372  * "r", "m.u", "E", "GeV", "Time", "ns")
373  enddo
374  endif
375 
376 
377  if(tkrezspec) then
378  do i = 1, hnsites
379  do j = 1, 2
380  call kwhisti3( rezspec(j, i),
381  * 0.1, 0.2, 15, b'00011',
382  * 500.e-6, 0.25, 14, b'01001',
383  * 0.0, 1.0, 20, b'10000')
384  call kwhiststep3(rezspec(j, i), 3,2)
385  call kwhistai3(rezspec(j, i),
386  * "cos zenith dist. of "//ptcl(j)//" with diff. r&E",
387  * "rez", "ptcls", .true., 0.,
388  * "r", "m.u", "E", "GeV", "cosz", " ")
389  enddo
390  call kwhisti3(rezspec(3, i),
391  * 0.1, 0.2, 15, b'00011',
392  * 100.e-3, 0.25, 10, b'01011',
393  * 0.0, 1.0, 20, b'10000')
394  call kwhiststep3(rezspec(3, i), 3,2)
395  call kwhistai3(rezspec(3, i),
396  * "cos zenith dist. of "//ptcl(3)//" with diff. r&E",
397  * "rez", "ptcls", .true., 0.,
398  * "r", "m.u", "E", "GeV", "cosz", " ")
399  enddo
400  endif
401 
402 
403 
404  if(tkrzfspec) then
405  do i = 1, hnsites
406  do j = 1, 3
407  call kwhisti3(rzfspec(j, i),
408  * 0.1, 0.2, 15, b'00011',
409  * 0.0, 1.0, 10, b'10000',
410  * -1.0, 1.0, 20, b'10000')
411  call kwhiststep3(rzfspec(j, i), 3,2)
412  call kwhistai3(rzfspec(j, i),
413  * "f spectrum of "//ptcl(j)//" with diff r&cosz",
414  * "rzf", "ptcls", .true., 0.,
415  * "r", "m.u", "cosz", " ", "f", " ")
416  enddo
417  enddo
418  endif
419 
420  if(tkrefspec) then
421  do i = 1, hnsites
422  do j = 1, 2
423  call kwhisti3( refspec(j, i),
424  * 0.1, 0.2, 15, b'00011',
425  * 500.e-6, 0.25, 16, b'01001',
426  * -1.0, 1.0, 20, b'10000')
427  call kwhiststep3(refspec(j, i), 3,3)
428  call kwhistai3(refspec(j, i),
429  * "f spectrum of "//ptcl(j)//" with diff. r&E",
430  * "ref", "ptcls", .true., 0.,
431  * "r", "m.u", "E", "GeV", "f", " ")
432  enddo
433  call kwhisti3(refspec(3, i),
434  * 0.1, 0.2, 15, b'00011',
435  * 100.e-3, 0.25, 10, b'01001',
436  * -1.0, 1.0, 20, b'10000')
437  call kwhiststep3(refspec(3, i), 3,2)
438  call kwhistai3(refspec(3, i),
439  * "f spectrum of "//ptcl(3)//" with diff. r&E",
440  * "ref", "ptcls", .true., 0.,
441  * "r", "m.u", "E", "GeV", "f", " ")
442  enddo
443  endif
444 
445  return
446 ! *********************************** hook for Beginning of 1 event
447 ! * All system-level initialization for 1 event generation has been
448 ! * eneded at this moment.
449 ! * After this is executed, event generation starts.
450 ! *
451  entry xbgevent
452 
453 
454  call cqincident(inci, angle)
455 
456 
457 #if FNODATDEF > 0
458  bufc=0
459 #endif
460 
461  cosz = -angle.r(3)
462  e0 = inci.p.fm.p(4)
463  if(inci.p.code .eq. 9) then
464  nn= inci.p.subcode
465  elseif(inci.p.code .eq. 1) then
466  nn=0
467  else
468  nn=1
469  endif
470 ! next is only available for parallel job. For normal job,
471 ! fisrt col.depth is not yet fixed.
472  firstz= zfirst.pos.depth*0.1
473  write(0,'(a,1pE11.3,a,E11.3,a,E11.3,a)')
474  * ' E0=',e0, ' cosz=',cosz, ' firstz=',
475  * firstz, ' g/cm2'
476 ! here we use enhanced limit
477  call howmuch(limite, e0, nn, cosz)
478 
479 !
480  do i = 1, noofassites
481  sumeloss(i) = 0.
482  ng(i) = 0.
483  ne(i) = 0.
484  nmu(i) = 0.
485  nhad(i) = 0.
486  enddo
487 
488  if(tklat) then
489  do i = 1, hnsites
490  do j = 1, 4
491  call kwhistc(rspec(j,i))
492  enddo
493  enddo
494  endif
495 
496  if(tkelosslat) then
497  do i = 1, hnsites
498  do j = 1, 3
499  call kwhistc(lossrspec(j,i))
500  enddo
501  enddo
502  endif
503 
504 
505  if(tkarspec) then
506  do i = 1, hnsites
507  do j = 1, 4
508  call kwhistc2(arspec(j,i))
509  enddo
510  enddo
511  endif
512 
513  if(tkrespec) then
514  do i = 1, hnsites
515  do j = 1, 3
516  call kwhistc2(respec(j,i))
517  enddo
518  enddo
519  endif
520 
521 
522  if(tkrzspec) then
523  do i = 1, hnsites
524  do j = 1, 3
525  call kwhistc2(rzspec(j,i))
526  enddo
527  enddo
528  endif
529 
530  if(tkrfspec) then
531  do i = 1, hnsites
532  do j = 1, 3
533  call kwhistc2(rfspec(j,i))
534  enddo
535  enddo
536  endif
537 
538 
539 
540  if(tkefspec) then
541  do i = 1, hnsites
542  do j = 1, 3
543  call kwhistc2(efspec(j,i))
544  enddo
545  enddo
546  endif
547 
548 
549 
550 
551  if(tkrtspec) then
552  do i = 1, hnsites
553  do j = 1, 3
554  call kwhistc2(rtspec1(j,i))
555  call kwhistc2(rtspec2(j,i))
556  enddo
557  enddo
558  endif
559 
560 
561  if(tkretspec) then
562  do i = 1, hnsites
563  do j = 1, 3
564  call kwhistc3(retspec1(j,i))
565  call kwhistc3(retspec2(j,i))
566  enddo
567  enddo
568  endif
569 
570 
571 
572  if(tkrezspec) then
573  do i = 1, hnsites
574  do j = 1, 3
575  call kwhistc3(rezspec(j,i))
576  enddo
577  enddo
578  endif
579 
580 
581 
582  if(tkrzfspec) then
583  do i = 1, hnsites
584  do j = 1, 3
585  call kwhistc3(rzfspec(j,i))
586  enddo
587  enddo
588  endif
589 
590 
591  if(tkrefspec) then
592  do i = 1, hnsites
593  do j = 1, 3
594  call kwhistc3(refspec(j,i))
595  enddo
596  enddo
597  endif
598 
599  do i = 1, hnsites
600  do j= 1, 4
601  do k = 1, nfai
602  do l = 1, nrbin
603  nrfairec(l, k, j, i) = 0.
604  nrfaiall(l, k, j, i) = 0.
605  derfai(l, k, i) = 0. ! no j (for all charged ptcls)
606  enddo
607  enddo
608  enddo
609  enddo
610 
611  obstimes = 0.
612 
613  return
614 ! ***************
615  entry xobs(atrack, id)
616 !
617 ! For id =2, you need not output the z value, because it is always
618 ! 0 (within the computational accuracy).
619 !
620 ! **************************
621  obstimes = obstimes + 1.d0
622  if(mod(obstimes, 500000.d0) .eq. 0. ) then
623  dosort=.false.
624  do i = 1, min(4,stack_pos)
625  if(stack(i).p.fm.p(4) .ne. savederg(i)) then
626  savederg(i)=stack(i).p.fm.p(4)
627  dosort=.true.
628  endif
629  enddo
630  if(dosort) then
631  call csortstack
632  endif
633  write(0, *) ' obstimes=', obstimes, ' ptclE=',atrack.p.fm.p(4)
634  do i = 1, min(4,stack_pos)
635  write(0,*)' stack tops=', stack(i).p.fm.p(4)
636  enddo
637  endif
638 ! ***************
639  code = atrack.p.code
640  if(id .eq. 2 .and. code .le. 6) then
641  codex=min(code, 4)
642  wz = atrack.vec.w.r(3) ! downgoing < 0
643  if(wz .gt. 0) return
644  wz = -wz
645  ldep = atrack.where
646 
647  if(code .eq. kphoton) then
648  ng(ldep) = ng(ldep) + enhance
649  elseif(code .eq. kelec) then
650  ne(ldep) = ne(ldep) + enhance
651  elseif(code .eq. kmuon) then
652  nmu(ldep) = nmu(ldep) + enhance
653  else
654  nhad(ldep)=nhad(ldep) + enhance
655  endif
656 
657  r = sqrt( atrack.pos.xyz.x**2 +
658  * atrack.pos.xyz.y**2 )
659 
660  ek = atrack.p.fm.p(4) -atrack.p.mass
661 
662 ! ---------- compute energy loss rate
663  if(atrack.p.charge .ne. 0 ) then
664  rho = cvh2den(atrack.pos.height)
665 ! get energy loss when aTrack goes 1 g/cm2 along the
666 ! primary direction. Gramage the particle can run is
667 ! 1/cos where cos is the cos of angle relative to the
668 ! primary angle . 1g/cm^2 = 10-3kg/10-4 m^2 =10 kg/m^2.
669 ! To travel 1 g/cm^2 along shower axis, the ptcl must
670 ! run dist kg/m^2
671  if(abs(wz) .gt. 1.d-2) then
672  dist =10./wz ! in kg/m2/(g/cm2)
673  else
674 ! for safety
675  dist =1000.
676  endif
677  if( heobs ) then
678 ! the ptcls is the one obsrved at skeleton making time
679 ! we must compute dedt here
680  call cdedxinair(atrack.p, rho, dedt) ! dedt; GeV/(kg/m2)
681  if(atrack.p.code .eq. kmuon ) then
682 ! dE/dx due to muon pair, brem, nuc.i
683  call cmudedx(muni, mubr, mupr, atrack.p.fm.p(4),
684  * dedxmu) ! dedxmu in GeV/(g/cm2)
685  dedxmu = dedxmu /10. ! GeV/(kg/m2)
686  dedt = dedt + dedxmu
687  endif
688  else
689 ! we can use already computed one
690  call cqelossrate(dedt) ! loss rate GeV/(kg/m^2)
691  endif
692 ! energy loss rate
693  eloss = dedt*dist*enhance ! GeV/(g/cm2)
694 
695  sumeloss(ldep)=sumeloss(ldep) + eloss
696  else
697  eloss=0.
698  endif
699 
700  j = w2il(ldep)
701  if(j .eq. 0) goto 100
702 
703  molu =obssites(ldep).mu
704  rinmu =r/molu
705  sr = rinmu
706  if(rinmu .gt. rminm) then
707  ridx= log10(rinmu/rminm)/bin + 1
708  else
709  ridx =0
710  endif
711 ! fai is in -15 to 345 (for dfai=30.)
712  fai=atan2(atrack.pos.xyz.y, atrack.pos.xyz.x)*
713  * 57.29577951308230d0
714  fai= mod(fai + 360.d0, 360.d0)
715  if(fai .gt. (360.d0-dfai/2.0d0)) fai= fai-360.d0
716  faiidx=(fai+dfai/2.0d0) /dfai + 1
717 
718  if(ridx .gt. 0 .and. ridx .le. nrbin ) then
719  call rndc(u)
720  accept = u .le. recprob(ridx, codex, j)
721  if(accept) then
722  nrfairec(ridx, faiidx, codex, j) =
723  * nrfairec(ridx, faiidx, codex, j) + 1.0
724  endif
725  nrfaiall(ridx, faiidx, codex, j) =
726  * nrfaiall(ridx, faiidx, codex, j) + 1.0
727  derfai(ridx, faiidx, j) =
728  * derfai(ridx, faiidx, j) + eloss
729  else
730  accept = .false.
731  endif
732  if(accept) then
733 #if FNODATDEF > 0
734  if(bufc .lt. bufsize) then
735  bufc = bufc + 1
736  buf(bufc).ldep=ldep
737  buf(bufc).code=code
738  buf(bufc).subcode=atrack.p.subcode
739  buf(bufc).charge = atrack.p.charge
740  buf(bufc).ridx=ridx
741  buf(bufc).faiidx= faiidx
742  buf(bufc).rinmu = rinmu
743  buf(bufc).fai= fai
744  buf(bufc).ek = ek
745  buf(bufc).t = atrack.t
746  buf(bufc).wx=-atrack.vec.w.r(1)
747  buf(bufc).wy=-atrack.vec.w.r(2)
748  buf(bufc).wz=wz
749  else
750  write(fnodat) bufc, buf
751  bufc= 0
752  endif
753 #else
754  if(wz .lt. 0.999) then
755  write(*,
756  * '(6i3, 1pE11.3, 0p,f6.1,1p2E11.3,0p, 2f8.4,f10.6)')
757  * ldep, code, atrack.p.subcode,
758  * atrack.p.charge, ridx, faiidx,
759  * rinmu, fai,
760  * ek, atrack.t,
761  * -atrack.vec.w.r(1), -atrack.vec.w.r(2), wz
762  else
763  write(*,'(6i3,1pE11.3, 0p, f6.1, 1p2E11.3, a)')
764  * ldep, code, atrack.p.subcode,
765  * atrack.p.charge, ridx, faiidx,
766  * rinmu, fai,
767  * ek, atrack.t, " 0 0 1"
768  endif
769 #endif
770  endif
771  endif
772  100 continue
773 ! ------------------
774  i=w2hl(ldep)
775 
776  if(id .eq. 2 .and. code .le. 6 .and. i .gt. 0) then
777  if( tklat ) then
778 ! lateral
779  call kwhist(rspec(codex, i), sr, enhance )
780  endif
781  if( tkelosslat
782  * .and. atrack.p.charge .ne. 0 ) then
783 !c de = Eloss*enhance
784  de = eloss
785  call kwhist(lossrspec(3,i), sr, de)
786  if( code .eq. kelec ) then
787 ! by electrons
788  call kwhist(lossrspec(1,i), sr, de)
789  else
790 ! by other charged particles
791  call kwhist(lossrspec(2,i), sr, de)
792  endif
793 ! end of eloss late
794  endif
795 
796 
797  if(tkarspec) then
798 ! arspec
799  call kwhist2( arspec(codex, i), sngl(fai), sr, enhance)
800  endif
801  endif
802  if(id .eq. 2 .and. code .le. 3 .and. i .gt. 0) then
803  if( tkrespec ) then
804 ! re spectrum
805  call kwhist2( respec(code, i), sr, ek, enhance)
806  endif
807 
808  if( tkrzspec ) then
809  call kwhist2( rzspec(code, i), sr,
810  * sngl(wz), enhance )
811  endif
812 
813  if( wz .lt. 1.0d0 ) then
814  temp = 1.d0 - wz**2
815  if(temp .gt. 0.) then
816  temp = sqrt(temp)
817  f = (
818  * atrack.pos.xyz.x* atrack.vec.w.r(1) +
819  * atrack.pos.xyz.y* atrack.vec.w.r(2) ) /r
820  * /temp
821  else
822  f = 1.0
823  endif
824  else
825  f = 1.0
826  endif
827  if( tkzfspec .and. f .lt. 1.0 ) then
828  call kwhist2( zfspec(code, i),
829  * sngl(wz), f, enhance )
830  endif
831 
832 
833  if( tkefspec .and. f .lt. 1.0 ) then
834  call kwhist2( efspec(code, i),
835  * ek, f, enhance )
836  endif
837 
838  if( tkrfspec .and. f .lt. 1.0 ) then
839  call kwhist2( rfspec(code, i),
840  * sr, f, enhance )
841  endif
842 
843  if( tkrtspec ) then
844  call kwhist2( rtspec1(code, i), sr,
845  * sngl( atrack.t ), enhance )
846  call kwhist2( rtspec2(code, i), sr,
847  * sngl( atrack.t ), enhance )
848  endif
849 
850  if(tkretspec ) then
851  call kwhist3( retspec1(code, i), sr, ek,
852  * sngl( atrack.t ), enhance)
853  call kwhist3( retspec2(code, i), sr, ek,
854  * sngl( atrack.t ), enhance)
855  endif
856 
857  if( tkrezspec ) then
858  call kwhist3( rezspec(code,i), sr, ek,
859  * sngl(wz), enhance )
860  endif
861 
862  if( tkrzfspec ) then
863  call kwhist3( rzfspec(code, i), sr,
864  * sngl(wz), f, enhance )
865  endif
866 
867  if( tkrefspec ) then
868  call kwhist3( refspec(code, i), sr,
869  * ek, f, enhance)
870  endif
871  endif
872  return
873 ! **************
874  entry xenevent
875 ! **************
876 #if FNODATDEF > 0
877 ! if(fnodat .gt. 0) then
878  if(bufc .gt. 0) then
879  write(fnodat) bufc, buf
880  bufc=0
881  endif
882 ! endif
883 #endif
884 
885 
886 ! call cqFirstID(firstz)
887  firstz= zfirst.pos.depth*0.1
888 
889  if(observeas) then
890  cog = 0.
891  sumne = 0.
892 
893  do i = 1, noofassites
894 
895  asobssites(i).esize = asobssites(i).esize* enhance
896 
897  if(i .gt. 1 .and. i .lt. noofassites ) then
898  dd =(asdepthlist(i+1) - asdepthlist(i-1))/2.0
899  elseif(i .eq. 1) then
900  dd =(asdepthlist(2) - asdepthlist(1))
901  else
902  dd =(asdepthlist(noofassites) -
903  * asdepthlist(noofassites-1))
904  endif
905  cog = cog + asobssites(i).esize*dd*asdepthlist(i)
906  sumne= sumne +asobssites(i).esize*dd
907  enddo
908 ! 0.1 is for g/cm2
909  cog = cog*0.1/sumne
910 
911  cog2 = 0.
912  sumne = 0.
913  do i = 1, noofassites
914  if( asobssites(i).age .gt.
915  * (2.0-asobssites(noofassites).age)) then
916  if(i .gt. 1 .and. i .lt. noofassites ) then
917  dd =( asdepthlist(i+1) - asdepthlist(i-1))/2.0
918  elseif(i .eq. 1) then
919  dd =(asdepthlist(2) - asdepthlist(1))
920  else
921  dd =(asdepthlist(noofassites) -
922  * asdepthlist(noofassites-1))
923  endif
924  dd = dd
925  cog2 = cog2 + asobssites(i).esize*asdepthlist(i)*dd
926  sumne= sumne +asobssites(i).esize*dd
927  endif
928  enddo
929  if(sumne .gt. 0.) then
930  cog2 = cog2*0.1/sumne
931  else
932 ! too deep penetration
933  cog2 = asdepthlist(noofassites)*0.1
934  endif
935 
936 
937  if(fnob .ge. 0 ) then
938  write(fnob,
939  * '("h ", i4, 3i3, 1pE11.3, 0p 3f11.7, f7.2, 2f7.0)')
940  * eventno, inci.p.code,
941  * inci.p.subcode, inci.p.charge,
942  * inci.p.fm.e, -angle.r(1), -angle.r(2), -angle.r(3),
943  * firstz, cog, cog2
944  else
945  write(*,
946  * '("h ", i4, 3i3, 1pE11.3, 0p 3f11.7, f7.2, 2f7.0)')
947  * eventno, inci.p.code,
948  * inci.p.subcode, inci.p.charge,
949  * inci.p.fm.e, -angle.r(1), -angle.r(2), -angle.r(3),
950  * firstz, cog, cog2
951  endif
952 
953  do i = 1, noofassites
954  if(fnob .ge. 0) then
955  write(fnob, '("t ", i3, 2f7.1, 2f6.3,
956  * 1p6E11.3)')
957  * i,
958  * asdepthlist(i)*0.1, asobssites(i).mu,
959  * asobssites(i).age, asdepthlist(i)*0.1/cog2,
960  * ng(i), ne(i), nmu(i), nhad(i),
961  * asobssites(i).esize, sumeloss(i)
962  else
963  write(*, '("t ", i3, 2f7.1, 2f6.3,
964  * 1p6E11.3)')
965  * i,
966  * asdepthlist(i)*0.1, asobssites(i).mu,
967  * asobssites(i).age, asdepthlist(i)*0.1/cog2,
968  * ng(i), ne(i), nmu(i), nhad(i),
969  * asobssites(i).esize, sumeloss(i)
970  endif
971  enddo
972  if(fnob .gt. 0 ) then
973  write(fnob,*)
974  else
975  write(*,*)
976  endif
977  endif
978 
979 
980 ! here we write limit; at reduce process, non enhanced value
981 ! is needed
982  write(fnonrfai,'(i2,1pE11.3, 0p,i3, f8.4, 1pE11.3,3i4)')
983  * eventno, e0, nn, cosz, limit, nrbin, nfai, ansites
984  do i = 1, ansites
985  do j = 1, 4
986  do k = 1, nfai
987  l = indivdep(i)
988  write(fnonrfai, '("rec",f7.1, 4i4)' )
989  * asdepthlist(l)*0.1, l, i, j, k
990  write(fnonrfai, '(1p10E11.3)')
991  * ( nrfairec(l,k,j,i), l=1,nrbin )
992  enddo
993  enddo
994  enddo
995 
996  do i = 1, ansites
997  do j = 1, 4
998  do k = 1, nfai
999  l = indivdep(i)
1000  write(fnonrfai, '("all",f7.1, 4i4)' )
1001  * asdepthlist(l)*0.1, l, i, j, k
1002  write(fnonrfai, '(1p10E11.3)')
1003  * ( nrfaiall(l,k,j,i)*enhance, l=1,nrbin )
1004  enddo
1005  enddo
1006  enddo
1007 ! dE/dx lateral
1008  do i = 1, ansites
1009  do k = 1, nfai
1010  l = indivdep(i)
1011  write(fnonrfai, '("dE/dx",f7.1, 3i4)' )
1012  * asdepthlist(l)*0.1, l, i, k
1013  write(fnonrfai, '(1p10E11.3)')
1014  * ( derfai(l,k,i)*enhance, l=1,nrbin )
1015  enddo
1016  enddo
1017 
1018  write(fnonrfai, *)
1019 
1020 
1021  do i = 1, hnsites
1022  j=histdep(i)
1023  write(evid(i),
1024  * '(i3, f7.1, f6.3, f6.3,
1025  * i5, i4)')
1026  * histdep(i), asdepthlist(j)*0.1,
1027  * asobssites(j).age, asdepthlist(j)*0.1/cog,
1028  * int(asobssites(j).mu), int(cog)
1029  enddo
1030 
1031 
1032  if( tklat ) then
1033  do j = 1, 4
1034  do i = 1, hnsites
1035  k = histdep(i)
1036  call kwhists( rspec(j,i), 0. ) ! 0. means area norm.
1037  call kwhistev(rspec(j,i), eventno)
1038  call kwhistid(rspec(j,i), evid(i))
1039  call kwhistdir(rspec(j,i), ptcl(j)//"/")
1040  call kwhistp( rspec(j,i), fno)
1041  enddo
1042  enddo
1043  endif
1044 
1045  if(tkelosslat) then
1046  do j = 1, 3
1047  do i = 1, hnsites
1048  call kwhists( lossrspec(j,i), 0. ) ! area norm
1049  call kwhistev(lossrspec(j,i), eventno)
1050  call kwhistid( lossrspec(j, i), evid(i))
1051  call kwhistdir( lossrspec(j, i), ptcl2(j)//"/")
1052  call kwhistp( lossrspec(j, i), fno)
1053  enddo
1054  enddo
1055 
1056  endif
1057  120 continue
1058 
1059  if(tkarspec) then
1060  do i = 1, hnsites
1061  k=histdep(i)
1062  do j = 1, 4
1063  call kwhists2( arspec(j, i), 0. )
1064  call kwhistev2( arspec(j,i), eventno)
1065  call kwhistid2( arspec(j, i), evid(i))
1066  write(dirstr,'(a,"/d",i4, "/")')
1067  * ptcl(j), int( asdepthlist(k)*0.1 )
1068  call kseblk(dirstr, "|", nstr)
1069  call kwhistdir2(arspec(j, i), dirstr)
1070  call kwhistp2( arspec(j, i), fno)
1071  enddo
1072  enddo
1073  endif
1074 
1075  if(tkrespec) then
1076  do j = 1, 3
1077  do i = 1, hnsites
1078  call kwhists2( respec(j, i), 0. )
1079  call kwhistev2( respec(j,i), eventno)
1080  call kwhistid2( respec(j, i), evid(i))
1081  k=histdep(i)
1082  write(dirstr,'(a,"/d",i4, "/")')
1083  * ptcl(j), int( asdepthlist(k)*0.1 )
1084  call kseblk(dirstr, "|", nstr)
1085  call kwhistdir2(respec(j, i), dirstr)
1086  call kwhistp2( respec(j, i), fno)
1087  enddo
1088  enddo
1089  endif
1090 
1091 
1092  if( tkrzspec ) then
1093  do j = 1, 3
1094  do i = 1, hnsites
1095  call kwhists2( rzspec(j, i), 0.)
1096  call kwhistid2( rzspec(j, i), evid(i))
1097  k=histdep(i)
1098  write(dirstr,'(a,"/d",i4, "/")')
1099  * ptcl(j), int( asdepthlist(k)*0.1 )
1100  call kseblk(dirstr, "|", nstr)
1101  call kwhistdir2( rzspec(j, i), dirstr)
1102  call kwhistp2( rzspec(j, i), fno)
1103  enddo
1104  enddo
1105  endif
1106 
1107 
1108  if( tkzfspec ) then
1109  do j = 1, 3
1110  do i = 1, hnsites
1111  call kwhists2( zfspec(j, i), 0.)
1112  call kwhistid2( zfspec(j, i), evid(i) )
1113  k=histdep(i)
1114  write(dirstr,'(a,"/d",i4, "/")')
1115  * ptcl(j), int( asdepthlist(k)*0.1 )
1116  call kseblk(dirstr, "|", nstr)
1117  call kwhistdir2( zfspec(j, i), dirstr)
1118  call kwhistp2( zfspec(j, i), fno)
1119  enddo
1120  enddo
1121  endif
1122 
1123 
1124 
1125  if( tkrfspec ) then
1126  do j = 1, 3
1127  do i = 1, hnsites
1128  call kwhists2( rfspec(j, i), 0.)
1129  call kwhistid2( rfspec(j, i), evid(i))
1130  k=histdep(i)
1131  write(dirstr,'(a,"/d",i4, "/")')
1132  * ptcl(j), int( asdepthlist(k)*0.1 )
1133  call kseblk(dirstr, "|", nstr)
1134  call kwhistdir2(rfspec(j, i), dirstr)
1135  call kwhistp2( rfspec(j, i), fno)
1136  enddo
1137  enddo
1138  endif
1139 
1140 
1141  if( tkefspec ) then
1142  do j = 1, 3
1143  do i = 1, hnsites
1144  call kwhists2( efspec(j, i), 0.)
1145  call kwhistev2(efspec(j,i), eventno)
1146  call kwhistid2( efspec(j, i), evid(i))
1147  k=histdep(i)
1148  write(dirstr,'(a,"/d",i4, "/")')
1149  * ptcl(j), int( asdepthlist(k)*0.1 )
1150  call kseblk(dirstr, "|", nstr)
1151  call kwhistdir2(efspec(j, i), dirstr)
1152  call kwhistp2( efspec(j, i), fno)
1153  enddo
1154  enddo
1155  endif
1156 
1157 
1158  if( tkrtspec ) then
1159  do j = 1, 3
1160  do i = 1, hnsites
1161  call kwhists2( rtspec1(j,i), 0.)
1162  call kwhistev2( rtspec1(j,i), eventno)
1163  call kwhistid2( rtspec1(j,i), evid(i))
1164  k=histdep(i)
1165  write(dirstr,'(a,"/d",i4, "/")')
1166  * ptcl(j), int( asdepthlist(k)*0.1 )
1167  call kseblk(dirstr, "|", nstr)
1168  call kwhistdir2(rtspec1(j,i), dirstr)
1169  call kwhistp2( rtspec1(j,i), fno)
1170  enddo
1171  enddo
1172  do j = 1, 3
1173  do i = 1, hnsites
1174  call kwhists2( rtspec2(j,i), 0.)
1175  call kwhistev2(rtspec2(j,i), eventno)
1176  call kwhistid2( rtspec2(j,i), evid(i))
1177  k=histdep(i)
1178  write(dirstr,'(a,"/d",i4, "/")')
1179  * ptcl(j), int( asdepthlist(k)*0.1 )
1180  call kseblk(dirstr, "|", nstr)
1181  call kwhistdir2(rtspec2(j,i), dirstr)
1182  call kwhistp2( rtspec2(j,i), fno)
1183  enddo
1184  enddo
1185  endif
1186 
1187 
1188  if(tkretspec ) then
1189  do j = 1, 3
1190  do i = 1, hnsites
1191  call kwhists3( retspec1(j, i), 0.)
1192  call kwhistev3(retspec1(j,i), eventno)
1193  call kwhistid3( retspec1(j, i), evid(i))
1194  k=histdep(i)
1195  write(dirstr,'(a,"/d",i4,"/")' )
1196  * ptcl(j), int( asdepthlist(k)*0.1 )
1197  call kseblk(dirstr, "|", nstr)
1198  call kwhistdir3( retspec1(j, i), dirstr)
1199  call kwhistp3( retspec1(j, i), fno)
1200  enddo
1201  enddo
1202  do j = 1, 3
1203  do i = 1, hnsites
1204  call kwhists3( retspec2(j, i), 0.)
1205  call kwhistev3(retspec2(j,i), eventno)
1206  call kwhistid3( retspec2(j, i), evid(i))
1207  k=histdep(i)
1208  write(dirstr,'(a,"/d",i4, "/")' )
1209  * ptcl(j), int( asdepthlist(k)*0.1 )
1210  call kseblk(dirstr, "|", nstr)
1211  call kwhistdir3(retspec2(j, i), dirstr)
1212  call kwhistp3( retspec2(j, i), fno)
1213  enddo
1214  enddo
1215  endif
1216 
1217 
1218  if( tkrezspec ) then
1219  do j = 1, 3
1220  do i = 1, hnsites
1221  call kwhists3( rezspec(j, i), 0.)
1222  call kwhistev3(rezspec(j,i), eventno)
1223  call kwhistid3( rezspec(j, i), evid(i))
1224  k=histdep(i)
1225  write(dirstr,'(a,"/d",i4,"/")' )
1226  * ptcl(j), int( asdepthlist(k)*0.1 )
1227  call kseblk(dirstr, "|", nstr)
1228  call kwhistdir3(rezspec(j, i), dirstr)
1229  call kwhistp3( rezspec(j, i), fno)
1230  enddo
1231  enddo
1232  endif
1233 
1234 
1235  if( tkrzfspec ) then
1236  do j = 1, 3
1237  do i = 1, hnsites
1238  call kwhists3( rzfspec(j, i), 0.)
1239  call kwhistev3(rzfspec(j,i), eventno)
1240  call kwhistid3( rzfspec(j, i), evid(i))
1241  k=histdep(i)
1242  write(dirstr,'(a,"/d",i4,"/")' )
1243  * ptcl(j), int( asdepthlist(k)*0.1 )
1244  call kseblk(dirstr, "|", nstr)
1245  call kwhistdir3(rzfspec(j, i), dirstr)
1246  call kwhistp3( rzfspec(j, i), fno)
1247  enddo
1248  enddo
1249  endif
1250 
1251 
1252  if( tkrefspec ) then
1253  do j = 1, 3
1254  do i = 1, hnsites
1255  call kwhists3( refspec(j, i), 0.)
1256  call kwhistev3(refspec(j,i), eventno)
1257  call kwhistid3( refspec(j, i), evid(i))
1258  k=histdep(i)
1259  write(dirstr,'(a,"/d",i4, "/")' )
1260  * ptcl(j), int( asdepthlist(k)*0.1 )
1261  call kseblk(dirstr, "|", nstr)
1262  call kwhistdir3( refspec(j, i), dirstr)
1263  call kwhistp3( refspec(j, i), fno)
1264  enddo
1265  enddo
1266  endif
1267  end
*Zfirst p fm *Zfirst p Zfirst p Zfirst p *Zfirst p *Zfirst pos xyz Zfirst pos xyz *Zfirst pos *Zfirst pos depth
Definition: ZavoidUnionMap.h:1
void kwhistid(struct histogram1 *h, char *id)
subroutine howmuch(limit, E0, NN, cosz)
Definition: howmuch.f:2
subroutine cqincident(incident, AngleAtObs)
Definition: cmkIncident.f:486
integer lengid integer lengdir character *dir integer kgetenv2 character *numb character *execid character *msg logical takehist save do nsites histdep(i)=0indivdep(i)=0enddo leng
integer nsites ! max real bin
Definition: Zprivate0.h:2
integer function kgetenv2(envname, envresult)
Definition: cgetLoginN.f:77
dE dx *! Nuc Int sampling table e
Definition: cblkMuInt.h:130
others if is ng
Definition: cblkManager.h:9
const int kphoton
Definition: Zcode.h:6
latitude latitude this system is used *****************************************************************! type coord sequence union map real z z in m endmap xyz map real ! latitude in deg is to the north ! longitude in deg is to the east *h ! height in m endmap llh map real ! polar angle ! azimuthal angle *radius ! radial distance endmap sph endunion character *sys ! which system xyz
Definition: Zcoord.h:25
void kwhistp(struct histogram1 *h, FILE *fno)
Definition: Ztrack.h:44
int nmu[nl][nth]
Definition: Zprivate.h:12
subroutine cmudedx(flagN, flagBr, flagPr, Emu, dEdx, dEdxF)
Definition: cmudEdx.f:2
void kwhisti(struct histogram1 *h, float ixmin, float ibinORxmax, int inbin, int itklg)
max ptcl codes in the kelec
Definition: Zcode.h:2
block data cblkElemag data *AnihiE ! Eposi< 1 TeV, anihilation considered *X0/365.667/, ! radiation length of air in kg/m2 *Ecrit/81.e-3/, ! critical energy of air in GeV *MaxComptonE/1./, ! compton is considered below 1 GeV *MaxPhotoE/1.e-3/, ! above this, PhotoElectric effect neg. *MinPhotoProdE/153.e-3/, ! below 153 MeV, no gp --> hadrons ! scattering const not MeV *Knockon true
Definition: cblkElemag.h:7
void kwhistev(struct histogram1 *h, int evno)
logical KeepWeight ! see setupenv sh logical tkarspec ! get lateral histo in a web fai bin logical tkrtspec ! get time histo in each web bin logical SeeLowdE common Zprivatec2 tkrtspec
Definition: Zprivate4.h:7
void kwhist(struct histogram1 *h, float x, float w)
void kwhistdir(struct histogram1 *h, char *dir)
int ne[nl][nth]
Definition: Zprivate.h:11
subroutine rndc(u)
Definition: rnd.f:91
! timing nrbin
Definition: Zprivate2.h:12
********************block data cblkHeavy ********************integer j data *HeavyG2symbol p
Definition: cblkHeavy.h:7
void kwhistso(int binw)
subroutine cdedxinair(aPtcl, rhoin, dedt, dedtF)
Definition: cdedxInAir.f:49
logical KeepWeight ! see setupenv sh logical tkarspec ! get lateral histo in a web fai bin logical tkrtspec ! get time histo in each web bin logical SeeLowdE common Zprivatec2 tkarspec
Definition: Zprivate4.h:7
integer w2hl(nsites)
block data cblkEvhnp ! currently usable models data RegMdls ad *special data *Cekaon d0
Definition: cblkEvhnp.h:5
*Zfirst p fm *Zfirst p Zfirst p Zfirst p *Zfirst p *Zfirst pos xyz Zfirst pos xyz *Zfirst pos *Zfirst pos Zfirst pos height
Definition: ZavoidUnionMap.h:1
latitude latitude this system is used *****************************************************************! type coord sequence union map real y
Definition: Zcoord.h:25
dE dx *! Nuc Int sampling table d
Definition: cblkMuInt.h:130
void kwhistc(struct histogram1 *h)
dE dx *! Nuc Int sampling table b
Definition: cblkMuInt.h:130
*Zfirst p fm *Zfirst p Zfirst p Zfirst p *Zfirst p charge
Definition: ZavoidUnionMap.h:1
block data cblkElemag data *AnihiE ! Eposi< 1 TeV, anihilation considered *X0/365.667/, ! radiation length of air in kg/m2 *Ecrit/81.e-3/, ! critical energy of air in GeV *MaxComptonE/1./, ! compton is considered below 1 GeV *MaxPhotoE/1.e-3/, ! above this, PhotoElectric effect neg. *MinPhotoProdE/153.e-3/, ! below 153 MeV, no gp --> hadrons ! scattering const not MeV *Knockon ! knockon is considered Obsolete *PhotoProd false
Definition: cblkElemag.h:7
subroutine csortstack
Definition: cstack.f:102
nodes t
void kwhistai(struct histogram1 *h, char *title, char *categ, char *dNunit, int logv, float pw, char *label, char *unit)
void kwhists(struct histogram1 *h, float inorm)
Definition: Zptcl.h:75
Definition: Zcoord.h:43
*Zfirst p fm *Zfirst p mass
Definition: ZavoidUnionMap.h:1
subroutine xbgrun
Definition: interface.f:10
max ptcl codes in the kmuon
Definition: Zcode.h:2
subroutine kseblk(text, c, lc)
Definition: kseblk.f:18
! timing nfai
Definition: Zprivate2.h:12
! structure defining a particle at production ! Basic idea of what is to be contained in ! the particle structue is that dynamical ones should be included those derivable from the particle code ! is not included ******************************************************type fmom momentum sequence union map real e endmap map real * x
Definition: Zptcl.h:21
*Zfirst p fm *Zfirst p Zfirst p Zfirst p subcode
Definition: ZavoidUnionMap.h:1