      SUBROUTINE MRQCOF(X,Y,SIG,NDATA,A,MA,LISTA,MFIT,ALPHA,BETA,NALP,
     *CHISQ,FUNCS)
	implicit real*8(a-h,o-z)
      PARAMETER (MMAX=20)
      DIMENSION X(NDATA),Y(NDATA),SIG(NDATA),ALPHA(NALP,NALP),BETA(MA),
     *    DYDA(MMAX),LISTA(MFIT),A(MA)
	common		elong,as,be, 		! extra parameters used by make_model
     &                  zs,ze,xs,xe,ccc,PI,au,au2,anu,alamb,freq,twopi,qxf,
     &			L1,L2,L3,L4,L5,L6,L7,L8,L9,A1,A2,A3,A4,A5,A6,A7,A8,A9,
     &			z,z0,R,p,Vx,Xnorm,XMnorm,ahpow,aanoise,
     &			er1,er2,er3,er4,er5,er6,er7,er8,er9,er10,
     &			ijsing,iFs(500),iFt(500),iFd(500),aispecq(500),ispecq1,iq1(500),ispecq2,iq2(500)
      DO 12 J=1,MFIT
        DO 11 K=1,J
          ALPHA(J,K)=0.0d0
11      CONTINUE
        BETA(J)=0.0d0
12    CONTINUE
      CHISQ=0.0d0
C
      DO 15 I=1,NDATA
C
C Begin B. Jackson modification 11/2013
C
	if(I.lt.500) then		
	  iq1(I)     = I
	  iq2(I)     = I
	  aispecq(I) = X(I)
	end if
	if(X(I) .lt. xs .or. X(I) .gt. xe) go to 15
C
C	write(*,'(A,4E15.7)') 'before funcs',A(1), A(2), A(3), A(4)
C
C End B. Jackson modification 11/2013
C
        CALL FUNCS(X(I),A,YMOD,DYDA,MA)
C	print *, I, X(I), A(1), A(2), Y(I), YMOD
        SIG2I=1.0d0/(SIG(I)*SIG(I))
        DY=Y(I)-YMOD					! Original (not weighted by power) BVJ 26 November 2013
C        DY=(Y(I)-YMOD)/Y(I)
        DO 14 J=1,MFIT
          WT=DYDA(LISTA(J))*SIG2I
          DO 13 K=1,J
            ALPHA(J,K)=ALPHA(J,K)+WT*DYDA(LISTA(K))
13        CONTINUE
          BETA(J)=BETA(J)+DY*WT
14      CONTINUE
        CHISQ=CHISQ+DY*DY*SIG2I
	print *, ' '
	write(*,'(A,I4,F10.4,F10.5,2F10.6,E15.7)'), 'In mrqcof', I, A(1)/1.0d3,A(2),A(3)/4.84813681d-6,A(4),CHISQ
	write(*,'(A,F12.7,A,F12.7,A,F8.5)') 'Data Power input =', Y(I),', Model Power = ',YMOD,', Frequency =',X(I)
	if(I.ne.NDATA) print *, ' '
15    CONTINUE
      DO 17 J=2,MFIT
        DO 16 K=1,J-1
          ALPHA(K,J)=ALPHA(J,K)
16      CONTINUE
17    CONTINUE
      RETURN
      END
