diff options
Diffstat (limited to 'math/bevington/gridls.f')
| -rw-r--r-- | math/bevington/gridls.f | 102 |
1 files changed, 102 insertions, 0 deletions
diff --git a/math/bevington/gridls.f b/math/bevington/gridls.f new file mode 100644 index 00000000..e62d3219 --- /dev/null +++ b/math/bevington/gridls.f @@ -0,0 +1,102 @@ +C SUBROUTINE GRIDLS.F +C +C SOURCE +C BEVINGTON, PAGES 212-213. +C +C PURPOSE +C MAKE A GRID-SEARCH LEAST-SQUARES FIT TO DATA WITH A SPECIFIED +C FUNCTION WHICH IS NOT LINEAR IN COEFFICIENTS +C +C USAGE +C CALL GRIDLS (X, Y, SIGMAY, NPTS, NTERMS, MODE, A, DELTAA, +C SIGMAA, 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 SIGMAA - ARRAY OF STANDARD DEVIATIONS 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 DELTAA VALUES ARE MODIFIED BY THE PROGRAM +C + SUBROUTINE GRIDLS (X,Y,SIGMAY,NPTS,NTERMS,MODE,A,DELTAA, + *SIGMAA,YFIT,CHISQR) + DIMENSION X(1),Y(1),SIGMAY(1),A(1),DELTAA(1),SIGMAA(1), + *YFIT(1) + REAL FUNCTN + EXTERNAL FUNCTN + +C +11 NFREE=NPTS-NTERMS + FREE=NFREE + CHISQR=0. + IF (NFREE) 100,100,20 +20 DO 90 J=1,NTERMS +C +C EVALUATE CHI SQUARE AT FIRST TWO SEARCH POINTS +C +21 DO 22 I=1,NPTS +22 YFIT(I)=FUNCTN (X,I,A) +23 CHISQ1=FCHISQ (Y,SIGMAY,NPTS,NFREE,MODE,YFIT) + FN=0. + DELTA=DELTAA(J) +41 A(J)=A(J)+DELTA + DO 43 I=1,NPTS +43 YFIT(I)=FUNCTN (X,I,A) +44 CHISQ2=FCHISQ (Y,SIGMAY,NPTS,NFREE,MODE,YFIT) +45 IF (CHISQ1-CHISQ2) 51,41,61 +C +C REVERSE DIRECTION OF SEARCH IF CHI SQUARE IS INCREASING +C +51 DELTA=-DELTA + A(J)=A(J)+DELTA + DO 54 I=1,NPTS +54 YFIT(I)=FUNCTN (X,I,A) + SAVE=CHISQ1 + CHISQ1=CHISQ2 +57 CHISQ2=SAVE +C +C INCREMENT A(J) UNTIL CHI SQUARE INCREASES +C +61 FN=FN+1. + A(J)=A(J)+DELTA + DO 64 I=1,NPTS +64 YFIT(I)=FUNCTN (X,I,A) + CHISQ3=FCHISQ (Y,SIGMAY,NPTS,NFREE,MODE,YFIT) +66 IF (CHISQ3-CHISQ2) 71,81,81 +71 CHISQ1=CHISQ2 + CHISQ2=CHISQ3 + GOTO 61 +C +C FIND MINIMUM OF PARABOLA DEFINED BY LAST THREE POINTS +C +81 DELTA=DELTA*(1./(1.+(CHISQ1-CHISQ2)/(CHISQ3-CHISQ2))+0.5) +82 A(J)=A(J)-DELTA +83 SIGMAA(J)=DELTAA(J)*SQRT(2./(FREE*(CHISQ3-2.*CHISQ2+CHISQ1))) +84 DELTAA(J)=DELTAA(J)*FN/3. +90 CONTINUE +C +C EVALUATE FIT AND CHI SQUARE FOR FINAL PARAMETERS +C +91 DO 92 I=1,NPTS +92 YFIT(I)=FUNCTN (X,I,A) +93 CHISQR=FCHISQ (Y,SIGMAY,NPTS,NFREE,MODE,YFIT) +100 RETURN + END |