
function xrt_align, iindex, x_fudge=x_fudge, y_fudge=y_fudge, x_off=x_off, y_off=y_off, $
  pnt_str=pnt_str, refresh=refresh, quiet=quiet, info=info, verbose=verbose, qstop=qstop

;+
;   Name: XRT_ALIGN
;   Purpose:
;      Correct XCEN and YCEN tags of XRT index records using Shimizu-generated Hinode
;      alignment results.
;   Input Parameters:
;      IINDEX
;   Output Paramters:
;      OINDEX
;   Keyword Parameters:
;   Calling Sequence:
;      IDL> oindex = xrt_align(iindex)
;   History:
;      18-oct-2010 - G.L.Slater
;-

; Define prognam, progver variables:
  prognam = 'XRT_ALIGN.PRO'
  progver = 'V1.0' ; 2011-08-04 (GLS)

if keyword_set(info) then begin
  print, "ec_fw1_ : 'Al_poly', 'Open'"
  print, "ec_fw1  :        0 ,     1"
  print, "ec_fw2_ : 'Al_mesh', 'Gband', 'Open', 'Ti_poly'"
  print, "ec_fw2  :        0 ,     1        2          3"
  return, ''
endif

oindex = iindex

if not exist(sign_xrt_ew) then sign_xrt_ew = 1
if not exist(sign_xrt_ns) then sign_xrt_ns = 1

; Define temporary variables to contain reset xcen and ycen values:
xcen_reset = iindex.xcen
ycen_reset = iindex.ycen

; Back out existing x and y offsets (either from header or reformatter defaults):

if tag_exist(iindex, 'x_offset') then begin

  xcen_reset = xcen_reset - iindex.x_offset
  ycen_reset = ycen_reset - iindex.y_offset

endif else begin

; Define defaults:
  if not exist(x_fudge) then x_fudge = 0 ; -2
  if not exist(y_fudge) then y_fudge = 0

; Define starting times of interpolation regions for Shimizu offset data (currently unused):
  sttim_region = ['01-oct-2006','23-apr-2007','12-aug-2007','06-may-2008','06-aug-2008','01-jan-2020']

; Define existing estimated offsets introduced by reformatter:
  x_off_old_gband  = -23.94       ; G band, ec_fw2 = 3
  y_off_old_gband  = -30.05

  x_off_old_xray   = -40.0
  y_off_old_xray   = -22.3

; Define fixed XRT gband to X-ray offsets:
  x_off_gband2xray = -16.06
  y_off_gband2xray =  07.25

; Back out appropriate reformatter defaults and put in gband offsets for gband headers:
  ss_gband = where(iindex.ec_fw2 eq 3, n_gband)
  if n_gband gt 0 then begin
    xcen_reset[ss_gband] = xcen_reset[ss_gband] - x_off_old_gband + x_off_gband2xray
    ycen_reset[ss_gband] = ycen_reset[ss_gband] - y_off_old_gband + y_off_gband2xray
  endif

  ss_xray = where(iindex.ec_fw2 ne 3, n_xray)
  if n_xray gt 0 then begin
    xcen_reset[ss_xray] = xcen_reset[ss_xray] - x_off_old_xray + x_fudge
    ycen_reset[ss_xray] = ycen_reset[ss_xray] - y_off_old_xray + y_fudge
  endif

endelse

; Now add in interpolated offsets from alignment data base:
oindex.xcen = xcen_reset
oindex.ycen = ycen_reset

oindex.crval1 = xcen_reset
oindex.crval2 = ycen_reset

; Convert XRT image times to seconds and 'internal' format:
t_sec_xrt = anytim(iindex.date_obs)

; Read in Shimizu pointing data if not passed ('sc-boresight-to-SOT' and 'SOT-to-XRT' offsets):
if ( (not exist(pnt_str)) or (keyword_set(refresh)) ) then $
  rd_shimizu_data, pnt_str=pnt_str, refresh=refresh, /do_corrections

; Interpolation methods used:
; sot2sc_ns, sot2sc_ew - linear interpolation (not piecewise)
; sot2xrt_ns, sot2xrt_ew - linear interpolation (not piecewise) on both a and b, and then averaging the two

sot2sc_ns_img = interpol(pnt_str.sot2sc_ns, pnt_str.t_sec_sot2sc_ns, t_sec_xrt)
sot2sc_ew_img = interpol(pnt_str.sot2sc_ew, pnt_str.t_sec_sot2sc_ew, t_sec_xrt)

sot2xrt_ns_a_img = interpol(pnt_str.sot2xrt_ns_a, pnt_str.t_sec_sot2xrt_ns_a, t_sec_xrt)
sot2xrt_ns_b_img = interpol(pnt_str.sot2xrt_ns_b, pnt_str.t_sec_sot2xrt_ns_b, t_sec_xrt)
sot2xrt_ns_img = (sot2xrt_ns_a_img + sot2xrt_ns_b_img)/2.

sot2xrt_ew_a_img = interpol(pnt_str.sot2xrt_ew_a, pnt_str.t_sec_sot2xrt_ew_a, t_sec_xrt)
sot2xrt_ew_b_img = interpol(pnt_str.sot2xrt_ew_b, pnt_str.t_sec_sot2xrt_ew_b, t_sec_xrt)
sot2xrt_ew_img = (sot2xrt_ew_a_img + sot2xrt_ew_b_img)/2.

; Calculate x (ew) and y (ns) offsets - depending on the signs of the sot2xrt components:
x_off =   sot2sc_ew_img + sign_xrt_ew * sot2xrt_ew_img
y_off =   sot2sc_ns_img + sign_xrt_ns * sot2xrt_ns_img

ss_inf = where(finite(x_off) eq 0, n_inf)
if n_inf gt 0 then begin
  print, 'Interpolated x_off gives Inf or NaN for the following times:'
  prstr, oindex[ss_inf].date_obs
  print, 'Forcing x_off = 0 for these times.'
  x_off[ss_inf] = 0
endif

ss_inf = where(finite(y_off) eq 0, n_inf)
if n_inf gt 0 then begin
  print, 'Interpolated y_off gives Inf or NaN for the following times:'
  prstr, oindex[ss_inf].date_obs
  print, 'Forcing y_off = 0 for these times.'
  y_off[ss_inf] = 0
endif

; Make corrections to XCEN, YCEN:
oindex.xcen = oindex.xcen + x_off
oindex.ycen = oindex.ycen + y_off

; Make corrections to CRVAL1, CRVAL2:
oindex.crval1 = oindex.crval1 + x_off
oindex.crval2 = oindex.crval2 + y_off

; Modify or add tags for x and y pointing offsets:
if tag_exist(oindex, 'x_offset') then $
  oindex.x_offset = x_off else $
  oindex = add_tag(oindex, x_off, 'x_offset')
if tag_exist(oindex, 'y_offset') then $
  oindex.y_offset = y_off else $
  oindex = add_tag(oindex, y_off, 'y_offset')

; Add history record giving pointing modification done:
  history_pnt = anytim(!stime, /ccsds, /date) + ' POINTING UPDATED: ' + $
    ' x_off = ' + string(x_off, format='(f8.2)') + ' ' + $
    ' y_off = ' + string(y_off, format='(f8.2)')
  update_history, oindex, history_pnt, caller=prognam, version=progver, mode=1

if keyword_set(qstop) then stop, ' Stopping on request.'

return, oindex

end
