diff options
Diffstat (limited to 'math/bevington/gradls.f')
| -rw-r--r-- | math/bevington/gradls.f | 113 |
1 files changed, 113 insertions, 0 deletions
diff --git a/math/bevington/gradls.f b/math/bevington/gradls.f new file mode 100644 index 00000000..f6b88d68 --- /dev/null +++ b/math/bevington/gradls.f @@ -0,0 +1,113 @@ +C SUBROUTINE GRADLS.F +C +C SOURCE +C BEVINGTON, PAGES 220-222. +C +C PURPOSE +C MAKE A GRADIENT-SEARCH LEAST-SQUARES FIT TO DATA WITH A +C SPECIFIED FUNCTION WHICH IS NOT LINEAR IN COEFFICIENTS +C +C USAGE +C CALL GRADLS (X, Y, SIGMAY, NPTS, NTERMS, MODE, A, DELTAA, +C YFIT, CHISQR) +C +C DESCRIPTION OF PARAMETERS +C X - ARRAY OF DATA POINTS FOR INDEPENDENT VARIABLE +C Y - ARRAY OF DATA POINTS FOR DEPENDENT VARIABLE +C SIGMAY - ARRAY OF STANDARD DEVIATIONS FOR Y DATA POINTS +C NPTS - NUMBER OF PAIRS OF DATA POINTS +C NTERMS - NUMBER OF PARAMETERS +C MODE - DETERMINES METHOD OF WEIGHTING LEAST-SQUARES FIT +C +1 (INSTRUMENTAL) WEIGHT(I) = 1./SIGMAY(I)**2 +C 0 (NO WEIGHTING) WEIGHT(I) = 1. +C -1 (STATISTICAL) WEIGHT(I) = 1./Y(I) +C A - ARRAY OF PARAMETERS +C DELTAA - ARRAY OF INCREMENTS FOR PARAMETERS A +C YFIT - ARRAY OF CALCULATED VALUES OF Y +C CHISQR - REDUCED CHI SQUARE FOR FIT +C +C SUBROUTINES AND FUNCTION SUBPROGRAMS REQUIRED +C FUNCTN (X, I, A) +C EVALUATES THE FITTING FUNCTION FOR THE ITH TERM +C FCHISQ (Y, SIGMAY, NPTS, NFREE, MODE, YFIT) +C EVALUATES REDUCED CHI SQUARED FOR FIT TO DATA +C +C COMMENTS +C DIMENSION STATEMENT VALID FOR NTERMS UP TO 10 +C + SUBROUTINE GRADLS (X,Y,SIGMAY,NPTS,NTERMS,MODE,A,DELTAA, + *YFIT,CHISQR) + DIMENSION X(1),Y(1),SIGMAY(1),A(1),DELTAA(1),YFIT(1) + DIMENSION GRAD(10) + REAL FUNCTN + EXTERNAL FUNCTN + +C +C EVALUATE CHI SQUARE AT BEGINNING +C +11 NFREE=NPTS-NTERMS + IF (NFREE) 13,13,21 +13 CHISQR=0. + GOTO 110 +21 DO 22 I=1,NPTS +22 YFIT(I)=FUNCTN (X,I,A) + CHISQ1=FCHISQ (Y,SIGMAY,NPTS,NFREE,MODE,YFIT) +C +C EVALUATE GRADIENT OF CHI SQUARE +C +31 SUM=0. +32 DO 39 J=1,NTERMS + DELTA=0.1*DELTAA(J) + A(J)=A(J)+DELTA + DO 36 I=1,NPTS +36 YFIT(I)=FUNCTN (X,I,A) + A(J)=A(J)-DELTA + GRAD(J)=CHISQ1-FCHISQ (Y,SIGMAY,NPTS,NFREE,MODE,YFIT) +39 SUM=SUM+GRAD(J)**2 +41 DO 42 J=1,NTERMS +42 GRAD(J)=DELTAA(J)*GRAD(J)/SQRT(SUM) +C +C EVALUATE CHI SQUARE AT NEW POINT +C +51 DO 52 J=1,NTERMS +52 A(J)=A(J)+GRAD(J) +53 DO 54 I=1,NPTS +54 YFIT(I)=FUNCTN (X,I,A) + CHISQ2=FCHISQ (Y,SIGMAY,NPTS,NFREE,MODE,YFIT) +C +C MAKE SURE CHI SQUARE DECREASES +C +61 IF (CHISQ1-CHISQ2) 62,62,71 +62 DO 64 J=1,NTERMS + A(J)=A(J)-GRAD(J) +64 GRAD(J)=GRAD(J)/2. + GOTO 51 +C +C INCREMENT PARAMETERS UNTIL CHI SQUARE STARTS TO INCREASE +C +71 DO 72 J=1,NTERMS +72 A(J)=A(J)+GRAD(J) + DO 74 I=1,NPTS +74 YFIT(I)=FUNCTN (X,I,A) +75 CHISQ3=FCHISQ (Y,SIGMAY,NPTS,NFREE,MODE,YFIT) +76 IF (CHISQ3-CHISQ2) 81,91,91 +81 CHISQ1=CHISQ2 +82 CHISQ2=CHISQ3 + GOTO 71 +C +C FIND MINIMUM OF PARABOLA DEFINED BY LAST THREE POINTS +C +91 DELTA=1./(1.+(CHISQ1-CHISQ2)/(CHISQ3-CHISQ2))+0.5 + DO 93 J=1,NTERMS +93 A(J)=A(J)-DELTA*GRAD(J) + DO 95 I=1,NPTS +95 YFIT(I)=FUNCTN (X,I,A) + CHISQR=FCHISQ (Y,SIGMAY,NPTS,NFREE,MODE,YFIT) +101 IF (CHISQ2-CHISQR) 102,110,110 +102 DO 103 J=1,NTERMS +103 A(J)=A(J)+(DELTA-1.)*GRAD(J) +104 DO 105 I=1,NPTS +105 YFIT(I)=FUNCTN (X,I,A) +106 CHISQR=CHISQ2 +110 RETURN + END |