pro new_corr, index_sot, data_sot, indexp_sot, datap_sot, file_sot=file_sot, $
  xcen_orig=xcen_orig, ycen_orig=ycen_orig, $
  sot_noscale=sot_noscale, $
  wave_aia=wave_aia, index_sdo, indexp_sdo, file_sdo=file_sdo, sdo_noscale=sdo_noscale, $
  dt_max=dt_max, use_maps=use_maps, $
  xquery=xquery, $
  xdisp_ref_arr=xdisp_ref_arr, ydisp_ref_arr=ydisp_ref_arr, $
  use_shimizu=use_shimizu, norollcorrect=norollcorrect, $
  do_prep_sot=do_prep_sot, $
  ss_axis3=ss_axis3, ss_axis4=ss_axis4, $
  do_restore_test=do_restore_test, $
  nocrop2disk=nocrop2disk, do_crop_nocrop_blink=do_crop_nocrop_blink, $
  xdisp_rot=xdisp_rot, ydisp_rot=ydisp_rot, $
  peak_str=peak_str, $
  last_good_xcen=last_good_xcen, last_good_ycen=last_good_ycen, $
  do_img_plot=do_img_plot, do_hist_plot=do_hist_plot, $
  do_plot_map=do_plot_map, do_raw_xstep=do_raw_xstep, do_corr_xstep=do_corr_xstep, $
  do_blink=do_blink, do_raw_blink=do_raw_blink, do_corr_blink=do_corr_blink, $
  do_new_corr_debug=do_new_corr_debug, do_pause=do_pause, $
  lil_reg_cube=lil_reg_cube, $
  add_waves=add_waves, cutoutp_add_wave_cube=cutoutp_add_wave_cube, $
  error_corr=error_corr, verbose=verbose, debug=debug, dir_top=dir_top, _extra=_extra

;+
; NAME:
;	sdo_corr
; PURPOSE:
;	Cross correlate input image with fov-and-time-matched SDO cutout image
;       of appropriate wavelength, to determine translational offsets.
; CATEGORY:
;	Image registration
; CALL:
;	sdo_corr, index_sot, data_sot
; OPTIONAL KEYWORDS INPUTS:
;	verbose
;	_extra
; OUTPUTS:
;	None
; OPTIONAL KEYWORDS OUTPUTS:
;	  date_obs_sot
;	  date_obs_sdo
;	  file_sdo
;	  x_off
;	  y_off
;	  wave_sot
;	  wave_sdo
;	  raster_xsize
;	  raster_ysize
;	  peak_str
; HISTORY:
;       2013-11-21 - JPW : - added logic to correct for AIA CROTA2
;                          - removed loop around interpolate
;                          - changed logic to rotate xdisp with PC matrix
;                            and to properly add it to the old CRVALi
;       2014-01-13 - JPW : - now uses an on-disk subfield for images at the limb
;                            turn off with /nocrop2disk
;                          - also fixed xcen,ycen and removed scale stuff and
;                            some other unused code
; TODO:
;       - Change instrument-specific parameter names into generic ones:
;       - Handle vectors of images
;       - Caching of reusables (fit coefficients, etc)
;       - Auto-scale SDO images using same method as for SOT
;-

t0_exec_sec = anytim(!stime)

if get_logenv('SOT_CORR_DB') eq '' then STOP
;  set_logenv, 'SOT_CORR_DB', concat_dir(concat_dir(get_logenv('SSWDB'), 'hinode/gen'), $
;                                        'hinode_sdo_alignment')
print, 'SOT_CORR_DB: ' + get_logenv('SOT_CORR_DB')
dir_top = get_logenv('SOT_CORR_DB')

if get_logenv('MINI_ARCHIVES') eq '' then STOP
;  set_logenv, 'MINI_ARCHIVES', '/archive1/hinode/sot/mini_archives'
parent_mini_archives = get_logenv('MINI_ARCHIVES')

if keyword_set(do_restore_test) then restore, './sdo_corr_20100628.sav'
dir_db = concat_dir(dir_top, 'database')
file_log_exec_times = concat_dir(dir_top, 'sot_l1_corr_exec_times.log')
if not exist(max_t_diff_sec) then max_t_diff_sec = 60
if not exist(n_smooth) then n_smooth = 5
if not exist(peak_str_min) then peak_str_min = 3

if not exist(wing_frac) then wing_frac = 0.1
if not exist(rebin) then rebin = 1
if not exist(do_smooth) then do_smooth = 1
if not exist(do_cubic) then do_cubic = 1
if not exist(rescale) then rescale = 1
if not exist(do_plots) then do_plots = 1
if not exist(do_corr) then do_corr = 1

if not exist(wait_glitch) then wait_glitch = 0.4
if not exist(wait_sec) then wait_sec = 0 ; 0.4
if not exist(n_blink) then n_blink = 5 ; 10
if not exist(wait_blink) then wait_blink = 0.5 ; 0.1 ; 0

if not exist(dt_max) then dt_max = 12 ; 60 ; seconds

t_sot = index_sot.date_obs
sat_rot = index_sot.sat_rot
sunpb0r = pb0r(index_sot.date_obs)
rsun = sunpb0r[2]*60.

; Determine which image to use in s 3+ dimensional image cube:
; (Use temporary kluge table to determine when to use first image in cube)
; TODO: Make this logic less klugey.) 

; Totally Ad Hoc Kluge 1:
mode_sot = index_sot.wave + ' ' + index_sot.obs_type
if mode_sot eq 'TF Na I 5896 FG shutterless I and V with 0.2s intervals 1' then ss_axis3 = 0

if tag_exist(index_sot, 'naxis3') then begin
  if not exist(ss_axis3) then ss_axis3 = index_sot.naxis3-1
endif
if tag_exist(index_sot, 'naxis4') then begin
  if not exist(ss_axis4) then ss_axis4 = index_sot.naxis4-1
endif

; Reduce the cube to a single 2D image:
if ( tag_exist(index_sot, 'naxis3') and (not tag_exist(index_sot, 'naxis4')) ) then $
  data_sot = reform(data_sot[*,*,ss_axis3])
if tag_exist(index_sot, 'naxis4') then $
  data_sot = reform(data_sot[*,*,ss_axis3,ss_axis4])

if mode_sot eq 'TF Na I 5896 FG shutterless I and V with 0.2s intervals 1' then begin
  ss_neg = where(data_sot lt 0, n_neg)
  data_sot[ss_neg] = -data_sot[ss_neg]
endif

; Optionally prep SOT image:
if keyword_set(do_prep_sot) then begin
  fg_prep, index_sot, data_sot, indexp_sot, datap_sot, /no_flat, /no_dark
endif else begin
  indexp_sot = index_sot
  datap_sot = data_sot
endelse

if keyword_set(use_shimizu) then begin
  STOP, 'NEW_CORR: Stopping on request at call to GET_SHIMIZU.'
; In case of no fg_prep, then update pointing keywords.
; For now, use Shimizu table values, eventually to be replaced by
;   continuously and iteratively improved table from the correlation results:
  sot_offsets = get_shimizu(t_sot)
PRINT, SOT_OFFSETS
  crval1_offset = sot_offsets[0]
  crval2_offset = sot_offsets[1]
  indexp_sot.crval1 = indexp_sot.crval1 + crval1_offset
  indexp_sot.crval2 = indexp_sot.crval2 + crval2_offset
  indexp_sot.xcen = indexp_sot.crval1
  indexp_sot.ycen = indexp_sot.crval2
endif

; Add CROTA2 tag to IRIS header if not present:
;if tag_exist(indexp_sot, 'crota2') then indexp_sot.crota2 = sat_rot 

; Update rotation angle from fit to measured data:
;if not keyword_set(norollcorrect) then begin
; restore, '/net/topaz/Users/slater/soft/idl/idl_startup/sot/coeff_fit_sot_delta_rotang.sav'
  restore, concat_dir(dir_db, 'coeff_fit_sot_delta_rotang.sav')
; indexp_sot.crota1 = indexp_sot.crota1 + (10*sine_func(anytim(t_sot), coeff_fit_da))
; indexp_sot.crota1 = indexp_sot.crota1 + (05*sine_func(anytim(t_sot), coeff_fit_da))
  indexp_sot.crota1 = indexp_sot.crota1 - 0.8
  indexp_sot.crota2 = indexp_sot.crota1
; indexp_sot.crota1 = indexp_sot.crota1 + 5
; indexp_sot.crota2 = indexp_sot.crota2 + 5
;endif

; Update platescale from fit to measured data:
;restore, '/net/topaz/Users/slater/soft/idl/idl_startup/sot/coeff_fit_sot_delta_platescale.sav'
restore, concat_dir(dir_db, 'coeff_fit_sot_delta_platescale.sav')
;indexp_sot.cdelt1 = indexp_sot.cdelt1 / (1.005*poly_func(anytim(t_sot), reform(coeff_fit_ds)))
;indexp_sot.cdelt2 = indexp_sot.cdelt1
;indexp_sot.cdelt1 = indexp_sot.cdelt1 / 1.02
 indexp_sot.cdelt1 = indexp_sot.cdelt1 / 1.015
 indexp_sot.cdelt2 = indexp_sot.cdelt1

; Generate regulation WCS structure from SSW index structure for SOT image:
wcsp_sot = fitshead2wcs(indexp_sot)

; Define size and FOV of SOT image (CCD readout region):
tsc1 = 1
tec1 = indexp_sot.naxis1
tsr1 = 1
ter1 = indexp_sot.naxis2
iNx = tec1-tsc1+1
iNy = ter1-tsr1+1

; Define arrays to contain x and y Helioprojective (arcsec) coordinates of all image pixels:
solx = fltarr(iNx,iNy)
soly = fltarr(iNx,iNy)

for i=0,iNx-1 do begin
	for j=0,iNy-1 do begin
		solx[i,j] = (tsc1 + i - indexp_sot.crpix1)*wcsp_sot.pc[0,0] + $
					(tsr1 + j - indexp_sot.crpix2)*wcsp_sot.pc[0,1]
		soly[i,j] = (tsc1 + i - indexp_sot.crpix1)*wcsp_sot.pc[1,0] + $
					(tsr1 + j - indexp_sot.crpix2)*wcsp_sot.pc[1,1]
	endfor
endfor

solx *= indexp_sot.cdelt1
soly *= indexp_sot.cdelt2
solx += indexp_sot.crval1
soly += indexp_sot.crval2
;solx += indexp_sot.xcen
;soly += indexp_sot.ycen

; Optionally check if limb crosses IRIS cutout field.  If so, crop to on-disk portion:
crop_flag = 0
;if ~keyword_set(nocrop2disk) then begin
if exist(hell_freezes_over) then begin
; fovmar =  5.0
; fovmin = 20.0
  fovmar =  1.0
  fovmin = 10.0
  ; 5 arcsec margin:
  mx_solr = rsun - fovmar
  solr = sqrt(solx*solx+soly*soly)
  ; corners of IRIS image :
  solr_corn = solr[[[0,iNx-1],[0,iNx-1]],[[0,0],[iNy-1,iNy-1]]]
  if min(solr_corn) lt mx_solr and max(solr_corn) gt mx_solr then begin
    print, 'The SOT FOV crosses the limb. Attempting to crop it to larges on disk rectangle.'
    ; corner closest to suncenter "cctsc" :
    wcorn = where(solr_corn eq min(solr_corn))
    cornx = wcorn[0] mod 2
    corny = wcorn[0] / 2
    ; find points where solar limb (minus margin) crosses edges of image
    ; edges adjacent to cctsc (only need nxe1,nye1 = # of edge pixels inside limb)
    dmy = where(solr[*,corny*(iNy-1)] lt mx_solr,nxe1)
    dmy = where(solr[cornx*(iNx-1),*] lt mx_solr,nye1)
    ; edges not adjacent to cctsc
    dmy = where(solr[*,(1-corny)*(iNy-1)] lt mx_solr,nxe2)
    dmy = where(solr[(1-cornx)*(iNx-1),*] lt mx_solr,nye2)
    ; Use limb crossing pts to determine optimum subfield size.
    ; Assuming straight limb for this calculation, which provides
    ; a somewhat conservative result (draw 4 cases to derive eq's)
    case 1 of                                  ; 4 cases
      (nxe1 eq inx) and (nye1 eq iny): begin   ; crossing only non-adjacent edges
        inx1 = nxe2 > fix((float(inx)+float(nye2)*float(inx-nxe2)/float(iny-nye2))/2.0) < inx
        iny1 = nye2 > fix((float(iny)+float(nxe2)*float(iny-nye2)/float(inx-nxe2))/2.0) < iny
        end
      (nxe1 eq inx) and (nye1 lt iny): begin   ; crossing both y edges
        inx1 = fix(float(nye1)*float(inx)/float(2*(nye1-nye2) > nye1))
        iny1 = nye1/2 > nye2
        end
      (nxe1 lt inx) and (nye1 eq iny): begin   ; crossing both x edges
        inx1 = nxe1/2 > nxe2
        iny1 = fix(float(nxe1)*float(iny)/float(2*(nxe1-nxe2) > nxe1))
        end
      else: begin                              ; crossing only adjacent edges
        inx1 = nxe1/2
        iny1 = nye1/2
        end
    endcase
    ; update IRIS image cutout if new subimage large enough (i.e., at least 20 " on the side)
    iFOVx1 = inx1*indexp_sot.cdelt1
    iFOVy1 = iny1*indexp_sot.cdelt2
    if iFOVx1 ge fovmin and iFOVy1 ge fovmin then begin
      otsc1 = cornx * (iNx-inx1)
      otsr1 = corny * (iNy-iny1)
      tsc1 = tsc1+otsc1
      tec1 = tsc1+inx1-1
      tsr1 = tsr1+otsr1
      ter1 = tsr1+iny1-1
      iNx = inx1
      iNy = iny1
      solx = solx[otsc1:otsc1+inx1-1,otsr1:otsr1+iny1-1]
      soly = soly[otsc1:otsc1+inx1-1,otsr1:otsr1+iny1-1]
      crop_flag = 1
    endif
  endif
endif

; Define SOT cutout:
;cutoutp_sot = datap_sot[tsc1-1:tec1-1,tsr1-1:ter1-1] / indexp_sot.exptime
 cutoutp_sot = datap_sot[tsc1-1:tec1-1,tsr1-1:ter1-1]

; Optionally blink movie the uncropped and cropped FOV's :
if ( (crop_flag eq 1) and keyword_set(do_crop_nocrop_blink) ) then begin
  nocrop_img = datap_sot
    crop_img = datap_sot*0
  crop_img[tsc1-1,tsr1-1] = cutoutp_sot
  xstepper, [[[nocrop_img]],[[crop_img]]]
STOP
endif

; Handle NaN's:
ss_fin = where(finite(cutoutp_sot) eq 1, n_fin, comp=ss_inf, ncomp=n_inf)
if n_inf gt 0 then cutoutp_sot[ss_inf] = min(cutoutp_sot[ss_fin])

; Find SDO wave to match, as well as scaling information for both SOT and SDO images:
if not exist(wave_aia) then $
  wave_aia = aia_wave_match(indexp_sot, mode=mode, ss_axis3=ss_axis3, ss_axis4=ss_axis4, $
                            loval_sot=loval_sot, hival_sot=hival_sot, $
                            do_invert_sot=do_invert_sot, nsig_sot=nsig_sot_tab, $
                            loval_sdo=loval_sdo, hival_sdo=hival_sdo, $
                            do_invert_sdo=do_invert_sdo, nsig_sdo=nsig_sdo_tab, $
                            _extra=_extra)

; Return if no match found:
if wave_aia eq '' then begin
  print, 'No AIA wave match found for SOT wave.  Returning.'
; stop,  'No AIA wave match found for SOT wave.  Stopping.'
  peak_str = 0
  error_corr = 1
  return
endif

; Optionally Clip and scale SOT image:
;if not keyword_set(sot_noscale) then begin
if exist(hell_freezes_over) then begin

  if not exist(nsig_sot) then nsig_sot = nsig_sot_tab

  if ( (loval_sot eq 0) and (hival_sot eq 0) ) then begin
    print, 'No limit values in database.'
    cutoutp_sot_scaled = bytscl(cutoutp_sot)
  endif else begin
    if (do_invert_sot or keyword_set(force_invert_sot)) then $
      cutoutp_sot_scaled = 255b-bytscl(loval_sot>cutoutp_sot<hival_sot) else $
      cutoutp_sot_scaled =      bytscl(loval_sot>cutoutp_sot<hival_sot)
  endelse

  if keyword_set(do_img_plot) then begin
    wdef,0,indexp_sot.naxis1,indexp_sot.naxis2
    tvscl, cutoutp_sot_scaled
  endif

  hist = histogram(cutoutp_sot, loc=xloc)
  ss_nz = where(hist ne 0, n_nz)

  if n_nz gt 0 then begin
    hist_nz = hist[ss_nz]
    xloc_nz = xloc[ss_nz]
; For now, at least, set hist, xloc to nz versions:
    hist = hist_nz
    xloc = xloc_nz

    if n_elements(hist) gt n_smooth then begin
      hist_integral_smooth = smooth(total(hist, /cum), n_smooth) / total(hist)
    endif else begin
      stop, 'Problem with SOT image histogram scaling. Stopping.'
      peak_str = 0
      error_corr = 1
      return
    endelse

;   hist_integral_smooth = total(hist, /cum) / total(hist)
;   deriv_hist_integral_smooth =  hist_integral_smooth[1:*] - hist_integral_smooth
    deriv_hist_integral_smooth = [hist_integral_smooth[1:*] - hist_integral_smooth,0]
; Try to fit the dominant guassian to automatically determine lo and hi cutoff values:
    n_pts = n_elements(hist)
    n0_crop = fix(0.1*n_pts)
    x_crop = xloc[n0_crop:*]
;   hist_crop = hist[n0_crop:*]
    hist_crop = deriv_hist_integral_smooth[n0_crop:*]
    y_crop = hist_crop
    yfit_crop = gaussfit(x_crop, y_crop, afit, chisq=chisqr, nterms=3, sigma=sigma, yerror=yerror)
    loval_sot_auto = afit[1] - nsig_sot*afit[2]
    hival_sot_auto = afit[1] + nsig_sot*afit[2]

    if (do_invert_sot or keyword_set(force_invert_sot)) then $
      cutoutp_sot_scaled_auto = 255b-bytscl(loval_sot_auto>cutoutp_sot<hival_sot_auto) else $
      cutoutp_sot_scaled_auto =      bytscl(loval_sot_auto>cutoutp_sot<hival_sot_auto)

    if keyword_set(do_hist_plot) then begin
      wdef,1,1024,768
      !p.multi = [0,1,3]
      plot, xloc, hist, charsize=2
      plot, xloc, hist_integral_smooth, charsize=2
      plot, xloc, deriv_hist_integral_smooth, charsize=2

      oplot, x_crop, yfit_crop, psym=2

      oplot, [loval_sot_auto, loval_sot_auto], [0, max(yfit_crop)], linestyle=4, thick=2
      oplot, [hival_sot_auto, hival_sot_auto], [0, max(yfit_crop)], linestyle=4, thick=2
    endif

  endif else begin
    stop, ' No non-zero values of histogram. Stopping.'
  endelse

  if keyword_set(do_img_plot) then begin
    wset,0
    tvscl, cutoutp_sot_scaled_auto
    if keyword_set(do_pause) then stop, 'NEW_CORR: Stopping on request after image display.'
  endif

endif else begin
  cutoutp_sot_scaled = cutoutp_sot
endelse


; Get an appropriate SDO image:

if exist(index_sdo) then STOP, 'index_sdo has been passed in. Stopping.'

; Construct xquery argument for ssw_jsoc_time2data call:

xquery = ['quality = 0']

;if not exist(xquery) then xquery = ['quality = 0']
;if not exist(xquery) then xquery = ['quality != 1075998720', 'quality != 1075990528']
;if not exist(xquery) then xquery = ['quality != 397312']
;if not exist(xquery) then xquery = ['quality != 397312', 'aiftsid>49152']
;if not exist(xquery) then xquery = [ $
;  'AIFTSID  >  49151', $ ; Cal mode image
;  'AIFCPS   >     80', $ ; Focus not nominal
;  'AIFCPS   <      0', $ ; Focus not nominal
;  'AIFWEN   =     74', $ ; Filterwheel in open position
;  'AIFWEN   =     75', $ ; Filterwheel in open position
;  'AIMGSHCE <    100', $ ; Very short exposure
;  'AIMGSHCE >   7000', $ ; Very long exposure (likely cal mode image)
;  'IMG_TYPE = 'LIGHT', $

; ********
TRY_AGAIN:
; ********

expand_interval = 1
n_tries_max = 04
n_tries = 0l
n_iter = 0l
re_read = 0
t_buff_sec = 25 ; seconds (normal UV wave sample time + 1)
index_sdo = -1
t0_exec_jsoc_sec = anytim(!stime)

while ( (expand_interval eq 1) and (t_buff_sec le 64024l) ) do begin

  delvarx, file_sdo0, file_sdo, index_sdo0, index_sdo, data_sdo, indexp_sdo, datap_sdo, $
           mapp_sdo, rmapp_sdo

  if ( (not exist(file_sdo)) or (size(index_sdo, /type) ne 8) ) then begin
    t_query = anytim(!stime, /ccsds)
    ssw_jsoc_time2data, anytim(anytim(t_sot)-t_buff_sec, /ccsds), $
                        anytim(anytim(t_sot)+t_buff_sec, /ccsds), $
    wave=wave_aia, index_sdo, file_sdo, /files_only, xquery=xquery
;   wave=171, index_sdo, file_sdo, /files_only, xquery=xquery
;   ssw_jsoc_time2data, t_sot, t_sot, wave=wave_aia, index_sdo, file_sdo, /files_only
    t_post_query = anytim(!stime, /ccsds)
    buff = t_query + '  ' + t_post_query
    if exist(index_sdo) then begin
      if size(index_sdo, /type) eq 8 then begin
        if tag_exist(index_sdo, 'runtime') then begin
          jsoc_runtime = index_sdo[0].runtime
          buff = buff + '  JSOC query execution time: ' + $
                 string(jsoc_runtime, '$(f9.4)') + 'seconds'
	  buff = buff + '  re_read: ' + string(re_read, '$(i1.1)') + $
			'  expand_interval: ' + string(expand_interval, '$(i1.1)') + $
			'  n_tries: ' + string(n_tries, '$(i2.2)') + $
			'  n_n_iter: ' + string(n_iter, '$(i2.2)')
        endif
      endif
    endif
    file_append, '/sanhome/slater/logs/jsoc_query_log', buff
  endif  
   
  re_read = 0
; If no valid file found, try one more time, in case there was a system glitch:
  if ( (not exist(file_sdo)) or (size(index_sdo, /type) ne 8) ) then begin
    re_read = 1
  endif else begin
    if ( (file_sdo[0] eq '') or (file_sdo[0] eq 'Remote listing failed.') ) then $
      re_read = 1
  endelse
  
  if re_read eq 1 then begin

    delvarx, file_sdo0, file_sdo, index_sdo0, index_sdo, data_sdo, indexp_sdo, datap_sdo, $
             mapp_sdo, rmapp_sdo

    wait, wait_glitch
    t_query = anytim(!stime, /ccsds)
    ssw_jsoc_time2data, anytim(anytim(t_sot)-t_buff_sec, /ccsds), $
                        anytim(anytim(t_sot)+t_buff_sec, /ccsds), $
    wave=wave_aia, index_sdo, file_sdo, /files_only, xquery=xquery
;   ssw_jsoc_time2data, t_sot, t_sot, wave=wave_aia, index_sdo, file_sdo, /files_only
    t_post_query = anytim(!stime, /ccsds)
    buff = t_query + '  ' + t_post_query
    if exist(index_sdo) then begin
      if size(index_sdo, /type) eq 8 then begin
        if tag_exist(index_sdo, 'runtime') then begin
          jsoc_runtime = index_sdo[0].runtime
          buff = buff + '  JSOC query execution time: ' + $
                 string(jsoc_runtime, '$(f9.4)') + 'seconds'
          buff = buff + '  re_read: ' + string(re_read, '$(i1.1)') + $
                        '  expand_interval: ' + string(expand_interval, '$(i1.1)') + $
                        '  n_tries: ' + string(n_tries, '$(i2.2)') + $
                        '  n_n_iter: ' + string(n_iter, '$(i2.2)')
        endif
      endif
    endif
    file_append, '/sanhome/slater/logs/jsoc_query_log', buff

  endif

HELP, INDEX_SDO,EXIST(FILE_SDO),N_ITER,T_BUFF_SEC,RE_READ
;STOP
; QUESTION FOR SAM - WHY IS INDEX_SDO = -1 WHEN FILE_SDO EXISTS?

  expand_interval = 0
  if ( (not exist(file_sdo)) or (size(index_sdo, /type) ne 8) ) then begin
    expand_interval = 1
  endif else begin
    if ( (file_sdo[0] eq '') or (file_sdo[0] eq 'Remote listing failed.') ) then $
      expand_interval = 1
  endelse

  if expand_interval eq 1 then begin
    n_iter = n_iter+1 
    t_buff_sec = 2l^(8+n_iter) ; 512, 1024, 2048, 4096, 8192 seconds
  endif

endwhile
t1_exec_jsoc_sec = anytim(!stime)

give_up = 0
if ( (not exist(file_sdo)) or (size(index_sdo, /type) ne 8) ) then begin
  give_up = 1
endif else begin
  if ( (file_sdo[0] eq '') or (file_sdo[0] eq 'Remote listing failed.') ) then $
    give_up = 1
endelse

if give_up eq 1 then begin
  print, 'No usable close in time SDO images found.  Giving up and returning.'
  if keyword_set(do_new_corr_debug) then $
    STOP, 'NEW_CORR: No usable close in time SDO images found.  Giving up.'
  peak_str = 0
  error_corr = 1
  return
endif
  
; If one or more good SDO files have been found, then attempt to use the closest
;   in time:

ss_exist_arr = file_exist(file_sdo)
ss_exist = where(ss_exist_arr ne -1, n_exist)

if n_exist eq 0 then begin
  print, 'None of the files returned by SSW_JSOC_TIME2DATA seem to exist. Returning.'
  if keyword_set(do_new_corr_debug) then $
    STOP, 'NEW_CORR: None of the files returned by SSW_JSOC_TIME2DATA seem to exist. Stopping.'
  peak_str = 0
  error_corr = 1
  return
endif else begin
  file_sdo = file_sdo[ss_exist]
  index_sdo = index_sdo[ss_exist]
endelse

min_t_off_sec_sdo = min(abs(anytim(index_sdo.t_obs) - $
                            anytim(indexp_sot.date_obs)), ss_min_t_off)

; Temporary KLUGE for as yet not understaood error resulting from ss_min_t_off equal to 0:
;if ss_min_t_off eq 0 then begin
;  error_message = 'NEW_CORR: Problem finding a suitable SDO file. ss_min_t_off equal to 0. Skipping this SOT image.'
;  print, 'Problem finding a suitable SDO file. Returning.'
;  if keyword_set(do_new_corr_debug) then $
;    STOP, 'NEW_CORR: Problem finding a suitable SDO file. ss_min_t_off equal to 0. Stopping.'
;  goto, failure
;  return
;endif

file_sdo0 = file_sdo[ss_min_t_off]
index_sdo0 = index_sdo[ss_min_t_off]
help, ss_min_t_off, file_sdo, file_sdo0, index_sdo, index_sdo0
help, index_sdo0.date_obs
print, index_sdo0.date_obs

; *********************************************************
if strlen(index_sdo0.date_obs) gt 23 then begin
   n_tries = n_tries + 1
   if n_tries le n_tries_max then begin
      print, 'FAILED on N_TRIES = ' + strtrim(n_tries,2) + ' due to munged header (DATE_OBS).  Trying again.'
      GOTO, TRY_AGAIN
   endif else begin
      print, 'Quitting this SOT image (' + indexp_sot.date_obs + ') after ' + strtrim(n_tries,2) + ' tries.'
      GOTO, failure
   endelse
endif
; *********************************************************

; Calculate cutout coordinates for read_sdo:
x_pix_sot_llur_fov = [0, indexp_sot.naxis1-1]
y_pix_sot_llur_fov = [0, indexp_sot.naxis2-1]
pix_sot_llur_fov = transpose([[x_pix_sot_llur_fov], [y_pix_sot_llur_fov]])

help, pix_sot_llur_fov, indexp_sot, ss_min_t_off, index_sdo0.date_obs
pix_sdo_llur_fov = wcs_pix2pix_2(pix_sot_llur_fov, indexp_sot, index_sdo0)

;xllp_sdo =      min([pix_sdo_llur_fov[0,0], pix_sdo_llur_fov[0,1]])
;yllp_sdo =      min([pix_sdo_llur_fov[1,0], pix_sdo_llur_fov[1,1]])
xllp_sdo = round(min([pix_sdo_llur_fov[0,0], pix_sdo_llur_fov[0,1]]))
yllp_sdo = round(min([pix_sdo_llur_fov[1,0], pix_sdo_llur_fov[1,1]]))
xurp_sdo = max([pix_sdo_llur_fov[0,0], pix_sdo_llur_fov[0,1]])
yurp_sdo = max([pix_sdo_llur_fov[1,0], pix_sdo_llur_fov[1,1]])
;nxp_sdo = ceil(xurp_sdo) - floor(xllp_sdo) + 1
;nyp_sdo = ceil(yurp_sdo) - floor(yllp_sdo) + 1
nxp_sdo = ceil(indexp_sot.naxis1*(indexp_sot.cdelt1 / index_sdo0.cdelt1) ) + 1
nyp_sdo = ceil(indexp_sot.naxis2*(indexp_sot.cdelt2 / index_sdo0.cdelt2) ) + 1

if ( (wave_aia eq 'blos') or (wave_aia eq 'cont') ) then begin
  uncomp_delete = 1
  use_index = 1
  use_shared = 0 ; NOTA BENA - 'use_shared' flag apparently incompatible with 'use_index' flag.
  noshell = 1
endif else begin
  uncomp_delete = 1
  use_index = 0 ;1
  use_shared = 0
  noshell = 1
endelse

t0_exec_read_sdo_sec = anytim(!stime)

if file_exist(file_sdo0) eq 0 then begin
  print, 'SDO file does not exist. Returning.'
; stop,  'SDO file does not exist. Stopping.' 
  peak_str = 0
  error_corr = 1
  return
endif

if ( (xllp_sdo le 0) or (yllp_sdo le 0) or $
     ((xllp_sdo + nxp_sdo) ge 4096)     or $
     ((yllp_sdo + nyp_sdo) ge 4096) )   then begin
  help, xllp_sdo, yllp_sdo, nxp_sdo, nyp_sdo
  print, 'Attempt to extract roi outside of SDO image.  Returning.'
  if keyword_set(do_new_corr_debug) then $
    STOP, 'NEW_CORR: Attempt to extract roi outside of SDO image.  Stopping.'
PRINT, 'TEST PRINT'
  peak_str = 0
  error_corr = 1
  return
endif

read_sdo, file_sdo0, index_sdo0, data_sdo, $
          xllp_sdo, yllp_sdo, nxp_sdo, nyp_sdo, $
          uncomp_delete=uncomp_delete, use_index=use_index, use_shared=use_shared, noshell=noshell

; If read_sdo failed, try once more:
if not exist(index_sdo) then begin
  print, 'READ_SDO failed.  Waiting 100 seconds and trying again.'
  wait, wait_glitch

  read_sdo, file_sdo0, index_sdo0, data_sdo, $
            xllp_sdo, yllp_sdo, nxp_sdo, nyp_sdo, $
            uncomp_delete=uncomp_delete, use_index=use_index, use_shared=use_shared, noshell=noshell
endif
t1_exec_read_sdo_sec = anytim(!stime)

if ( (not exist(index_sdo)) or (not exist(data_sdo)) ) then begin
  print, 'Either INDEX_SDO or DATA_SDO does not exist. Returning.'
  if keyword_set(do_new_corr_debug) then STOP, ' Either INDEX_SDO or DATA_SDO does not exist.  Stopping.' 
  peak_str = 0
  error_corr = 1
  return
endif

if not tag_exist(index_sdo, 'crota2') then begin
  print, " SDO image header has no 'CROTA2' tag. Returning."
  if keyword_set(do_new_corr_debug) then STOP, " SDO image header has no 'CROTA2' tag.  Stopping."
  peak_str = 0
  error_corr = 1
  return
endif

file_sdo  = file_sdo0
index_sdo = index_sdo0

; Optionally prep SDO image:
;do_prep_sdo=1
if keyword_set(do_prep_sdo) then begin
  aia_prep, index_sdo, data_sdo, indexp_sdo, datap_sdo
endif else begin
  indexp_sdo = index_sdo
  datap_sdo = data_sdo
endelse

print, 'De-rotate SDO image.'
index2map, indexp_sot, datap_sot, mapp_sot
index2map, indexp_sdo, datap_sdo, mapp_sdo

; KLUGE - Attempt to correct for apparent indexing error in READ_SDO coordinate reference:
;mapp_sdo.xc = mapp_sdo.xc - indexp_sdo.cdelt1
;mapp_sdo.yc = mapp_sdo.yc - indexp_sdo.cdelt2

dt_sdo = anytim(indexp_sdo.date_obs) - anytim(indexp_sot.date_obs)
if abs(dt_sdo) le dt_max then no_drotate = 1 else no_drotate = 0

file_log = './sot_corr_info.log'
print, 'indexp_sot.date_obs = ' + indexp_sot.date_obs
print, 'indexp_sdo.date_obs = ' + indexp_sdo.date_obs
print, 'DT_SDO = ' + strtrim(dt_sdo,2)
print, 'DT_MAX = ' + strtrim(dt_max,2)
print, 'NO_DROTATE = ' + strtrim(no_drotate, 2)
file_append, file_log, 'indexp_sot.date_obs = ' + indexp_sot.date_obs
file_append, file_log, 'indexp_sdo.date_obs = ' + indexp_sdo.date_obs
file_append, file_log, 'DT_SDO = ' + strtrim(dt_sdo,2)
file_append, file_log, 'DT_MAX = ' + strtrim(dt_max,2)
file_append, file_log, 'NO_DROTATE = ' + strtrim(no_drotate, 2)
;wait, 2
;rmapp_sdo = drot_map(mapp_sdo, ref_map=mapp_sot, no_drotate=0, $
rmapp_sdo = drot_map(mapp_sdo, ref_map=mapp_sot, no_drotate=no_drotate, /keep_limb, $
                     rigid=rigid, no_project=no_project, _extra=_extra)
cutoutp_sdo = rmapp_sdo.data
if crop_flag eq 1 then cutoutp_sdo = cutoutp_sdo[tsc1-1:tec1-1,tsr1-1:ter1-1]

t_sot = indexp_sot.date_obs
t_sdo = indexp_sdo.date_obs
t_delt_sec = anytim(t_sdo) - anytim(t_sot)

; Optionally Clip and scale SDO image:
;if not keyword_set(sdo_noscale) then begin
if exist(hell_freezes_over) then begin

  if not exist(nsig_sdo) then nsig_sdo = nsig_sdo_tab

  if ( (loval_sdo eq 0) and (hival_sdo eq 0) ) then begin
    print, 'No limit values for SDO image in database.'
    cutoutp_sdo_scaled = bytscl(cutoutp_sdo)
  endif else begin
    if (do_invert_sdo or keyword_set(force_invert_sdo)) then $
      cutoutp_sdo_scaled = 255b-bytscl(loval_sdo>cutoutp_sdo<hival_sdo) else $
      cutoutp_sdo_scaled =      bytscl(loval_sdo>cutoutp_sdo<hival_sdo)
  endelse

  if keyword_set(do_img_plot) then begin
    wset,0
    tvscl, cutoutp_sdo_scaled
  endif

  hist = histogram(cutoutp_sdo, loc=xloc)
  ss_nz = where(hist ne 0, n_nz)

  if n_nz gt 0 then begin
    hist_nz = hist[ss_nz]
    xloc_nz = xloc[ss_nz]
; For now, at least, set hist, xloc to nz versions:
    hist = hist_nz
    xloc = xloc_nz

    if n_elements(hist) gt n_smooth then begin
      hist_integral_smooth = smooth(total(hist, /cum), n_smooth) / total(hist)
    endif else begin
      stop, 'Problem with SDO image histogram scaling. Stopping.'
      peak_str = 0
      error_corr = 1
      return
    endelse

;   hist_integral_smooth = total(hist, /cum) / total(hist)
;   deriv_hist_integral_smooth =  hist_integral_smooth[1:*] - hist_integral_smooth
    deriv_hist_integral_smooth = [hist_integral_smooth[1:*] - hist_integral_smooth,0]

; Try to fit the dominant guassian to automatically determine lo and hi cutoff values:
    n_pts = n_elements(hist)
    n0_crop = fix(0.00*n_pts)
    x_crop = xloc[n0_crop:*]
;   hist_crop = hist[n0_crop:*]
    hist_crop = deriv_hist_integral_smooth[n0_crop:*]
    y_crop = hist_crop
    yfit_crop = gaussfit(x_crop, y_crop, afit, chisq=chisqr, nterms=3, sigma=sigma, yerror=yerror)

    loval_sdo_auto = afit[1] - nsig_sdo*afit[2]
    hival_sdo_auto = afit[1] + nsig_sdo*afit[2]

    if (do_invert_sdo or keyword_set(force_invert_sdo)) then $
      cutoutp_sdo_scaled_auto = 255b-bytscl(loval_sdo_auto>cutoutp_sdo<hival_sdo_auto) else $
      cutoutp_sdo_scaled_auto =      bytscl(loval_sdo_auto>cutoutp_sdo<hival_sdo_auto)

    if keyword_set(do_hist_plot) then begin
      wdef,1,1024,768
      !p.multi = [0,1,3]
      plot, xloc, hist, charsize=2
      plot, xloc, hist_integral_smooth, charsize=2
      plot, xloc, deriv_hist_integral_smooth, charsize=2

      oplot, x_crop, yfit_crop, psym=2

      oplot, [loval_sdo_auto, loval_sdo_auto], [0, max(yfit_crop)], linestyle=4, thick=2
      oplot, [hival_sdo_auto, hival_sdo_auto], [0, max(yfit_crop)], linestyle=4, thick=2
    endif

  endif else begin
    stop, ' No non-zero values of histogram. Stopping.'
  endelse

endif else begin
  cutoutp_sdo_scaled = cutoutp_sdo
endelse

; Create an index structure for the SDO cutout:
map2index, rmapp_sdo, index_sdo_cutout, data_sdo_cutout

; Create SDO/SOT cutout cube and cross correlate to get offsets:
cutout_cube = [[[cutoutp_sdo_scaled]],[[cutoutp_sot_scaled]]]
cutout_cube_unshifted = cutout_cube
if keyword_set(do_new_corr_debug) then STOP, 'NEW_CORR: Stopping on request prior to TR_GET_DISP.'

if keyword_set(do_raw_xstep) then xstepper, cutout_cube_unshifted
if keyword_set(do_raw_blink) then begin
  xs_blink = n_elements(cutout_cube_unshifted[*,0,0])
  ys_blink = n_elements(cutout_cube_unshifted[0,*,0])
  wdef,0,xs_blink,ys_blink
  for ii=0, n_blink-1 do begin
    tvscl, cutout_cube_unshifted[*,*,0] & wait, wait_blink
    tvscl, cutout_cube_unshifted[*,*,1] & wait, wait_blink
  endfor
endif

; ************************************************************************************************
disp = tr_get_disp_2d_corr(cutout_cube, /shift, do_surf=do_surf, peak_str=peak_str, _extra=_extra)
; ************************************************************************************************

; Save reduced-size versions of cutouts:
t_cutout = time2file(indexp_sot.date_obs, /sec)
file_cutout = concat_dir('/sanhome/slater/missions/sot/corr/cutout_db_synop', $
                         'sot_sdo_cutout_pair_' + t_cutout + '.sav')
cutout_cube_save = [[[cutout_cube_unshifted[*,*,0]]],[[cutout_cube_unshifted[*,*,1]]], $
                    [[cutout_cube [*,*,1]]]]
size_cube = size(cutout_cube_save)
xsiz = size_cube[1] & ysiz = size_cube[2]
cutout_cube_save = congrid(cutout_cube_save, xsiz/4., ysiz/4., 3)
save, cutout_cube_save, file=file_cutout

if keyword_set(do_corr_xstep) then xstepper, cutout_cube
if keyword_set(do_corr_blink) then begin
  xs_blink = n_elements(cutout_cube[*,0,0])
  ys_blink = n_elements(cutout_cube[0,*,0])
  wdef,0,xs_blink,ys_blink
  for ii=0, n_blink-1 do begin
    tvscl, cutout_cube[*,*,0] & wait, wait_blink
    tvscl, cutout_cube[*,*,1] & wait, wait_blink
  endfor
  wait, 0.3
endif

print, ''
print, ' PEAK_STR = ' + string(peak_str, '(f5.1)')
print, ''
;xstepper, [[cutout_cube_unshifted],[cutout_cube]]

if ~ finite(peak_str) then begin
  print, 'PEAK_STR is not finite. Returning.'
  peak_str = 0
  error_corr = 1
  return
endif

xdisp = disp[0,1]*indexp_sot.cdelt1
ydisp = disp[1,1]*indexp_sot.cdelt2

; Rotate the X and Y offsets from the correlation to the roll angle of the SOT image:
if keyword_set(use_maps) then begin
  xdisp_rot = xdisp*wcsp_sot.pc[0,0] + ydisp*wcsp_sot.pc[0,1]
  ydisp_rot = xdisp*wcsp_sot.pc[1,0] + ydisp*wcsp_sot.pc[1,1]
endif else begin
  xdisp_rot = xdisp
  ydisp_rot = ydisp
endelse

xcen_sot_cor = indexp_sot.xcen + xdisp_rot ;+ total(xdisp_ref_arr)
ycen_sot_cor = indexp_sot.ycen + ydisp_rot ;+ total(ydisp_ref_arr)
;xcen_sot_cor = indexp_sot.xcen + total(xdisp_ref_arr)
;ycen_sot_cor = indexp_sot.ycen + total(ydisp_ref_arr)

xdisp_tot = xdisp_rot + total(xdisp_ref_arr)
ydisp_tot = ydisp_rot + total(ydisp_ref_arr)

last_good_wcsp_sot = wcsp_sot
last_good_xcen_sot = xcen_sot_cor
last_good_ycen_sot = ycen_sot_cor

crval1_sot_cor = $
  comp_fits_crval(xcen_sot_cor, indexp_sot.cdelt1, indexp_sot.naxis1, indexp_sot.crpix1)
crval2_sot_cor = $
  comp_fits_crval(ycen_sot_cor, indexp_sot.cdelt2, indexp_sot.naxis2, indexp_sot.crpix2)

crpix1_sot_cor = $
  comp_fits_crpix(xcen_sot_cor, indexp_sot.cdelt1, indexp_sot.naxis1, crval1_sot_cor)
crpix2_sot_cor = $
  comp_fits_crpix(ycen_sot_cor, indexp_sot.cdelt2, indexp_sot.naxis2, crval2_sot_cor)

; Backup the original indexp_sot before updating pointing, etc tags:
indexp_sot_orig = indexp_sot

; Update SOT X,YCEN, CRVAL1,2, CRPIX1,2 tags (add X,YCEN tags if they don't exist):
if not tag_exist(indexp_sot, 'xcen') then $
  indexp_sot = add_tag(indexp_sot, xcen_sot_cor, 'xcen') else $
  indexp_sot.xcen = xcen_sot_cor
if not tag_exist(indexp_sot, 'ycen') then $
  indexp_sot = add_tag(indexp_sot, ycen_sot_cor, 'ycen') else $
  indexp_sot.ycen = ycen_sot_cor
PRINT, 'x offsets: ', xdisp_ref_arr
PRINT, 'y offsets: ', ydisp_ref_arr

indexp_sot.crval1 = crval1_sot_cor
indexp_sot.crval2 = crval2_sot_cor
indexp_sot.crpix1 = crpix1_sot_cor
indexp_sot.crpix2 = crpix2_sot_cor

indexp_sot = add_tag(indexp_sot, xdisp,     'xdisp')
indexp_sot = add_tag(indexp_sot, ydisp,     'ydisp')
indexp_sot = add_tag(indexp_sot, xdisp_rot, 'xdisp_rot')
indexp_sot = add_tag(indexp_sot, ydisp_rot, 'ydisp_rot')
indexp_sot = add_tag(indexp_sot, xdisp_tot, 'xdisp_tot')
indexp_sot = add_tag(indexp_sot, ydisp_tot, 'ydisp_tot')
indexp_sot = add_tag(indexp_sot, peak_str,  'peak_str')
indexp_sot = add_tag(indexp_sot, file_sot,  'file_sot')
indexp_sot = add_tag(indexp_sot, file_sdo,  'file_sdo')
indexp_sot = add_tag(indexp_sot, indexp_sot.date_obs, 't_sot')
indexp_sot = add_tag(indexp_sot, indexp_sdo.date_obs, 't_sdo')
indexp_sot = add_tag(indexp_sot, dt_sdo, 'dt_sdo')
indexp_sot = add_tag(indexp_sot, indexp_sdo.quality, 'quality_sdo')

; Create reduced-scale images and index recs of sot, sdo for archiving:
sdata_sot_unshift = confac(reform(cutout_cube_unshifted[*,*,1]), 0.25)
index_sot_unshift = indexp_sot
mreadfits_fixup, index_sot_unshift, sdata_sot_unshift
sdata_sot_shifted = confac(reform(cutout_cube[*,*,1]), 0.25)
index_sot_shifted = indexp_sot_orig
mreadfits_fixup, index_sot_shifted, sdata_sot_shifted
sdata_mini_sdo = confac(reform(cutout_cube[*,*,0]), 0.25)
index_mini_sdo = indexp_sdo
mreadfits_fixup, index_mini_sdo, sdata_mini_sdo

parent_sot_unshift = concat_dir(parent_mini_archives, 'mini_sot_unshift')
outdir_sot_unshift = ssw_time2paths(parent=parent_sot_unshift, times=index_sot_unshift.date_obs) & mk_dir, outdir_sot_unshift
filnam = str_replace(strcompress(get_infox(index_sot_unshift,'wave,naxis1,naxis2')),' ','_') + '_' + $
         time2file(index_sot_unshift.date_obs) + '.fits'
outfile = concat_dir(outdir_sot_unshift, filnam)
;mwritefits, index_sot_unshift, sdata_sot_unshift, outdir=outdir_sot_unshift, prefix=fprefix 
 mwritefits, index_sot_unshift, sdata_sot_unshift, outfile=outfile

parent_sot_shifted = concat_dir(parent_mini_archives, 'mini_sot_shifted')
outdir_sot_shifted = ssw_time2paths(parent=parent_sot_shifted, times=index_sot_shifted.date_obs) & mk_dir, outdir_sot_shifted
filnam = str_replace(strcompress(get_infox(index_sot_shifted,'wave,naxis1,naxis2')),' ','_') + '_' + $
         time2file(index_sot_shifted.date_obs) + '.fits'
outfile = concat_dir(outdir_sot_shifted, filnam)
;mwritefits, index_sot_shifted, sdata_sot_shifted, outdir=outdir_sot_shifted, prefix=fprefix 
 mwritefits, index_sot_shifted, sdata_sot_shifted, outfile=outfile

parent_mini_sdo = concat_dir(parent_mini_archives, 'mini_sdo')
outdir_mini_sdo = ssw_time2paths(parent=parent_mini_sdo, times=index_sot_shifted.date_obs) & mk_dir, outdir_mini_sdo
filnam = 'sdo_' + str_replace(strtrim(strcompress(get_infox(index_mini_sdo,'wavelnth,naxis1,naxis2')),2),' ','_') + '_' + $
         time2file(index_sot_shifted.date_obs) + '.fits'
outfile = concat_dir(outdir_mini_sdo, filnam)
;mwritefits, index_mini_sdo, sdata_mini_sdo, outdir=outdir_mini_sdo, prefix=fprefix 
 mwritefits, index_mini_sdo, sdata_mini_sdo, outfile=outfile

t1_exec_sec = anytim(!stime)
;t1_exec_sec = systime(/sec)

jsoc_time_sec = t1_exec_jsoc_sec - t0_exec_jsoc_sec
sbuff = 'JSOC time for this image       = ' + strtrim(jsoc_time_sec,2) + ' seconds.'
print, sbuff
file_append, file_log_exec_times, sbuff

read_sdo_time_sec = t1_exec_read_sdo_sec - t0_exec_read_sdo_sec
sbuff = 'READ_SDO time for this image   = ' + strtrim(read_sdo_time_sec,2) + ' seconds.'
print, sbuff
file_append, file_log_exec_times, sbuff

processing_time_sec = t1_exec_sec - t0_exec_sec
sbuff = 'Processing time for this image = ' + strtrim(processing_time_sec,2) + ' seconds.'
print, sbuff
file_append, file_log_exec_times, sbuff

indexp_sot = add_tag(indexp_sot, processing_time_sec, 'processing_time_sec')

return

failure: begin
;STOP, 'NEW_CORR: Stopping at beginning of failure GOTO block.'
  peak_str = 0
  error_corr = 1
  if exist(error_message) then file_append, file_error_log, error_message
; If available create reduced-scale images and index recs of the failed SOT image,
; as well as its sdo match (if one was found), for archiving:
  if exist(cutoutp_sot) then begin
  
     sdata_sot_failure = confac(reform(cutoutp_sot), 0.25)
     index_sot_failure = indexp_sot
     mreadfits_fixup, index_sot_failure, sdata_sot_failure

     parent_sot_failure = concat_dir(parent_mini_archives, 'mini_sot_failure')
     outdir_sot_failure = ssw_time2paths(parent=parent_sot_failure, times=indexp_sot.date_obs) & mk_dir, outdir_sot_failure
     filnam = str_replace(strcompress(get_infox(index_sot_failure,'wave,naxis1,naxis2')),' ','_') + '_' + $
              time2file(indexp_sot.date_obs) + '.fits'
     outfile = concat_dir(outdir_sot_failure, filnam)
;    mwritefits, index_sot_unshift, sdata_sot_unshift, outdir=outdir_sot_unshift, prefix=fprefix 
     mwritefits, index_sot_failure, sdata_sot_failure, outfile=outfile

     return
  endif
end

;sdata_sot_shifted = confac(reform(cutout_cube[*,*,1]), 0.25)
;index_sot_shifted = indexp_sot
;mreadfits_fixup, index_sot_shifted, sdata_sot_shifted

;parent_sot_shifted = concat_dir(parent_mini_archives, 'mini_sot_shifted')
;outdir_sot_shifted = ssw_time2paths(parent=parent_sot_shifted, times=index_sot_shifted.date_obs) & mk_dir, outdir_sot_shifted
;filnam = str_replace(strcompress(get_infox(index_sot_shifted,'wave,naxis1,naxis2')),' ','_') + '_' + $
;         time2file(index_sot_shifted.date_obs) + '.fits'
;outfile = concat_dir(outdir_sot_shifted, filnam)
;;mwritefits, index_sot_shifted, sdata_sot_shifted, outdir=outdir_sot_shifted, prefix=fprefix 
; mwritefits, index_sot_shifted, sdata_sot_shifted, outfile=outfile

;sdata_mini_sdo = confac(reform(cutout_cube[*,*,0]), 0.25)
;index_mini_sdo = indexp_sdo
;mreadfits_fixup, index_mini_sdo, sdata_mini_sdo

;parent_mini_sdo = concat_dir(parent_mini_archives, 'mini_sdo')
;outdir_mini_sdo = ssw_time2paths(parent=parent_mini_sdo, times=index_sot_shifted.date_obs) & mk_dir, outdir_mini_sdo
;filnam = 'sdo_' + str_replace(strtrim(strcompress(get_infox(index_mini_sdo,'wavelnth,naxis1,naxis2')),2),' ','_') + '_' + $
;         time2file(index_sot_shifted.date_obs) + '.fits'
;outfile = concat_dir(outdir_mini_sdo, filnam)
;;mwritefits, index_mini_sdo, sdata_mini_sdo, outdir=outdir_mini_sdo, prefix=fprefix 
; mwritefits, index_mini_sdo, sdata_mini_sdo, outfile=outfile

end