64 common/medparam/centrmin,centrmax,breal,centr,rau,nf
66 DOUBLE PRECISION centrmin,centrmax,breal,centr,rau
67 common/medparamint/taui,ti,
tc,d3,
zeta3,d,
68 &
n0,sigmann,a,woodssaxon
69 DOUBLE PRECISION taui,ti,
tc,
alpha,
beta,
gamma,d3,
zeta3,d,
n0,
74 common/rapmax2/etamax2
75 double precision etamax2
80 common/mdfac/mdfactor,mdscalefac
81 DOUBLE PRECISION mdfactor,mdscalefac
83 COMMON /thickfnc/ rmax,ta(100,2)
84 DOUBLE PRECISION rmax,ta
86 COMMON /crosssec/ impmax,cross(200,3)
87 DOUBLE PRECISION impmax,cross
98 CHARACTER*100 buffer,label,tempbuf
129 write(logfid,*)
'Reading medium parameters from ',
file
130 OPEN(unit=lun,
file=
file,status=
'old',err=10)
132 READ(lun,
'(A)', iostat=ios)
buffer
133 if (ios.ne.0) goto 30
135 if (firstchar.eq.
'#') goto 20
139 IF (label==
"TAUI")
THEN
140 READ(
buffer,*,iostat=ios) taui
141 ELSE IF (label==
"TI")
THEN
142 READ(
buffer,*,iostat=ios) ti
143 ELSE IF (label==
"TC")
THEN
145 ELSE IF (label==
"WOODSSAXON")
THEN
146 READ(
buffer,*,iostat=ios) woodssaxon
147 ELSE IF (label==
"CENTRMIN")
THEN
148 READ(
buffer,*,iostat=ios) centrmin
149 ELSE IF (label==
"CENTRMAX")
THEN
150 READ(
buffer,*,iostat=ios) centrmax
151 ELSE IF (label==
"NF")
THEN
152 READ(
buffer,*,iostat=ios) nf
153 ELSE IF (label==
"N0")
THEN
155 ELSE IF (label==
"D")
THEN
156 READ(
buffer,*,iostat=ios) d
157 ELSE IF (label==
"SIGMANN")
THEN
158 READ(
buffer,*,iostat=ios) sigmann
159 ELSE IF (label==
"MDFACTOR")
THEN
160 READ(
buffer,*,iostat=ios) mdfactor
161 ELSE IF (label==
"MDSCALEFAC")
THEN
162 READ(
buffer,*,iostat=ios) mdscalefac
164 write(logfid,*)
'unknown label ',label
168 30
close(lun,status=
'keep')
169 write(logfid,*)
'...done'
172 10
write(logfid,*)
'Could not open medium parameter file, '//
173 &
'will run with default settings.'
176 write(logfid,*)
'No medium parameter file found, '//
177 &
'will run with default settings.'
180 40
write(logfid,*)
'using parameters:'
181 write(logfid,*)
'TAUI =',taui
182 write(logfid,*)
'TI =',ti
183 write(logfid,*)
'TC =',
tc
184 write(logfid,*)
'WOODSSAXON =',woodssaxon
185 write(logfid,*)
'CENTRMIN =',centrmin
186 write(logfid,*)
'CENTRMAX =',centrmax
187 write(logfid,*)
'NF =',nf
188 write(logfid,*)
'A =',a
189 write(logfid,*)
'N0 =',
n0
190 write(logfid,*)
'D =',d
191 write(logfid,*)
'SIGMANN =',sigmann
192 write(logfid,*)
'MDFACTOR =',mdfactor
193 write(logfid,*)
'MDSCALEFAC =',mdscalefac
210 common/medparam/centrmin,centrmax,breal,centr,rau,nf
212 DOUBLE PRECISION centrmin,centrmax,breal,centr,rau
213 common/medparamint/taui,ti,
tc,d3,
zeta3,d,
214 &
n0,sigmann,a,woodssaxon
215 DOUBLE PRECISION taui,ti,
tc,
alpha,
beta,
gamma,d3,
zeta3,d,
n0,
220 COMMON /crosssec/ impmax,cross(200,3)
221 DOUBLE PRECISION impmax,cross
227 r=(
pyr(0)*(centrmax-centrmin)+centrmin)/100.
230 if ((
r-cross(
j,3)/cross(200,3)).ge.0.)
then
239 breal = (b2*(cross(
i,3)/cross(200,3)-
r)
240 & +
b1*(
r-cross(
i+1,3)/cross(200,3)))/
241 & (cross(
i,3)/cross(200,3)-cross(
i+1,3)/cross(200,3))
247 common/medparam/centrmin,centrmax,breal,centr,rau,nf
249 DOUBLE PRECISION centrmin,centrmax,breal,centr,rau
259 common/medparam/centrmin,centrmax,breal,centr,rau,nf
261 DOUBLE PRECISION centrmin,centrmax,breal,centr,rau
263 DOUBLE PRECISION x1,x2,y1,y2,
z,xval,yval,zval,
nthick,
pyr
267 y1=-sqrt(4*rau**2-breal**2)/2.
268 y2=sqrt(4*rau**2-breal**2)/2.
269 131 xval=
pyr(0)*(x2-x1)+x1
270 yval=
pyr(0)*(y2-y1)+y1
271 IF((
nthick(xval-breal/2.,yval).EQ.0.d0).OR.
272 &
nthick(xval+breal/2.,yval).EQ.0.d0) goto 131
275 IF(zval.GT.
z) goto 131
283 common/medparam/centrmin,centrmax,breal,centr,rau,nf
285 DOUBLE PRECISION centrmin,centrmax,breal,centr,rau
286 DOUBLE PRECISION bval
295 common/medparam/centrmin,centrmax,breal,centr,rau,nf
297 DOUBLE PRECISION centrmin,centrmax,breal,centr,rau
299 common/medparamint/taui,ti,
tc,d3,
zeta3,d,
300 &
n0,sigmann,a,woodssaxon
301 DOUBLE PRECISION taui,ti,
tc,
alpha,
beta,
gamma,d3,
zeta3,d,
n0,
306 common/boostmed/
boost
311 common/logfile/logfid
314 DOUBLE PRECISION x,
y,
z,
t,ms,
px,
py,
pz,
e,md,temp
316 DOUBLE PRECISION r,
pyr,pmax,wt,
tau,theta,
phi,
pi,
p,ys,pz2,e2
317 DATA pi/3.141592653589793d0/
320 IF(
r.LT.(2.*12.*nf*d3/3.)/(2.*12.*nf*d3/3.+3.*16.*
zeta3/2.))
THEN
336 IF(temp.LT.1.d-2)
THEN
337 write(logfid,*)
'asking for a scattering centre without medium:'
338 write(logfid,*)
'at (x,y,z,t)=',
x,
y,
z,
t
339 write(logfid,*)
'making one up to continue but '//
340 &
'something is wrong!'
345 ms=
getms(0.d0,0.d0,0.d0,0.d0)
346 md=
getmd(0.d0,0.d0,0.d0,0.d0)
347 e=sqrt(
px**2+
py**2+
pz**2+ms**2)
351 10
p =
pyr(0)**0.3333333*pmax
352 e2 = sqrt(
p**2+ms**2)
354 wt = (exp(ms/temp)-1.)/(exp(e2/temp)-1.)
356 wt = (exp(ms/temp)+1.)/(exp(e2/temp)+1.)
358 if (wt.gt.1.)
write(logfid,*)
'Error in getscatterer: weight = ',wt
359 if (wt.lt.0.)
write(logfid,*)
'Error in getscatterer: weight = ',wt
360 if (
pyr(0).gt.wt) goto 10
362 theta = -acos(2.*
pyr(0)-1.)+
pi
363 px =
p*sin(theta)*cos(
phi)
364 py =
p*sin(theta)*sin(
phi)
366 e = cosh(ys)*e2 + sinh(ys)*pz2
367 pz = sinh(ys)*e2 + cosh(ys)*pz2
374 common/boostmed/
boost
377 common/rapmax2/etamax2
378 double precision etamax2
380 double precision x,
y,
z,
t,
px,
py,
pz,
e,
getms,m,ys
383 if ((
z.eq.0.d0).and.(
t.eq.0.d0)) ys =0.d0
384 if (ys.gt.etamax2) ys=etamax2
385 if (ys.lt.-etamax2) ys=-etamax2
400 common/boostmed/
boost
403 common/rapmax2/etamax2
404 double precision etamax2
421 DOUBLE PRECISION FUNCTION getmd(X1,Y1,Z1,T1)
424 common/mdfac/mdfactor,mdscalefac
425 DOUBLE PRECISION mdfactor,mdscalefac
433 DOUBLE PRECISION FUNCTION getms(X2,Y2,Z2,T2)
435 DOUBLE PRECISION x2,y2,z2,t2,
getmd
441 DOUBLE PRECISION FUNCTION getneff(X3,Y3,Z3,T3)
443 common/medparam/centrmin,centrmax,breal,centr,rau,nf
445 DOUBLE PRECISION centrmin,centrmax,breal,centr,rau
446 common/medparamint/taui,ti,
tc,d3,
zeta3,d,
447 &
n0,sigmann,a,woodssaxon
448 DOUBLE PRECISION taui,ti,
tc,
alpha,
beta,
gamma,d3,
zeta3,d,
n0,
454 DATA pi/3.141592653589793d0/
455 tau = sqrt(t3**2-z3**2)
464 DOUBLE PRECISION FUNCTION gettemp(X4,Y4,Z4,T4)
467 common/medparam/centrmin,centrmax,breal,centr,rau,nf
469 DOUBLE PRECISION centrmin,centrmax,breal,centr,rau
470 common/medparamint/taui,ti,
tc,d3,
zeta3,d,
471 &
n0,sigmann,a,woodssaxon
472 DOUBLE PRECISION taui,ti,
tc,
alpha,
beta,
gamma,d3,
zeta3,d,
n0,
477 common/rapmax2/etamax2
478 double precision etamax2
480 DOUBLE PRECISION x4,y4,z4,t4,
tau,
npart,eps0,epsin,tempin,
pi,
482 DATA pi/3.141592653589793d0/
486 IF(abs(z4).GT.t4)
RETURN
488 tau=sqrt(t4**2-z4**2)
490 IF((
nthick(x4-breal/2.,y4).EQ.0.d0).OR.
491 &
nthick(x4+breal/2.,y4).EQ.0.d0)
RETURN
493 ys = 0.5*
log((t4+z4)/(t4-z4))
494 if (abs(ys).gt.etamax2)
return
497 eps0=(16.*8.+7.*2.*6.*nf)*
pi**2*ti**4/240.
498 epsin=eps0*
npart(x4-breal/2.,y4,x4+breal/2.,y4)
500 tempin=(epsin*240./(
pi**2*(16.*8.+7.*2.*6.*nf)))**0.25
519 common/medparam/centrmin,centrmax,breal,centr,rau,nf
521 DOUBLE PRECISION centrmin,centrmax,breal,centr,rau
522 common/medparamint/taui,ti,
tc,d3,
zeta3,d,
523 &
n0,sigmann,a,woodssaxon
524 DOUBLE PRECISION taui,ti,
tc,
alpha,
beta,
gamma,d3,
zeta3,d,
n0,
538 common/mdfac/mdfactor,mdscalefac
539 DOUBLE PRECISION mdfactor,mdscalefac
550 common/medparam/centrmin,centrmax,breal,centr,rau,nf
552 DOUBLE PRECISION centrmin,centrmax,breal,centr,rau
553 common/medparamint/taui,ti,
tc,d3,
zeta3,d,
554 &
n0,sigmann,a,woodssaxon
555 DOUBLE PRECISION taui,ti,
tc,
alpha,
beta,
gamma,d3,
zeta3,d,
n0,
560 common/mdfac/mdfactor,mdscalefac
561 DOUBLE PRECISION mdfactor,mdscalefac
580 common/medparam/centrmin,centrmax,breal,centr,rau,nf
582 DOUBLE PRECISION centrmin,centrmax,breal,centr,rau
583 common/medparamint/taui,ti,
tc,d3,
zeta3,d,
584 &
n0,sigmann,a,woodssaxon
585 DOUBLE PRECISION taui,ti,
tc,
alpha,
beta,
gamma,d3,
zeta3,d,
n0,
590 common/rapmax2/etamax2
591 double precision etamax2
593 common/mdfac/mdfactor,mdscalefac
594 DOUBLE PRECISION mdfactor,mdscalefac,
pi
595 DATA pi/3.141592653589793d0/
608 common/medparam/centrmin,centrmax,breal,centr,rau,nf
610 DOUBLE PRECISION centrmin,centrmax,breal,centr,rau
611 common/medparamint/taui,ti,
tc,d3,
zeta3,d,
612 &
n0,sigmann,a,woodssaxon
613 DOUBLE PRECISION taui,ti,
tc,
alpha,
beta,
gamma,d3,
zeta3,d,
n0,
618 common/rapmax2/etamax2
619 double precision etamax2
629 common/medparam/centrmin,centrmax,breal,centr,rau,nf
631 DOUBLE PRECISION centrmin,centrmax,breal,centr,rau
632 common/medparamint/taui,ti,
tc,d3,
zeta3,d,
633 &
n0,sigmann,a,woodssaxon
634 DOUBLE PRECISION taui,ti,
tc,
alpha,
beta,
gamma,d3,
zeta3,d,
n0,
639 common/rapmax2/etamax2
640 double precision etamax2
643 DATA pi/3.141592653589793d0/
650 DOUBLE PRECISION FUNCTION npart(XX1,YY1,XX2,YY2)
652 common/medparamint/taui,ti,
tc,d3,
zeta3,d,
653 &
n0,sigmann,a,woodssaxon
654 DOUBLE PRECISION taui,ti,
tc,
alpha,
beta,
gamma,d3,
zeta3,d,
n0,
659 DOUBLE PRECISION xx1,yy1,xx2,yy2,
nthick
670 common/medparam/centrmin,centrmax,breal,centr,rau,nf
672 DOUBLE PRECISION centrmin,centrmax,breal,centr,rau
673 common/medparamint/taui,ti,
tc,d3,
zeta3,d,
674 &
n0,sigmann,a,woodssaxon
675 DOUBLE PRECISION taui,ti,
tc,
alpha,
beta,
gamma,d3,
zeta3,d,
n0,
680 common/logfile/logfid
683 COMMON /thickfnc/ rmax,ta(100,2)
684 DOUBLE PRECISION rmax,ta
685 INTEGER line,lmin,lmax,
i
686 DOUBLE PRECISION x1,y1,xa(4),ya(4),
y,
dy,
r,
c,b,delta
689 IF(
r.GT.ta(100,1))
THEN
692 line=int(
r*99.d0/ta(100,1)+1)
695 IF((
r.LT.ta(lmin,1)).OR.(
r.GT.ta(lmin+1,1)))
696 &
write(logfid,*)
line,lmin,
r,ta(lmin,1),ta(lmin+1,1)
701 c=(ya(2)-ya(1))/(xa(2)-xa(1))
712 common/medparam/centrmin,centrmax,breal,centr,rau,nf
714 DOUBLE PRECISION centrmin,centrmax,breal,centr,rau
715 common/medparamint/taui,ti,
tc,d3,
zeta3,d,
716 &
n0,sigmann,a,woodssaxon
717 DOUBLE PRECISION taui,ti,
tc,
alpha,
beta,
gamma,d3,
zeta3,d,
n0,
722 COMMON /thickfnc/ rmax,ta(100,2)
723 DOUBLE PRECISION rmax,ta
729 DOUBLE PRECISION eps,hfirst,
y
735 r=1.12*a**(0.33333)-0.86*a**(-0.33333)
740 b=(
i-1)*2.d0*
r/nsteps
754 common/medparam/centrmin,centrmax,breal,centr,rau,nf
756 DOUBLE PRECISION centrmin,centrmax,breal,centr,rau
757 common/medparamint/taui,ti,
tc,d3,
zeta3,d,
758 &
n0,sigmann,a,woodssaxon
759 DOUBLE PRECISION taui,ti,
tc,
alpha,
beta,
gamma,d3,
zeta3,d,
n0,
764 COMMON /crosssec/ impmax,cross(200,3)
765 DOUBLE PRECISION impmax,cross
782 p=(1.d0-prod)*8.8d0/14.d0*b
788 cross(ib,3)=cross(ib-1,3)+(
p+pprev)/2.*0.1
799 DOUBLE PRECISION xval
802 common/medparamint/taui,ti,
tc,d3,
zeta3,d,
803 &
n0,sigmann,a,woodssaxon
804 DOUBLE PRECISION taui,ti,
tc,
alpha,
beta,
gamma,d3,
zeta3,d,
n0,