# Written, prompt="Error Code" by Phil Massey, Lowell Obs.  March 21, 2002
# Added support for sround
# Added support for finding the "best" fwhm for fitpsf
# some bug!
#Changed how we did mode April 27, 2005
#Changed how we did min
#Same as "character6" except that the gain is included in the call.
#Included erorr checking ierr=0 is OK; ierr>1 is not OK
#Necessitated by txdump...
        procedure character6 (img,gain,rdnoise,saturate,fname,timestamp)
        string img {prompt="Image"}
        real gain {3., prompt="gain"}
        real rdnoise {6., prompt="readnoise"}
        real saturate {30000., prompt="saturation"}
	string fname {"filter", prompt="filter key word"}
	string timestamp {"time-obs", prompt="time key word"}

	string *inlist

        begin

        string s1, s2, s3, s4, s5, s6, s7
        int ierr


real fwhmv,gainv,rdnoisev,saturatev,f1,fd1,fd2,fd3,rms,bot,top,round,best,val,num

# Get parameters
	s1=img
        gainv=gain
	rdnoisev=rdnoise
	saturatev=saturate
	s6=fname
	s7=timestamp
	delete("errorfile")
        hedit(s1,"gain",gainv,add+,verify-,show+)
	hedit(s1,"rdnoise",rdnoisev,add+,verify-,show+)
	hedit(s1,"saturate",saturatev,add+,verify-,show+)
	hedit(s1,"fname",s6,add+,verify-,show+)
	hedit(s1,"timestamp",s7,add+,verify-,show+)
# Now find mode:
	delete("temp")
#REPLACED MODAL DEFINTION
#	imstat(s1,fields="mode",lower=-500,upper=20000,format-,>"temp")
#
	imhist(s1,z1=-500,z2=20000,binwidth=1,nbins=INDEF,listout=yes,>"temp")
	delete("temp2")
	sort("temp",col=2,num+,rev+,>"temp2")
	inlist="temp2"
	i=fscan(inlist,f1)
	hedit(s1,"mode",f1,add+,verify-,show+)
	print "mode=",(f1)
	rms=sqrt(f1*gainv+rdnoisev**2)/gainv
	hedit(s1,"sigma",rms,add+,verify-,show+)
	print "rms=",(rms)
	val=f1-5*rms
        print "Final value for DATAMI=",(val)
	bot=val
	hedit(s1,"DATAMI",bot,add+,veri-,show+)
	top=saturatev-1500.
	hedit(s1,"DATAMA",top,add+,veri-,show+)
	delete("tempxy")
	print("running daofind in character")
daofind(s1,output="tempxy",sigma=rms,fwhm=4.0,thresh=30.0,datamin=bot,datamax=top,interact-,verify-,verbose-,filter=s6,obstime=s7)
flpr
	delete("tempround")
	if(!access("tempxy")) {
        ierr=1
        goto err_exit
        }
	txdump("tempxy","sround","yes",>"tempround")
	delete("temproundans")
	pickround
	inlist="temproundans"
	i=fscan(inlist,round)
	hedit(s1,"sround",round,add+,verify-,show+,update+)
	print("running fitpsf in character:",(rms),(bot),(top))
# 4.0
	print("PSF=4.0")
	delete("temprad4.0")
fitpsf(s1,12,sigma=rms,coords="tempxy",output="temprad4.0",datamin=bot,datamax=top,interact-,verify-,verbose-,cache+,filter=s6,obstime=s7,fwhmpsf=4.0)
flpr
	print("done running fitpsf")
	delete("rad4.0")
	if(!access("temprad4.0")){
	ierr=2
	goto err_exit
        }
	txdump("temprad4.0",fields="rsigma",expr="rsigma!=INDEF",>"rad4.0")

# 5.0:
	print("PSF=5.0")
	delete("temprad5.0")
fitpsf(s1,15,sigma=rms,coords="tempxy",output="temprad5.0",datamin=bot,datamax=top,interact-,verify-,verbose-,cache+,filter=s6,obstime=s7,fwhmpsf=5.0)
flpr
	print("done running fitpsf")
	delete("rad5.0")
	if(!access("temprad5.0")){
	ierr=3
	goto err_exit
        }
	txdump("temprad5.0",fields="rsigma",expr="rsigma!=INDEF",>"rad5.0")
# 6.0:
        print("PSF=6.0")
	delete("temprad6.0")
fitpsf(s1,18,sigma=rms,coords="tempxy",output="temprad6.0",datamin=bot,datamax=top,interact-,verify-,verbose-,cache+,filter=s6,obstime=s7,fwhmpsf=6.0)
flpr
        print("done running fitpsf")
        delete("rad6.0")
	if(!access("temprad6.0")){
	ierr=4
	goto err_exit
        }
        txdump("temprad6.0",fields="rsigma",expr="rsigma!=INDEF",>"rad6.0")

# 7.0:
	print("PSF=7.0")
	delete("temprad7.0")
fitpsf(s1,21,sigma=rms,coords="tempxy",output="temprad7.0",datamin=bot,datamax=top,interact-,verify-,verbose-,cache+,filter=s6,obstime=s7,fwhmpsf=7.0)
flpr
	print("done running fitpsf")
	delete("rad7.0")
	if(!access("temprad7.0")){
	ierr=5
        goto err_exit
        }
	txdump("temprad7.0",fields="rsigma",expr="rsigma!=INDEF",>"rad7.0")
# 8.0:
	print("PSF=8.0")
	delete("temprad8.0")
fitpsf(s1,24,sigma=rms,coords="tempxy",output="temprad8.0",datamin=bot,datamax=top,interact-,verify-,verbose-,cache+,filter=s6,obstime=s7,fwhmpsf=8.0)
	print("done running fitpsf")
flpr
	delete("rad8.0")
	if(!access("temprad8.0")){
	ierr=6
        goto err_exit
        }
	txdump("temprad8.0",fields="rsigma",expr="rsigma!=INDEF",>"rad8.0")
# 9.0:
	print("PSF=9.0")
	delete("rad9.0")
	delete("temprad9.0")
fitpsf(s1,27,sigma=rms,coords="tempxy",output="temprad9.0",datamin=bot,datamax=top,interact-,verify-,verbose-,cache+,filter=s6,obstime=s7,fwhmpsf=9.0)
flpr
	print("done running fitpsf")
	if(!access("temprad9.0")){
	ierr=7
        goto err_exit
        }
	txdump("temprad9.0",fields="rsigma",expr="rsigma!=INDEF",>"rad9.0")
#10.0
       print("PSF=10.")
       delete("temprad10X")
fitpsf(s1,30,sigma=rms,coords="tempxy",output="temprad10X",datamin=bot,datamax=top,interact-,verify-,verbose-,cache+,filter=s6,obstime=s7,fwhmpsf=10.0)
flpr        
print("done running fitpsf")        
delete("rad10X")      
if(!access("temprad10X")) {
	ierr=8
        goto err_exit
        }
txdump("temprad10X",fields="rsigma",expr="rsigma!=INDEF",>"rad10X")
#Now pick one
	delete("tempone")
	pickbig
	inlist="tempone"
	i=fscan(inlist,best)
#Do it again...
	print("BEST PSF=",(best))
	delete("temprad")
fitpsf(s1,15,sigma=rms,coords="tempxy",output="temprad",datamin=bot,datamax=top,interact-,verify-,verbose-,filter=s6,obstime=s7,fwhmpsf=best)
flpr
	print("done running fitpsf")
	delete("allrad")
	if(!access("temprad")){
	ierr=9
	goto err_exit
         }
	txdump("temprad",fields="rsigma",expr="rsigma!=INDEF",>"allrad")
	delete("tempfwhm")
        radbig
flpr
	inlist="tempfwhm"
	i=fscan(inlist,fwhmv)
#For consistency do it all again...
        print("One more time from the top...")
	delete("tempxy")
	print("running daofind in character...again...with"//(fwhmv)//" as a value")
daofind(s1,output="tempxy",sigma=rms,fwhm=fwhmv,thresh=30.0,datamin=bot,datamax=top,interact-,verify-,verbose-,filter=s6,obstime=s7)
flpr
	delete("tempround")
	if(!access("tempxy")) {
	ierr=10
	goto err_exit
	}
	txdump("tempxy","sround","yes",>"tempround")
	delete("temproundans")
	pickround
	inlist="temproundans"
	i=fscan(inlist,round)
	hedit(s1,"sround",round,add+,verify-,show+,update+)
# I'm sorry. This is dumb.  But apparently this is what we have to do in
# order to get a good value.
	print("running fitpsf in character:",(rms),(bot),(top))
# 4.0
	print("PSF=4.0")
	delete("temprad4.0")
fitpsf(s1,12,sigma=rms,coords="tempxy",output="temprad4.0",datamin=bot,datamax=top,interact-,verify-,verbose-,cache+,filter=s6,obstime=s7,fwhmpsf=4.0)
flpr
	print("done running fitpsf")
	delete("rad4.0")
	if(!access("temprad4.0")){
	ierr=11
	goto err_exit
	}
	txdump("temprad4.0",fields="rsigma",expr="rsigma!=INDEF",>"rad4.0")

# 5.0:
	print("PSF=5.0")
	delete("temprad5.0")
fitpsf(s1,15,sigma=rms,coords="tempxy",output="temprad5.0",datamin=bot,datamax=top,interact-,verify-,verbose-,cache+,filter=s6,obstime=s7,fwhmpsf=5.0)
flpr
	print("done running fitpsf")
	delete("rad5.0")
	if(!access("temprad5.0")){
	ierr=12
	goto err_exit
        }
	txdump("temprad5.0",fields="rsigma",expr="rsigma!=INDEF",>"rad5.0")
# 6.0:
        print("PSF=6.0")
	delete("temprad6.0")
fitpsf(s1,18,sigma=rms,coords="tempxy",output="temprad6.0",datamin=bot,datamax=top,interact-,verify-,verbose-,cache+,filter=s6,obstime=s7,fwhmpsf=6.0)
flpr
        print("done running fitpsf")
        delete("rad6.0")
	if(!access("temprad6.0")){
	ierr=13
	goto err_exit
        }
        txdump("temprad6.0",fields="rsigma",expr="rsigma!=INDEF",>"rad6.0")

# 7.0:
	print("PSF=7.0")
	delete("temprad7.0")
fitpsf(s1,21,sigma=rms,coords="tempxy",output="temprad7.0",datamin=bot,datamax=top,interact-,verify-,verbose-,cache+,filter=s6,obstime=s7,fwhmpsf=7.0)
flpr
	print("done running fitpsf")
	delete("rad7.0")
	if(!access("temprad7.0")){
	ierr=14
	goto err_exit
	}
	txdump("temprad7.0",fields="rsigma",expr="rsigma!=INDEF",>"rad7.0")
# 8.0:
	print("PSF=8.0")
	delete("temprad8.0")
fitpsf(s1,24,sigma=rms,coords="tempxy",output="temprad8.0",datamin=bot,datamax=top,interact-,verify-,verbose-,cache+,filter=s6,obstime=s7,fwhmpsf=8.0)
	print("done running fitpsf")
flpr
	delete("rad8.0")
	if(!access("temprad8.0")){
	ierr=15
	goto err_exit
         }
	txdump("temprad8.0",fields="rsigma",expr="rsigma!=INDEF",>"rad8.0")
# 9.0:
	print("PSF=9.0")
	delete("temprad9.0")
fitpsf(s1,27,sigma=rms,coords="tempxy",output="temprad9.0",datamin=bot,datamax=top,interact-,verify-,verbose-,cache+,filter=s6,obstime=s7,fwhmpsf=9.0)
flpr
	print("done running fitpsf")
	delete("rad9.0")
	if(!access("temprad9.0")){
	ierr=16
	goto err_exit
	}
	txdump("temprad9.0",fields="rsigma",expr="rsigma!=INDEF",>"rad9.0")

# 10.0:
        print("PSF=10.0")
        delete("temprad10X")
fitpsf(s1,30,sigma=rms,coords="tempxy",output="temprad10X",datamin=bot,datamax=top,interact-,verify-,verbose-,cache+,filter=s6,obstime=s7,fwhmpsf=10.0)
flpr
        print("done running fitpsf")
        delete("rad10X")
	if(!access("temprad10X")){
	ierr=17
	goto err_exit
	}
        txdump("temprad10X",fields="rsigma",expr="rsigma!=INDEF",>"rad10X")
#Now pick one
	delete("tempone")
	pickbig
	inlist="tempone"
	i=fscan(inlist,best)
#Now do it one last time for real:
	print("BEST PSF=",(best))
	delete("temprad")
fitpsf(s1,20,sigma=rms,coords="tempxy",output="temprad",datamin=bot,datamax=top,interact-,verify-,verbose-,filter=s6,obstime=s7,fwhmpsf=best)
flpr
	print("done running fitpsf")
	delete("allrad")
	if(!access("temprad")){
	ierr=18
	goto err_exit
	}
	txdump("temprad",fields="rsigma",expr="rsigma!=INDEF",>"allrad")
	delete("tempfwhm")
        radbig
flpr
	inlist="tempfwhm"
	i=fscan(inlist,fwhmv)
	hedit(s1,"fwhm",fwhmv,add+,ver-,show+)
	bye
err_exit: print(s1//" Error code "//ierr//" in character---ABORTING", >> "../ERRORLOG")
        print(s1//" Error code "//ierr//" in character---ABORTING")
	print(ierr,>>"errorfile")
end