pro interp3d, infil_pair, index_in=index_in, data_in=data_in, index_out=index_out, data_out=data_out, $
  t_delt=t_delt, display=display, disp_xsiz=disp_xsiz, write_file=write_file, outdir=outdir, $
  verbose=verbose

;+
; NAME:
;   interp3d
; PURPOSE:
;   Interpolate intermediate images between a pair input images.
;   Optionally write them out, preserving the naming convention.
; CALLING SEQUENCE:
;   IDL> interp3d, ['file1.fts','file2.fts'], /display, /write_file, t_delt=600
;   IDL> interp3d, infil_pair, index_out=index_out, data_out=data_out, /display
;   IDL> interp3d, index_in=index_in, data_in=data_in, index_out=index_out, data_out=data_out, /display
; INPUTS:
;   infil_pair - Two element string vector of FITS file names to be interpolated
; OUTPUTS:
;   index_out  - Vector of index records for interpolated images
;   data_out   - Interpolated image cube
; OPTIONAL KEYWORD PARAMETERS:
;   t_delt     - Cadence for interpolated images (seconds, default = 10 sec)
;   display    - Set to 1 to display interpolated image cube
;   disp_xsiz  - X dimension of display version of interpolated image cube (Y dim derived)
;   write_file - Set to 1 to write interpolated FITS files (file names are consistent with
;                input file naming convention)
;   outdir     - Path specification for output files (default is current dir)
; MODIFICATION HISTORY:
;   2009-08-17 - GLS - Written

; Set defaults for optional inputs:
if not exist(t_delt) then t_delt = 10.
if not exist(disp_xsiz) then disp_xsiz = 512

; Read in FITS file pair if infil_pair supplied (otherwise index_in and data_in are assumed):
if exist (infil_pair) then begin
  mreadfits, infil_pair, index_pair, data_pair
endif else begin
  index_pair = index_in
  data_pair = data_in
endelse

infil0 = infil_pair[0]
index0 = index_pair[0]

t_pair = index_pair.date_d$obs
tsec_pair = anytim(t_pair)
tdiff_pair = tsec_pair[1] - tsec_pair[0]

if (tdiff_pair gt t_delt) then begin

; Create time grid for interpolated images and index records for those images:
  t_grid = timegrid(t_pair[0], t_pair[1], seconds=tdiff_pair)
  n_samp = fix(tdiff_pair/t_delt)
  tsec_samp = tsec_pair[0] + (indgen(n_samp)+1)*t_delt
  t_samp = anytim(tsec_samp, /ccsds)
  index_samp = replicate(index0, n_samp)
  index_samp.date_d$obs = anytim(t_samp, /ccsds)

; Do interpolation (CONGRID with linear interpolation):
  naxis1 = index0.naxis1
  naxis2 = index0.naxis2
  buff = congrid(data_pair, naxis1, naxis2, n_samp+2)
  data_samp = buff[*,*,1:n_samp]

; Optionally display movie of interpolated files interactively:
  if keyword_set(display) then begin
    if not exist(disp_xsiz) then disp_xsiz=512
    disp_ysiz = disp_xsiz*(float(naxis2)/naxis1)
    lil_data_samp = congrid(data_samp, disp_xsiz, disp_ysiz, n_samp, /interp)
    xstepper, safe_log10(lil_data_samp)
  endif

; Optionally write out interpolated FITS files:
  if keyword_set(write_file) then begin
    break_file, infil0, disk0, dir0, filnam0, ext0
    file_spec0 = strmid(filnam0, 15, 6) + '_' + strtrim(index_samp.wavelnth,2)
    outfil_samp = time2file(index_samp.date_d$obs, /seconds) + file_spec0 + ext0
    if exist(outdir) then outfil_samp = concat_dir(outdir, outfil_samp)
    mwritefits, index_samp, data_samp, outfile=outfil_samp
  endif

endif else begin

  if keyword_set(verbose) then $
    print, 'Requested sample interval exceeds ime difference between image pair.  No files created.'

endelse

end