!ORDER=1
print,'yr,mo,da,cam,filt,exp'
read,yr,mo,da,cam,filt,exp
img = INTARR(512,512,3)
rimg = img & srimg=img & ptimg = img
;g1 = img & g2 = img & g3 = img & g4 = img
stim = INTARR(512,512,3)
stimg = INTARR(512,512) & finimg = stimg
gstar = fltarr(512,512) & stimg2 = stimg
;rdmeandk,FIX(cam),FIX(exp),dk
FOR i=0,2 DO BEGIN
	path =STRING(cam,filt,da,cam,filt,exp,2*i+24,$
	   FORMAT='("~steele/data/aur/c",I1,"/f",I1,"/a",I2.2,2I1,I3.3,".",I3.3)')
	rdkimg,path,hb,ximg
	kh = gethd(hb)
	help,/st,kh.misc.tm
	;img(*,*,i) = ximg(*,*)-dk
	rimg(*,*,i) = ximg(*,*)
	IF i EQ 1 THEN head = hb
ENDFOR
stop
imgg=rimg
FOR i=0,2 DO BEGIN
	tpt=0 & pt=0 & pt2=0 & bpt= 0 & x=0

	ns=1.5
	ng=2.
	dr=4.

	dim=SIZE(rimg(*,*,i))
	g1=ABS(rimg(*,*,i)-SHIFT(rimg(*,*,i),dr,0))
	g1(dr-1,*)=0
	sd = STDEV(g1,m) & gt1=m+ns*sd
	g1big=WHERE(g1 GT gt1)
	IF (g1big(0) NE -1) THEN big1=g1(g1big) ELSE big1=-1

	g2=ABS(rimg(*,*,i)-SHIFT(rimg(*,*,i),0,dr))
	g2(*,dim(2)-dr)=0
	sd = STDEV(g2,m) & gt2=m+ns*sd
	g2big=WHERE(g2 GT gt2)
	IF (g2big(0) NE -1) THEN big2=g2(g2big) ELSE big2=-1

	g3=ABS(rimg(*,*,i)-SHIFT(rimg(*,*,i),-dr,0))
	g3(dim(1)-dr,*)=0
	sd = STDEV(g3,m) & gt3=m+ns*sd
	g3big=WHERE(g3 GT gt3)
	IF (g3big(0) NE -1) THEN big3=g3(g3big) ELSE big3=-1

	g4=ABS(rimg(*,*,i)-SHIFT(rimg(*,*,i),0,-dr))
	g4(*,dr-1)=0
	sd = STDEV(g4,m) & gt4=m+ns*sd
	g4big=WHERE(g4 GT gt4)
	IF (g4big(0) NE -1) THEN big4=g4(g4big) ELSE big4=-1

	gradthresh=MAX([gt1,gt2,gt3,gt4])      ; select the largest gradient
	
	vbig=WHERE(big1 GT gradthresh)            ; restrict big gradients to
	IF (vbig(0) NE -1) THEN g1big=g1big(vbig) ELSE g1big=-1;those greater than gradthresh
	vbig=WHERE(big2 GT gradthresh)
	IF (vbig(0) NE -1) THEN g2big=g2big(vbig) ELSE g2big=-1
	vbig=WHERE(big3 GT gradthresh)
	IF (vbig(0) NE -1) THEN g3big=g3big(vbig) ELSE g3big=-1
	vbig=WHERE(big4 GT gradthresh)
	IF (vbig(0) NE -1) THEN g4big=g4big(vbig) ELSE g4big=-1
	nbig=INTARR(dim(1),dim(2))    ; no. of big gradients at each pixel
	IF (g1big(0) NE -1) THEN nbig(g1big)=nbig(g1big)+1
	IF (g2big(0) NE -1) THEN nbig(g2big)=nbig(g2big)+1
	IF (g3big(0) NE -1) THEN nbig(g3big)=nbig(g3big)+1
	IF (g4big(0) NE -1) THEN nbig(g4big)=nbig(g4big)+1
	toobig=WHERE(nbig GE ng)
	r=toobig/dim(1)
	c=toobig MOD dim(1) 
	tmp=WHERE(((dr-1 LT c) AND (c LT dim(1)-dr)) AND $
	  ((dr-1 LT r) AND (r LT dim(2)-dr)))
	toobig=toobig(tmp)
	r=r(tmp) & c=c(tmp)
	stimg = INTARR(512,512)
	stimg(toobig) = rimg(toobig +i*512.^2.)

;here, we have either stars or hot spots
	pt = where((stimg GT 0) AND stimg GE SHIFT(stimg,1,0) $
	 AND stimg GE SHIFT(stimg,0,1) $
	 AND stimg GE SHIFT(stimg,-1,0) $
	 AND stimg GE SHIFT(stimg,0,-1) $
	 AND stimg GE SHIFT(stimg,0,3) $
	 AND stimg GE SHIFT(stimg,3,0) $
	 AND stimg GE SHIFT(stimg,0,-3) $
	 AND stimg GE SHIFT(stimg,-3,0) $
	 AND stimg GE SHIFT(stimg,1,1) $
	 AND stimg GE SHIFT(stimg,1,-1) $
	 AND stimg GE SHIFT(stimg,-1,1) $
	 AND stimg GE SHIFT(stimg,-1,-1))

	stimg2 = INTARR(512,512)
	stimg2(pt) = stimg(pt)
ptimg(*,*,i)=stimg2 gt 0

;	 s1 = SHIFT(stimg2 ,1,0)
;	 s3 = SHIFT(stimg2 ,0,1)
;	 s5 = SHIFT(stimg2,1,1)
;	 s2 = SHIFT(stimg2 ,-1,0)
;	 s4 = SHIFT(stimg2 ,0,-1)
;	 s6 = SHIFT(stimg2,1,-1)
;	 s7 = SHIFT(stimg2,-1,-1)
;;	 s8 = SHIFT(stimg2,-1,1)
;;
;	pt2 = WHERE(stimg2 AND NOT s1 AND NOT s2 AND NOT s3 AND NOT s4 $
;	 AND NOT s5 AND NOT s6 AND NOT s7 AND NOT s8)
;	IF pt2(0) NE -1 THEN tpt =  [tpt,pt2]
;
;	x = [(WHERE(stimg2 AND s1)) , (WHERE(stimg2 AND s3)) , $
;	    (WHERE(stimg2 AND s5))] 
;	
;	REPEAT BEGIN
;
;		bpt=x(WHERE(x NE -1))
;		IF bpt(0) NE -1 THEN BEGIN
;	   	     timg=INTARR(512,512)
;		     timg(bpt) = stimg2(bpt)
; 		     stimg2 = timg
;
;		     s1 = SHIFT(stimg2 ,1,0)
;		     s3 = SHIFT(stimg2 ,0,1)
;		     s5 = SHIFT(stimg2,1,1)
;		     s2 = SHIFT(stimg2 ,-1,0)
;		     s4 = SHIFT(stimg2 ,0,-1)
;		     s6 = SHIFT(stimg2,1,-1)
;	             s7 = SHIFT(stimg2,-1,-1)
;		     s8 = SHIFT(stimg2,-1,1)
;
;		     pt2 = WHERE(stimg2 AND NOT s1 AND NOT s2 AND NOT s3 AND NOT s4 $
;;;		      AND NOT s5 AND NOT s6 AND NOT s7 AND NOT s8)
;		     IF pt2(0) NE -1 THEN tpt =  [tpt,pt2]
;		ENDIF
;
;		x = [(WHERE(stimg2 AND s1)) , (WHERE(stimg2 AND s3)) , $
;		    (WHERE(stimg2 AND s5))]
;
;	ENDREP UNTIL((x(0) EQ -1) AND (x(1) EQ -1) AND (x(2) EQ -1))
;
;	tpt = tpt(1:*)
;	ptimg(tpt+i*65536) = stimg(tpt)
;	ptimg(*,*,i) = ptimg(*,*,i) GT 0
;;tvscl, stimg,i
ENDFOR

surr=INTARR(9,9)+1
;surr(0,0)=0 & surr(2,0)=0 & surr(0,2)=0 & surr(2,2)=0
d0img = DILATE(ptimg(*,*,0),surr)
d2img = DILATE(ptimg(*,*,2),surr)
finimg = d0img AND ptimg(*,*,1) AND d2img
;tvscl,ptimg(*,*,1),2
;tvscl,finimg,0
im = ptimg(*,*,1)-finimg
;tvscl,im,3
;tvscl,rimg(*,*,1)<6000,1
starpts = WHERE(d0img AND ptimg(*,*,1) AND d2img)

;window,2
;;tvscl,rimg(*,*,1)<500

xy,starpts,512,c2,r2
k=2.00
;plot,k*c2,k*(-r2+512),clip=[0,0,511,511],/device,xmargin=[0,0],ymargin=[0,0],/noerase,$
;symsize=2,position=[0,0,511,511],psym=4,xrange=[0,511],yrange=[0,511],xstyle=5, ystyle=5
plotmod,c2,r2,head,cam,filt,exp
;splotbig,c2,r2,head,cam,filt,exp
;window,2
;;tvscl,rimg(*,*,1)<500
;plot,c2,-r2+512,psym=1,/device,xmargin=[0,0],ymargin=[0,0],position=[0,0,512,512],/noerase,$
;xrange=[0,512],yrange=[0,512]

END