!********************************
!***Basic Analysis of Data Set***
!********************************
	SUBROUTINE PARMESTIMATE()
	USE WINTERACTER
	USE RESOURCE
	USE IONTOXMODEL_INTERFACEVARS
	USE IONTOXMODEL_MODELVARS
	INTEGER RESETERROR,IPARM,JPARM,ILINE,IDIR,IHIST,IPT,I,IHISTSAVE,K,IERR,IFLAG,IDAT
	REAL*8 MODELCALC,FUNCCALC,X,DEV,SUM,SUM2,PI/3.14159265D0/
    EXTERNAL GETDATA,GRAPHDATA,FUNCCALC,COMPLEXOPTIMIZE,MODELCALC,GRAPHLINE
!***
!***Update data from Data Set grid and plot data if this is basic estimation phase (not QA phase)
!***
	CALL GETDATA()
	CALL GRAPHDATA()
!***
!***Extract parameter guesses and limits from parameter grid
!***	
	TYPE *,'MODOPT',MODOPT
	DO IPARM=1,NPARM
		RESETERROR=InfoError(1)
		CALL WGridGetCellDouble(IDF_MAIN_PARAMETERS_GRID,1,IPARM,PARMGUESS(IPARM))
		CALL WGridGetCellDouble(IDF_MAIN_PARAMETERS_GRID,2,IPARM,PARMDELTA(IPARM))
		CALL WGridGetCellDouble(IDF_MAIN_PARAMETERS_GRID,3,IPARM,PARMMIN(IPARM))
		CALL WGridGetCellDouble(IDF_MAIN_PARAMETERS_GRID,4,IPARM,PARMMAX(IPARM))
	ENDDO
!***
!***Calculate AIC for null hypothesis (constant YDAT)
!***
	DATAMEAN=0.0
	SUM=0.0
	SUM2=0.0
	DO IDAT=1,NDAT
		SUM=SUM+YDAT(IDAT)
		SUM2=SUM2+YDAT(IDAT)**2
	ENDDO
	DATAMEAN=SUM/NDAT
	DATASTDV=SQRT((SUM2-SUM**2/NDAT)/NDAT)
	TYPE *,'MEAN,STDV',DATAMEAN,DATASTDV
	LIKENULL=0.0
	DO IDAT=1,NDAT
		DEV=YDAT(IDAT)-DATAMEAN
		LIKENULL=LIKENULL-LOG(EXP(-DEV**2/(2*DATASTDV**2))/(SQRT(2*PI)*DATASTDV))
	END DO
	AICNULL=2+2*LIKENULL
     CALL WDialogPutDouble(IDF_MAIN_AIC_DOUBLE2,AICNULL)
!***
!***Do optimization for best parameter fit using initial guesses and limits
!***
	NPOINT=10
	PRINTFLAG=.FALSE.
	TYPE "('OPTIMIZATION; INITIAL FUNC:',F15.8)",FUNCCALC(PARMGUESS)
	CALL COMPLEXOPTIMIZE(NPOINT,NPARM,PARMGUESS,PARMDELTA,PARMMIN,PARMMAX,PARMSOLUT,FUNCSOLUT,MODELERROR,NITER,NHIST,HISTITER,HISTFUNC,HISTPARM,PRINTFLAG)
	TYPE "('Optimization Errors',7I2,F15.8)",(MODELERROR(I),I=1,7),FUNCSOLUT
	DO IPARM=1,NPARM
		IF(MODELERROR(IPARM+1).EQ.0) THEN
			TYPE "(I3,I5,4F12.8)",IPARM,NITER,PARMSOLUT(IPARM),PARMGUESS(IPARM),PARMMIN(IPARM),PARMMAX(IPARM)
		ELSE
			TYPE "(I3,I5,4F12.8,' ***Convergence Error')",IPARM,NITER,PARMSOLUT(IPARM),PARMGUESS(IPARM),PARMMIN(IPARM),PARMMAX(IPARM)
		ENDIF
		CALL WGridPutCellDouble(IDF_MAIN_PARAMETERS_GRID,5,IPARM,PARMSOLUT(IPARM))
		CALL WGridPutCellInteger(IDF_MAIN_PARAMETERS_GRID,6,IPARM,MODELERROR(IPARM+1))
	ENDDO
	CALL WDialogPutInteger(IDF_MAIN_CONVERGE_INTEGER,MODELERROR(1))
	CALL WDialogPutDouble(IDF_MAIN_MAXLIKE_DOUBLE,FUNCSOLUT)
	NITERSAVE1=NITER
     AIC=2*NPARM+2*FUNCSOLUT
     CALL WDialogPutDouble(IDF_MAIN_AIC_DOUBLE,AIC)
!***
!***Optimization history for main parameter optimization saved in window grid
!***
	DO IHIST=1,NHIST
		CALL WGridPutCellInteger(IDF_MAIN_OPTIMIZATION_GRID,1,IHIST,HISTITER(IHIST))
		CALL WGridPutCellDouble(IDF_MAIN_OPTIMIZATION_GRID,10,IHIST,HISTFUNC(IHIST,1))
		CALL WGridPutCellDouble(IDF_MAIN_OPTIMIZATION_GRID,11,IHIST,HISTFUNC(IHIST,2))
		DO IPARM=1,NPARM
			CALL WGridPutCellDouble(IDF_MAIN_OPTIMIZATION_GRID,2*IPARM,IHIST,HISTPARM(IHIST,IPARM,1))
			CALL WGridPutCellDouble(IDF_MAIN_OPTIMIZATION_GRID,2*IPARM+1,IHIST,HISTPARM(IHIST,IPARM,2))
		ENDDO		
	ENDDO
!***
!***Save selected history values for later output 
!***
	DO IHISTSAVE=1,5
		IF(IHISTSAVE.EQ.1) THEN
			IHIST=1
		ELSEIF(IHISTSAVE.EQ.2) THEN
			IHIST=NHIST*0.25
		ELSEIF(IHISTSAVE.EQ.3) THEN
			IHIST=NHIST*0.50
		ELSEIF(IHISTSAVE.EQ.4) THEN
			IHIST=NHIST*0.75
		ELSE
			IHIST=NHIST
		ENDIF
		HISTITERSAVE(IHISTSAVE)=HISTITER(IHIST)
		HISTFUNCSAVE(IHISTSAVE,1)=HISTFUNC(IHIST,1)
		HISTFUNCSAVE(IHISTSAVE,2)=HISTFUNC(IHIST,2)
		DO IPARM=1,NPARM
			HISTPARMSAVE(IHISTSAVE,IPARM,1)=HISTPARM(IHIST,IPARM,1)
			HISTPARMSAVE(IHISTSAVE,IPARM,2)=HISTPARM(IHIST,IPARM,2)
		ENDDO
	ENDDO
 
!***
!***Compute and plot curve and enter points into grid
!***
	DO IPARM=1,NPARM
		PARMTEMP(IPARM)=PARMSOLUT(IPARM)
	ENDDO
	DO ILINE=1,100
		IF(ILINE.LE.90) THEN
    			XLINE(ILINE)=-2.0+0.98*(PARMSOLUT(1)+1.5)*(ILINE-1)/89.0
		ELSE
			XLINE(ILINE)=-2.0+0.98*(PARMSOLUT(1)+1.5)+0.0199999*(PARMSOLUT(1)+1.5)*(ILINE-90)/10
		ENDIF
	    	YLINE(ILINE)=MODELCALC(DBLE(XLINE(ILINE)))
		CALL WGridPutCellReal(IDF_MAIN_PREDLINE_GRID,1,ILINE,XLINE(ILINE))
		CALL WGridPutCellReal(IDF_MAIN_PREDLINE_GRID,2,ILINE,YLINE(ILINE))
	END DO
	CALL GRAPHLINE()
!***
!***Document whether solution is at local optimum
!***
	DO IPARM=1,NPARM
		DO JPARM=1,NPARM
			PARMTEMP(JPARM)=PARMSOLUT(JPARM)
		ENDDO
		DO I=-5,5
			PARMTEMP(IPARM)=PARMSOLUT(IPARM)+(I/20D0)*PARMDELTA(IPARM)*(PARMMAX(IPARM)-PARMMIN(IPARM))
			X=FUNCCALC(PARMTEMP)
			TYPE "(2I3,<NPARM>F11.6,F14.8)",IPARM,I,(PARMTEMP(JPARM),JPARM=1,NPARM),ABS(X)
		ENDDO
	ENDDO
	RETURN
!***
!***CONFIDENCE LIMITS PHASE
!***
	ENTRY CONFLIMITS()
!***
!***Extract parameter guesses and limits from parameter grid in case user changed them
!***
	DO IPARM=1,4
		CALL WGridGetCellDouble(IDF_MAIN_PARAMETERS_GRID,1,IPARM,CLPARMGUESS(IPARM))
		CALL WGridGetCellDouble(IDF_MAIN_PARAMETERS_GRID,2,IPARM,CLPARMDELTA(IPARM))
		CALL WGridGetCellDouble(IDF_MAIN_PARAMETERS_GRID,3,IPARM,CLPARMMIN(IPARM))
		CALL WGridGetCellDouble(IDF_MAIN_PARAMETERS_GRID,4,IPARM,CLPARMMAX(IPARM))
	ENDDO
!***
!***Initiate Loops and Set Index Variables 
!***
	IF(MODOPT.EQ.1) THEN
		NPARMCL=2
	ELSE
		NPARMCL=3
	ENDIF
	NPOINT=10
	DO IPARM=1,6
	DO IERR=0,4
		CLFLAGLCL(IPARM,IERR)=0
		CLFLAGUCL(IPARM,IERR)=0
	ENDDO
	ENDDO
	CLLIKETARGET=FUNCSOLUT+1.92
	TYPE "('***CLCALC*** START, CLLIKETARGET=',F10.3)",CLLIKETARGET
!***
!***Loop through parameters and set search range for CL based on parameter solutions and parameter limits
!***
	DO FIXPARM=1,NPARM
        	IF(FIXPARM.EQ.1) THEN
        		FITPARM1=2
        		FITPARM2=3
        		FITPARM3=4
        	ELSEIF(FIXPARM.EQ.2) THEN
        		FITPARM1=1
        		FITPARM2=3
        		FITPARM3=4
        	ELSEIF(FIXPARM.EQ.3) THEN
        		FITPARM1=1
        		FITPARM2=2
        		FITPARM3=4
        	ELSEIF(FIXPARM.EQ.4) THEN
        		FITPARM1=1
        		FITPARM2=2
        		FITPARM3=3
        	ENDIF
		PARMLCL(FIXPARM)=-1D6
		PARMUCL(FIXPARM)=+1D6
		FIXPARMINIT=PARMSOLUT(FIXPARM)
		FIXPARMMIN=CLPARMMIN(FIXPARM)
		FIXPARMMAX=CLPARMMAX(FIXPARM)
!***
!***Loop between lower and upper confidence limits
!***
		DO IDIR=-1,1,2
!***
!***If basic parameter solutions are at min or max, set confidence limits to limits without further calculation
!*** 
        		IF(IDIR.EQ.-1.AND.ABS((PARMSOLUT(FIXPARM)-PARMMIN(FIXPARM))/PARMMIN(FIXPARM)).LT.0.001) THEN
        			PARMLCL(FIXPARM)=PARMMIN(FIXPARM)
        			CLFLAGLCL(FIXPARM,0)=-1
				CLFLAGLCL(FIXPARM,FIXPARM)=-1
        			FUNCLCL(FIXPARM)=FUNCSOLUT 
        			FITPARMVALatLCL(FIXPARM,1)=PARMSOLUT(FITPARM1)
        			FITPARMVALatLCL(FIXPARM,2)=PARMSOLUT(FITPARM2)
        			FITPARMVALatLCL(FIXPARM,3)=PARMSOLUT(FITPARM3)
				TYPE "('***CLCALC*** OFFRAMP PARMSOLUT AT MIN',5I5)",fixparm,idir
        			CYCLE
			ELSEIF(IDIR.EQ.1.AND.ABS((PARMSOLUT(FIXPARM)-PARMMAX(FIXPARM))/PARMMAX(FIXPARM)).LT.0.001) THEN
        			PARMUCL(FIXPARM)=PARMMAX(FIXPARM)
        			CLFLAGUCL(FIXPARM,0)=-1
				CLFLAGUCL(FIXPARM,FIXPARM)=1
        			FUNCUCL(FIXPARM)=FUNCSOLUT
        			FITPARMVALatUCL(FIXPARM,1)=PARMSOLUT(FITPARM1)
        			FITPARMVALatUCL(FIXPARM,2)=PARMSOLUT(FITPARM2)
        			FITPARMVALatUCL(FIXPARM,3)=PARMSOLUT(FITPARM3)
				TYPE "('***CLCALC*** OFFRAMP PARMSOLUT AT MAX',5I5)",fixparm,idir
        			CYCLE
        		ENDIF

!***
!***Want to identify two values for FIXPARM that bracket target likelihood for confidence limits
!***First point is FIXPARMINIT, so that FITPARM values and likelihood is maximum likelihood solutions
!***For second point, incrementally increase deviation of FIXPARM from PARMSOLUT, reoptimizing  
!***FITPARMs each time, until the likelihood at FIXPARM exceeds CLLIKETARGET
!***FIXPARM incremented by 1/10 distance from FIXPARMINIT to FIXPARMMIN or FIXPARMMAX    
!***
			IF(IDIR.EQ.-1) THEN
				FIXPARMDELTA=(FIXPARMMIN-FIXPARMINIT)/10
			ELSE
				FIXPARMDELTA=(FIXPARMMAX-FIXPARMINIT)/10
			ENDIF
			FITPARMSOLUT(1)=PARMSOLUT(FITPARM1)
			FITPARMSOLUT(2)=PARMSOLUT(FITPARM2)
			FITPARMSOLUT(3)=PARMSOLUT(FITPARM3)
		 	FITPARMMIN(1)=CLPARMMIN(FITPARM1)
 			FITPARMMIN(2)=CLPARMMIN(FITPARM2)
 			FITPARMMIN(3)=CLPARMMIN(FITPARM3)
	 		FITPARMMAX(1)=CLPARMMAX(FITPARM1)
			FITPARMMAX(2)=CLPARMMAX(FITPARM2)
			FITPARMMAX(3)=CLPARMMAX(FITPARM3)
	 		FITPARMDELTA(1)=CLPARMDELTA(FITPARM1)
	 		FITPARMDELTA(2)=CLPARMDELTA(FITPARM2)
	 		FITPARMDELTA(3)=CLPARMDELTA(FITPARM3)
			FIXPARMVALOLD=PARMSOLUT(FIXPARM)
			LIKEVALUEOLD=FUNCSOLUT
			DO IPT=1,10
				FIXPARMVAL=FIXPARMINIT+IPT*FIXPARMDELTA
       			FITPARMGUESS(1)=FITPARMSOLUT(1)
        			FITPARMGUESS(2)=FITPARMSOLUT(2)
        			FITPARMGUESS(3)=FITPARMSOLUT(3)
	 			CALL COMPLEXOPTIMIZE(NPOINT,NPARMCL,FITPARMGUESS,FITPARMDELTA,FITPARMMIN,FITPARMMAX,FITPARMSOLUT,LIKEVALUE,CLERROR,NITER,NHISTCL,HISTITERCL,HISTFUNCCL,HISTPARMCL,PRINTFLAG)
				IF(IPHASE.EQ.2) TYPE "('***CLCALC CONFLIMIT BRACKET***',3I3,F12.5,2X,3F10.5,2X,3F10.5,2X,F10.3,2X,4I3)",IPARM,IDIR,IPT,FIXPARMVAL,(FITPARMGUESS(K),K=1,3),(FITPARMSOLUT(K),K=1,3),LIKEVALUE,(CLERROR(K),K=1,4)
				IF(LIKEVALUE.GT.CLLIKETARGET) EXIT
				FIXPARMVALOLD=FIXPARMVAL
				LIKEVALUEOLD=LIKEVALUE
			ENDDO
!***
!***If loop reaches its limit, routine terminates with LCL or UCL set to min/max and CL flags set to indicate this
!***
			IF(IPT.GT.10) THEN
            		IF(IDIR.EQ.-1) THEN
            			PARMLCL(FIXPARM)=FIXPARMVAL
           			FUNCLCL(FIXPARM)=LIKEVALUE
            			FITPARMVALatLCL(FIXPARM,1)=FITPARMSOLUT(1)
            			FITPARMVALatLCL(FIXPARM,2)=FITPARMSOLUT(2)
            			FITPARMVALatLCL(FIXPARM,3)=FITPARMSOLUT(3)
            			CLFLAGLCL(FIXPARM,0)=-1
    					CLFLAGLCL(FIXPARM,FIXPARM)=-1
    					IF(IPHASE.EQ.2) TYPE "('***CLCALC*** LCL MIN OFFRAMP',5I5)",IDIR,FIXPARM
            			CYCLE
            		ELSE
            			PARMUCL(FIXPARM)=FIXPARMVAL
            			FUNCUCL(FIXPARM)=LIKEVALUE
            			FITPARMVALatUCL(FIXPARM,1)=FITPARMSOLUT(1)
            			FITPARMVALatUCL(FIXPARM,2)=FITPARMSOLUT(2)
            			FITPARMVALatUCL(FIXPARM,3)=FITPARMSOLUT(3)
            			CLFLAGUCL(FIXPARM,0)=-1
    					CLFLAGUCL(FIXPARM,FIXPARM)=1
    					IF(IPHASE.EQ.2) TYPE "('***CLCALC*** UCL MAX OFFRAMP',5I5)",IDIR,FIXPARM
            			CYCLE
            		ENDIF
			ENDIF
!***
!***Have bracketed points for confidence limit, so set binary search points 1 and 3 to these bracketing points, search point 2 to their midpoint
!***(FITPARMGUESS is set to most recent solution rather than midpoint of most recent solution and previous one) 
!***
        		FIXPARMSEARCH(1)=FIXPARMVALOLD
        		FIXPARMSEARCH(2)=(FIXPARMVALOLD+FIXPARMVAL)/2
			FIXPARMSEARCH(3)=FIXPARMVAL
        		LIKESEARCH(1)=LIKEVALUEOLD
			LIKESEARCH(3)=LIKEVALUE
        		FIXPARMVAL=FIXPARMSEARCH(2)
       		FITPARMGUESS(1)=FITPARMSOLUT(1)
        		FITPARMGUESS(2)=FITPARMSOLUT(2)
        		FITPARMGUESS(3)=FITPARMSOLUT(3)
        	 	CALL COMPLEXOPTIMIZE(NPOINT,NPARMCL,FITPARMGUESS,FITPARMDELTA,FITPARMMIN,FITPARMMAX,FITPARMSOLUT,LIKEVALUE,CLERROR,NITER,NHISTCL,HISTITERCL,HISTFUNCCL,HISTPARMCL,PRINTFLAG)
        		LIKESEARCH(2)=LIKEVALUE
        		IF(IPHASE.EQ.2) TYPE "('***CLCALC*** CONF LIMIT BISECT:',I3,3F10.5,3F10.3,3F10.5,2X,7I4)",0,(FIXPARMSEARCH(IPT),IPT=1,3),(LIKESEARCH(IPT),IPT=1,3),(FITPARMSOLUT(IPARM),IPARM=1,3),(CLERROR(I),I=1,4)
        		DO K=1,20
        			IF(LIKESEARCH(2).GT.CLLIKETARGET) THEN
        				FIXPARMSEARCH(3)=FIXPARMSEARCH(2)
        				LIKESEARCH(3)=LIKESEARCH(2)
        				FIXPARMSEARCH(2)=0.5*(FIXPARMSEARCH(1)+FIXPARMSEARCH(2))
        				FITPARMGUESS(1)=FITPARMSOLUT(1)
        				FITPARMGUESS(2)=FITPARMSOLUT(2)
        				FITPARMGUESS(3)=FITPARMSOLUT(3)
        				FIXPARMVAL=FIXPARMSEARCH(2)
        	 			CALL COMPLEXOPTIMIZE(NPOINT,NPARM,FITPARMGUESS,FITPARMDELTA,FITPARMMIN,FITPARMMAX,FITPARMSOLUT,LIKEVALUE,CLERROR,NITER,NHISTCL,HISTITERCL,HISTFUNCCL,HISTPARMCL,PRINTFLAG)
        				LIKESEARCH(2)=LIKEVALUE
        				IF(IPHASE.EQ.2) THEN
						TYPE "('***CLCALC*** CONF LIMIT BISECT:',I3,3F10.5,3F10.3,3F10.5,2X,7I4)",K,(FIXPARMSEARCH(IPT),IPT=1,3),(LIKESEARCH(IPT),IPT=1,3),(FITPARMSOLUT(IPARM),IPARM=1,3),(CLERROR(I),I=1,4)
					ENDIF
        			ELSE
        				FIXPARMSEARCH(1)=FIXPARMSEARCH(2)
        				LIKESEARCH(1)=LIKESEARCH(2)
        				FIXPARMSEARCH(2)=0.5*(FIXPARMSEARCH(3)+FIXPARMSEARCH(2))
        				FITPARMGUESS(1)=FITPARMSOLUT(1)
        				FITPARMGUESS(2)=FITPARMSOLUT(2)
        				FITPARMGUESS(3)=FITPARMSOLUT(3)
        				FIXPARMVAL=FIXPARMSEARCH(2)
        	 			CALL COMPLEXOPTIMIZE(NPOINT,NPARM,FITPARMGUESS,FITPARMDELTA,FITPARMMIN,FITPARMMAX,FITPARMSOLUT,LIKEVALUE,CLERROR,NITER,NHISTCL,HISTITERCL,HISTFUNCCL,HISTPARMCL,PRINTFLAG)
        				LIKESEARCH(2)=LIKEVALUE
        				IF(IPHASE.EQ.2) THEN
						TYPE "('***CLCALC*** CONF LIMIT BISECT:',I3,3F10.5,3F10.3,3F10.5,2X,7I4)",K,(FIXPARMSEARCH(IPT),IPT=1,3),(LIKESEARCH(IPT),IPT=1,3),(FITPARMSOLUT(IPARM),IPARM=1,3),(CLERROR(I),I=1,4)
					ENDIF
        			ENDIF
        			IF(ABS(FIXPARMSEARCH(3)-FIXPARMSEARCH(1)).LT.1D-6) EXIT
				!ACCEPT *,IGO
        		ENDDO
			PRINTFLAG=.FALSE.
!***
!***At end of binary search, store calculated confidence limit value, the likelihood statistic at the CL, 
!***the values of FITPARMs at Cl, and error flags (0 if converged at LIKETARGET, other values for reaching min/max or nonconvergence  
!***
        		IF(IDIR.EQ.-1) THEN
        			PARMLCL(FIXPARM)=FIXPARMSEARCH(2)
        			FUNCLCL(FIXPARM)=LIKESEARCH(2)
        			FITPARMVALatLCL(FIXPARM,1)=FITPARMSOLUT(1)
        			FITPARMVALatLCL(FIXPARM,2)=FITPARMSOLUT(2)
        			FITPARMVALatLCL(FIXPARM,3)=FITPARMSOLUT(3)
				CLFLAGLCL(FIXPARM,0)=CLERROR(1)
				IF(ABS(CLLIKETARGET-FUNCLCL(FIXPARM)).GT.0.01) CLFLAGLCL(FIXPARM,0)=-1
				CLFLAGLCL(FIXPARM,FITPARM1)=CLERROR(2)
				CLFLAGLCL(FIXPARM,FITPARM2)=CLERROR(3)
				CLFLAGLCL(FIXPARM,FITPARM3)=CLERROR(4)
        		ELSE
        			PARMUCL(FIXPARM)=FIXPARMSEARCH(2)
        			FUNCUCL(FIXPARM)=LIKESEARCH(2)
        			FITPARMVALatUCL(FIXPARM,1)=FITPARMSOLUT(1)
        			FITPARMVALatUCL(FIXPARM,2)=FITPARMSOLUT(2)
        			FITPARMVALatUCL(FIXPARM,3)=FITPARMSOLUT(3)
				CLFLAGUCL(FIXPARM,0)=CLERROR(1)
				IF(ABS(CLLIKETARGET-FUNCUCL(FIXPARM)).GT.0.01) CLFLAGUCL(FIXPARM,0)=1
				CLFLAGUCL(FIXPARM,FITPARM1)=CLERROR(2)
				CLFLAGUCL(FIXPARM,FITPARM2)=CLERROR(3)
				CLFLAGUCL(FIXPARM,FITPARM3)=CLERROR(4)
        		ENDIF
        	ENDDO
	ENDDO
!***
!***Display CLs on parameter grid, with error flags
!***
	DO FIXPARM=1,NPARM
		CALL WGridPutCellDouble(IDF_MAIN_PARAMETERS_GRID,7,FIXPARM,PARMLCL(FIXPARM))
		DO IFLAG=0,4
			CALL WGridPutCellInteger(IDF_MAIN_PARAMETERS_GRID,8+IFLAG,FIXPARM,CLFLAGLCL(FIXPARM,IFLAG))
		ENDDO
		CALL WGridPutCellDouble(IDF_MAIN_PARAMETERS_GRID,13,FIXPARM,PARMUCL(FIXPARM))
		DO IFLAG=0,4
			CALL WGridPutCellInteger(IDF_MAIN_PARAMETERS_GRID,14+IFLAG,FIXPARM,CLFLAGUCL(FIXPARM,IFLAG))
		ENDDO
	ENDDO
	ANALYSISDONE=.TRUE.
	RETURN
	END SUBROUTINE PARMESTIMATE 



!******************************************
!***Routine for Function Being Minimized***
!******************************************
	REAL*8 FUNCTION FUNCCALC(PARMCURR)
	USE IONTOXMODEL_MODELVARS
	REAL*8 PARMCURR(6),FUNCCALCFUL,FUNCCALCRED
    EXTERNAL FUNCCALCFUL,FUNCCALCRED,MODELCALC
	IF(IPHASE.EQ.1) THEN
		FUNCCALC=FUNCCALCFUL(PARMCURR)
	ELSE
		FUNCCALC=FUNCCALCRED(PARMCURR)
	ENDIF
	RETURN
	END
	
	REAL*8 FUNCTION FUNCCALCFUL(PARMCURR)
	USE IONTOXMODEL_MODELVARS
	REAL*8 MODELCALC,YO,YP
	REAL*8 PARMCURR(6),SUM,DEV,XCALL,STDDEV
	REAL*8 PI/3.14159265D0/
	INTEGER IDAT,IPARM
    EXTERNAL MODELCALC
	DO IPARM=1,NPARM
		PARMTEMP(IPARM)=PARMCURR(IPARM)
	ENDDO
	STDDEV=PARMTEMP(NPARM)
	SUM=0.0
	DO IDAT=1,NDAT
		KDAT=IDAT
		YO=YDAT(IDAT)
		XCALL=MIN(0.9999D0*PARMTEMP(1),XDAT(IDAT))
		YP=MODELCALC(XCALL)
		DEV=YO-YP
		SUM=SUM-LOG(EXP(-DEV**2/(2*STDDEV**2))/(SQRT(2*PI)*STDDEV))
	END DO
	FUNCCALCFUL=SUM
	RETURN
	END

	REAL*8 FUNCTION FUNCCALCRED(PARMCURR)
	USE IONTOXMODEL_MODELVARS
	REAL*8 MODELCALC,YO,YP
	REAL*8 PARMCURR(6),SUM,DEV,XCALL
	REAL*8 PI/3.14159265D0/
	INTEGER IDAT,IPARM
    EXTERNAL MODELCALC
	DO IPARM=1,FIXPARM-1
		PARMTEMP(IPARM)=PARMCURR(IPARM)
	ENDDO
	PARMTEMP(FIXPARM)=FIXPARMVAL
	DO IPARM=FIXPARM+1,NPARM
		PARMTEMP(IPARM)=PARMCURR(IPARM-1)
	ENDDO
	SUM=0.0
	DO IDAT=1,NDAT
		KDAT=IDAT
		YO=YDAT(IDAT)
		XCALL=MIN(0.9999D0*PARMTEMP(1),XDAT(IDAT))
		YP=MODELCALC(XCALL)
		DEV=YO-YP
		SUM=SUM-LOG(EXP(-DEV**2/(2*PARMTEMP(NPARM)**2))/(SQRT(2*PI)*PARMTEMP(NPARM)))
	END DO
	IF(PRINTFLAG) TYPE "('FUNC',I4,F12.6,5F10.6)",NDAT,SUM,(PARMTEMP(IPARM),IPARM=1,4),YP
	FUNCCALCRED=SUM
	RETURN
	END

!********************************************
!***Routine for Predicted Effects Variable
!***Equations Documented in "ModelEquations.docx"
!********************************************
	REAL*8 FUNCTION MODELCALC(LOGCA)
	USE IFPORT
	USE IONTOXMODEL_MODELVARS
	REAL*8 SLOPE,CAT,SHARP,LOGCAMIN
	REAL*8 LOGCA,CA
	REAL*8 LCXONLY,LCXONLY1
	REAL*8 LOGLCXONLY,LOGLCXONLY1
	IF(MODOPT.EQ.1) THEN      !NoCaTox, ModelA-LinearLogXvsLogCa
		LOGLCXONLY1=PARMTEMP(1)
		SLOPE=PARMTEMP(2)
		LOGLCXONLY=LOGLCXONLY1+SLOPE*LOGCA
		MODELCALC=LOGLCXONLY
	ELSEIF(MODOPT.EQ.2) THEN     !NoCaTox, ModelB-Linear@LoCa,ExponentialAsymptote@HiCa
		CA=10D0**LOGCA
		LCXONLY1=10D0**PARMTEMP(1)
		SLOPE=PARMTEMP(2)
		CAT=10D0**PARMTEMP(3)
		LCXONLY=LCXONLY1*(1D0-10D0**(-(CA/CAT)**SLOPE))
		LOGLCXONLY=LOG10(LCXONLY) 
		MODELCALC=LOGLCXONLY
	ELSEIF(MODOPT.EQ.3) THEN     !NoCaTox, ModelC-Hyperbola(onLogScale)
		LOGLCXONLY1=PARMTEMP(1)
		LOGCAMIN=PARMTEMP(2)
		SHARP=10D0**PARMTEMP(3)
		LOGLCXONLY=LOGLCXONLY1-SHARP/(LOGCA-LOGCAMIN)
		MODELCALC=LOGLCXONLY
	ELSEIF(MODOPT.EQ.4) THEN     !NoCaTox, ModelD-Rationale expression raised to power
		CA=10D0**LOGCA
		LCXONLY1=10D0**PARMTEMP(1)
		SLOPE=PARMTEMP(2)
		CAT=10D0**PARMTEMP(3)
		LCXONLY=LCXONLY1*(CA/(CA+CAT))**SLOPE
		LOGLCXONLY=LOG10(LCXONLY) 
		MODELCALC=LOGLCXONLY
	ELSEIF(MODOPT.EQ.5) THEN     !NoCaTox, ModelF-Hockey stick
		LOGLCXONLY1=PARMTEMP(3)
		SLOPE=PARMTEMP(2)
		LOGLCXONLY=MIN(PARMTEMP(1),LOGLCXONLY1+SLOPE*LOGCA)
		MODELCALC=LOGLCXONLY
	ENDIF
	RETURN
	END FUNCTION MODELCALC

