;+
; NAME: ircal_regrid
; PURPOSE: resample IRCAL images to compensate for anamorphic magnification
; USAGE:
; 		IRCAL is anamorphic, so images have different plate scales in
; 		X and Y. Congrid the data in order to fix that. 
;
; 		This routine doesn't do a perfect job of it since it's constrained
; 		to deal with integer pixels...
;
; INPUTS:
; KEYWORDS:
; OUTPUTS:
;
; HISTORY:
; 	Began 2005-07-26 02:57:29 by Marshall Perrin 
; 	2006-05-09	Some updates
;-

PRO ircal_regrid,image,header,platescale=platescale

	@instruments.include
	scalex = ircal.platescale
	scaley = ircal.platescale_y

	; what should be the new (regridded) platescale?
	;if ~(keyword_set(platescale)) then platescale=scalex
	if ~(keyword_set(platescale)) then platescale=0.040 ; Nyquist at J for Lick.
							; TODO would this be better at twice Nyquist?
	
	; The following is not perfect, since we assume that we can regrid the 
	; data into an integer number of pixels. 
	; This is imperfect. But the error is less than about 1 part in 1000
	; (for typical values)
	; and we don't know the IRCAL platescale to a better accuracy than that!
	;
	; ---
	;
	; An alternative algorithm uses interpolate with the /grid option; see
	; ircal_dewarp for that one. That's probably more correct ultimately
	; 
	sz = size(image)
	sx = sz[1]
	sy = sz[2]
	nsx = round(sz[1]*scalex/platescale)
	nsy = round(sz[2]*scaley/platescale)

	; Pixel area conversion factor, to preserve flux
	pixelarearatio = platescale^2 / (scalex*scaley)


	; should I set /MINUS_ONE when calling congrid?
	; I don't really understand it.
	;image_scl = congrid(image,sx,nsy,cubic=-0.5)
	;
	; Actually, since we want the extra to be evenly spaced on both
	; sides, the most appropriate thing is to use CMCONGRID with
	; /HALF_HALF. MDP 2006-05-10 00:49:26

	image_scl = cmcongrid(image,sx,nsy,cubic=-0.5,/half_half)


	

	; TODO do we need to adjust the image header at all here?
	; Well, there's this:
	sxaddpar,header,"NAXIS2",nsy

	; perhaps we should update the WCS keywords to reflect the fact that the
	; Y axis won't be perfectly correct here, since we rounded nsy to be an
	; integer.  TBD later.
	image=image_scl2

end