COSMOS v7.655  COSMOSv7655
(AirShowerMC)
mkNrfaiFromDat.f
Go to the documentation of this file.
1  implicit none
2 ! input: h0 .nrfai file with invalid "rec" but
3 ! correct "all" and "dE/dx" part
4 ! give the name by env. NRFAIFILE
5 ! f0 .dat file. combined main output file
6 ! stdin
7 ! output:h1 .nrfai with correct "rec", "all" and
8 ! "dE/dx"'
9 ! stdout
10 ! this creates "rec" part of .nrfai data from
11 ! .dat file
12 !
13 
14 !
15 #include "Zmaxdef.h"
16 #include "Zobs.h"
17 #include "../FleshHist/Zprivate0.h"
18 
19  integer ndepth
20  parameter(ndepth= nsites)
21  real nrfaiRec(nrbin, nfai, 4, ndepth)
22  real nrfaiAll(nrbin, nfai, 4, ndepth)
23  real dErfai(nrbin, nfai, ndepth)
24 
25  integer EvNo,EvNo0
26  integer fn0
27  integer kgetenv2
28  integer klena
29  real intdep(ndepth)
30  real E0, cosz, limit(4)
31  integer newfmt, nr, maxsites
32  real Emin
33  integer NN
34  integer nrbina, nfaia, ansites0
35  integer leng, i, j, k, l
36  integer i0, j0, k0
37  integer l0(ndepth)
38  integer icon
39  character*128 nrfai
40  character*128 input
41  integer code, subcode, charge, ldep
42  integer w2il(ndepth)
43  character*20 field(15)
44  character*128 input0
45  integer ridx, faiidx
46  real rinmu, fai, Ek, time, wx, wy, wz
47  fn0 = 11
48  leng = kgetenv2("NRFAIFILE", nrfai)
49  call copenfw2(fn0, nrfai, 1, icon)
50  if(icon .ne. 1) then
51  write(0,*) nrfai(1:leng)
52  if( icon .eq. 0) then
53  write(0,*) 'not exists'
54  else
55  write(0,*) ' cannot be opened '
56  endif
57  write(0,*) ' icon=',icon
58  stop 9999
59  else
60  write(0,*) nrfai(1:leng), ' opened'
61  endif
62 
63  do while(.true.)
64  input = ' '
65  read( fn0, '(a)', end=1000 ) input
66  if(input .ne. " ") then
67 
68  call ksplit(input0, 20, 15, field, nr)
69  if(nr .eq. 8) then
70  read(input0(1:klena(input0)), *)
71  * evno0, e0,nn, cosz, limit(1), nrbina, nfaia, ansites0
72  maxsites = ansites0
73  emin=500.d-6
74  newfmt = 0
75  elseif(nr .eq. 9 ) then
76  read(input0(1:klena(input0)), *)
77  * evno0, e0,nn, cosz, limit(1), nrbina, nfaia, ansites0,
78  * maxsites
79  emin=500.d-6
80  newfmt = 1
81  elseif( nr .eq. 13) then
82  read(input0(1:klena(input0)), *)
83  * evno0, e0,nn, cosz, emin, nrbina, nfaia, ansites0,
84  * maxsites, limit
85  newfmt = 2
86  endif
87  if(nrbina .ne. nrbin .or. nfaia .ne. nfai) then
88  write(0,*)' nrbina=',nrbina, 'or nfaia=',nfaia,
89  * ' differ from the def. in this prog'
90  stop 5555
91  endif
92 ! ********
93  do i = 1, ansites0
94  do j = 1, 4
95  do k= 1, nfai
96  read(fn0, '(3x, f7.1, 4i4)' )
97  * intdep(i), l0(i), i0, j0, k0
98  w2il( l0(i) ) = i
99  if(i0 .ne. i .or. j .ne. j0 .or. k .ne. k0) then
100  write(0,*) ' intdep, i0,j0,k0=',
101  * intdep(i), i0, j0, k0, ' strange'
102  stop 8888
103  endif
104  read(fn0, *)
105  * ( nrfairec(l,k,j,i), l=1,nrbin )
106  do l = 1, nrbin
107  nrfairec(l,k,j,i) = 0.
108  enddo
109  enddo
110  enddo
111  enddo
112 ! ************
113  do i = 1, ansites0
114  do j = 1, 4
115  do k = 1, nfai
116  read(fn0, '(3x,f7.1, 4i4)' )
117  * intdep(i), l0(i), i0, j0, k0
118  if(i0 .ne. i .or. j .ne. j0 .or. k .ne. k0) then
119  write(0,*) ' intdep, i0,j0,k0=',
120  * intdep(i), i0, j0, k0, ' strange'
121  stop 9876
122  endif
123  read(fn0, *)
124  * ( nrfaiall(l,k,j,i), l=1,nrbin )
125  enddo
126  enddo
127  enddo
128  do i = 1, ansites0
129  do k = 1, nfai
130  read(fn0, '(5x,f7.1, 3i4)' )
131  * intdep(i), l0(i), i0, k0
132  if(i0 .ne. i .or. k .ne. k0) then
133  write(0,*) ' intdep, i0,k0=',
134  * intdep(i), i0, k0, ' strange'
135  stop 9875
136  endif
137  read(fn0, *)
138  * ( derfai(l,k,i), l=1,nrbin )
139  enddo
140  enddo
141  endif
142  enddo
143  1000 continue
144 
145 ! nrfai for 1 event read
146 ! read .dat and fill "rec"
147  do while(.true.)
148  read(*, *, end=100)
149  * ldep, code, subcode,
150  * charge, ridx, faiidx,
151  * rinmu, fai, ek, time, wx, wy, wz
152  if(code .gt. 4 ) code = 4
153  i = w2il(ldep)
154  nrfairec(ridx, faiidx, code, i) =
155  * nrfairec(ridx, faiidx, code, i) + 1.0
156  enddo
157  100 continue
158 
159 ! outpute .nrfai
160  if(newfmt .eq. 0 ) then
161  0 write(*,
162  * '(i2,1pE11.3, 0p,i3, f8.4, 1pE11.3,0p,3i4)')
163  * evno, e0, nn, cosz, limit(1),
164  * nrbina, nfaia, ansites0
165  elseif( newfmt .eq. 1) then
166  0 write(*,
167  * '(i2,1pE11.3, 0p,i3, f8.4, 1pE11.3,0p, 4i4)')
168  * evno, e0, nn, cosz, limit(1),
169  * nrbina, nfaia, ansites0, maxsites
170  elseif( newfmt .eq. 2) then
171  0 write(*,
172  * '(i2,1pE11.3, 0p,i3, f8.4, 1pE11.3,0p, 4i4,4f10.0)')
173  * evno, e0, nn, cosz, emin,
174  * nrbina, nfaia, ansites0, maxsites, limit
175  endif
176  do i = 1, ansites0
177  do j = 1, 4
178  do k = 1, nfaia
179  write(*, '("rec",f7.1, 4i4)' )
180  * intdep(i), l0(i), i, j, k
181  write(*, '(1p10E11.3)')
182  * ( nrfairec(l,k,j,i), l=1,nrbin )
183  enddo
184  enddo
185  enddo
186  do i = 1, ansites0
187  do j = 1, 4
188  do k = 1, nfaia
189  write(*, '("all",f7.1, 4i4)' )
190  * intdep(i), l0(i), i, j, k
191  write(*, '(1p10E11.3)')
192  * ( nrfaiall(l,k,j,i), l=1,nrbin )
193  enddo
194  enddo
195  enddo
196 
197  do i = 1, ansites0
198  do k = 1, nfaia
199  write(*, '("dE/dx",f7.1, 3i4)' )
200  * intdep(i), l0(i), i, k
201  write(*, '(1p10E11.3)')
202  * ( derfai(l,k,i), l=1,nrbin )
203  enddo
204  enddo
205  write(*, *)
206 
207  write(0,*) 'all data in the event processed '
208  end
parameter
integer npitbl real *nx parameter(n=101, npitbl=46, nx=n-1) real *8 uconst
ksplit
subroutine ksplit(a, m, n, b, nr)
Definition: ksplit.f:2
kgetenv2
integer function kgetenv2(envname, envresult)
Definition: cgetLoginN.f:77
i
nodes i
Definition: cstdatmos0.fNew2.h:93
time
subroutine time(xxx)
Definition: chook.f:5
true
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
leng
integer leng
Definition: interface2.h:1
nrbin
! timing nrbin
Definition: Zprivate2.h:12
j
********************block data cblkHeavy ********************integer j data *HeavyG2symbol *data *HeavyG2code kiron data j
Definition: cblkHeavy.h:36
copenfw2
subroutine copenfw2(io, fnin, form, icon)
Definition: copenf.f:205
code
*Zfirst p fm *Zfirst p Zfirst p code
Definition: ZavoidUnionMap.h:1
d
dE dx *! Nuc Int sampling table d
Definition: cblkMuInt.h:130
charge
*Zfirst p fm *Zfirst p Zfirst p Zfirst p *Zfirst p charge
Definition: ZavoidUnionMap.h:1
klena
integer function klena(cha)
Definition: klena.f:20
nfai
! timing nfai
Definition: Zprivate2.h:12
subcode
*Zfirst p fm *Zfirst p Zfirst p Zfirst p subcode
Definition: ZavoidUnionMap.h:1