function fg_make, index, data, files=files, t0=t0, t1=t1, obs_type=obs_type, $
                  do_display=do_display, qstop=qstop

;+
;   Name: fg_make
;
;   Purpose: Create composite images from a set of SOT FG FITS files and return
;            a multi-dimensional array with dimensions:
;            
;   Input Parameters:
;   
;   Keyword Paramters:
;
;   Calling Sequence:
;     IDL> 
;
;   Calling Examples:
;     IDL>
;
;   Restrictions:
;
;   Nota Bena:
;
;   Suggestions:
;
;   TODO:
;
;   History:
;     18-nov-2014 - G L Slater - Contributed.
;-

if not keyword_set(level0) then level0 = 1

; Aribitrarily (for now), define time window
if not exist(t_off_sec) then t_off_sec = 120l ; 600l ; 3600l

if not exist(index) or not exist(data) then begin
  if exist(files) then begin
    mreadfits_header, files, index_sot
    ss_good = where(index_sot.fgmode eq 'shutterless', n_good)
    if n_good gt 0 then begin
      files_sot = files[ss_good]
      index_sot = index_sot[ss_good]
      t_sec_sot = anytim(index_sot)
      t_sec_sot = t_sec_sot[sort(t_sec_sot)]
      t0_sec = t_sec[0]
      t1_sec = t_sec_sot[n_elements(t_sec_sot)-1]
      sttim_cat = anytim(t0 - t_off_sec, /ccsds)
      entim_cat = anytim(t1 + t_off_sec, /ccsds)
      sot_cat, sttim_cat, entim_cat, search_arr=['obs_type=*shutterless*'], $
               cat_sot, files_sot, count=count, level=level0
      if count gt 0 then begin
        mreadfits_header, files_sot, index_sot
      endif else begin
        print, ' No shutterless files found in time interval.  Returning.'
        return, ''
      endelse
    endif else begin
      print, ' None of the files are for shutterless images.  Returning.'
      return, ''
    endelse
  endif else begin
    if exist(t0) then begin
      t0_sec = anytim(t0)
      if not exist(t1) then t1 = t0
      t1_sec = anytim(t1)
      sttim_cat = anytim(t0_sec - t_off_sec, /ccsds)
      entim_cat = anytim(t1_sec + t_off_sec, /ccsds)
      sot_cat, sttim_cat, entim_cat, search_arr=['obs_type=*shutterless*'], $
               cat_sot, files_sot, count=count, level=level0
      if count gt 0 then begin
        mreadfits_header, files_sot, index_sot
      endif else begin
        print, ' No shutterless files found in time interval.  Returning.'
        return, ''
      endelse
    endif else begin
      print, ' None of index, files, or time range passed.  Returning.'
      return, ''
    endelse
  endelse
endif else begin
  index_sot = index
  ss_good = where(index_sot.fgmode eq 'shutterless', n_good)
  if n_good gt 0 then begin
    index_sot = index_sot[ss_good]
        t_sec_sot = anytim(index_sot.date_obs)
    t_sec_sot = t_sec[sort(t_sec)]
    t0_sec = t_sec_sot[0]
    t1_sec = t_sec_sot[n_elements(t_sec_sot)-1]
    sttim_cat = anytim(t0_sec - t_off_sec, /ccsds)
    entim_cat = anytim(t1_sec + t_off_sec, /ccsds)
    sot_cat, sttim_cat, entim_cat, search_arr=['obs_type=*shutterless*'], $
             cat_sot, files_sot, count=count, level=level0
    if count gt 0 then begin
      mreadfits_header, files_sot, index_sot
    endif else begin
      print, ' No shutterless files found in time interval.  Returning.'
      return, ''
    endelse
  endif else begin
    print, ' None of the images are shutterless.  Returning.'
    return, ''
  endelse
endelse

; Select mode:
mode_tags = 'naxis1,naxis2,wave,waveid,obs_type'
ss_mode = ssw_uniq_modes(index_sot, mode_tags, /interactive)
files_sot = files_sot[ss_mode]
read_sot, files_sot, index_sot, data_sot

n_files_sot = n_elements(files_sot)

t_sot = index_sot.date_obs
t_sec_sot = anytim(t_sot)
ss_t0_closest = min(where(t_sec_sot ge t0_sec, n_match_t0))
ss_t1_closest = max(where(t_sec_sot le t1_sec, n_match_t1))

roiloop = index_sot.roiloop
waveoff = index_sot.waveoff
max_roiloop = max(roiloop)
n_roi = max_roiloop > 1
n_waveoff = n_elements(all_vals(waveoff))
siz_data_sot = size(data_sot)
n_polar = siz_data_sot[3]

n_img_per_cycle = n_roi*n_waveoff

tdelt_sot = t_sec_sot[1:*] - t_sec_sot
tdelt_fix = fix(tdelt_sot*10)/10

main_cnt = index_sot.main_cnt
delta_main_cnt = main_cnt[1:*] - main_cnt
ss_change_main = where(delta_main_cnt ne 0, n_change_main)
if n_change_main gt 0 then begin
  length_vals_main = ss_change_main[1:*] - ss_change_main
  length_mode_main = mode_val(length_vals_main, freq_mode_main)
  n_mode_main = n_elements(freq_mode_main)
endif else begin
  n_mode_main = 0
endelse

subr_cnt = index_sot.subr_cnt
delta_subr_cnt = subr_cnt[1:*] - subr_cnt
ss_change_subr = where(delta_subr_cnt ne 0, n_change_subr)
if n_change_subr gt 0 then begin
  length_vals_subr = ss_change_subr[1:*] - ss_change_subr
  length_mode_subr = mode_val(length_vals_subr, freq_mode_subr)
  n_mode_subr = n_elements(freq_mode_subr)
endif else begin
  n_mode_subr = 0
endelse

seqn_cnt = index_sot.seqn_cnt
delta_seqn_cnt = seqn_cnt[1:*] - seqn_cnt
ss_change_seqn = where(delta_seqn_cnt ne 0, n_change_seqn)
if n_change_seqn gt 0 then begin
  length_vals_seqn = ss_change_seqn[1:*] - ss_change_seqn
  length_mode_seqn = mode_val(length_vals_seqn, freq_mode_seqn)
  n_mode_seqn = n_elements(freq_mode_seqn)
endif else begin
  n_mode_seqn = 0
endelse

max_n_mode = max([n_mode_main, n_mode_subr, n_mode_seqn], ss_max_mode)
case ss_max_mode of
  0: begin
       delta_cnt = delta_main_cnt
       ss_change = ss_change_main
       length_vals = length_vals_main
       length_mode = length_mode_main
       freq_mode = freq_mode_main
       n_mode = n_mode_main
       end
  1: begin
       delta_cnt = delta_subr_cnt
       ss_change = ss_change_subr
       length_vals = length_vals_subr
       length_mode = length_mode_subr
       freq_mode = freq_mode_subr
       n_mode = n_mode_subr
       end
  2: begin
       delta_cnt = delta_seqn_cnt
       ss_change = ss_change_seqn
       length_vals = length_vals_seqn
       length_mode = length_mode_seqn
       freq_mode = freq_mode_seqn
       n_mode = n_mode_seqn
       end
  else: begin
       print, ' No max freq mode found.'
       end
endcase

if length_mode ne n_img_per_cycle then begin
  stop, 'length_mode ne n_img_per_cycle. Stopping.'
endif

ss_length_mode = where(length_vals eq length_mode, n_length_mode)

;for i=0, n_length_mode-1 do begin
;  if not exist(ss_good_arr) then $
;    ss_good_arr = ss_length_mode[i] + indgen(length_mode)
;  endif else begin
;    ss_good_arr = [ss_good_arr, ss_length_mode[i] + indgen(length_mode)]
;  endif
;endfor

ss_good_cycle_start = ss_change[ss_length_mode] + 1
n_good_cycle_start = n_length_mode
;STOP, 'STOPPING AFTER CALCULATION OF SS_GOOD_CYCLE_START.'
xsiz_comp = siz_data_sot[1]*n_roi
comp_arr = fltarr(siz_data_sot[1]*n_roi, siz_data_sot[2], n_polar, n_waveoff, n_good_cycle_start)

for i=0, n_good_cycle_start-1 do begin
  for j=0, n_img_per_cycle-1 do begin
    data0 = reform(data_sot[*,*,*,ss_good_cycle_start[i]+j])
    for k=0, n_polar-1 do begin
      polar = k
      woff  = j  /  n_roi
      roi   = j mod n_roi

if roi gt 0 then $
  subimg0 = shift(data0[0:(siz_data_sot[1]/2),*,polar], 1) else $
  subimg0 =       data0[0:(siz_data_sot[1]/2),*,polar]
      subimg1 = data0[(siz_data_sot[1]/2):*,*,polar]
      if k eq 0 then begin
        subimg0 = uint(subimg0)
        subimg1 = uint(subimg1)
      endif
      comp_arr[(siz_data_sot[1]/2d)*(n_roi-(roi+1)), 0, polar, woff, i] = subimg0
      comp_arr[(siz_data_sot[1]/2d)*(n_roi+(roi  )), 0, polar, woff, i] = subimg1
    endfor
  endfor
endfor    

if keyword_set(do_display) then begin
; if not exist(lo_cut_non_mag) then lo_cut_non_mag =  16000
; if not exist(hi_cut_non_mag) then hi_cut_non_mag =  33000
  if not exist(lo_cut_non_mag) then lo_cut_non_mag =  16000
  if not exist(hi_cut_non_mag) then hi_cut_non_mag =  33000
  if not exist(lo_cut_mag)     then lo_cut_mag     = - 1000
  if not exist(hi_cut_mag)     then hi_cut_mag     =   1000

  ans = ''
  help, comp_arr
  print, ' Enter (polarization, waveoff) pair for movie (e.g. "1,0"),'
  read,  ' Return to quit, default is quit: ', ans
; while strlowcase(strmid(ans,0,1)) ne '' do begin
  while ans ne '' do begin
    comp_ind_arr = str2arr(ans)
    polar_ind = comp_ind_arr[0]
    wave_ind = comp_ind_arr[1]
    comp_arr0 = reform(comp_arr[*,*,polar_ind,wave_ind,*])
    if polar_ind eq 0 then begin
;     lo_cut = lo_cut_non_mag
;     hi_cut = hi_cut_non_mag
      lo_cut = min(comp_arr0)
      hi_cut = max(comp_arr0)
    endif else begin
      lo_cut = lo_cut_mag
      hi_cut = hi_cut_mag
    endelse

;   xstepper, sot_filt(lo_cut>reform(comp_arr[*,*,polar_ind,wave_ind,*])<hi_cut)
    xstepper, lo_cut>comp_arr0<hi_cut

    help, comp_arr
    print, ' Enter (polarization,waveoff) pair to movie (e.g. "1,0"),'
    read,  ' Return to quit, default is quit: ', ans
  endwhile
endif

help, comp_arr, n_roi, n_polar, n_waveoff

if keyword_set(qstop) then stop, 'FG_MAKE: Stopping on request before return.'
return, comp_arr

end