;
; Correlates the (calibrated) regional reconstructions with the instrumental
; regional series in sliding windows, hi and low freq.
;
showband=0        ; 0=do Hug , 1=do Age
dolongeuro=0      ; 0=don't include long euro instrumental series, 1=do
docol=1           ; 0=B&W 1=color
dopub=2           ; 0=draft 1=publication quality 2=just ALL
;
if showband eq 1 then begin
  print,'THE ALL REGION FOR AGE-BANDED DATA IS NOT THE PCR VERSION!!!!'
endif
;
if showband eq 0 then begin
  maints='reg'
endif else begin
  maints='band'
endelse
;
hiwin=20.
lowin=50.
thalf=10.
;
; The following are uncertainty ranges for a 20-yr sample, given the mean
; correlation as an estimate of the population correlation coefficient.
;
rrmin=[0.66,0.35,0.36,0.24,0.08,-0.1,0.60,0.54,0.36,0.46]
rrmax=[0.92,0.83,0.83,0.78,0.69,0.60,0.91,0.89,0.83,0.86]
;
if dopub eq 0 then begin
  xran=[1725,1994]
endif else begin
  xran=[1861,1994]
endelse
;
; Read in long instrumental records from Europe
;
restore,filename='/cru/u2/f055/data/obs/station/surface/longeuro/longeuro_seas.idlsave'
;  Gets: nyr,timey,nseas,seasname,nstat,statname,nreg,regname,seasts,regseas
;
liyear=timey
lineur=reform(regseas(0,*,5))
liseur=reform(regseas(2,*,5))
; Add 2001 NaNs to records
lineur=[lineur,!values.f_nan]
liseur=[liseur,!values.f_nan]
filter_cru,thalf,/nan,tsin=lineur,tslow=neurlow,tshigh=neurhi
filter_cru,thalf,/nan,tsin=liseur,tslow=seurlow,tshigh=seurhi
;ml=where(liyear lt 1750.)
;neurlow(ml)=!values.f_nan
;seurlow(ml)=!values.f_nan
;
; Read in the calibrated series
;
restore,filename=maints+'temp_calibrated.allversion.idlsave'
; Gets: nyr,nreg,calregts,regname,timey,tempregts,tempnyr,temptimey
;
; Prepare for plotting
;
loadct,39
case dopub of
  0: multi_plot,nrow=9,layout='large'
  1: multi_plot,nrow=5,ncol=2,layout='large'
  2: multi_plot,nrow=4,layout='large'
endcase
if !d.name eq 'X' then begin
  window,ysize=850
  !p.font=-1
endif else begin
  !p.font=0
  device,/helvetica,/bold,font_size=12
endelse
!y.margin=[0,0]
!y.omargin=[4,4]
if dopub eq 0 then begin
  !x.margin=[10,4]
  !x.omargin=[0,0]
endif else begin
  !x.margin=[0,0]
  !x.omargin=[10,4]
endelse
if docol eq 1 then begin
  def_1color,20,color='lsand'
  def_1color,21,color='lyellow'
  def_1color,22,color='red'
  def_1color,23,color='mblue'
  def_1color,24,color='lgreen'
endif else begin
  def_1color,20,color='mgrey'
  def_1color,21,color='palegrey'
  def_1color,22,color='black'
  def_1color,23,color='black'
  def_1color,24,color='black'
endelse
;
; Repeat for each region
;
partno='('+['a','b','c','d','e','f','g','h','i','j']+') '
if dopub eq 2 then streg=nreg-1 else streg=0
for ireg = streg , nreg-1 do begin
  ;
  if dopub ne 1 then begin
    if ((ireg mod 2) eq 0) then ytit='' $
                           else ytit='Sliding correlations'
    ytf=''
  endif else begin
    if ((ireg mod 2) eq 0) then begin
      ytit='Sliding correlations'
      ytf=''
    endif else begin
      ytit=''
      ytf='nolabels'
    endelse
  endelse
  ;ypos=-0.4
  ypos=-0.25
  ;
  pause
  plot,[0,1],/nodata,$
    /xstyle,xrange=xran,xtickformat='nolabels',$
    /ystyle,yrange=[-0.49,0.99],ytickformat=ytf,$
    ytitle=ytit
  ;
  tem1=reform(tempregts(*,ireg))
  filter_cru,thalf,/nan,tsin=tem1,tslow=temlow,tshigh=temhi
  tyr=temptimey
  ;
  cal1=reform(calregts(*,ireg))
  filter_cru,thalf,/nan,tsin=cal1,tslow=callow,tshigh=calhi
  ;
  ; Now compute the running correlations
  ;
  minyr=tyr(0)
  maxyr=timey(nyr-1)
  elecal=where(timey eq minyr)
  elecal=elecal(0)
  nover=maxyr-minyr+1
  ;
  ; High frequency first
  ;
  nwin=nover-hiwin+1
  xhi=findgen(nwin)+minyr+float(hiwin-1)/2.
  rhi=fltarr(nwin)
  rhi(*)=!values.f_nan
  for iwin = 0 , nwin-1 do begin
    onex=calhi(elecal+iwin:elecal+iwin+hiwin-1)
    oney=temhi(iwin:iwin+hiwin-1)
    rhi(iwin)=correlate(onex,oney)
  endfor
  ;
  ; Low frequency first
  ;
if dopub eq 0 then begin
  nwin=nover-lowin+1
  xlo=findgen(nwin)+minyr+float(lowin-1)/2.
  rlo=fltarr(nwin)
  rlo(*)=!values.f_nan
  for iwin = 0 , nwin-1 do begin
    onex=callow(elecal+iwin:elecal+iwin+lowin-1)
    oney=temlow(iwin:iwin+lowin-1)
    rlo(iwin)=correlate(onex,oney)
  endfor
endif
  rhimean=total(rhi,/nan)
  rhinum=total(finite(rhi))
  rhimean=rhimean/float(rhinum)
  print,rhimean
  polyfill,[!x.crange,reverse(!x.crange)],$
    [replicate(rrmin(ireg),2),replicate(rrmax(ireg),2)],color=20
  axis,yaxis=0,/ystyle,ytickformat='nolabels'
  axis,yaxis=1,/ystyle,ytickformat='nolabels'
  ;
  if dopub ne 2 then begin
    if !d.name eq 'PS' then device,font_size=9
    xyouts,1910+25*dopub,ypos,partno(ireg)+regname(ireg)
    if !d.name eq 'PS' then device,font_size=12
  endif
  ;
  if ireg le dopub then axis,xaxis=1,xtitle='Year',/xstyle
  if ireg ge nreg-1-dopub then axis,xaxis=0,xtitle='Year',/xstyle
  ;
  oplot,!x.crange,[0.,0.],linestyle=1
  ;
  oplot,xhi,rhi,thick=2+2*dopub
  oplot,!x.crange,[rhimean,rhimean]
  oplot,!x.crange,replicate(rrmin(ireg),2)
  oplot,!x.crange,replicate(rrmax(ireg),2)
  if dopub eq 0 then oplot,xlo,rlo,thick=5
  ;
  ; Now correlate with the long European instrumental series if required
  ;
  if (dolongeuro ne 0) and (ireg le 1) then begin
    if ireg eq 0 then begin
      temlow=neurlow
      temhi=neurhi
    endif else begin
      temlow=seurlow
      temhi=seurhi
    endelse
    tyr=liyear
    ; Now compute the running correlations
    minyr=tyr(0)
    maxyr=timey(nyr-1)
    elecal=where(timey eq minyr)
    elecal=elecal(0)
    nover=maxyr-minyr+1
    ; High frequency first
    nwin=nover-hiwin+1
    xhi=findgen(nwin)+minyr+float(hiwin-1)/2.
    rhi=fltarr(nwin)
    rhi(*)=!values.f_nan
    for iwin = 0 , nwin-1 do begin
      onex=calhi(elecal+iwin:elecal+iwin+hiwin-1)
      oney=temhi(iwin:iwin+hiwin-1)
      rhi(iwin)=correlate(onex,oney)
    endfor
    ; Low frequency first
    nwin=nover-lowin+1
    xlo=findgen(nwin)+minyr+float(lowin-1)/2.
    rlo=fltarr(nwin)
    rlo(*)=!values.f_nan
    for iwin = 0 , nwin-1 do begin
      onex=callow(elecal+iwin:elecal+iwin+lowin-1)
      oney=temlow(iwin:iwin+lowin-1)
      rlo(iwin)=correlate(onex,oney)
    endfor
    oplot,xhi,rhi,thick=4,color=20
    oplot,xlo,rlo,thick=8,color=20
  endif
  ;
endfor
;
end