;
; Plots low-frequency timeseries and their differences
;
;plot,[0,1]
loadct,0
def_1color,r,g,b,10,color='lgrey'
multi_plot,nrow=11,ncol=4,layout='large'
if !d.name eq 'X' then begin
  window,ysize=900
endif else begin
  device,xsize=18,ysize=18
endelse
!x.omargin=[6.5,2.5]
!x.ticklen=0.06
!y.omargin=[4,2]
!y.margin=[0,0]
!y.ticklen=0.04
!y.minor=2
!p.charsize=!p.charsize
yr=[-2.99,2.99]
;
thalf=10.
pp=[1881,1940]
refp=[ [1881,1975] , [pp] , [pp] , [pp] , [pp] , [pp] , [pp] , [pp] , [pp] , $
  [pp] , [pp] ]
regorder=[0,1,2,3,4,5,6,7,9,10,8]
;
; Get region info
;
restore,filename='regboxes.idlsave'
for jreg = 0 , nreg-1 do begin
  ireg=regorder(jreg)
  refperiod=reform(refp(*,ireg))
  restore,filename='instradj_'+regname(ireg)+'.idlsave'
  ;
  openw,1,'reginstr_aprsep_'+regname[ireg]+'.txt'
  printf,1,'Regional-mean CRUTEM1 instrumental temperature for Apr-Sep for region '+regname[ireg]
  for iyr = 0 , n_elements(x)-1 do printf,1,x[iyr],instradj[iyr],format='(I4,F10.3)'
  close,1
  ;
  wdens=instradj
  mknormal,wdens,x,refperiod=refperiod,refmean=refmean,refsd=refsd
  filter_cru,thalf,tsin=wdens,tslow=instrlow,tshigh=instrhi,/nan
  ;
  restore,filename='jjaadj_'+regname(ireg)+'.idlsave'
  ;
  openw,1,'reginstr_junaug_'+regname[ireg]+'.txt'
  printf,1,'Regional-mean CRUTEM1 instrumental temperature for Jun-Aug for region '+regname[ireg]
  for iyr = 0 , n_elements(x)-1 do printf,1,x[iyr],jjaadj[iyr],format='(I4,F10.3)'
  close,1
  ;
  wdens=jjaadj
  mknormal,wdens,x,refperiod=refperiod,refmean=refmean,refsd=refsd
  filter_cru,thalf,tsin=wdens,tslow=jjalow,tshigh=jjahi,/nan
  ;
  restore,filename='densadj_'+regname(ireg)+'.idlsave'
  mknormal,densadj,x,refperiod=refperiod,refmean=refmean,refsd=refsd
  filter_cru,thalf,tsin=densadj,tslow=denslow,tshigh=denshi,/nan
  ;
  restore,filename='rwidadj_'+regname(ireg)+'.idlsave'
  mknormal,rwidadj,x,refperiod=refperiod,refmean=refmean,refsd=refsd
  filter_cru,thalf,tsin=rwidadj,tslow=rwidlow,tshigh=rwidhi,/nan
  ;
  if regname(ireg) eq 'NWNA' then begin
    dellist=where(x le 1890.)
    instrlow(dellist)=!values.f_nan
    jjalow(dellist)=!values.f_nan
  endif
  ;
  if !p.multi(0) eq 4 then begin
    xt=' '
    xtf=''
  endif else begin
    xt=' '
    xtf='nolabels'
  endelse
  yt=' '
  ytf=''
  if !p.multi(0) eq 24 then $
    yt='!5Smoothed anomalies (standard deviation units)!3'
  ;
  !x.margin=[0.5,0]
  plot,x,/nodata,$
    xrange=[1881,1992],xstyle=1,xtitle=xt,xtickformat=xtf,$
    yrange=yr,ystyle=1,ytickformat=ytf
  oplot,x,instrlow,thick=4
  oplot,x,denslow
  oplot,!x.crange,[0.,0.],thick=0.5
  ;
  if yt ne ' ' then xyouts,1860,0,yt,orient=90.,align=0.5,charsize=0.8*!p.charsize
  ;
  txt=regname(ireg)
  if ireg eq 0 then txt='ESIB'
  if ireg eq 1 then txt='CSIB'
  if ireg eq 2 then txt='WSIB'
  xyouts,1885,1.6,'!5'+txt+'!3',charsize=!p.charsize
  ;
  if jreg eq 0 then xyouts,1992,!y.crange(1)+0.5,'!5Ring density!3',align=0.5,$
    charsize=!p.charsize*0.9
  if jreg eq nreg-1 then xyouts,1992,!y.crange(0)-2.4,'!5Year!3',align=0.5,$
    charsize=!p.charsize*0.8
  ;
  yt=' '
  ytf='nolabels'
  ;
  !x.margin=[0,0.5]
  plot,x,/nodata,$
    xrange=[1881,1992],xstyle=1,xtitle=xt,xtickformat=xtf,$
    yrange=yr,ystyle=1,ytitle=yt,ytickformat=ytf
  ;
  aa=denslow-instrlow
  xx=x
  keeplist=where(xx ge 1881.)
  aa=aa(keeplist)
  xx=xx(keeplist)
  misslist=where(finite(aa) eq 0,nmiss)
  if nmiss gt 0 then aa(misslist)=0.
  poslist=where(aa gt 0.,npos)
  if npos gt 0 then aa(poslist)=0.
  nend=n_elements(aa)-1
  xx=[xx(0),xx,xx(nend),xx(0)]
  aa=[0.,aa,0.,0.]
  polyfill,xx,aa,color=10
  oplot,x,denslow-instrlow
  oplot,!x.crange,[0.,0.],thick=0.5
  ;
  axis,yaxis=1,ystyle=1,ytickformat='nolabels'
  ;
  if regname(ireg) eq 'ALL' then begin
    xlen=!x.crange(1)-!x.crange(0)
    xrig=!x.crange(1)
    xbox=xrig-[xlen*2.,0.,0.,xlen*2.,xlen*2.]
    plots,xbox,$
      [!y.crange(0),!y.crange,reverse(!y.crange)],thick=6
  endif
  ;
  !x.margin=[0.5,0]
  plot,x,/nodata,$
    xrange=[1881,1992],xstyle=1,xtitle=xt,xtickformat=xtf,$
    yrange=yr,ystyle=1,ytitle=yt,ytickformat=ytf
  oplot,x,jjalow,thick=4
  oplot,x,rwidlow
  oplot,!x.crange,[0.,0.],thick=0.5
  ;
  if jreg eq 0 then xyouts,1992,!y.crange(1)+0.5,'!5Ring width!3',align=0.5,$
    charsize=!p.charsize*0.9
  if jreg eq nreg-1 then xyouts,1992,!y.crange(0)-2.4,'!5Year!3',align=0.5,$
    charsize=!p.charsize*0.8
  ;
  !x.margin=[0,0.5]
  plot,x,/nodata,$
    xrange=[1881,1992],xstyle=1,xtitle=xt,xtickformat=xtf,$
    yrange=yr,ystyle=1,ytitle=yt,ytickformat=ytf
  ;
  aa=rwidlow-jjalow
  xx=x
  keeplist=where(xx ge 1881.)
  aa=aa(keeplist)
  xx=xx(keeplist)
  misslist=where(finite(aa) eq 0,nmiss)
  if nmiss gt 0 then aa(misslist)=0.
  poslist=where(aa gt 0.,npos)
  if npos gt 0 then aa(poslist)=0.
  nend=n_elements(aa)-1
  xx=[xx(0),xx,xx(nend),xx(0)]
  aa=[0.,aa,0.,0.]
  polyfill,xx,aa,color=10
  oplot,x,rwidlow-jjalow
  oplot,!x.crange,[0.,0.],thick=0.5
  ;
  axis,yaxis=1,ystyle=1,yticks=4,ytickv=findgen(5)-2,$
    ytickname=['-2','-1',' 0',' 1',' 2']
  ;
  if regname(ireg) eq 'ALL' then begin
    xlen=!x.crange(1)-!x.crange(0)
    xrig=!x.crange(1)
    xbox=xrig-[xlen*2.,0.,0.,xlen*2.,xlen*2.]
    plots,xbox,$
      [!y.crange(0),!y.crange,reverse(!y.crange)],thick=6
  endif
  ;
endfor
;
end