COSMOS v7.655  COSMOSv7655
(AirShowerMC)
cinipipx.f File Reference

Go to the source code of this file.

Functions/Subroutines

subroutine cinipipx
 
real *8 function cintdndx2 (x)
 

Function/Subroutine Documentation

◆ cinipipx()

subroutine cinipipx ( )

Definition at line 2 of file cinipipx.f.

References cerrormsg(), cintdndx2(), d, d0, dx, false, intendndx2, kbinchop(), ktrpzintt2(), n, ndndx, pipsx, uconst, and xval.

Referenced by cs2lp().

2  implicit none
3 ! initialize leading particle x sampling for pi+p-->pi leading
4 !
5 #include "Zcinippxc.h"
6 
7  external cintdndx2
8 
9  integer i, icon
10  logical first/.true./
11  save first
12 
13  real*8 s
14  real*8 x, x1, eps, ans
15  real*8 xa(npitbl)
16  real*8 dndx(npitbl) ! unnormalized dn/dx of leading particle at
17  data xa/0,
18  * .01, .02, .03, .04, .05, .06, .07, .08, .09, .1,
19  * .11, .15, .20, .25, .30, .35, .40, .45, .50, .55,
20  * .60, .65, .70, .725, .75, .77, .78, .79, .80,.82,
21  * .84, .85, .87, .88, .90, .91, .92, .93, .94, .95,
22  * .96, .97, .98, .99, 1.00/
23  data dndx/0,
24  * 31., 38., 43., 48., 52., 53.5, 55.5, 58., 59., 60.,
25  * 58.5, 58., 52, 44, 36, 30, 25, 21, 17, 14.5,
26  * 12.5, 11, 9.6, 9.0, 8.8, 8.6, 8.5, 8.4, 8.4, 8.4,
27  * 8.6, 8.8, 9.2, 9.5, 10.1, 11.5, 12.2, 14.2, 16., 21.,
28  * 26.0, 33., 50, 80, 130./
29 
30  eps = 1.d-4
31  if(first) then
32 ! move xa into xval;
33  do i = 1, npitbl
34  xval(i) = xa(i)
35  enddo
36  first = .false.
37  endif
38 
39 ! integral of dn/dx from 0 to 1.
40  call ktrpzintt2(dndx, 1, npitbl, xa, 1, 0.d0, 1.d0, s)
41 !
42 ! make normalized table
43 !
44  do i = 1, npitbl
45  ndndx(i) = dndx(i)/s
46  enddo
47 ! make table of intgral(0:x) of ndndx. for x = 0 to 1.0 step
48 ! 0.01
49  x = 0.
50  dx = 0.01d0
51  do i=2, n-1
52  x = x + dx
53  call ktrpzintt2(ndndx, 1, npitbl, xa, 1, 0.d0,
54  * x, intendndx2(i) )
55  enddo
56  intendndx2(1) = 0.
57  intendndx2(n) = 1.
58 ! solve inte(0:x) of ndndx = u for u = 0 to 1.0 step 0.01
59  x = dx
60  uconst = 0.
61  x1 = 0.
62  do i = 2, n-1
63  uconst = uconst + dx
64  call kbinchop(cintdndx2, x1, 1.d0, x, eps, ans, icon)
65  if(icon .ne. 0) then
66  call cerrormsg(
67  * 'failed in making pi-leading particle sampling table', 0)
68  endif
69  pipsx(i) = ans
70  x = ans
71  x1 = x
72  enddo
73  pipsx(1) = 0.
74  pipsx(n) = 1.
subroutine cerrormsg(msg, needrtn)
Definition: cerrorMsg.f:4
block data include Zlatfit h c fitting region data x1(1)/0.03/
real *8 function cintdndx2(x)
Definition: cinipipx.f:77
atmos%rho(atmos%nodes) **exp(-(z-atmos%z(atmos%nodes))/Hinf) elseif(z .lt. atmos%z(1)) then ans=atmos%rho(1) **exp((atmos%z(1) -z)/atmos%H(1)) else call kdwhereis(z, atmos%nodes, atmos%z, 1, i) a=atmos%a(i) if(a .ne. 0.d0) then ans=atmos%rho(i) **(1+a *(z-atmos%z(i))/atmos%H(i)) **(-1.0d0-1.d0/a) else ans=*atmos%rho(i) *exp(-(z-atmos%z(i))/atmos%H(i)) endif endif ! zsave=z ! endif cvh2den=ans end ! ---------------------------------- real *8 function cvh2temp(z) implicit none ! vettical height to temperatur(Kelvin) real *8 z ! input. vertical height in m ! output is temperature of the atmospher in Kelvin real *8 ans integer i if(z .gt. atmos%z(atmos%nodes)) then ans=atmos%T(atmos%nodes) elseif(z .lt. atmos%z(1)) then ans=atmos%T(1)+atmos%b(1) *(z - atmos%z(1)) else call kdwhereis(z, atmos%nodes, atmos%z, 1, i) ans=atmos%T(i)+atmos%b(i) *(z-atmos%z(i)) endif cvh2temp=ans end !--------------------------------------------- real *8 function cthick2h(t) implicit none real *8 t ! input. air thickness in kg/m^2 real *8 logt, ans integer i real *8 dod0, fd, a logt=log(t) if(t .ge. atmos%cumd(1)) then ans=atmos%z(1) - *(logt - atmos%logcumd(1)) *atmos%H(1) elseif(t .le. atmos%cumd(atmos%nodes)) then ans=atmos%z(atmos%nodes) - *Hinf *log(t/atmos%cumd(atmos%nodes)) else call kdwhereis(t, atmos%nodes, atmos%cumd, 1, i) ! i is such that X(i) > x >=x(i+1) ans
logical, save first
Definition: cNRLAtmos.f:8
nodes i
integer npitbl real *nx dx real dx
Definition: Zcinippxc.h:10
integer npitbl real *nx dx real ndndx
Definition: Zcinippxc.h:10
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
integer npitbl real *nx dx real xval
Definition: Zcinippxc.h:10
subroutine kbinchop(f, x1, x2, x, eps, ans, icon)
Definition: kbinChop.f:20
integer npitbl real *nx dx real uconst
Definition: Zcinippxc.h:10
subroutine ktrpzintt2(t, intv, n, xt, intvx, a, b, ans)
Definition: ktrpzIntT2.f:5
block data cblkEvhnp ! currently usable models data RegMdls ad *special data *Cekaon d0
Definition: cblkEvhnp.h:5
dE dx *! Nuc Int sampling table d
Definition: cblkMuInt.h:130
integer npitbl real *nx dx real * intendndx2
Definition: Zcinippxc.h:10
integer npitbl real *nx dx real pipsx
Definition: Zcinippxc.h:10
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
integer n
Definition: Zcinippxc.h:1
! 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
Here is the call graph for this function:
Here is the caller graph for this function:

◆ cintdndx2()

real*8 function cintdndx2 ( real*8  x)

Definition at line 77 of file cinipipx.f.

References d0, ktrpzintt2(), ndndx, uconst, and xval.

Referenced by cinipipx().

77  implicit none
78  real*8 x
79 
80 #include "Zcinippxc.h"
81 
82  real*8 ans
83 
84  call ktrpzintt2(ndndx, 1, npitbl, xval, 1, 0.d0, x, ans)
85  cintdndx2 = ans - uconst
real *8 function cintdndx2(x)
Definition: cinipipx.f:77
atmos%rho(atmos%nodes) **exp(-(z-atmos%z(atmos%nodes))/Hinf) elseif(z .lt. atmos%z(1)) then ans=atmos%rho(1) **exp((atmos%z(1) -z)/atmos%H(1)) else call kdwhereis(z, atmos%nodes, atmos%z, 1, i) a=atmos%a(i) if(a .ne. 0.d0) then ans=atmos%rho(i) **(1+a *(z-atmos%z(i))/atmos%H(i)) **(-1.0d0-1.d0/a) else ans=*atmos%rho(i) *exp(-(z-atmos%z(i))/atmos%H(i)) endif endif ! zsave=z ! endif cvh2den=ans end ! ---------------------------------- real *8 function cvh2temp(z) implicit none ! vettical height to temperatur(Kelvin) real *8 z ! input. vertical height in m ! output is temperature of the atmospher in Kelvin real *8 ans integer i if(z .gt. atmos%z(atmos%nodes)) then ans=atmos%T(atmos%nodes) elseif(z .lt. atmos%z(1)) then ans=atmos%T(1)+atmos%b(1) *(z - atmos%z(1)) else call kdwhereis(z, atmos%nodes, atmos%z, 1, i) ans=atmos%T(i)+atmos%b(i) *(z-atmos%z(i)) endif cvh2temp=ans end !--------------------------------------------- real *8 function cthick2h(t) implicit none real *8 t ! input. air thickness in kg/m^2 real *8 logt, ans integer i real *8 dod0, fd, a logt=log(t) if(t .ge. atmos%cumd(1)) then ans=atmos%z(1) - *(logt - atmos%logcumd(1)) *atmos%H(1) elseif(t .le. atmos%cumd(atmos%nodes)) then ans=atmos%z(atmos%nodes) - *Hinf *log(t/atmos%cumd(atmos%nodes)) else call kdwhereis(t, atmos%nodes, atmos%cumd, 1, i) ! i is such that X(i) > x >=x(i+1) ans
integer npitbl real *nx dx real ndndx
Definition: Zcinippxc.h:10
integer npitbl real *nx dx real xval
Definition: Zcinippxc.h:10
integer npitbl real *nx dx real uconst
Definition: Zcinippxc.h:10
subroutine ktrpzintt2(t, intv, n, xt, intvx, a, b, ans)
Definition: ktrpzIntT2.f:5
block data cblkEvhnp ! currently usable models data RegMdls ad *special data *Cekaon d0
Definition: cblkEvhnp.h:5
! 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
Here is the call graph for this function:
Here is the caller graph for this function: