;
; Compares PCR-reconstructed series using separate segments or using
; just the overall regression model with zero's filled in for missing
; predictors.
;
regname='ALL'
thalf=10.
;
; Read in the different PCR models
;
nver=10
for iver = 0 , nver-1 do begin
  ;
  fn='bandtemp'+regname+string(iver,format='(I1)')+'_fixed1950_calpcr.idlsave'
  print,fn
  restore,filename=fn
  ; Gets: nyr,nhtit,yrmxd,prednh,fullnh
  if nhtit ne regname then message,'Re-calibrate the age-banded!'
  ;
  if iver eq 0 then begin
    allnyr=nyr
    alltimey=yrmxd
    allts=fltarr(nyr,nver)
    cmpts=fullnh
  endif
  allts(*,iver)=prednh(*)
  ;
endfor
;
; Make a combined reconstruction using different models for different periods
; Some models are used for the early period, some for the late period
;
kearly=where(alltimey le 1960,nearly)
klate=where(alltimey ge 1960,nlate)
segts=allts(*,0)
for iver = 1 , 6 do begin
  filllist=where(finite(segts(kearly)) eq 0,nfill)
  gotlist=where(finite(segts(kearly)),ngot)
  print,iver-1,nearly,ngot,nfill
  if nfill gt 0 then segts(kearly(filllist))=allts(kearly(filllist),iver)
  if ngot gt 0 then allts(kearly(gotlist),iver)=!values.f_nan
  allts(klate,iver)=!values.f_nan
endfor
for iver = 7 , 9 do begin
  filllist=where(finite(segts(klate)) eq 0,nfill)
  gotlist=where(finite(segts(klate)),ngot)
  print,iver-1,nlate,ngot,nfill
  if nfill gt 0 then segts(klate(filllist))=allts(klate(filllist),iver)
  if ngot gt 0 then allts(klate(gotlist),iver)=!values.f_nan
  allts(kearly,iver)=!values.f_nan
endfor
;
loadct,39
multi_plot,nrow=3,layout='large'
if !d.name eq 'X' then begin
  window,ysize=900
  !p.font=-1
endif else begin
  !p.font=0
  device,/helvetica,bold=0,font_size=14
endelse
def_1color,10,color='black'
def_1color,11,color='red'
;
plot,alltimey,cmpts,$
  title=nhtit,$
  /xstyle,xrange=[1400,2000],xtitle='Year',$
  /ystyle,yrange=[-1.75,0.6],ytitle='Temperature anomaly  (!Uo!NC wrt 1961-90)'
oplot,!x.crange,[0.,0.],linestyle=1
;oplot,alltimey,segts,color=15
filter_cru,thalf,/nan,tsin=cmpts,tslow=ylow
oplot,alltimey,ylow,thick=3
filter_cru,thalf,/nan,tsin=segts,tslow=ylow
oplot,alltimey,ylow,color=11,thick=3
;
plot,alltimey,segts,$
  title=nhtit,$
  /xstyle,xrange=[1400,2000],xtitle='Year',$
  /ystyle,yrange=[-1.75,0.6],ytitle='Temperature anomaly  (!Uo!NC wrt 1961-90)'
oplot,!x.crange,[0.,0.],linestyle=1
filter_cru,thalf,/nan,tsin=segts,tslow=ylow
oplot,alltimey,ylow,thick=3,color=11
;
plot,alltimey,cmpts,/nodata,$
  title=nhtit,$
  /xstyle,xrange=[1400,2000],xtitle='Year',$
  /ystyle,yrange=[-1.75,0.3],ytitle='Temperature anomaly  (!Uo!NC wrt 1961-90)'
oplot,!x.crange,[0.,0.],linestyle=1
for iver = 0 , nver-1 do begin
  oplot,alltimey,allts(*,iver),color=(iver mod 2)+10
endfor
;
; Save the composite, multi-model reconstruction for later use
;
;yrmxd=alltimey
;prednh=segts
;nhtit=regname
;nyr=allnyr
;save,filename='bandtemp'+regname+'_calmultipcr.idlsave',$
;  yrmxd,prednh,nhtit,nyr
;
; Also output the series for MATLAB spectral analysis
; Stop in 1994 cos the Hugershoff doesn't have 1995 in it
; Start in 1402 cos 1400-1 are missing in age-banded multi-PCR
;
;if regname eq 'ALL' then begin
;  openw,1,'ts2.matlab_mxdband_all'
;  printf,1,prednh(2:nyr-2),format='(F7.2)'
;  close,1
;endif
;
;
; Output the final time series to an ASCII file too
;
;openw,1,regname+'temp_agebandbriffa.dat'
;if regname eq 'ALL' then begin
;  printf,1,'Calibrated against observed Apr-Sep temperature over 1881-1960'
;  printf,1,'averaged over all grid boxes containing tree-ring chronologies'
;  printf,1,'(i.e., co-located temperatures).'
;endif else begin
;  printf,1,'Calibrated against observed Apr-Sep temperature over 1881-1960'
;  printf,1,'averaged over all land grid boxes north of 20N'
;  printf,1
;endelse
;printf,1
;printf,1,'Year  Reconstructed temperature anomaly (degrees C wrt 1961-90)
;for i = 0 , nyr-1 do printf,1,yrmxd(i),prednh(i),format='(I4,F10.3)'
;close,1
;
end
