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,
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
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
144 READ(
buffer,*,iostat=ios) tc
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
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
267 y1=-
sqrt(4*rau**2-breal**2)/2.
268 y2=
sqrt(4*rau**2-breal**2)/2.
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,
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!'
346 md=
getmd(0.d0,0.d0,0.d0,0.d0)
351 10
p =
pyr(0)**0.3333333*pmax
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
366 e = cosh(ys)*
e2 + sinh(ys)*pz2
367 pz = sinh(ys)*
e2 + cosh(ys)*pz2
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
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)
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/
457 getneff=(2.*6.*nf*d3*2./3. + 16.*zeta3*3./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,
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
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.
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,
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,
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,
640 double precision etamax2
643 DATA pi/3.141592653589793d0/
644 getneffmax=(2.*6.*nf*d3*2./3. + 16.*zeta3*3./2.)
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,
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
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,