pro sotsp_l21struct2data, l21struct, index, data, $
   scale=scale,autoscale=autoscale,scan_info=scan_info

;+
;   Name: sotsp_l21struct2data
;
;   Purpose: SOT/SP Level2.1 struct (from binext FITS) -> "index,data"
;
;   Input Parameters:
;      l21struct - structure, from Level2.1 BinExt FITS (ME0 disambiguation)
;
;   Output:
;      returns index,data  
;
;   Keyword Parameters:
;     scale (switch) - if set, apply PI teams suggested scale (def=floating)
;     scan_info (input) ; structure with 1D tags TIMES, SLITPOS, SCATTPRO
;          (scan time, slitpos, scattered light profile), usually from Level2
;          file called before (see 2020 version of read_sotsp.pro)
;
;   History:
;     9 Jul 2020 - M.DeRosa - created, based on SLF's sotsp_l2struct2data.pro
;     2 Sep 2020 - M.DeRosa - changed array indexing from () to [] to
;                  for image[ssinf] to avoid namespace collision w/new
;                  IDL function 'image()'
;    23-Dec-2021 - M.DeRosa - (version 1.1) fixed two issues when /scale
;                  switch is set: (1) routine now removes non-finite values
;                  (Inf or NaN) before performing the scaling, and (2) routine
;                  now ensures positive values when doing log scaling for
;                  cases where the minimum image value is 0
;
;   Restrictions:
;      Not for casual users; 'read_sotsp.pro' is recommended for them
;      Currently includes All data extensions in "data" - plan for user
;      defined subset which is easy but not done as of now 
;
;-
version=1.1

if ~required_tags(l21struct,'br,bt,bp,field_azimuth_disamb') then begin 
  box_message,'Expect SOTSP Level2.1 FITS structure'
  return
endif
if ~keyword_set(scan_info) then begin
  box_message,'Required keyword scan_info not present, returning.'
  return
endif

time_window,scan_info.times,date_obs,date_end

autoscale=keyword_set(autoscale)
scale=keyword_set(scale) or autoscale

fitstags=['BP','BT','BR','FIELD_AZIMUTH_DISAMB']  ;  fits tag names
fieldlabs=['Magnetic Field Zonal Component (positive west)',$
           'Magnetic Field Meridional Component (positive north)', $
           'Magnetic Field Radial Component (positive up)', $
           'Magnetic Field Azimuth (disambiguated)']
ranges=[[-1e3,1e3],[-1e3,1e3],[-1e3,1e3],[0,360]]
logs=[0,0,0,0]
units=['Gauss','Gauss','Gauss','Degrees']

nout=4

index=replicate(fitshead2struct(l21struct.header),nout)
mag=gt_tagval(l21struct,/br)
data=make_array(data_chk(mag,/nx),data_chk(mag,/ny), nout, /float)

for i=0,nout-1 do begin
  tind=tag_index(l21struct,fitstags(i))
  if tind gt 0 then data(0,0,i)=l21struct.(tind) 
endfor

;  loop through the various output quantities and scale
if keyword_set(scale) then begin
  sdata=make_array(data_chk(data,/nx),data_chk(data,/ny),nout,/byte) 
  low=reform(ranges(0,*))
  high=reform(ranges(1,*))
  for i=0,nout-1 do begin  ; optional scaling loop

    image=data(*,*,i)

    ;  remove bad scaling values
    ssinf=where(1-finite(image),icnt)
    if icnt gt 0 then begin 
      box_message,'Inf fixup.'
      image[ssinf]=0.
    endif
    
    ;  re-scale, if desired/needed
    if (low(i) eq high(i)) or autoscale then begin 
      top=fix(get_logenv('sotsp_l2top'))
      hstats=ssw_image_stats(index(i),image,percent=[3,top],status=status) ; 
      if not status then begin 
        box_message,'Problem with Hist, reducing percentile..'
        hstats=ssw_image_stats(index(i), image,percent=[3,top-6], status=status)
      endif
      if not status then hstats={low:min(image),high:max(image)}
      low(i)=hstats.(0)
      high(i)=hstats.(1)
    endif

    ;  scale the image
    if logs(i) then begin  
      if low(i) eq 0 then low(i)=min(image+(image eq 0)*high(i))
      sdata(0,0,i) = bytscl(alog10(image>low(i)), $
                            min=alog10(low(i)),max=alog10(high(i)))
    endif else sdata(0,0,i)=bytscl(image,min=low(i),max=high(i))
  endfor
  delvarx,data & data=temporary(sdata)
endif

info={FITSTAGS:'',CUNIT1:'',CUNIT2:'',FTYPE:'',INSTRUME:'',TELESCOP:'', $
      cdelt1:0.0, cdelt2:0.0, crota1:0, crota2:0, date_obs:'', date_end:'', $
      BUNIT:''}
info=replicate(info,nout)
info.fitstags=fitstags
info.cunit1='arcsec'
info.cunit2='arcsec'
info.ftype=fieldlabs
info.telescop='HINODE'
info.INSTRUME='SOT/SP ' + fieldlabs
info.date_obs=date_obs
info.date_end=date_end
info.cdelt1=gt_tagval(l21struct.header,/xscale,missing=0.29)
info.cdelt2=gt_tagval(l21struct.header,/yscale,missing=0.32)
info.bunit=units
index=join_struct(temporary(index),info)

update_history,index,version=version,/caller

return
end