C SUBROUTINE POLFIT.F 
C
C SOURCE
C   BEVINGTON, PAGES 140-142.
C
C PURPOSE
C   MAKE A LEAST-SQUARES FIT TO DATA WITH A POLYNOMIAL CURVE
C      Y = A(1) + A(2)*X + A(3)*X**2 + A(3)*X**3 + . . .
C
C USAGE 
C   CALL POLFIT (X, Y, SIGMAY, NPTS, NTERMS, MODE, A, 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 COEFFICIENTS (DEGREE OF POLYNOMIAL + 1)
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 COEFFICIENTS OF POLYNOMIAL
C   CHISQR - REDUCED CHI SQUARE FOR FIT 
C
C SUBROUTINES AND FUNCTION SUBPROGRAMS REQUIRED 
C   DETERM (ARRAY, NORDER)
C      EVALUATES THE DETERMINANT OF A SYMETRICAL
C      TWO-DIMENSIONAL MATRIX OF ORDER NORDER
C
C COMMENTS
C   DIMENSION STATEMENT VALID FOR NTERMS UP TO 10
C
	SUBROUTINE POLFIT (X,Y,SIGMAY,NPTS,NTERMS,MODE,A,CHISQR)
	DOUBLE PRECISION SUMX,SUMY,XTERM,YTERM,ARRAY,CHISQ
	DIMENSION X(1),Y(1),SIGMAY(1),A(1)
	DIMENSION SUMX(19),SUMY(10),ARRAY(10,10)
C
C ACCUMULATE WEIGHTED SUMS
C
11	NMAX=2*NTERMS-1 
	DO 13 N=1,NMAX
13	SUMX(N)=0.
	DO 15 J=1,NTERMS
15	SUMY(J)=0.
	CHISQ=0.
21	DO 50 I=1,NPTS
	XI=X(I) 
	YI=Y(I) 
31	IF (MODE) 32,37,39
32	IF (YI) 35,37,33
33	WEIGHT=1./YI
	GOTO 41 
35	WEIGHT=1./(-YI) 
	GOTO 41 
37	WEIGHT=1.
	GOTO 41 
39	WEIGHT=1./SIGMAY(I)**2
41	XTERM=WEIGHT
	DO 44 N=1,NMAX
	SUMX(N)=SUMX(N)+XTERM
44	XTERM=XTERM*XI
45	YTERM=WEIGHT*YI 
	DO 48 N=1,NTERMS
	SUMY(N)=SUMY(N)+YTERM
48	YTERM=YTERM*XI
49	CHISQ=CHISQ+WEIGHT*YI**2
50	CONTINUE
C
C CONSTRUCT MATRICES AND CALCULATE COEFFICIENTS 
C
51	DO 54 J=1,NTERMS
	DO 54 K=1,NTERMS
	N=J+K-1 
54	ARRAY(J,K)=SUMX(N)
	DELTA=DETERM (ARRAY,NTERMS)
	IF (DELTA) 61,57,61
57	CHISQR=0.
	DO 59 J=1,NTERMS
59	A(J)=0. 
	GOTO 80 
61	DO 70 L=1,NTERMS
62	DO 66 J=1,NTERMS
	DO 65 K=1,NTERMS
	N=J+K-1 
65	ARRAY(J,K)=SUMX(N)
66	ARRAY(J,L)=SUMY(J)
70	A(L)=DETERM (ARRAY,NTERMS) /DELTA
C
C CALCULATE CHI SQUARE
C
71	DO 75 J=1,NTERMS
	CHISQ=CHISQ-2.*A(J)*SUMY(J)
	DO 75 K=1,NTERMS
	N=J+K-1 
75	CHISQ=CHISQ+A(J)*A(K)*SUMX(N)
76	FREE=NPTS-NTERMS
77	CHISQR=CHISQ/FREE
80	RETURN
	END