| 1 | SUBROUTINE MMOL4(Y,X,VAL,ARG,EPS,IER)
|
|---|
| 2 |
|
|---|
| 3 | C-----------------------------------------------------------------------
|
|---|
| 4 | C M(UON) MOL(IERE SCATTERING) 4 (POINT CONTINUED FRACT. INTERPOLATION)
|
|---|
| 5 | C
|
|---|
| 6 | C ROUTINE TAKEN FROM IBM SCIENTIFIC SUBROUTINE PACKAGE
|
|---|
| 7 | C ROUTINE TAKEN FROM GEANT321 (CERN)
|
|---|
| 8 | C 4 POINT CONTINUED FRACTION INTERPOLATION
|
|---|
| 9 | C THIS SUBROUTINE IS CALLED FROM MMOLIE
|
|---|
| 10 | C ARGUMENTS:
|
|---|
| 11 | C Y = INTERPOLATED VALUE FOR THE ARGUMENT X
|
|---|
| 12 | C X = ARGUMENT FOR Y
|
|---|
| 13 | C VAL = VALUE ARRAY
|
|---|
| 14 | C ARG = ARGUMENT ARRAY
|
|---|
| 15 | C EPS = DESIRED ACCURACY
|
|---|
| 16 | C IER = OUTPUT ERROR PARAMETER
|
|---|
| 17 | C 0 ACCURACY O.K.
|
|---|
| 18 | C 1 ACCURACY CAN NOT BE TESTED IN 4TH ORDER INTERPOLATION
|
|---|
| 19 | C 2 TWO IDENTICAL ELEMENTS IN THE ARGUMENT ARRAY
|
|---|
| 20 | C-----------------------------------------------------------------------
|
|---|
| 21 |
|
|---|
| 22 | IMPLICIT NONE
|
|---|
| 23 | REAL ARG(4),AUX,DELT,EPS,H,P1,P2,P3,Q1,Q2,Q3,VAL(4),X,Y,Z
|
|---|
| 24 | INTEGER I,II,III,IER,J,JEND
|
|---|
| 25 | C-----------------------------------------------------------------------
|
|---|
| 26 |
|
|---|
| 27 | IER = 1
|
|---|
| 28 | Y = VAL(1)
|
|---|
| 29 | P2 = 1.
|
|---|
| 30 | P3 = Y
|
|---|
| 31 | Q2 = 0.
|
|---|
| 32 | Q3 = 1.
|
|---|
| 33 | DO 16 I = 2,4
|
|---|
| 34 | II = 0
|
|---|
| 35 | P1 = P2
|
|---|
| 36 | P2 = P3
|
|---|
| 37 | Q1 = Q2
|
|---|
| 38 | Q2 = Q3
|
|---|
| 39 | Z = Y
|
|---|
| 40 | JEND = I - 1
|
|---|
| 41 | 3 AUX = VAL(I)
|
|---|
| 42 | DO 10 J = 1,JEND
|
|---|
| 43 | H = VAL(I) - VAL(J)
|
|---|
| 44 | IF ( ABS(H) .GT. 1.E-6*ABS(VAL(I)) ) GOTO 9
|
|---|
| 45 | IF ( ARG(I) .EQ. ARG(J) ) GOTO 17
|
|---|
| 46 | IF ( J .LT. JEND ) GOTO 8
|
|---|
| 47 | II = II + 1
|
|---|
| 48 | III = I + II
|
|---|
| 49 | IF ( III .GT. 4 ) GOTO 19
|
|---|
| 50 | VAL(I) = VAL(III)
|
|---|
| 51 | VAL(III) = AUX
|
|---|
| 52 | AUX = ARG(I)
|
|---|
| 53 | ARG(I) = ARG(III)
|
|---|
| 54 | ARG(III) = AUX
|
|---|
| 55 | GO TO 3
|
|---|
| 56 | 8 VAL(I) = 1.E36
|
|---|
| 57 | GO TO 10
|
|---|
| 58 | 9 VAL(I) = ( ARG(I)-ARG(J) ) / H
|
|---|
| 59 | 10 CONTINUE
|
|---|
| 60 | P3 = VAL(I) * P2 + ( X - ARG(I-1) ) * P1
|
|---|
| 61 | Q3 = VAL(I) * Q2 + ( X - ARG(I-1) ) * Q1
|
|---|
| 62 | IF ( Q3. NE. 0. ) THEN
|
|---|
| 63 | Y = P3 / Q3
|
|---|
| 64 | ELSE
|
|---|
| 65 | Y = 1.E36
|
|---|
| 66 | ENDIF
|
|---|
| 67 | DELT = ABS(Z-Y)
|
|---|
| 68 | IF ( DELT .LE. EPS ) GOTO 19
|
|---|
| 69 | 16 CONTINUE
|
|---|
| 70 | RETURN
|
|---|
| 71 | 17 IER = 2
|
|---|
| 72 | RETURN
|
|---|
| 73 | 19 IER = 0
|
|---|
| 74 | RETURN
|
|---|
| 75 | END
|
|---|
