;
; Plots Apr-Sep means of NOAA index timeseries, and saves them too.
;
loadct,39
multi_plot,nrow=5,ncol=2,layout='large'
if !d.name eq 'X' then begin
  window,ysize=850
  !p.font=-1
endif else begin
  !p.font=0
  device,/helvetica,/bold,font_size=9
endelse
;
; Read in all indices
;
openr,1,'noaa_indices.txt'
headst=strarr(21)
readf,1,headst
nyr=1997-1950+1
nmon=12
nt=nyr*nmon
rawdat=fltarr(16,nt)
readf,1,rawdat,format='(I4,I3,1X,14F5.1)'
close,1
timey=findgen(nyr)+1950.
;
; Deal with 'missing' values.  They are not really missing, it is just that
; some modes do not have any power in some months.  We want an Apr-Sep mean,
; so we could ignore months in this period that have no values, but it is
; probably better to treat them as zero.  First they're set to NaN, but
; treated as zero later.  Need 3 out of 6 months to get a mean
; for the season!
;
ml=where(rawdat eq -9.9,nmiss)
if nmiss gt 0 then rawdat(ml)=!values.f_nan
;
ind=['NAO','EA','EA-JET','WP','EP','NP','PNA','EA/WR','SCA','TNH','POL',$
  'PT','SZ','ASU']
nind=n_elements(ind)
;
nsummer=10
summerind=strarr(nsummer)
summerts=fltarr(nsummer,nyr)
summertslow=fltarr(nsummer,nyr)
isummer=0
;
for i = 0 , nind-1 do begin
  monts=reform(rawdat(i+2,*),nmon,nyr)
  monts=monts(3:8,*)
  seastot=total(monts,1,/nan)
  seasnum=total(finite(monts),1)
  seasts=seastot/6.           ; always divide by six regardless of missing!
  ml=where(seasnum lt 3.,nmiss)
  if nmiss gt 0 then seasts(ml)=!values.f_nan
  ;
  kl=where(finite(seasts),nkeep)
  if nkeep gt 20 then begin        ; only plot if at least 20 yrs
    filter_cru,10.,/nan,tsin=seasts,tslow=tslow
    pause
    plot,timey,seasts,$
      /xstyle,xtitle='Year',$
      ytitle='Apr-Sep mean index value',$
      title=ind(i)
    oplot,!x.crange,[0.,0.],linestyle=1
    oplot,timey,tslow,thick=3
    summerind(isummer)=ind(i)
    summerts(isummer,*)=seasts(*)
    summertslow(isummer,*)=tslow(*)
    isummer=isummer+1
  endif
  ;
endfor
;
save,filename='noaa_summer.idlsave',$
  nsummer,summerind,summerts,summertslow,timey
;
end