COSMOS v7.655  COSMOSv7655
(AirShowerMC)
ctotx.f
Go to the documentation of this file.
1 ! ctotx: x A collsion total cross section
2 !
3  module modcosxsec
4  logical,save:: nisdiff=.false.
5  end module modcosxsec
6  subroutine ctotx(pjin, Ain, xs)
7  implicit none
8 #include "Zptcl.h"
9 #include "Zcode.h"
10 #include "Zmass.h"
11 #include "Zevhnp.h"
12  type(ptcl)::pjin !input projectile particle
13  real(8),intent(in):: Ain ! input. effective target mass number
14  real(8),intent(out):: xs ! output. in mb. total xsection.
15 ! ctotx0 gives small Xs so we renormalize by
16 ! using cPDGsigmaTotpA(A) at @200 GeV
17  real(8),intent(in):: Zin
18  real(8):: cPDGsigmaTotpA
19 
20  type(ptcl):: pj
21  type(ptcl):: proton
22  save proton
23  real(8),save:: xst
24  logical,save::first=.true.
25  real(8),save::Asave=-1.
26  real(8),save::ratio
27  real(8):: Z, A
28  real(8),external:: cA2Z
29 
30  z = ca2z(ain)
31  goto 20
32 
33 
34  entry ctotx2(pjin, ain, zin, xs)
35 
36 
37  z = zin
38  20 continue
39  pj = pjin
40  a = ain
41  if( pjin%code == kgnuc .and. ain == 1.0) then
42  call cmkptc(knuc, -1, int(z), pj)
43  pj%fm%p(4) = pjin%fm%p(4)/pjin%subcode
44  a = pjin%subcode
45  z = pjin%charge
46  endif
47 
48 
49  call ctotx0(pj, a, z, xs)
50  if( a /= 1.0d0 ) then
51  if(first) then
52  call cmkptc(knuc, -1, 1, proton)
53  proton%fm%p(4) = 200. + masp
54  first = .false.
55  endif
56  if( a /= asave ) then
57  call ctotx0(proton, a, z, xst)
58  ratio = cpdgsigmatotpa(a)/xst
59  asave = a
60  endif
61  xs = xs *ratio
62  endif
63  end
64 
65  subroutine ctotx0(pj, A, Z, xs)
66  use modcosxsec
67  implicit none
68 #include "Zptcl.h"
69 #include "Zcode.h"
70 #include "Zevhnp.h"
71  type(ptcl)::pj !input projectile particle
72  real(8),intent(in):: A ! effective target mass number
73  real(8),intent(in):: Z ! effective target charge
74  real(8),intent(out)::xs ! in mb. total xsection.
75 
76  real*8 p
77 
78  real(8)::shp, shn
79  if( pj%fm%p(4) .le. pj%mass) then
80  if(pj%code .eq. knuc .and. pj%subcode .eq. antip) then
81  xs = largexs
82  elseif(pj%code .eq. kelec .and. pj%charge .eq. 1) then
83  xs = largexs
84  else
85  xs = smallxs
86  endif
87  else
88  p = sqrt(pj%fm%p(4)**2 - pj%mass**2)
89 ! p = max(p, 0.1d0)
90 ! if( p .ge. 20.) then
91 ! call cerrorMsg('Momentum is >20 GeV/c for ctotx',0)
92 ! endif
93 
94  if(pj%code .eq. knuc) then
95  if(pj%charge .eq. 1) then
96 ! proton
97 ! shp = ctotpp1(p)
98  if( a == 1.0 ) then
99  if( z == 1.0) then ! pp
100  call cpptotxs(p, xs)
101  elseif( z== 0.0 ) then ! pn
102  call cnptotxs(p, xs)
103  else
104  write(0,*) 'pj=proton target A,Z=',a,z
105  write(0,*) ' stragne for ctotx/ctotx2'
106  stop
107  endif
108  elseif(a > 1.) then
109  call cpptotxs(p, shp)
110  if( nisdiff ) then
111  call cnptotxs(p, shn)
112 ! use average of pp,pn
113  shp = (shp + shn)/2.0
114  endif
115  call cxp2xaxsec(a, shp, xs)
116  else
117  write(0,*) 'pj: p, target A,Z=',a,z
118  write(0,*) 'strange for ctotx/ctotx2'
119  stop
120  endif
121  elseif(pj%charge .eq. -1) then
122 ! pbar no n target
123  call cpbarptotxs(p, shp)
124  if( a .gt. 1.) then
125  call cxp2xaxsec(a, shp, xs)
126  else
127  xs = shp
128  endif
129  elseif(pj%subcode .eq. antip) then
130 ! anti-neutron; assume the same one as pbar
131 ! shp = ctotpbarp1(p)
132  call cnbarptotxs(p, shp)
133  if( a .gt. 1.) then
134  call cxp2xaxsec(a, shp, xs)
135  else
136  xs = shp
137  endif
138  else
139 ! neutron
140 ! shp = ctotnp1(p)
141  call cnptotxs(p, shp)
142  if(a .gt. 1.) then
143  call cxp2xaxsec(a, shp, xs)
144  else
145  xs = shp
146  endif
147  endif
148  elseif(pj%code .eq. kpion) then
149  if(pj%charge .eq. 1) then
150  call cpipptotxs(p, shp)
151  elseif(pj%charge .eq. -1) then
152  call cpimptotxs(p, shp)
153  else
154 ! at low energy, pi0 would not interact. any would be o.k
155 ! shp = ctotpiMp1(p)
156  call cpimptotxs(p, shp)
157  endif
158  if(a .ne. 1.0) then
159  call cxp2xaxsec(a, shp, xs)
160  else
161  xs = shp
162  endif
163  elseif(pj%code .eq. kkaon) then
164  if(pj%charge .eq. 1) then
165  if( a == 1.0 ) then
166  if( z == 1.0 ) then
167  call ckpptotxs(p, xs)
168  elseif( z == 0. ) then
169  call ckpntotxs(p, xs)
170  else
171  write(0,*) ' A,Z=',a,z, 'for K+ pj'
172  write(0,*) ' strange to ctotx/ctotx2 '
173  stop
174  endif
175  elseif(a > 1.) then
176  call ckpptotxs(p, shp)
177  if( nisdiff ) then
178  call ckpntotxs(p, shn)
179  shp = (shp + shn)/2.0
180  endif
181  call cxp2xaxsec(a, shp, xs)
182  else
183  write(0,*) ' A,Z=',a,z, 'for K+ p'
184  write(0,*) ' strange to ctotx/ctotx2 '
185  stop
186  endif
187  elseif(pj%charge .eq. -1) then ! K-
188  call ckmptotxs(p, shp) ! kmp = kmn at present
189  if(a .gt. 1.) then
190  if( nisdiff ) then
191  call ckmntotxs(p, shn)
192  shp = (shp + shn)/2.
193  endif
194  call cxp2xaxsec(a, shp, xs)
195  else
196  xs = shp
197  endif
198  else
199 ! k0; don't worry so much
200  call ckmptotxs(p, shp)
201  call ckmntotxs(p, shn)
202  shp = (shp+ shn)/2.
203  if(a .gt. 1.) then
204  call cxp2xaxsec(a, shp, xs)
205  else
206  xs = shp
207  endif
208  endif
209  elseif(pj%code .eq. kgzai .or. pj%code .eq. ksigma .or.
210  * pj%code .eq. kbomega .or. pj%code .eq. klambda .or.
211  * pj%code .eq. klambdac) then
212 ! don't worry, not used almost at all; use proton
213 ! shp = ctotpp1(p)
214  call cpptotxs(p, shp)
215  call cxp2xaxsec(a, shp, xs)
216  elseif( pj%code == kgnuc ) then
217  write(0,*) 'Sorry: ctotx is not usable for heavy ions'
218  write(0,*) ' only inelastic xs is used for heavy ions'
219  write(0,*) ' so that you may use cinela or cAAXsec2'
220  stop
221  else
222 ! use pion
223 ! shp =ctotpiPp1(p)
224  call cpipptotxs(p, shp)
225  call cxp2xaxsec(a, shp, xs)
226  endif
227  endif
228  end
229 
230  function ca2z(A) result(Z)
231  implicit none
232 ! very rough charge assignment for mass # A nucleus
233  real(8),intent(in):: A ! mass #
234 
235  real(8)::Z
236  z = int( a/2.15 + 0.7 )
237  end function ca2z
238 
239 
kgzai
max ptcl codes in the kgzai
Definition: Zcode.h:2
cpbarptotxs
subroutine cpbarptotxs(p, xs)
Definition: cpbarpTotXs.f:9
cpptotxs
subroutine cpptotxs(p, xs)
Definition: cppTotXs.f:3
kgnuc
max ptcl codes in the kgnuc
Definition: Zcode.h:2
cnptotxs
subroutine cnptotxs(p, xs)
Definition: cnpTotXs.f:3
klambdac
max ptcl codes in the klambdac
Definition: Zcode.h:2
ckpptotxs
subroutine ckpptotxs(p, xs)
Definition: ckppTotXs.f:3
kkaon
max ptcl codes in the kkaon
Definition: Zcode.h:2
kelec
max ptcl codes in the kelec
Definition: Zcode.h:2
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
cxp2xaxsec
subroutine cxp2xaxsec(a, xsxp, xsxa)
Definition: cxp2xAXsec.f:30
ckpntotxs
subroutine ckpntotxs(p, xs)
Definition: ckpnTotXs.f:3
modcosxsec::nisdiff
logical, save nisdiff
Definition: ctotx.f:4
ctotx
subroutine ctotx(pjin, Ain, xs)
Definition: ctotx.f:7
ckmptotxs
subroutine ckmptotxs(p, xs)
Definition: ckmpTotXs.f:3
d0
block data cblkEvhnp ! currently usable models data RegMdls ad *special data *Cekaon d0
Definition: cblkEvhnp.h:5
klambda
max ptcl codes in the klambda
Definition: Zcode.h:2
knuc
********************block data cblkHeavy ********************integer j data *HeavyG2symbol *data *HeavyG2code knuc
Definition: cblkHeavy.h:7
ca2z
real(8) function ca2z(A)
Definition: ctotx.f:231
ckmntotxs
subroutine ckmntotxs(p, xs)
Definition: ckmnTotXs.f:2
cpipptotxs
subroutine cpipptotxs(p, xs)
Definition: cpippTotXs.f:3
false
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
cmkptc
subroutine cmkptc(code, subcode, charge, p)
Definition: cmkptc.f:15
masp
masp
Definition: Zmass.h:5
ptcl
Definition: Zptcl.h:75
antip
max ptcl codes in the kseethru ! subcode integer antip
Definition: Zcode.h:2
kpion
max ptcl codes in the kpion
Definition: Zcode.h:2
ksigma
max ptcl codes in the ksigma
Definition: Zcode.h:2
kbomega
max ptcl codes in the kbomega
Definition: Zcode.h:2
cpimptotxs
subroutine cpimptotxs(p, xs)
Definition: cpimpTotXs.f:3
ctotx0
subroutine ctotx0(pj, A, Z, xs)
Definition: ctotx.f:66
cnbarptotxs
subroutine cnbarptotxs(p, xs)
Definition: cnbarpTotXs.f:3