pro ipc_op_dir

	cd,'\\hawk.cis.rit.edu\RAID1\IPC\Op_Dir
	im=float(readfits('IPC_Diff_Image2.fits',h))
	y0=216
	x0=163.
	num=450
	ipcs=fltarr(5,5,15,450)
	sumipcs=fltarr(15,450)
	offset=0.

	; loop over the 15 repeated patterns across the array
	for j=0,14 do begin
		; loop over 450 central pixels going down the array
		for i=0,num-1 do begin
			x=x0+256.*j
			y=y0+8.*i
			ipcs(0:2,*,j,i)=(im(x-2:x,y-2:y+2)-offset)/(im(x,y)-offset)
			x=x0+256.*j-71
			ipcs(3:4,*,j,i)=(im(x+1:x+2,y-2:y+2)-offset)/(im(x,y)-offset)
			sumipcs(j,i)=round((total(ipcs(*,*,j,i))-1.)*10000.)/100.
		endfor
	endfor

	; collapse all 5x5 maps into a single median map
	medipc=median(median(ipcs,dimension=4),dimension=3)
	medipc_percent=round(medipc*10000.)/100.
	medsumipcs=median(sumipcs)
	print,medipc_percent

	openw,1,'ipc.txt'
	printf,1,medipc_percent
	close,1

	set_plot,'win'
	surface,medipc_percent,charsize=3,/lego,ztitle='percent'
	set_plot,'ps'
	device, filename = 'ipc.ps',/landscape
	!p.font = -1    ; use vector fonts for PS
	!P.MULTI = 0    ; Set up for one plot per page
	cd,current=curdir
	fsub = curdir+'\ipc.ps'   ; set text to print at bottom of plot

	; set some size values
	asize = 1.3    ; charsize of annotations
	athick = 4.0    ; charthick of annotations
	lthick = 8.0    ; line thickness
	!x.thick=lthick
	!y.thick=lthick
	!z.thick=lthick

	surface,medipc_percent,charsize=asize*2,/lego,ztitle='percent',thick=lthick,charthick=athick,xtitle='x pixel',ytitle='y pixel'
	; Add labels & fitted line
	; write the file name outside the margin of the plot
	black=0
	xyouts, 1.0, 40, fsub, charsize = asize*0.75, charthick = athick, color = black, /device

	device,/close

	; work on sumipcs
	outliers=where(sumipcs gt median(sumipcs)+stddev(sumipcs,/nan)*2. or sumipcs lt median(sumipcs)-stddev(sumipcs,/nan)*2.,noutliers )
	if (noutliers ne 0) then sumipcs(outliers)=!VALUES.F_NAN
	outliers=where(sumipcs gt median(sumipcs)+stddev(sumipcs,/nan)*2. or sumipcs lt median(sumipcs)-stddev(sumipcs,/nan)*2.,noutliers )
	if (noutliers ne 0) then sumipcs(outliers)=!VALUES.F_NAN
	device, filename = 'sumipc.ps',/landscape
	fsub = curdir+'\sumipc.ps'   ; set text to print at bottom of plot
	surface,sumipcs,[indgen(15)*256.+x0],[indgen(450)*8.+y0],charsize=asize*2,/lego,ztitle='percent',thick=lthick,charthick=athick,xtitle='x pixel',ytitle='y pixel'
	xyouts, 1.0, 40, fsub, charsize = asize*0.75, charthick = athick, color = black, /device
	device,/close
stop
end

