]> git.uio.no Git - u/mrichter/AliRoot.git/blame - ISAJET/isasusy/ssmhn.F
Buffer Size can be defined
[u/mrichter/AliRoot.git] / ISAJET / isasusy / ssmhn.F
CommitLineData
0795afa3 1#include "isajet/pilot.h"
2 SUBROUTINE SSMHN(MHLNEG)
3C-----------------------------------------------------------------------
4C
5C Calculate HL, HH masses and ALFAH
6C (scalar Higgs mixing angle) using radiative
7C corrections calculated by M. Bisset
8C and save results in /SSPAR/.
9C
10C Both top and bottom couplings are now
11C included. Non-degenerate mixed squark
12C masses and A-terms are also included.
13C The D-terms from the squark mass matrix
14C (terms prop. to g**2 * Yukawa coupling)
15C are included as an option:
16C INRAD = 1 ==> D-TERMS ON
17C INRAD = 2 ==> D-TERMS OFF .
18C
19C 10/18/93 D-terms are now turned on.
20C INRAD = 1
21C
22C There is an arbitrary mass scale that must
23C chosen to avoid dimensionful logarithms.
24C The choice does not matter if D-terms are
25C not included, but it does matter if D-terms
26C are included.
27C
28C Arbitrary mass scale updated to
29C QQQ = HIGFRZ = SQRT(AMTLSS*AMTRSS)
30C with running masses to include dominant 2-loop
31C effects. 12/10/96 H. Baer
32C
33C It is assumed that the A-terms are real.
34C Complex A-terms are allowed
35C (unless RTT=0 or RBB=0 --see below) in
36C this subroutine, but the imaginary parts
37C are now set to zero.
38C
39C-----------------------------------------------------------------------
40#if defined(CERNLIB_IMPNONE)
41 IMPLICIT NONE
42#endif
43#include "isajet/sslun.inc"
44#include "isajet/sssm.inc"
45#include "isajet/sspar.inc"
46C
47 REAL PI,PI2,SR2,G2,GP2,GGP,GG1,GG2
48 REAL TANB,COTB,COSB,SINB,BE
49 REAL SINB2,COSB2,COS2B,SIN2B
50 REAL V2,VP2,V,VP,VVP,VPVM,VVPP
51 REAL MT2,MB2,FT2,FB2,MW2,ZAP,QQQ2
52 REAL EP,EP2,RR,MHP2
53 REAL ATI,ABI,ATR,ABR,AT2,AB2
54 REAL TLRM,BLRM,TLRP,BLRP
55 REAL MST1SQ,MST2SQ,MSB1SQ,MSB2SQ
56
57 REAL RTT,TTT1,TEMPT,TM1BT
58 REAL TEMPS,T1RD,T2RD,T1RPD,T2RPD
59 REAL CT1,A1,A2,T1RR,T2RR
60 REAL CT5,A5,A6,T1RPRP,T2RPRP
61 REAL A9,T1RRP,T2RRP
62 REAL TEMPSQ,DT1,DT2,VRRT,VRPRPT,VRRPT
63 REAL ALPHAT,LAT
64C
65 REAL RBB,BBB1,TEMPB,TM1BB
66 REAL B1RD,B2RD,B1RPD,B2RPD
67 REAL CB3,A3,A4,B1RR,B2RR
68 REAL CB7,A7,A8,B1RPRP,B2RPRP
69 REAL A10,B1RRP,B2RRP
70 REAL DB1,DB2,VRRB,VRPRPB,VRRPB
71 REAL ALPHAB,LAB
72C
73 REAL DVRR,DVRPRP,DVRRP,TEMPH
74 REAL MHL2,MHH2,TRACEM,TPAL,TANAH
75 REAL ASMB,MBMB,MBQ,ASMT,MTMT,MTQ,SUALFS,HIGFRZ
76 DOUBLE PRECISION SSMQCD
77 INTEGER INRAD,MHLNEG
78C
79 MHLNEG=0
80 PI=4.*ATAN(1.)
81 PI2 = PI**2
82 SR2=SQRT(2.)
83 G2=4.*PI*ALFAEM/SN2THW
84 GP2=G2*SN2THW/(1.-SN2THW)
85 HIGFRZ=SQRT(AMTLSS*AMTRSS)
86 ASMB=SUALFS(AMBT**2,.36,AMTP,3)
87 MBMB=AMBT*(1.-4*ASMB/3./PI)
88 MBQ=SSMQCD(DBLE(MBMB),DBLE(HIGFRZ))
89 ASMT=SUALFS(AMTP**2,.36,AMTP,3)
90 MTMT=AMTP/(1.+4*ASMT/3./PI+(16.11-1.04*(5.-6.63/AMTP))*
91 $(ASMT/PI)**2)
92 MTQ=SSMQCD(DBLE(MTMT),DBLE(HIGFRZ))
93 MT2=MTQ**2
94 MB2=MBQ**2
95 MW2=AMW**2
96 EP=TWOM1
97 EP2=EP**2
98 RR=RV2V1
99 MHP2=AMHA**2
100 TANB=1.0/RR
101 COTB=RR
102 BE=ATAN(1./RV2V1)
103 SINB=SIN(BE)
104 COSB=COS(BE)
105 SINB2=SINB**2
106 COSB2=COSB**2
107 SIN2B=SIN(2.*BE)
108 COS2B=COS(2.*BE)
109 V2=2.0*MW2*SINB2/G2
110 VP2=2.0*MW2*COSB2/G2
111 V=SQRT(V2)
112 VP=SQRT(VP2)
113 VVP=SQRT(V2*VP2)
114 VPVM=VP2-V2
115 GGP=G2+GP2
116 GG1=G2-5.0*GP2/3.0
117 GG2=G2-GP2/3.0
118 VVPP=2.0*AMZ**2/GGP
119 FT2=MT2/V2
120 FB2=MB2/VP2
121C
122 TLRM=AMTLSS**2-AMTRSS**2
123 BLRM=AMBLSS**2-AMBRSS**2
124 TLRP=AMTLSS**2+AMTRSS**2
125 BLRP=AMBLSS**2+AMBRSS**2
126C
127C Higgs mass matrix
128C
129C (AAT and AAB are also assumed to be real)
130C
131 ATR=AAT
132 ABR=AAB
133 ATI=0.0
134 ABI=0.0
135 AT2=ATR**2+ATI**2
136 AB2=ABR**2+ABI**2
137C
138 MST1SQ=AMT1SS**2
139 MST2SQ=AMT2SS**2
140 MSB1SQ=AMB1SS**2
141 MSB2SQ=AMB2SS**2
142 INRAD=1
143 QQQ2=HIGFRZ**2
144C
145 ZAP = 1.0
146C
147C STOP TERMS
148C
149 RTT=(TLRM+VPVM*ZAP*GG1/4.0)**2
150 $ +4.0*MT2*(EP*COTB+ATR)**2+4.0*MT2*ATI**2
151 RTT=SQRT(RTT)
152C
153C calculate 2M1*B term
154C
155 TTT1=0.5*TLRP+MT2+VPVM*ZAP*GGP/8.0
156 IF(RTT.NE.0.0) THEN
157 TEMPT=4.0*EP*FT2*VVP*ATI**2/(RTT**2)
158 TM1BT=-2.0*FT2*(TEMPT+ATR)*TTT1
159 $ *LOG(MST2SQ/MST1SQ)/RTT
160 TM1BT=TM1BT-FT2*ATR
161 $ *LOG(MST1SQ*MST2SQ/QQQ2/QQQ2)
162 TM1BT=TM1BT+FT2*(2.0*TEMPT-ATR)
163 TM1BT=3.0*EP*TM1BT/32.0/PI2
164C
165C calculate first derivatives w.r.t H_R
166C divided by sqrt(2) * v
167C
168 TEMPS=-ZAP*GG1*(TLRM+ZAP*GG1*VPVM/4.0)/2.0
169 TEMPS=TEMPS+4.0*FT2*(AT2+EP*COTB*ATR)
170 TEMPS=TEMPS/RTT/4.0
171 T1RD=FT2-ZAP*GGP/8.0-TEMPS
172 T2RD=FT2-ZAP*GGP/8.0+TEMPS
173C
174C calculate first derivatives w.r.t H_R'
175C divided by sqrt(2) * v'
176C
177 TEMPS=ZAP*GG1*(TLRM+ZAP*GG1*VPVM/4.0)/2.0
178 TEMPS=TEMPS+4.0*FT2*EP*(EP+TANB*ATR)
179 TEMPS=TEMPS/RTT/4.0
180 T1RPD=ZAP*GGP/8.0-TEMPS
181 T2RPD=ZAP*GGP/8.0+TEMPS
182C
183C calculate second derivatives w.r.t. H_R
184C
185 CT1=-V*ZAP*GG1*(TLRM+ZAP*GG1*VPVM/4.0)/SR2
186 CT1=CT1+4.0*SR2*FT2*V*(EP*COTB*ATR+AT2)
187 A1=-CT1**2/(RTT**3)/8.0
188 A2=-ZAP*GG1*(TLRM+ZAP*GG1*VPVM/4.0)/2.0
189 A2=A2+V2*ZAP*GG1**2/4.0+4.0*FT2*AT2
190 A2=A2/RTT/4.0
191 T1RR=FT2-ZAP*GGP/8.0-A1-A2
192 T2RR=FT2-ZAP*GGP/8.0+A1+A2
193C
194C calculate second derivatives w.r.t. H_R'
195C
196 CT5=VP*ZAP*GG1*(TLRM+ZAP*GG1*VPVM/4.0)/SR2
197 CT5=CT5+4.0*SR2*FT2*VP*EP*(EP+TANB*ATR)
198 A5=-CT5**2/(RTT**3)/8.0
199 A6=ZAP*GG1*(TLRM+ZAP*GG1*VPVM/4.0)/2.0
200 A6=A6+VP2*ZAP*GG1**2/4.0+4.0*FT2*EP2
201 A6=A6/RTT/4.0
202 T1RPRP=ZAP*GGP/8.0-A5-A6
203 T2RPRP=ZAP*GGP/8.0+A5+A6
204C
205C calculate second derivatives w.r.t. H_R and H_R'
206C
207 A9=-VVP*ZAP*(GG1**2)/4.0+4.0*FT2*EP*ATR
208 A9=A9/RTT/4.0
209 T1RRP=CT1*CT5/(RTT**3)/8.0-A9
210 T2RRP=-CT1*CT5/(RTT**3)/8.0+A9
211C
212C calculate D^2 V / D^2 H_R
213C
214 TEMPSQ=MST1SQ*(T1RR-T1RD)
215 DT1=2.0*(2.0*V2*T1RD**2+TEMPSQ)*LOG(MST1SQ/QQQ2)
216 DT1=DT1+6.0*V2*T1RD**2+TEMPSQ
217 TEMPSQ=MST2SQ*(T2RR-T2RD)
218 DT2=2.0*(2.0*V2*T2RD**2+TEMPSQ)*LOG(MST2SQ/QQQ2)
219 DT2=DT2+6.0*V2*T2RD**2+TEMPSQ
220 VRRT=DT1+DT2-8.0*FT2*MT2*LOG(MT2/QQQ2)-12.0*FT2*MT2
221 VRRT=-TM1BT*COTB+3.0*VRRT/32.0/PI2
222C
223C calculate D^2 V / D^2 H'_R
224C
225 TEMPSQ=MST1SQ*(T1RPRP-T1RPD)
226 DT1=2.0*(2.0*VP2*T1RPD**2+TEMPSQ)*LOG(MST1SQ/QQQ2)
227 DT1=DT1+6.0*VP2*T1RPD**2+TEMPSQ
228 TEMPSQ=MST2SQ*(T2RPRP-T2RPD)
229 DT2=2.0*(2.0*VP2*T2RPD**2+TEMPSQ)*LOG(MST2SQ/QQQ2)
230 DT2=DT2+6.0*VP2*T2RPD**2+TEMPSQ
231 VRPRPT=-TM1BT*TANB+3.0*(DT1+DT2)/32.0/PI2
232C
233C calculate D^2 V / D^H_R D^H_R'
234C
235 DT1=2.0*VVP*T1RD*T1RPD+MST1SQ*T1RRP
236 DT1=2.0*DT1*LOG(MST1SQ/QQQ2)
237 DT1=DT1+6.0*VVP*T1RD*T1RPD+MST1SQ*T1RRP
238 DT2=2.0*VVP*T2RD*T2RPD+MST2SQ*T2RRP
239 DT2=2.0*DT2*LOG(MST2SQ/QQQ2)
240 DT2=DT2+6.0*VVP*T2RD*T2RPD+MST2SQ*T2RRP
241 VRRPT=TM1BT+3.0*(DT1+DT2)/32.0/PI2
242C
243 ELSE IF(RTT.EQ.0.0) THEN
244C
245 ALPHAT=TLRP/2.0+MT2+ZAP*GGP*VPVM/8.0
246 LAT=2.0*LOG(ALPHAT/QQQ2)+3.0
247C
248C calculate D^2 V / D^2 H_R
249C
250 VRRT=V2*(GGP**2+GG1**2)/16.0-MT2*GGP
251 VRRT=ZAP*VRRT*LAT+8.0*FT2*MT2*LOG(ALPHAT/MT2)
252 VRRT=3.0*VRRT/32.0/PI2
253C
254C calculate D^2 V / D^2 H_R'
255C
256 VRPRPT=ZAP*VP2*(GGP**2+GG1**2)/16.0
257 VRPRPT=3.0*(VRPRPT*LAT)/32.0/PI2
258C
259C calculate D^2 V / D^H_R D^H_R'
260C
261 VRRPT=FT2*GGP-(GGP**2+GG1**2)/8.0
262 VRRPT=ZAP*VVP*VRRPT*LAT/2.0
263 VRRPT=3.0*VRRPT/32.0/PI2
264C
265C
266 ENDIF
267C
268C SBOTTOM TERMS
269C
270 RBB=(BLRM-VPVM*ZAP*GG2/4.0)**2
271 $ +4.0*MB2*(EP*TANB+ABR)**2+4.0*MB2*ABI**2
272 RBB=SQRT(RBB)
273C IF(RBB.EQ.0.0.AND.ABI.NE.0.0) THEN
274C WRITE(6,*) 'RBB=0, ABI NOT 0'
275C WRITE(6,*) 'ERROR: THIS CASE NOT COVERED YET'
276C GO TO 1000
277C ENDIF
278C
279 IF(RBB.NE.0.0) THEN
280C
281C calculate 2M1*B term
282C
283 BBB1=0.5*BLRP+MB2-VPVM*ZAP*GGP/8.0
284 TEMPB=4.0*EP*FB2*VVP*ABI**2/(RBB**2)
285 TM1BB=-2.0*FB2*(TEMPB+ABR)*BBB1
286 $ *LOG(MSB2SQ/MSB1SQ)/RBB
287 TM1BB=TM1BB-FB2*ABR
288 $ *LOG(MSB1SQ*MSB2SQ/QQQ2/QQQ2)
289 TM1BB=TM1BB+FB2*(2.0*TEMPB-ABR)
290 TM1BB=3.0*EP*TM1BB/32.0/PI2
291C
292C calculate first derivatives w.r.t H_R
293C divided by sqrt(2) * v
294C
295 TEMPS=ZAP*GG2*(BLRM-ZAP*GG2*VPVM/4.0)/2.0
296 TEMPS=TEMPS+4.0*FB2*EP*(EP+COTB*ABR)
297 TEMPS=TEMPS/RBB/4.0
298 B1RD=ZAP*GGP/8.0-TEMPS
299 B2RD=ZAP*GGP/8.0+TEMPS
300
301C calculate first derivatives w.r.t H_R'
302C divided by sqrt(2) * v'
303C
304 TEMPS=-ZAP*GG2*(BLRM-ZAP*GG2*VPVM/4.0)/2.0
305 TEMPS=TEMPS+4.0*FB2*(AB2+EP*TANB*ABR)
306 TEMPS=TEMPS/RBB/4.0
307 B1RPD=FB2-ZAP*GGP/8.0-TEMPS
308 B2RPD=FB2-ZAP*GGP/8.0+TEMPS
309C
310C calculate second derivatives w.r.t. H_R
311C
312 CB3=V*ZAP*GG2*(BLRM-ZAP*GG2*VPVM/4.0)/SR2
313 CB3=CB3+4.0*SR2*FB2*V*EP*(EP+COTB*ABR)
314 A3=-CB3**2/(RBB**3)/8.0
315 A4=ZAP*GG2*(BLRM-ZAP*GG2*VPVM/4.0)/2.0
316 A4=A4+V2*ZAP*GG2**2/4.0+4.0*FB2*EP2
317 A4=A4/RBB/4.0
318 B1RR=ZAP*GGP/8.0-A3-A4
319 B2RR=ZAP*GGP/8.0+A3+A4
320C
321C calculate second derivatives w.r.t. H_R'
322C
323 CB7=-VP*ZAP*GG2*(BLRM-ZAP*GG2*VPVM/4.0)/SR2
324 CB7=CB7+4.0*SR2*FB2*VP*(AB2+EP*TANB*ABR)
325 A7=-CB7**2/(RBB**3)/8.0
326 A8=-ZAP*GG2*(BLRM-ZAP*GG2*VPVM/4.0)/2.0
327 A8=A8+VP2*ZAP*GG2**2/4.0+4.0*FB2*AB2
328 A8=A8/RBB/4.0
329 B1RPRP=FB2-ZAP*GGP/8.0-A7-A8
330 B2RPRP=FB2-ZAP*GGP/8.0+A7+A8
331C
332C calculate second derivatives w.r.t. H_R and H_R'
333C
334 A10=-VVP*ZAP*(GG2**2)/4.0+4.0*FB2*EP*ABR
335 A10=A10/RBB/4.0
336 B1RRP=CB3*CB7/(RBB**3)/8.0-A10
337 B2RRP=-CB3*CB7/(RBB**3)/8.0+A10
338C
339C calculate D^2 V / D^2 H_R
340C
341 TEMPSQ=MSB1SQ*(B1RR-B1RD)
342 DB1=2.0*(2.0*V2*B1RD**2+TEMPSQ)*LOG(MSB1SQ/QQQ2)
343 DB1=DB1+6.0*V2*B1RD**2+TEMPSQ
344 TEMPSQ=MSB2SQ*(B2RR-B2RD)
345 DB2=2.0*(2.0*V2*B2RD**2+TEMPSQ)*LOG(MSB2SQ/QQQ2)
346 DB2=DB2+6.0*V2*B2RD**2+TEMPSQ
347 VRRB=-TM1BB*COTB+3.0*(DB1+DB2)/32.0/PI2
348C
349C calculate D^2 V / D^2 H'_R
350C
351 TEMPSQ=MSB1SQ*(B1RPRP-B1RPD)
352 DB1=2.0*(2.0*VP2*B1RPD**2+TEMPSQ)*LOG(MSB1SQ/QQQ2)
353 DB1=DB1+6.0*VP2*B1RPD**2+TEMPSQ
354 TEMPSQ=MSB2SQ*(B2RPRP-B2RPD)
355 DB2=2.0*(2.0*VP2*B2RPD**2+TEMPSQ)*LOG(MSB2SQ/QQQ2)
356 DB2=DB2+6.0*VP2*B2RPD**2+TEMPSQ
357 VRPRPB=DB1+DB2
358 VRPRPB=DB1+DB2-8.0*FB2*MB2*LOG(MB2/QQQ2)-12.0*FB2*MB2
359 VRPRPB=-TM1BB*TANB+3.0*VRPRPB/32.0/PI2
360C
361C calculate D^2 V / D H_R D H'_R
362C
363 DB1=2.0*VVP*B1RD*B1RPD+MSB1SQ*B1RRP
364 DB1=2.0*DB1*LOG(MSB1SQ/QQQ2)
365 DB1=DB1+6.0*VVP*B1RD*B1RPD+MSB1SQ*B1RRP
366 DB2=2.0*VVP*B2RD*B2RPD+MSB2SQ*B2RRP
367 DB2=2.0*DB2*LOG(MSB2SQ/QQQ2)
368 DB2=DB2+6.0*VVP*B2RD*B2RPD+MSB2SQ*B2RRP
369 VRRPB=TM1BB+3.0*(DB1+DB2)/32.0/PI2
370
371 ELSE IF(RBB.EQ.0.0) THEN
372C
373 ALPHAB=BLRP/2.0+MB2-ZAP*GGP*VPVM/8.0
374 LAB=2.0*LOG(ALPHAB/QQQ2)+3.0
375C
376C calculate D^2 V / D^2 H_R
377C
378 VRRB=ZAP*V2*(GGP**2 + GG2**2)/16.0
379 VRRB=3.0*(VRRB*LAB)/32.0/PI2
380C
381C calculate D^2 V / D^2 H_R'
382C
383 VRPRPB=VP2*(GGP**2+GG2**2)/16.0-MB2*GGP
384 VRPRPB=ZAP*VRPRPB*LAB+8.0*FB2*MB2*LOG(ALPHAB/MB2)
385 VRPRPB=3.0*VRPRPB/32.0/PI2
386C
387C calculate D^2 V / D^H_R D^H_R'
388C
389 VRRPB=FB2*GGP-(GGP**2+GG2**2)/8.0
390 VRRPB=ZAP*VVP*VRRPB*LAB/2.0
391 VRRPB=3.0*VRRPB/32.0/PI2
392C
393 ENDIF
394C
395 DVRR=VRRT+VRRB+VP2*MHP2/VVPP + V2*GGP/2.0
396 DVRPRP=VRPRPT+VRPRPB+V2*MHP2/VVPP + VP2*GGP/2.0
397 DVRRP=VRRPT+VRRPB-VVP*MHP2/VVPP - VVP*GGP/2.0
398C TEMPH is always non-negative:
399 TEMPH=(DVRR-DVRPRP)**2+4*DVRRP**2
400 TEMPH=0.5*SQRT(TEMPH)
401 MHL2=0.5*(DVRR+DVRPRP)-TEMPH
402 MHH2=0.5*(DVRR+DVRPRP)+TEMPH
403 IF(MHL2.LT.0.0) THEN
404 MHLNEG=1
405C WRITE(LOUT,*) 'SSMHN: ERROR: MHL**2 < 0.0 FOR PARAMETERS:'
406C WRITE(LOUT,*) 'MHP =', AMHA, 'TANB =', 1.0/RR
407C WRITE(LOUT,*) 'MSTL=', AMTLSS, 'MSBL=', AMBLSS
408C WRITE(LOUT,*) 'MSTR=', AMTRSS, 'MSBR=', AMBRSS
409C WRITE(LOUT,*) 'AT=', AAT, 'AB=', AAB
410C WRITE(LOUT,*) 'MU=-2M1=', -EP
411C WRITE(LOUT,*) 'MT=', AMTP, 'MB=', AMBT
412C WRITE(LOUT,*) 'D-TERMS? 1=YES 2=NO :', INRAD
413C WRITE(LOUT,*) 'MASS SCALE (QQQ)=', SQRT(QQQ2)
414 AMHH=SQRT(MHH2)
415 AMHL=SQRT(ABS(MHL2))
416 GO TO 1000
417 ENDIF
418 AMHL=SQRT(MHL2)
419 AMHH=SQRT(MHH2)
420
421C
422C Now calculate mixing angle ALFAH
423C
424 TRACEM=DVRR-DVRPRP
425 TPAL=TRACEM**2 + 4.0*DVRRP**2
426 TANAH=TRACEM+SQRT(TPAL)
427 IF(DVRRP.EQ.0.0) THEN
428 WRITE(LOUT,*) 'SSMHN: OFF-DIAGONAL TERM OF SCALAR HIGGS',
429 $ ' MASS MATRIX IS ZERO '
430 IF(TANAH.NE.0.0) THEN
431 WRITE(LOUT,*) 'SSMHN: WARNING: TAN(ALFAH) FORMULA',
432 $ ' YIELDS INFINITY'
433 ELSE IF(TANAH.EQ.0.0) THEN
434 WRITE(LOUT,*) 'SSMHN: WARNING: TAN(ALFAH) FORMULA',
435 $ ' YIELDS 0/0 '
436 ENDIF
437 IF(DVRR.GT.DVRPRP) THEN
438 WRITE(LOUT,*) 'SSMHN: DVRR > DVRPRP ==> SET ALFAH=PI/2'
439 ALFAH = PI/2.0
440 ELSE IF (DVRR .LT. DVRPRP) THEN
441 WRITE(LOUT,*) 'SSMHN: DVRR < DVRPRP ==> SET ALFAH=0'
442 ALFAH = 0.0
443 ELSE IF (DVRR .EQ. DVRPRP) THEN
444 WRITE(LOUT,*) 'SSMHN: DVRR = DVRPRP ==> ALFAH INDETERMINANT'
445 WRITE(LOUT,*) 'SETTING SCALAR MIXING ANGLE ALPHA=PI/4'
446 ALFAH=PI/4.0
447 ENDIF
448 GO TO 1000
449 ENDIF
450 TANAH = -0.5*TANAH/DVRRP
451 ALFAH = ATAN(TANAH)
452C
4531000 RETURN
454 END