
function aia_fix_header, iindex, odata, normalize=normalize

; TODO: Handle case of differing sizes of image and reference image
; Update header tag values as needed:

;+
; NAME:
;   AIA_FIX_HEADER
; PURPOSE:
;   Update headers for registered images
; CATEGORY:
; SAMPLE CALLS:
; INPUTS:
; OUTPUT:
;   Updated header structure
; KEYWORDS:
;   normalize - If set, and image is AIA, then document in HISTORY tag.
; TODO:
; HISTORY:
;   GLS - 2013-04-01 - If exposure normailzation has been done, then
;                      add this to history tag
;   GLS - 2018-06-21 - Added updating of data statistics tags
;   GLS - 2018-07-03 - Added check for existence of header tag 'WAVELNTH'
;                      as prerequisite for updating statistics tags
;-

; Define prognam, progver variables
prognam = 'AIA_FIX_HEADER.PRO'
progver = 'V1.00' ; 2013-05-28 (GLS)

oindex = iindex

instr_prefix = strupcase(strmid(oindex.instrume,0,3))

; Update LVL_NUM tag:
if tag_exist(oindex, 'LVL_NUM') then oindex.lvl_num = 1.5 else $
   oindex = add_tag(oindex, 1.5, 'LVL_NUM')

; If exposure normalization has been done, then update history for this:
if keyword_set(normalize) then begin
   new_hist_rec = 'Exposure normalization preformed.'
   update_history, oindex, new_hist_rec
endif

; Update or add XCEN, YCEN tags:
xcen = comp_fits_cen(oindex.crpix1, oindex.cdelt1, oindex.naxis1, oindex.crval1)
if tag_exist(oindex, 'XCEN') then oindex.xcen = xcen else $
   oindex = add_tag(oindex, xcen, 'XCEN')
ycen = comp_fits_cen(oindex.crpix2, oindex.cdelt2, oindex.naxis2, oindex.crval2)
if tag_exist(oindex, 'YCEN') then oindex.ycen = ycen else $
   oindex = add_tag(oindex, ycen, 'YCEN')

; Update or add R_SUN tag, if RSUN_OBS tag exists:
if tag_exist(oindex, 'RSUN_OBS') then begin
   if tag_exist(oindex, 'R_SUN') then $
      oindex.r_sun = oindex.rsun_obs/oindex.cdelt1-1 else $
         oindex = add_tag(oindex, oindex.rsun_obs/oindex.cdelt1-1, 'R_SUN')
;     oindex.r_sun = oindex.rsun_obs/oindex.cdelt1-1 else $
;        oindex = add_tag(oindex, oindex.rsun_obs/oindex.cdelt1-1, 'R_SUN')
endif

; Update data statistics header tags:

; Create mask for on-disk pixels, in preparation for calculating statistics tags:
r_mask = oindex.r_sun ; 0.9*oindex.r_sun
dum = (dindgen(oindex.naxis1)-oindex.naxis1/2)
mask_disk = (((dum##(dblarr(oindex.naxis1)+1.))^2 + $
             ((dblarr(oindex.naxis1)+1.)##dum)^2) lt (r_mask^2))
ss_ondisk  = where(mask_disk ne 0, n_ondisk)
ss_offdisk = where(mask_disk eq 0, n_offdisk)

; Mask limb for magnetogram:
; r_mask = index_blos.rsun_obs / index_blos.cdelt1-1
; dum = (findgen(4096)-2048)
; mask_mag = (((dum##(dblarr(4096)+1.))^2 + ((dblarr(4096)+1.)##dum)^2) lt (r_mask^2))
; ss_offlimb_mag = where(mask_mag eq 0, n_offlimb_mag)
; datap_scaled[ss_offlimb_mag] = 0

; Update/add data statistics keywords:
if not tag_exist(oindex, 'WAVELNTH') then begin
   print, "Header tag 'WAVELNTH' does not exits for input header."
   print, "Not updating data statistics header tags."
endif else begin
   if oindex.wavelnth eq 6173 then begin
      if tag_exist(oindex, 'DATAMIN')  then oindex.datamin = min(odata[ss_ondisk])
      if tag_exist(oindex, 'DATAMAX')  then oindex.datamax = max(odata[ss_ondisk])
      if tag_exist(oindex, 'DATAMEDN') then oindex.datamedn = median(odata[ss_ondisk])

;     moments_odata_ondisk = moment(finite(odata[ss_ondisk]))
      moments_odata_ondisk = moment(odata[ss_ondisk])
      if tag_exist(oindex, 'DATAMEAN') then oindex.datamean = moments_odata_ondisk[0]
      if tag_exist(oindex, 'DATARMS')  then oindex.datarms  = moments_odata_ondisk[1]^2
      if tag_exist(oindex, 'DATASKEW') then oindex.dataskew = moments_odata_ondisk[2]
      if tag_exist(oindex, 'DATAKURT') then oindex.datakurt = moments_odata_ondisk[3]
   endif else begin
      if tag_exist(oindex, 'DATAMIN')  then oindex.datamin = min(odata)
      if tag_exist(oindex, 'DATAMAX')  then oindex.datamax = max(odata)
      if tag_exist(oindex, 'DATAMEDN') then oindex.datamedn = median(odata)

;     moments_odata = moment(finite(odata))
      moments_odata = moment(odata)
      if tag_exist(oindex, 'DATAMEAN') then oindex.datamean = moments_odata[0]
      if tag_exist(oindex, 'DATARMS')  then oindex.datarms  = moments_odata[1]^2
      if tag_exist(oindex, 'DATASKEW') then oindex.dataskew = moments_odata[2]
      if tag_exist(oindex, 'DATAKURT') then oindex.datakurt = moments_odata[3]
   endelse
endelse

; Update or add miscellaneous header keyword records:

; Add miscellaneous HISTORY keyword records:
update_history, oindex, version=progver, caller=prognam

return, oindex

end


