pro read_soon, files, dset, index, data, pindex=pindex, $
	       primary_only=primary_only, imgperfile=imgperfile, $
	       loud=loud
;+ 
;   Name: read_soon
;
;   Purpose: read one or more soon 2D/3D files
;
;   Input Parameters:
;
;       files - list of files to read
;       dset  - subscripts of images relative to <files>
;               [ -1 implies all ]
;       index - output vector of 'header-structures'
;       data  - 2D/3D data (image or cube)
;
;   Keyword Parameters:
;      pindex (OUTPUT) - the primary headers (one per input file)
;      loud (switch) - if set, make soonfits loud
;  
;   History:
;     2-march-1999 - S.L.Freeland - written
;                    Based on earlier email to Alan Kiplinger
;     3-march-1999 - S.L.Freeland - set ref CRVALx->arcsecs using lat/lon
;                                   CROTA -> P angle
;     7-Apr-1999   - S.L.Freeland - use solar Pangle (soon.PA was position angle..)
;
;     13-May-1999  - Zarro (SM&A/GSFC) - added initial call to struct2ssw to
;                    get standard time input for conv_h2a.
;                    Also fixed P-angle sign
;      4-apr-2006 - S.L.Freeland - fix date_d$obs 
;
;   Method:
;      Use mreadfits (SLF), soonfits(Randy Miesner) and SSW
;
;   Referenced routines:
;      mreadfits, soonfits, fitshead2struct, struct2ssw, $
;      gt_tagval, mxf_dset_arr, conv_h2a
;  
; Calling Example
;----------------------------------------------------
;   IDL> read_soon,files,-1,index                            ; all headers->vect
;   IDL> ss=where(index.xxx eq aaa and index.yyy eq bbb etc) ; SOME filter
;   IDL> read_soon,files,SS, index, data                     ; read subset
;----------------------------------------------------
;-

if n_params() lt 3 then begin
  box_message,['Expect at least filelist(in), dset(in) and index(out)',$
	       'IDL> read_soon, filelist, dset, index [,data ,pindex=pindex]']
  return
endif    

mreadfits, files, pindex                       ; read primary header
imgperfile=gt_tagval(pindex,/FILE_D$LEN)       ; images per file (vector)
totimages=total(imgperfile)                    ; total for this filelist

if keyword_set(primary_only) then begin
   box_message,'Only reading primary headers->index and returning....'
   index=pindex
   return
endif

; ------- set a couple of useful variables --------
if dset(0) eq -1 then ss=lindgen(totimages) else ss=dset
nimage=n_elements(ss)
 
; ---- read 1st extension, convert-> str and replicate --------
head=soonfits(files(0),ext=0,/header) ; get extension header template
tdata=soonfits(files(0),ext=0)        ; get image template
headstr=fitshead2struct(head,/add_standard)         ; structure version
allstr=replicate(headstr,nimage)      ; all index output vector

; --------- map -----------------
darr=mxf_dset_map(pindex,ss)          ; map file/dset -> dset map

for i=0,nimage-1 do begin 
   head = soonfits(files(darr(i).ifile) $               ; read extended header
	   ,exten=darr(i).dset,/header_only, loud=loud)
   headstr=fitshead2struct(head,headstr)                ; header->structure
   allstr(i)=headstr                                    ; insert in output
endfor

if n_params() gt 3 then begin                   ; data output is requested
   data=make_array(data_chk(tdata,/nx),$        ; preformat output
		   data_chk(tdata,/ny),$
		   nimage,type=data_chk(tdata,/type))
   for i=0,nimage-1 do begin 
      image = soonfits(files(darr(i).ifile) $   ; read this image extension
	  ,exten=darr(i).dset,loud=loud)
      data(0,0,i)=image                         ; insert in cube
   endfor
endif

; -------- convert structures -> SSW -------
;   slf, 3-march-1999 convert LAT/LON->CRVALx

;-- 11 May 1999, DMZ added to get time right for conv_h2a

index=struct2ssw(allstr,/nopointing,/nosolar)                                         

if required_tags(allstr(0),'lon,lat,crpix1,crpix2') then begin 
   if max(abs(gt_tagval(allstr,/crval1,missing=0))) eq 0 then begin 
      arcs=conv_h2a(transpose([[allstr.lon],[allstr.lat]]),index)
      ew=reform(arcs(0,*)) & ns=reform(arcs(1,*))
      if n_elements(ns) eq 1 then begin
	  ns=ns(0)
	  ew=ew(0)
      endif	  
      index.crval1=ew & index.crval2=ns 
   endif
endif   

index=struct2ssw(index,/notime)                   ; point standards

if tag_exist(index,'date_d$obs') then index=rem_tag(temporary(index),'date_d$obs')

sunstuff=get_sun(index,p=pangle)
index.crota=-pangle

return
end
