; ; PLOTS 'ALL' REGION MXD timeseries from age banded and from hugershoff ; standardised datasets. ; Reads Harry's regional timeseries and outputs the 1600-1992 portion ; with missing values set appropriately. Uses mxd, and just the ; "all band" timeseries ;****** APPLIES A VERY ARTIFICIAL CORRECTION FOR DECLINE********* ; yrloc=[1400,findgen(19)*5.+1904] valadj=[0.,0.,0.,0.,0.,-0.1,-0.25,-0.3,0.,-0.1,0.3,0.8,1.2,1.7,2.5,2.6,2.6,$ 2.6,2.6,2.6]*0.75 ; fudge factor if n_elements(yrloc) ne n_elements(valadj) then message,'Oooops!' ; loadct,39 def_1color,20,color='red' plot,[0,1] multi_plot,nrow=4,layout='large' if !d.name eq 'X' then begin window, ysize=800 !p.font=-1 endif else begin !p.font=0 device,/helvetica,/bold,font_size=18 endelse ; ; Get regional tree lists and rbar ; restore,filename='reglists.idlsave' harryfn=['nwcan','wnam','cecan','nweur','sweur','nsib','csib','tib',$ 'esib','allsites'] ; rawdat=fltarr(4,2000) for i = nreg-1 , nreg-1 do begin fn='mxd.'+harryfn(i)+'.pa.mean.dat' print,fn openr,1,fn readf,1,rawdat close,1 ; densadj=reform(rawdat(2:3,*)) ml=where(densadj eq -99.999,nmiss) densadj(ml)=!values.f_nan ; x=reform(rawdat(0,*)) kl=where((x ge 1400) and (x le 1992)) x=x(kl) densall=densadj(1,kl) ; all bands densadj=densadj(0,kl) ; 2-6 bands ; ; Now normalise w.r.t. 1881-1960 ; mknormal,densadj,x,refperiod=[1881,1960],refmean=refmean,refsd=refsd mknormal,densall,x,refperiod=[1881,1960],refmean=refmean,refsd=refsd ; ; APPLY ARTIFICIAL CORRECTION ; yearlyadj=interpol(valadj,yrloc,x) densall=densall+yearlyadj ; ; Now plot them ; filter_cru,20,tsin=densall,tslow=tslow,/nan cpl_barts,x,densall,title='Age-banded MXD from all sites',$ xrange=[1399.5,1994.5],xtitle='Year',/xstyle,$ zeroline=tslow,yrange=[-7,3] oplot,x,tslow,thick=3 oplot,!x.crange,[0.,0.],linestyle=1 ; endfor ; ; Restore the Hugershoff NHD1 (see Nature paper 2) ; xband=x restore,filename='../tree5/densadj_MEAN.idlsave' ; gets: x,densadj,n,neff ; ; Extract the post 1600 part ; kl=where(x ge 1400) x=x(kl) densadj=densadj(kl) ; ; APPLY ARTIFICIAL CORRECTION ; yearlyadj=interpol(valadj,yrloc,x) densadj=densadj+yearlyadj ; ; Now plot it too ; filter_cru,20,tsin=densadj,tslow=tshug,/nan cpl_barts,x,densadj,title='Hugershoff-standardised MXD from all sites',$ xrange=[1399.5,1994.5],xtitle='Year',/xstyle,$ zeroline=tshug,yrange=[-7,3],bar_color=20 oplot,x,tshug,thick=3,color=20 oplot,!x.crange,[0.,0.],linestyle=1 ; ; Now overplot their bidecadal components ; plot,xband,tslow,$ xrange=[1399.5,1994.5],xtitle='Year',/xstyle,$ yrange=[-6,2],thick=3,title='Low-pass (20-yr) filtered comparison' oplot,x,tshug,thick=3,color=20 oplot,!x.crange,[0.,0.],linestyle=1 ; ; Now overplot their 50-yr components ; filter_cru,50,tsin=densadj,tslow=tshug,/nan filter_cru,50,tsin=densall,tslow=tslow,/nan plot,xband,tslow,$ xrange=[1399.5,1994.5],xtitle='Year',/xstyle,$ yrange=[-6,2],thick=3,title='Low-pass (50-yr) filtered comparison' oplot,x,tshug,thick=3,color=20 oplot,!x.crange,[0.,0.],linestyle=1 ; ; Now compute the full, high and low pass correlations between the two ; series ; perst=1400. peren=1992. ; openw,1,'corr_age2hug.out' thalf=[10.,30.,50.,100.] ntry=n_elements(thalf) printf,1,'Correlations between timeseries' printf,1,'Age-banded vs. Hugershoff-standardised' printf,1,' Region Full <10 >10 >30 >50 >100' ; kla=where((xband ge perst) and (xband le peren)) klh=where((x ge perst) and (x le peren)) ts1=densadj(klh) ts2=densall(kla) ; r1=correlate(ts1,ts2) rall=fltarr(ntry) for i = 0 , ntry-1 do begin filter_cru,thalf(i),tsin=ts1,tslow=tslow1,tshigh=tshi1,/nan filter_cru,thalf(i),tsin=ts2,tslow=tslow2,tshigh=tshi2,/nan if i eq 0 then r2=correlate(tshi1,tshi2) rall(i)=correlate(tslow1,tslow2) endfor ; printf,1,'ALL SITES',r1,r2,rall,$ format='(A11,2X,6F6.2)' ; printf,1,' ' printf,1,'Correlations carried out over the period ',perst,peren ; close,1 ; end