]> git.uio.no Git - u/mrichter/AliRoot.git/blame - HBTAN/fsiini.F
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Typo corrected (HP)
[u/mrichter/AliRoot.git] / HBTAN / fsiini.F
CommitLineData
f5ab1a71 1
2 SUBROUTINE fsiini
3C---Note:
4C-- ICH= 0 (1) if the Coulomb interaction is absent (present);
5C-- ISPIN= JJ= 1,2,..,MSPIN denote increasing values of the pair
6C-- total spin S.
7C-- To calculate the CF of two particles (with masses m1, m2 and
8C-- charges C1, C2) the following information is required:
9C-- AM= twice the reduced mass= 2*m1*m2/(m1+m2) in GeV/c^2,
10C-- DM= (m1-m2)/(m1+m2), required if NS=2;
11C-- AC= Bohr radius= 2*137.036*0.1973/(C1*C2*AMH) in fm;
12C-- AC > 1.D9 if C1*C2= 0, AC < 0 if C1*C2 < 0;
13C-- MSPIN= MSPINH(LL)= number of the values of the total pair spin S;
14C-- FD= FDH(LL,JJ), RD= RDH(LL,JJ)= scattering length and effective
15C-- radius for each value of the total pair spin S, JJ= 1,..,MSPIN; ;
16C-- the corresponding square well parameters EB= EBH(LL,JJ), RB=
17C-- RBH(LL,JJ) (required if NS=1) may be calculated by sear.f;
18C-- if the effective range approximation is not valid (as is the case,
19C-- e.g., for two-pion system) a code for calculation of the scattering
20C-- amplitude should be supplemented;
21C-- RHO= RHOH(LL,JJ), SF= SFH(LL,JJ), SE= SEH(LL) are spin factors;
22C-- RHO= the probability that the spins j1 and j2 of the two particles
23C-- will combine in a total spin S;
24C-- RHO= (2*S+1)/[(2j1+1)*(2j2+1)] for unpolarized particles;
25C-- RHO= (1-P1*P2)/4 and (3+P1*P2)/4 correspond to S=0 and 1 in the
26C-- case of spin-1/2 particles with polarizations P1 and P2;
27C-----------------------------------------------------------------------
28 IMPLICIT REAL*8 (A-H,O-Z)
29 COMMON/FSI_POC/AMN,AM1,AM2,CN,C1,C2,AC1,AC2
30 COMMON/FSI_PRF/PPX,PPY,PPZ,AK,AKS,
31 1 X,Y,Z,T,RP,RPS
32 COMMON/FSI_SPIN/RHO(10)
33 COMMON/FSI_ACH/HPR,AC,ACH,ACHR,ISPIN,MSPIN
34 COMMON/FSI_NS/LL,NS,ICH,ISI,IQS,I3C,I3S
35 COMMON/FSI_FD/FD(10),RD(10)
36 COMMON/FSI_C/C(10),AM,AMS,DM
37 COMMON/FSI_CONS/PI,PI2,SPI,DR,W
38 COMPLEX*16 C
39 COMMON/FSI_AA/AA
40 COMMON/FSI_AAPI/AAPI(20,2)/FSI_AAND/AAND(20,4)
41 COMMON/FSI_AAPIN/AAPIN(20,2)
42 COMMON/FSI_SW/RB(10),EB(10),BK(10),CDK(10),SDK(10),
43 1 SBKRB(10),SDKK(10)
44 COMMON/LEDWEIGHT/WEIF,WEI,WEIN,ITEST,IRANPOS
45
46c Include 'common_fsi_poc.inc'
47c Include 'common_fsi_prf.inc'
48c Include 'common_fsi_spin.inc'
49c Include 'common_fsi_ach.inc'
50c Include 'common_fsi_ns.inc'
51c Include 'common_fsi_fd.inc'
52c-mlv Include 'common_fsi_c.inc'
53c Include 'common_fsi_cons.inc'
54c Include 'common_fsi_aa.inc'
55c Include 'common_fsi_aapi.inc'
56c Include 'common_fsi_aapin.inc'
57c Include 'common_fsi_sw.inc'
58c Include 'common_fsi_aand.inc'
59
60
61 DIMENSION FDH(30,10),RDH(30,10),EBH(30,10),RBH(30,10)
62 DIMENSION RHOH(30,10)
63 DIMENSION AM1H(30),AM2H(30),C1H(30),C2H(30),MSPINH(30)
64C============= declarations pour l'appel de READ_FILE()============
65 CHARACTER*10 KEY
66 CHARACTER*8 CH8
67 INTEGER*4 INT4
68 REAL*8 REAL8
69 INTEGER*4 IERR
70C
71C--- mass of the first and second particle
72 DATA AM1H/.93956563D0,.93827231D0,.93956563D0,3.72737978D0,
73 C .13957D0,.13498D0,.13957D0, .93956563D0, .93827231D0,
74 C 4*.13957D0,4*.493677D0,
75 C 2*1.87561339D0,2*2.80892165D0,2*.497672D0,
76 C 1.87561339D0,3*.93827231D0,.93956563D0, 2*0.D0/
77 DATA AM2H/.93956563D0,.93827231D0,.93827231D0,3.72737978D0,
78 C .13957D0,.13498D0,.13957D0, 2*1.87561339D0,
79 C 2*.493677D0,2*.93827231D0,
80 C 2*.493677D0,2*.93827231D0,
81 C 1.87561339D0,3.72737978D0,2.80892165D0,3.72737978D0,
82 C 2*.497672D0,2*2.80892165D0,3.72737978D0,
83 C 2*1.115684D0,2*0.D0/
84c--------|---------|---------|---------|---------|---------|---------|----------
85C--- charge of the first and second particle
86 DATA C1H/0.D0,1.D0,0.D0,2.D0, 1.D0,0.D0,1.D0,0.D0,1.D0,
87 C 3*1.D0,-1.D0,3*1.D0,-1.D0,
88 C 4*1.D0,2*0.D0,4*1.D0,0.D0, 2*0.D0/
89 DATA C2H/0.D0,1.D0,1.D0,2.D0,-1.D0,0.D0,3*1.D0,
90 C -1.D0,3*1.D0,-1.D0,3*1.D0,
91 C 1.D0,2.D0,1.D0,2.D0,2*0.D0,2*1.D0,2.D0,2*0.D0,2*0.D0/
92C---MSPIN vs (LL)
93 DATA MSPINH/3*2,4*1,2*2,8*1,3,1,2,1,2*1,2*2,1,2*2, 2*0/
94C---Spin factors RHO vs (LL,ISPIN)
95 DATA RHOH/3*.25D0, 4*1.D0, 2*.3333D0, 8*1.D0,
96 1 .1111D0,1.D0,.25D0,1.D0,2*1.D0,
97 1 .3333D0,.25D0,1.D0,2*.25D0, 2*0.D0,
98 2 3*.75D0, 4*0.D0, 2*.6667D0, 8*0.D0,
99 2 .3333D0,.0D0,.75D0,.0D0,2*0.D0,
100 2 .6667D0,.75D0,0.D0,2*.75D0, 2*0.D0,
101 3 17*.0D0,.5556D0,3*0.D0, 7*0.D0,2*0.D0,210*0.D0/
102C---Scattering length FD and effective radius RD in fm vs (LL,ISPIN)
103 DATA FDH/17.0D0,7.8D0,23.7D0,2230.1218D0,.225D0,.081D0,-.063D0,
104 1 -.65D0,-2.73D0,
105 1 .137D0,-.071D0,-.148D0,.112D0,2*1.D-6,-.360D0,
106 1 2*1.D-6,1.344D0,6*1.D-6,-5.628D0,2.18D0,2.40D0, 2*0.D0,
107cc 2 -10.8D0,2*-5.4D0,4*0.D0,-6.35D0,-11.88D0,8*0.D0,9*0.D0,
108 2 3*-5.4D0,4*0.D0,-6.35D0,-11.88D0,8*0.D0,9*0.D0,
109 2 1.93D0,1.84D0,2*0.D0,
110 3 240*0.D0/
111c--------|---------|---------|---------|---------|---------|---------|----------
112 DATA RDH/2.7D0,2.8D0,2.7D0,1.12139906D0,-44.36D0,64.0D0,784.9D0,
113 1 477.9D0, 2.27D0, 9*0.D0,-69.973D0, 6*0.D0,3.529D0,
114 1 3.19D0,3.15D0, 2*0.D0,
115 2 3*1.7D0,4*0.D0,2.0D0,2.63D0, 17*0.D0,3.35D0,3.37D0, 2*0.D0,
116 3 240*0.D0/
117C---Corresponding square well parameters RB (width in fm) and
118C-- EB =SQRT(-AM*U) (in GeV/c); U is the well height
119 DATA RBH/2.545739D0, 2.779789D0, 2.585795D0, 5.023544D0,
120 1 .124673D0, .3925180D0,.09D0, 2.D0, 4.058058D0, 17*0.D0,
121 1 2.252623D0, 2.278575D0, 2*0.D0,
122 2 3*2.003144D0,
123 2 4*0.D0, 2.D0, 4.132163D0, 17*0.D0,
124 2 2.272703D0, 2.256355D0, 2*0.D0,
125 3 240*0.D0/
126 DATA EBH/.1149517D0, .1046257D0, .1148757D0, .1186010D0,
127 1 .7947389D0,2.281208D0,8.7D0,.4D0,.1561219D0,17*0.D0,
128 1 .1013293D0, .1020966D0, 2*0.D0,
129 2 3*.1847221D0,
130 2 4*0.D0, .4D0, .1150687D0, 17*0.D0,
131 2 .09736083D0, .09708310D0, 2*0.D0,
132 3 240*0.D0/
133C=======< constants >========================
134 W=1/.1973D0 ! from fm to 1/GeV
135 PI=4*DATAN(1.D0)
136 PI2=2*PI
137 SPI=DSQRT(PI)
138 DR=180.D0/PI ! from radian to degree
139
140c WRITE(*,*)'from C++ to fortran W PI PI2 SPI DR',W,PI,PI2,SPI,DR
141
142 AC1=1.D10
143 AC2=1.D10
144C=======< condition de calculs >=============
145 NUNIT=11 ! for IBM or HP
146C NUNIT=4 ! for SUN in Prague
147c-mlv CALL readint4(NUNIT,'ITEST ',ITEST)
148c-mlv CALL readint4(NUNIT,'LL ',LL) ! Two-particle system
149c-mlv CALL readint4(NUNIT,'NS ',NS)
150c CALL READ_FILE(NUNIT,'ITEST ',CHAR,ITEST,REAL8,IERR)
151c CALL READ_FILE(NUNIT,'LL ',CHAR,LL,REAL8,IERR)
152c CALL READ_FILE(NUNIT,'NS ',CHAR,NS,REAL8,IERR)
153
154
155
156C---setting particle masses and charges
157 AM1=AM1H(LL)
158 AM2=AM2H(LL)
159 C1=C1H(LL)
160 C2=C2H(LL)
161
162C-- Switches:
163C ISI=1(0) the strong interaction between the two particles ON (OFF)
164C IQS=1(0) the quantum statistics ON (OFF);
165C should be OFF for nonidentical particles
166C I3C=1(0) the Coulomb interaction with the nucleus ON (OFF)
167C I3S=1(0) the strong interaction with the nucleus ON (OFF)
168C ICH=1(0) if C1*C2 is different from 0 (is equal to 0)
169C- to switch off the Coulomb force between the two particles
170C put ICH=0 and substitute the strong amplitude parameters by
171C the ones not affected by Coulomb interaction
172
173 IF(ITEST.EQ.0)THEN
174 ICH=0
175 IF(C1*C2.NE.0.D0) ICH=1
176 IQS=0
177 IF(C1+AM1.EQ.C2+AM2) IQS=1
178 I3S=0 ! only this option is available
179 ISI=1
180 I3C=0
181 IF(CN*ICH.NE.0.D0) I3C=1
182 ENDIF
183
184
185c IF(ITEST.EQ.1)THEN
186c NS=4 ! SPHER. WAVE
187c ICH=0
188c IQS=1
189c ISI=0
190c I3C=0
191
192
193c-mlv CALL readint4(NUNIT,'ICH ',ICH)
194c-mlv CALL readint4(NUNIT,'IQS ',IQS)
195c-mlv CALL readint4(NUNIT,'ISI ',ISI)
196c-mlv CALL readint4(NUNIT,'I3C ',I3C)
197c CALL READ_FILE(NUNIT,'ICH ',CHAR,ICH,REAL8,IERR)
198c CALL READ_FILE(NUNIT,'IQS ',CHAR,IQS,REAL8,IERR)
199c CALL READ_FILE(NUNIT,'ISI ',CHAR,ISI,REAL8,IERR)
200c CALL READ_FILE(NUNIT,'I3C ',CHAR,I3C,REAL8,IERR)
201c ENDIF
202
203 write(*,*)'====itest ll ich iqs isi i3c===',itest,ll,ich, iqs, isi, i3c
204
205C==================================================================
206C---fm to 1/GeV
207 DO 3 J1=1,30
208 DO 3 J2=1,10
209 FDH(J1,J2)=FDH(J1,J2)*W
210 RDH(J1,J2)=RDH(J1,J2)*W
211 3 RBH(J1,J2)=RBH(J1,J2)*W
212C---calcul. twice the reduced mass (AM), the relative mass difference
213C-- (DM) and the Bohr radius (AC)
214 AM=2*AM1*AM2/(AM1+AM2)
215 AMS=AM*AM
216 DM=(AM1-AM2)/(AM1+AM2)
217 AC=1.D10
218 C12=C1*C2
219 IF(C12.NE.0.D0)AC=2*137.036D0/(C12*AM)
220C---Setting spin factors
221 MSPIN=MSPINH(LL)
222 MSP=MSPIN
223 DO 91 ISPIN=1,10
224 91 RHO(ISPIN)=RHOH(LL,ISPIN)
225C---Integration limit AA in the spherical wave approximation
226 AA=0.D0
227cc IF(NS.EQ.2.OR.NS.EQ.4)AA=.5D0 !!in 1/GeV --> 0.1 fm
228 IF(NS.EQ.2.OR.NS.EQ.4)AA=6.D0 !!in 1/GeV --> 1.2 fm
229C---Setting scatt. length (FD), eff. radius (RD) and, if possible,
230C-- also the corresp. square well parameters (EB, RB)
231 DO 55 JJ=1,MSP
232 ISPIN=JJ
233 FD(JJ)=FDH(LL,JJ)
234 RD(JJ)=RDH(LL,JJ)
235 EB(JJ)=EBH(LL,JJ)
236 RB(JJ)=RBH(LL,JJ)
237C---Resets FD and RD for a nucleon-deuteron system (LL=8,9)
238 IF(LL.EQ.8.OR.LL.EQ.9)THEN
239 JH=LL-7+2*JJ-2
240 FD(JJ)=AAND(1,JH)
241 RD(JJ)=AAND(2,JH)-2*AAND(3,JH)/AAND(1,JH)
242 ENDIF
243C---Resets FD and RD for a pion-pion system (LL=5,6,7)
244 IF(LL.EQ.5.OR.LL.EQ.6.OR.LL.EQ.7)THEN
245 IF(LL.EQ.7)FD(JJ)=AAPI(1,2)/AM
246 IF(LL.EQ.5)FD(JJ)=(.6667D0*AAPI(1,1)+.3333D0*AAPI(1,2))/AM
247 IF(LL.EQ.6)FD(JJ)=(.3333D0*AAPI(1,1)+.6667D0*AAPI(1,2))/AM
248 AKS=0.D0
249 DAKS=1.D-5
250 AKSH=AKS+DAKS
251 AKH=DSQRT(AKSH)
252 GPI1H=GPIPI(AKSH,1)
253 GPI2H=GPIPI(AKSH,2)
254 H=1/FD(JJ)
255 IF(LL.EQ.7)C(JJ)=1/DCMPLX(GPI2H,-AKH)
256 IF(LL.EQ.5)
257 + C(JJ)=.6667D0/DCMPLX(GPI1H,-AKH)+.3333D0/DCMPLX(GPI2H,-AKH)
258 IF(LL.EQ.6)
259 + C(JJ)=.3333D0/DCMPLX(GPI1H,-AKH)+.6667D0/DCMPLX(GPI2H,-AKH)
260 HH=DREAL(1/C(JJ))
261 RD(JJ)=2*(HH-H)/DAKS
262 ENDIF
263C---Resets FD and RD for a pion-nucleon system (LL=12,13)
264 IF(LL.EQ.12.OR.LL.EQ.13)THEN
265 IF(LL.EQ.12)FD(JJ)=AAPIN(1,2)
266 IF(LL.EQ.13)FD(JJ)=(.6667D0*AAPIN(1,1)+.3333D0*AAPIN(1,2))
267 AKS=0.D0
268 DAKS=1.D-5
269 AKSH=AKS+DAKS
270 AKH=DSQRT(AKSH)
271 GPI1H=GPIN(AKSH,1)
272 GPI2H=GPIN(AKSH,2)
273 H=1/FD(JJ)
274 IF(LL.EQ.12)C(JJ)=1/DCMPLX(GPI2H,-AKH)
275 IF(LL.EQ.13)
276 + C(JJ)=.6667D0/DCMPLX(GPI1H,-AKH)+.3333D0/DCMPLX(GPI2H,-AKH)
277 HH=DREAL(1/C(JJ))
278 RD(JJ)=2*(HH-H)/DAKS
279 ENDIF
280C---Calculation continues for any system (any LL)
281 55 CONTINUE
282 RETURN
283 END
284
285
286
287 FUNCTION GPIPI(X,J)
288C--- GPIPI = k*COTG(DELTA), X=k^2
289C-- J=1(2) corresponds to isospin=0(2)
290 IMPLICIT REAL*8 (A-H,O-Z)
291c-- Include 'common_fsi_aapi.inc'
292c-- Include 'common_fsi_c.inc'
293 COMMON/FSI_AAPI/AAPI(20,2)
294 COMMON/FSI_C/HELP(20),AM,AMS,DM
295 OM=DSQRT(X+AMS)
296 XX=X/AMS
297 GPIPI=OM/AAPI(1,J)
298 GPIPI=GPIPI*(1+(AAPI(3,J)-AAPI(1,J)**2)*XX+AAPI(4,J)*XX*XX)
299 GPIPI=GPIPI/(1+(AAPI(3,J)+AAPI(2,J)/AAPI(1,J))*XX)
300
301 END
302
303 FUNCTION GPIN(X,J)
304C--- GPIN = k*COTG(DELTA), X=k^2
305C-- J=1(2) corresponds to piN isospin=1/2(3/2)
306 IMPLICIT REAL*8 (A-H,O-Z)
307c-- Include 'common_fsi_aapin.inc'
308 COMMON/FSI_AAPIN/AAPIN(20,2)
309 GPIN=1/AAPIN(1,J)+.5D0*AAPIN(2,J)*X
310
311 END