IAP GITLAB

Commit 3874cb68 authored by Tanguy Pierog's avatar Tanguy Pierog

update Sibyll to version 2.3


git-svn-id: https://devel-ik.fzk.de/svn/mc/crmc/trunk@5268 c7a5e08c-de06-0410-9364-b41cf42a0b17
parent 813c6f4b
......@@ -29,8 +29,8 @@ MESSAGE(STATUS "build ${CMAKE_BINARY_DIR}")
## The version number. Only nummeric PLEASE!
SET (CRMC_VERSION_MAJOR 1)
SET (CRMC_VERSION_MINOR 5)
SET (CRMC_VERSION_PATCH 6)
SET (CRMC_VERSION_MINOR 6)
SET (CRMC_VERSION_PATCH 0)
##Compiler options
# hotfix for lxbatch not sure why he selects wrong compiler
......
CRMC v1.5.7 Last modifications 2016/03/23
CRMC v1.6.0 Last modifications 2016/03/29
**********************************************************
The Program "crmc"
......@@ -57,7 +57,9 @@ QGSJET01 (-m2) : N.N. Kalmykov, S. Ostapchenko, and A.I. Pavlov,
QGSJETII-03 (-m11) : S. Ostapchenko,
Nucl.Phys.Proc.Suppl. 151 (2006) 143
SIBYLL2.1 (-m6) : R. Engel, T.K. Gaisser, P. Lipari, and T. Stanev,
SIBYLL2.3 (-m6) : ???
and
R. Engel, T.K. Gaisser, P. Lipari, and T. Stanev,
Proc. 26th Int. Cosmic Ray Conf., Salt Lake City (USA), 1 (1999) 415
and
E.-J. Ahn, R. Engel, T.K. Gaisser, P. Lipari, and T. Stanev,
......
......@@ -2020,7 +2020,7 @@ c-------------------------------------------------------------------------------
integer istatus
character*3 code1,code2
parameter (ncode=5,nidt=341)
parameter (ncode=5,nidt=351)
integer idt(ncode,nidt)
double precision drangen,dummy
......@@ -2066,10 +2066,10 @@ c nxs|pdg|qgs|cor|sib
* , 231,313,99,62,30 !k*0
* , -231,-313,99,65,31 !k*0b
* , 331,333,99,99,33 !phi
* , -140,421,8,116,99 !D0(1.864)
* , 140,-421,8,119,99 !D0b(1.864)
* , -240,411,7,117,99 !D(1.869)+
* , 240,-411,7,118,99 !Db(1.869)-
* , -140,421,8,116,71 !D0(1.864)
* , 140,-421,8,119,72 !D0b(1.864)
* , -240,411,7,117,59 !D(1.869)+
* , 240,-411,7,118,60 !Db(1.869)-
* , 1120,2212,2,14,13 !proton
* , 1220,2112,3,13,14 !neutron
* , 2130,3122,6,18,39 !lambda
......@@ -2088,33 +2088,33 @@ c nxs|pdg|qgs|cor|sib
* , 1331, 3324,99,99,47 !xi*0
* , 2331, 3314,99,99,48 !xi*-
* , 3331, 3334,99,24,49 !omega-
* , 2140, 4122,9,137,99 !LambdaC(2.285)+
* , 2140, 4122,9,137,89 !LambdaC(2.285)+
* ,17,1000010020,99,201,1002 ! Deuteron
* ,18,1000010030,99,301,1003 ! Triton
* ,19,1000020040,99,402,1004 ! Alpha
* ,0,0,99,0,0 ! Air
* ,99,99,99,99,99 / ! unknown
* ,99,99,99,99,100 / ! unknown
data ((idt(i,j),i=1,ncode),j= 69,89)/
$ -340,431,99,120,99 ! Ds+
$ ,340,-431,99,121,99 ! Ds-
$ ,-241,413,99,124,99 ! D*+
$ ,241,-413,99,125,99 ! D*-
$ ,-141,423,99,123,99 ! D*0
$ ,141,-423,99,126,99 ! D*0b
$ ,-341,433,99,127,99 ! Ds*+
$ ,341,-433,99,128,99 ! Ds*-
$ ,440,441,99,122,99 ! etac
$ ,441,443,99,130,99 ! J/psi
$ ,2240,4112,99,142,99 ! sigmac0
$ ,1240,4212,99,141,99 ! sigmac+
$ ,1140,4222,99,140,99 ! sigmac++
$ ,2241,4114,99,163,99 ! sigma*c0
$ ,1241,4214,99,162,99 ! sigma*c+
$ ,1141,4224,99,161,99 ! sigma*c++
$ ,3240,4132,99,139,99 ! Xic0
$ ,2340,4312,99,144,99 ! Xi'c0
$ ,3140,4232,99,138,99 ! Xic+
$ ,1340,4322,99,143,99 ! Xi'c+
$ -340,431,99,120,74 ! Ds+
$ ,340,-431,99,121,75 ! Ds-
$ ,-241,413,99,124,78 ! D*+
$ ,241,-413,99,125,79 ! D*-
$ ,-141,423,99,123,80 ! D*0
$ ,141,-423,99,126,81 ! D*0b
$ ,-341,433,99,127,76 ! Ds*+
$ ,341,-433,99,128,77 ! Ds*-
$ ,440,441,99,122,73 ! etac
$ ,441,443,99,130,83 ! J/psi
$ ,2240,4112,99,142,86 ! sigmac0
$ ,1240,4212,99,141,85 ! sigmac+
$ ,1140,4222,99,140,84 ! sigmac++
$ ,2241,4114,99,163,96 ! sigma*c0
$ ,1241,4214,99,162,95 ! sigma*c+
$ ,1141,4224,99,161,94 ! sigma*c++
$ ,3240,4132,99,139,88 ! Xic0
$ ,2340,4312,99,144,98 ! Xi'c0
$ ,3140,4232,99,138,87 ! Xic+
$ ,1340,4322,99,143,97 ! Xi'c+
$ ,3340,4332,99,145,99 / ! omegac0
data ((idt(i,j),i=1,ncode),j= 90,nidt)/
$ 1112,32224,99,99,99 ! Delta(1600)++
......@@ -2135,7 +2135,7 @@ c nxs|pdg|qgs|cor|sib
$ ,2224,21114,99,99,99 ! Delta(1920)-
$ ,2224,11116,99,99,99 ! Delta(1930)-
$ ,2224, 1118,99,99,99 ! Delta(1950)-
$ ,1122,12212,99,99,99 ! N(1440)+
$ ,1122,12212,99,99,51 ! N(1440)+
$ ,1123, 2124,99,99,99 ! N(1520)+
$ ,1123,22212,99,99,99 ! N(1535)+
$ ,1124,32214,99,99,99 ! Delta(1600)+
......@@ -2145,7 +2145,7 @@ c nxs|pdg|qgs|cor|sib
$ ,1125,12216,99,99,99 ! N(1680)+
$ ,1126,12214,99,99,99 ! Delta(1700)+
$ ,1127,22124,99,99,99 ! N(1700)+
$ ,1127,42212,99,99,99 ! N(1710)+
$ ,1127,42212,99,99,53 ! N(1710)+
$ ,1127,32124,99,99,99 ! N(1720)+
$ ,1128,12122,99,99,99 ! Delta(1900)+
$ ,1128, 2126,99,99,99 ! Delta(1905)+
......@@ -2153,7 +2153,7 @@ c nxs|pdg|qgs|cor|sib
$ ,1128,22214,99,99,99 ! Delta(1920)+
$ ,1128,12126,99,99,99 ! Delta(1930)+
$ ,1128, 2218,99,99,99 ! Delta(1950)+
$ ,1222,12112,99,99,99 ! N(1440)0
$ ,1222,12112,99,99,52 ! N(1440)0
$ ,1223, 1214,99,99,99 ! N(1520)0
$ ,1223,22112,99,99,99 ! N(1535)0
$ ,1224,32114,99,99,99 ! Delta(1600)0
......@@ -2163,7 +2163,7 @@ c nxs|pdg|qgs|cor|sib
$ ,1225,12116,99,99,99 ! N(1680)0
$ ,1226,12114,99,99,99 ! Delta(1700)0
$ ,1227,21214,99,99,99 ! N(1700)0
$ ,1227,42112,99,99,99 ! N(1710)0
$ ,1227,42112,99,99,54 ! N(1710)0
$ ,1227,31214,99,99,99 ! N(1720)0
$ ,1228,11212,99,99,99 ! Delta(1900)0
$ ,1228, 1216,99,99,99 ! Delta(1905)0
......@@ -2254,8 +2254,8 @@ c nxs|pdg|qgs|cor|sib
$ ,451,543,99,99,99 ! B*c+
$ ,550,551,99,99,99 ! etab
$ ,551,553,99,99,99 ! Upsilon
$ ,2341,4314,99,99,99 ! Xi*c0
$ ,1341,4324,99,99,99 ! Xi*c+
$ ,2341,4314,99,99,99 ! Xi*c0(2645)
$ ,1341,4324,99,99,99 ! Xi*c+(2645)
$ ,3341,4334,99,99,99 ! omega*c0
$ ,2440,4412,99,99,99 ! dcc
$ ,2441,4414,99,99,99 ! dcc*
......@@ -2299,24 +2299,21 @@ c nxs|pdg|qgs|cor|sib
$ ,3551,5534,99,99,99 ! sbb*
$ ,4551,5544,99,99,99 ! cbb*
$ ,5551,5554,99,99,99 ! bbb*
$ ,123,10213,99,99,99 ! b1
$ ,122,10211,99,99,99 ! a0+
$ ,233,10313,99,99,99 ! K0_1
$ ,232,10311,99,99,99 ! K*0_1
$ ,133,10323,99,99,99 ! K+_1
$ ,132,10321,99,99,99 ! K*+_1
$ ,123,10213,99,99,99 ! b_1+
$ ,122,10211,99,99,62 ! a_0+
$ ,-122,-10211,99,99,63 ! a_0-
$ ,232,10311,99,99,66 ! K*0_1
$ ,-232,-10311,99,99,67 ! K*0b_1
$ ,132,10321,99,99,64 ! K*+_1
$ ,-132,-10321,99,99,65 ! K*-_1
$ ,143,10423,99,99,99 ! D0_1
$ ,132,10421,99,99,99 ! D*0_1
$ ,142,10421,99,99,99 ! D*0_1
$ ,243,10413,99,99,99 ! D+_1
$ ,242,10411,99,99,99 ! D*+_1
$ ,343,10433,99,99,99 ! D+s_1
$ ,342,10431,99,99,99 ! D*0s+_1
$ ,223,10113,99,99,99 ! b_10
$ ,222,10111,99,99,99 ! a_00
$ ,113,10223,99,99,99 ! h_10
$ ,112,10221,99,99,99 ! f_00
$ ,333,10333,99,99,99 ! h'_10
$ ,332,10331,99,99,99 ! f'_00
$ ,113,10113,99,99,99 ! b_10
$ ,112,10111,99,99,61 ! a_00
$ ,443,10443,99,99,99 ! h_1c0
$ ,442,10441,99,99,99 ! Xi_0c0
$ ,444,10443,99,99,99 ! psi'
......@@ -2332,22 +2329,31 @@ c nxs|pdg|qgs|cor|sib
$ ,552,10551,99,99,99 ! Upsilon'*
$ ,124,20213,99,99,99 ! a_1+
$ ,125,215,99,99,99 ! a_2+
$ ,234,20313,99,99,99 ! K*0_1
$ ,235,315,99,99,99 ! K*0_2
$ ,134,20323,99,99,99 ! K*+_1
$ ,135,325,99,99,99 ! K*+_2
$ ,126,10215,99,99,99 ! pi_2+(1670)
$ ,127,217,99,99,99 ! rho_3+(1690)
$ ,232,20313,99,99,99 ! K*0_1(1400)
$ ,233,315,99,99,99 ! K*0_2(1430)
$ ,234,10311,99,99,99 ! K*0_0(1430)
$ ,132,20323,99,99,99 ! K*+_1(1400)
$ ,133,325,99,99,99 ! K*+_2(1430)
$ ,134,10321,99,99,99 ! K*+_0(1430)
$ ,144,20423,99,99,99 ! D*_10
$ ,135,425,99,99,99 ! D*_20
$ ,145,425,99,99,99 ! D*_20
$ ,244,20413,99,99,99 ! D*_1+
$ ,245,415,99,99,99 ! D*_2+
$ ,344,20433,99,99,99 ! D*_1s+
$ ,345,435,99,99,99 ! D*_2s+
$ ,224,20113,99,99,99 ! a_10
$ ,225,115,99,99,99 ! a_20
$ ,114,20223,99,99,99 ! f_10
$ ,115,225,99,99,99 ! f_20
$ ,334,20333,99,99,99 ! f'_10
$ ,335,335,99,99,99 ! f'_20
$ ,114,20113,99,99,99 ! a_10
$ ,115,115,99,99,99 ! a_20
$ ,116,10115,99,99,99 ! pi_20(1670)
$ ,117,117,99,99,99 ! rho_30(1690)
$ ,222,9010221,99,99,99 ! f_00(980)
$ ,223,10223,99,99,99 ! h_10(1170)
$ ,224,225,99,99,99 ! f_20(1270)
$ ,225,20223,99,99,99 ! f_10(1285)
$ ,226,9030221,99,99,99 ! f_00(1500)
$ ,332,20333,99,99,99 ! f_10(1420)
$ ,333,335,99,99,99 ! f'_20(1525)
$ ,444,20443,99,99,99 ! Xi_1c0
$ ,445,445,99,99,99 ! Xi_2c0
$ ,254,20513,99,99,99 ! db*_10
......@@ -2362,11 +2368,15 @@ c nxs|pdg|qgs|cor|sib
$ ,555,555,99,99,99 ! bb*_20
$ ,11099,9900110,99,99,99 ! diff pi0 state
$ ,12099,9900210,99,99,99 ! diff pi+ state
$ ,13099,9900320,99,99,99 ! diff K+ state
$ ,22099,9900220,99,99,99 ! diff omega state
$ ,2099,9900310,99,99,99 ! diff K0 state
$ ,-2099,9900130,99,99,99 ! diff pi+ state
$ ,33099,9900330,99,99,99 ! diff phi state
$ ,44099,9900440,99,99,99 ! diff J/psi state
$ ,112099,9902210,99,99,99 ! diff proton state
$ ,122099,9902110,99,99,99 ! diff neutron state
$ ,213099,9903120,99,99,99 ! diff lambda state
$ ,800000110,110,99,99,99 ! Reggeon
$ ,800000990,990,99,99,99 / ! Pomeron
......@@ -2415,6 +2425,7 @@ c print *,'idtrafo',' ',code1,' ',code2,idi
endif
if(i.eq.2) nidtmx=nidt
if(i.eq.4) nidtmx=89
if(i.eq.5) nidtmx=nidt ! maximal reach in conversion table
elseif(code2.eq.'pdg')then
j=2
ji=j
......@@ -2489,6 +2500,8 @@ c print *,'idtrafo',' ',code1,' ',code2,idi
else
mm=int(drangen(dummy)*dble(m))
endif
else !m=0 only one line, take care of sign
if(idt(i,n).lt.0)isi=-isi
endif
idtrafo=idt(j,n+mm)*isi
if(abs(idtrafo).eq.99)call utstop('New particle not allowed ')
......
......@@ -324,6 +324,10 @@
double precision GTOT,GPROD,GABS,GDD,GQEL,GCOH
double precision e0
#ifdef __SIBYLL_
double precision engy_d,ssig,slope,rho,sigt_d,sige_d,sigqel_d,
& sigsd_d,sigqsd_d
#endif
dimension dumdif(3)
if(model.eq.2)then
#ifndef __QGSJET01__
......@@ -369,9 +373,15 @@
elseif(iclpro.eq.3)then
K=3
endif
CALL SIB_SIGMA_HP(K,engy,SSIG,dum0,dum1,dumdif,SLOPE,RHO)
CALL GLAUBER(matarg,SSIG,SLOPE,RHO,sigt,sige,sigqel)
sigi=sigt-sigqel
engy_d = dble(engy)
CALL SIB_SIGMA_HP(K,engy_d,SSIG,dum0,dum1,dumdif,SLOPE,RHO)
CALL GLAUBER2
& (matarg,SSIG,SLOPE,RHO,sigt_d,sige_d,sigqel_d,
& sigsd_d,sigqsd_d)
sigt=sngl(sigt_d)
sige=sngl(sige_d)
sigqel=sngl(sigqel_d)
sigi=sngl(sigt_d-sigqel_d)
sigc=sigi
#endif
elseif(model.eq.7.or.model.eq.11)then
......@@ -436,7 +446,11 @@
end
subroutine m6SIGMA(icl,engy,stot,sela,sine,sdifr,slela,Rho)
dimension sdifr0(3)
#ifdef __SIBYLL__
double precision engy_d,stot_d,sela_d,sine_d,slela_d,rho_d
double precision sdifr0_d(3)
#endif
#ifndef __SIBYLL__
stop'please compile with requested model'
print *, icl,engy,stot,sela,sine,sdifr,slela,Rho,sdifr0 !get rid of unused warning
......@@ -448,8 +462,15 @@
else
L=3
endif
call SIB_SIGMA_HP(L,engy,stot,sela,sine,sdifr0,slela,Rho)
sdifr=sdifr0(1)+sdifr0(2)+sdifr0(3)
engy_d = dble(engy)
CALL SIB_SIGMA_HP
& (L,engy_d,stot_d,sela_d,sine_d,sdifr0_d,slela_d,Rho_d)
stot = sngl( stot_d )
sela = sngl( sela_d )
sine = sngl( sine_d )
slela = sngl( slela_d )
Rho = sngl( Rho_d )
sdifr=sngl(sdifr0_d(1)+sdifr0_d(2)+sdifr0_d(3))
#endif
end
......
This source diff could not be displayed because it is too large. You can view the blob instead.
-----------------------
code particle mass
-----------------------
1 gam 0.0000
2 e+ 0.0005
3 e- 0.0005
4 mu+ 0.1057
5 mu- 0.1057
6 pi0 0.1350
7 pi+ 0.1396
8 pi- 0.1396
9 k+ 0.4937
10 k- 0.4937
11 k0l 0.4977
12 k0s 0.4977
13 p 0.9383
14 n 0.9396
15 nue 0.0000
16 nueb 0.0000
17 num 0.0000
18 numb 0.0000
19 empty
20 empty
21 k0 0.4977
22 k0b 0.4977
23 eta 0.5488
24 etap 0.9576
25 rho+ 0.7714
26 rho- 0.7714
27 rho0 0.7717
28 k*+ 0.8921
29 k*- 0.8921
30 k*0 0.8965
31 k*0b 0.8965
32 omeg 0.7826
33 phi 1.1020
34 SIG+ 1.1894
35 SIG0 1.1925
36 SIG- 1.1973
37 XI0 1.3149
38 XI- 1.3213
39 LAM 1.1156
40 DELT++ 1.2300
41 DELT+ 1.2310
42 DELT0 1.2320
43 DELT- 1.2330
44 SIG*+ 1.3828
45 SIG*0 1.3837
46 SIG*- 1.3872
47 XI*0 1.5318
48 XI*- 1.5350
49 OME*- 1.6724
Antibaryons have negative codes (antineutrons = -14 for example)
Other particles are unknown in this version.
A special numbering scheme has been introduced for nuclei.
A Helium nucleus (A=4) is indicated with the code L=1004
an Iron nucleus (A=56) has L=1056.
When a particle decays it is NOT deleted from the
event record, but it is kept there (in case one want to understand
what happened). In this case 10000 is added to
the code (or -10000 is subtracted if the code is negative)
of the particle that has decayed. The daughter particles are
added starting at the end of the event record.
This source diff could not be displayed because it is too large. You can view the blob instead.
This source diff could not be displayed because it is too large. You can view the blob instead.
......@@ -8,15 +8,19 @@ c-----------------------------------------------------------------------
c Primary initialization for Sibyll
c-----------------------------------------------------------------------
include 'epos.inc'
COMMON /S_DEBUG/ Ncall, Ndebug
COMMON /S_CSYDEC/ CBR(102), KDEC(612), LBARP(49), IDB(49)
INTEGER NCALL, NDEBUG, LUN
COMMON /S_DEBUG/ NCALL, NDEBUG, LUN
DOUBLE PRECISION CBR
INTEGER KDEC,LBARP,IDB
COMMON /S_CSYDEC/ CBR(223+16+12+8), KDEC(1338+6*(16+12+8)),
& LBARP(99), IDB(99)
call utpri('inisib',ish,ishini,6)
write(ifmt,'(a,i6)')'initialize Sibyll ...'
Ndebug=0
if(ish.ge.3)Ndebug=ish-2
lun = 7
C... SIBYLL initialization
CALL SIBYLL_INI
......@@ -24,7 +28,7 @@ C...Cross sections for nucleus-nucleus and hadron nucleus
CALL NUC_NUC_INI
C...define all particles as unstable
do i=1,49
do i=1,99
IDB(i) = abs(IDB(i)) ! >0 means unstable
enddo
......@@ -35,6 +39,7 @@ C...define all particles as unstable
call utprix('inisib',ish,ishini,6)
end
c-----------------------------------------------------------------------
......@@ -43,10 +48,12 @@ c-----------------------------------------------------------------------
c Initialization for each type of event (for given proj, targ and egy)
c-----------------------------------------------------------------------
include 'epos.inc'
COMMON /S_CSYDEC/ CBR(102), KDEC(612), LBARP(49), IDB(49)
DOUBLE PRECISION CBR
INTEGER KDEC,LBARP,IDB
COMMON /S_CSYDEC/ CBR(223+16+12+8), KDEC(1338+6*(16+12+8)),
& LBARP(99), IDB(99)
common/geom/rmproj,rmtarg,bmax,bkmx
if(matarg.gt.18.or.maproj.gt.64)
& call utstop('Mass too big for Sibyll (Mtrg<18, Mprj<64) !&')
id=idtrafo('nxs','sib',idproj)
......@@ -107,7 +114,7 @@ c (part taken from epos-dky: hdecas)
else
C...define all particles as stable
do i=1,49
do i=1,99
IDB(i) = -abs(IDB(i)) ! <0 means stable
enddo
......@@ -128,28 +135,38 @@ c-----------------------------------------------------------------------
include 'epos.inc'
common/geom/rmproj,rmtarg,bmax,bkmx
C SIBYLL
COMMON /S_PLIST/ P(8000,5), LLIST(8000), NP
COMMON /S_PLNUC/ PA(5,40000), LLA(40000), NPA
DOUBLE PRECISION P
INTEGER NP,LLIST,NP_max
PARAMETER (NP_max=8000)
COMMON /S_PLIST/ P(NP_max,5), LLIST(NP_max), NP
INTEGER NCALL, NDEBUG, LUN
COMMON /S_DEBUG/ NCALL, NDEBUG, LUN
INTEGER NW_max
PARAMETER (NW_max = 20)
PARAMETER (NS_max = 20, NH_max = 50)
PARAMETER (NJ_max = (NS_max+NH_max)*NW_max)
COMMON /S_CHIST/ X1J(NJ_max),X2J(NJ_max),
& X1JSUM(NW_max),X2JSUM(NW_max),PTJET(NJ_max),PHIJET(NJ_max),
& NNPJET(NJ_max),NNPSTR(2*NW_max),NNSOF(NW_max),NNJET(NW_max),
& JDIF(NW_max),NW,NJET,NSOF
INTEGER NNSOF,NNJET,JDIF,NWD,NJET,NSOF
COMMON /S_CHIST/ NNSOF(NW_max),NNJET(NW_max),
& JDIF(NW_max),NWD,NJET,NSOF
COMMON /S_PLNUC/ PA(5,40000), LLA(40000), NPA
double precision engy_dbl,ptm,pa
COMMON /S_CLDIF/ LDIFF
INTEGER LDIFF
iret=0
b1=bminim
b2=min(bmax,bmaxim)
a=pi*(b2**2-b1**2)
if(a.gt.0..and.rangen().gt.sibincs/10./a)goto 1001 !no interaction
if(ish.ge.3)call alist('Determine Sibyll Production&',0,0)
nptl=0
NP=0
NPA=0
LDIFF=0 !all types of events
nevt=1
kolevt=-1
koievt=-1
......@@ -170,11 +187,12 @@ C SIBYLL
call conre
call conwr
engy_dbl = dble(engy)
itrg=matarg
if(idtargin.eq.0)itrg=0
if(maproj.eq.1)then !hadronic projectile
L0=idtrafo('nxs','sib',idproj)
CALL SIBYLL (L0, itrg, engy)
CALL SIBYLL (L0, itrg, engy_dbl)
CALL DECSIB
if(ish.ge.5)write(ifch,'(a,i5)')
$ ' number of particles from Sibyll :',NP
......@@ -210,11 +228,11 @@ c LLIST is the code of final particle, P - its 4-momentum and mass.
$ ' epos particle ',nptl,' id :',id,' after conversion'
pptl(1,nptl)=P(k,1) !P_x
pptl(2,nptl)=P(k,2) !P_y
pptl(3,nptl)=P(k,3) !P_z
pptl(4,nptl)=abs(P(k,4)) !E
pptl(5,nptl)=P(k,5) !mass
pptl(1,nptl)=sngl(P(k,1)) !P_x
pptl(2,nptl)=sngl(P(k,2)) !P_y
pptl(3,nptl)=sngl(P(k,3)) !P_z
pptl(4,nptl)=abs(sngl(P(k,4))) !E
pptl(5,nptl)=sngl(P(k,5)) !mass
ityptl(nptl)=0
iorptl(nptl)=1
jorptl(nptl)=maproj+matarg
......@@ -235,9 +253,9 @@ c LLIST is the code of final particle, P - its 4-momentum and mass.
enddo
if(NW.lt.matarg)then
ntgevt=NW
do is=maproj+1+NW,maproj+matarg !make the ns last target nucleon final
if(NWD.lt.matarg)then
ntgevt=NWD
do is=maproj+1+NWD,maproj+matarg !make the ns last target nucleon final
iorptl(is)=0
istptl(is)=0
enddo
......@@ -245,7 +263,7 @@ c LLIST is the code of final particle, P - its 4-momentum and mass.
else !for nucleus projectile
nbar=0
IAP = maproj
CALL SIBNUC (IAP, itrg, engy)
CALL SIBNUC (IAP, itrg, engy_dbl)
if(ish.ge.5)write(ifch,'(a,i5)')
$ ' number of particles from Sibyll :',NPA
do 100 k=1,NPA
......@@ -262,7 +280,7 @@ c LLIST is the code of final particle, P - its 4-momentum and mass.
if(ic.ge.1001) then !count spectators
nNuc=ic-1000
if(infragm.le.1
& .or.PA(4,k).lt.0.5*egymin)then !nuclear interaction only above min energy, otherwise : fragmentation
& .or.PA(4,k).lt.0.5d0*dble(egymin))then !nuclear interaction only above min energy, otherwise : fragmentation
ns=ns+nNuc
goto 100
elseif(ic.eq.1001)then
......@@ -291,11 +309,11 @@ c LLIST is the code of final particle, P - its 4-momentum and mass.
nbar=nbar+nNuc
if(abs(id).gt.1000.and.nNuc.eq.0)nbar=nbar+sign(1,id)
pptl(1,nptl)=PA(1,k) !P_x
pptl(2,nptl)=PA(2,k) !P_y
pptl(3,nptl)=PA(3,k) !P_z
pptl(4,nptl)=PA(4,k) !E
pptl(5,nptl)=PA(5,k) !mass
pptl(1,nptl)=sngl(PA(1,k)) !P_x
pptl(2,nptl)=sngl(PA(2,k)) !P_y
pptl(3,nptl)=sngl(PA(3,k)) !P_z
pptl(4,nptl)=sngl(PA(4,k)) !E
pptl(5,nptl)=sngl(PA(5,k)) !mass
istptl(nptl)=0
ityptl(nptl)=0
iorptl(nptl)=1
......@@ -419,6 +437,10 @@ c egy - center of mass energy
c------------------------------------------------------------------------------
include 'epos.inc'
dimension SDIF(3)
double precision egy_dbl,ST,SEL,SINEL,SDIF,SL,RHO,AL,
& SIGTA,SIGELAdum,SIGQEA,SIGSDA,SIGQSDA
egy_dbl = dble(egy)
if(iclpro.eq.1)then
L=2
......@@ -427,10 +449,11 @@ c------------------------------------------------------------------------------
else
L=3
endif
call SIB_SIGMA_HP(L,egy,ST,SEL,SINEL,SDIF,SL,RHO)
call SIB_SIGMA_HP(L,egy_dbl,ST,SEL,SINEL,SDIF,SL,RHO)
if(matar.gt.1)then
C calculate hadron-A(matar) cross section
CALL GLAUBER(matar,ST,SL,RHO,SIGTA,SIGELAdum,SIGQEA)
CALL GLAUBER2
& (matar,ST,SL,RHO,AL,SIGTA,SIGELAdum,SIGQEA,SIGSDA,SIGQSDA)
fsibcrse=SIGTA-SIGQEA
else
fsibcrse=SINEL
......@@ -452,7 +475,7 @@ c matarg - projec mass number
c id - proj id (sibyll code)
c------------------------------------------------------------------------------
include 'epos.inc'
double precision egy
double precision egy, sigbmdif, sibcr,SSIGNUC,alnuc,sqs_dbl
COMMON /CLENNN/ SSIGNUC(60), ALNUC(60)