; ; Creates MEAN timeseries from adjusted regional timeseries and adjusts the ; MEAN for changing sample size using a time-dependent rbar. ; Timeseries are weighted by their EPS when computing the mean. ; doland=0 ; 0=portrait, 1=landscape ; ploton,1,/color8 loadct,39 def_1color,50,color='red' def_1color,51,color='deepblue' def_1color,52,color='green' def_1color,53,color='purple' def_1color,54,color='orange' if doland eq 0 then multi_plot,nrow=3,layout='large' $ else multi_plot,nrow=1,/landscape if !d.name eq 'X' then begin window,ysize=900 !p.font=-1 initx endif else begin !p.font=0 if doland eq 0 then device,/helvetica,/bold,font_size=12,xsize=17.9 $ else device,/helvetica,/bold,font_size=12,ysize=12.,$ xoffset=6. endelse ; ; Get adjusted regional timeseries ; restore,filename='reglists.idlsave' nreg=nreg-2 ; not ARCTIC or EXTRA-ARCTIC ; for i = 3 , nreg-1 do begin ; not 0, 1 & 2 (ALL, NORTH & SOUTH) ; ; Get adjusted regional timseries ; restore,filename='densadj_'+regname(i)+'.idlsave' ; if i eq 3 then begin nyr=n_elements(x) ndata=nreg-3 ; not 0, 1 & 2 allmxd=fltarr(nyr,ndata) alleps=fltarr(nyr,ndata) endif allmxd(*,i-3)=densadj(*) alleps(*,i-3)=denseps(*) ; endfor ; ; Get grand correlation matrix between regions ; restore,filename='rbartd_dens_MEAN.idlsave' ; ; Make a time-dependent rbar next ; print,'Time-dependent rbar for ',nyr,' years' rbar_td=fltarr(nyr)*!values.f_nan for it = 0 , nyr-1 do begin print,it,format='($,I4)' mask = where( finite(allmxd(it,*)) , nsamp) if nsamp gt 0 then begin if nsamp eq 1 then begin rbar_td(it)=1. ; if n=1 then n'=1 regardless of rbar endif else begin ; get reduced grandr redr=grandr(mask,*) redr=redr(*,mask) ; Total it and average it, but only do i<>j. totr=0. totn=0. for j = 1 , nsamp-1 do begin for k = 0 , j-1 do begin totr=totr+redr(j,k) ; should have no missing values I hope totn=totn+1. if redr(j,k) eq 1. then message,'Ooops!' endfor endfor rbar_td(it)=totr/totn endelse endif endfor print ; ; Now compute a weighted mean of the available regional timeseries ; n=float(total(finite(allmxd),2)) ; compute timeseries of number of regions totmxd=total(allmxd*alleps,2,/nan) toteps=total(alleps,2,/nan) meanmxd=totmxd/toteps ; ; Now need to adjust for varying sample size ; xadj=mkeffective(meanmxd,n,rbar=rbar_td,neff=neff) meanmxd=xadj ; ; Now normalise w.r.t. 1881-1960 ; mknormal,meanmxd,x,refperiod=[1881,1960] ; ; Now plot them ; filter_cru,20,tsin=meanmxd,tslow=tslow,/nan cpl_barts,x,meanmxd,$ ;title='Normalised adjusted MEAN of all regional MXD timeseries',$ xtitle='Year',/xstyle,$ ytitle=' Normalised mean density anomaly (standard deviation units)',$ yminor=4,yticks=4,ytickv=findgen(5)*2.-6.,$ zeroline=tslow,yrange=[-9,3],ystyle=9,bar_color=[50,51] oplot,x,tslow,thick=2 oplot,!x.crange,[0.,0.] ; axis,xaxis=1,/xstyle ;,xtitle='Year' axis,xaxis=0,/xstyle,xtickformat='nolabels' ; ;xyouts,1420.,1.9,/data,'NHD1',size=1.5,color=54 ; yrlab=['1448','1453','1495','1544','1587','1601','1641/2/3','1666',$ '1675','1695','1698/9','1740','1783','1816/7/8','1836/7','1884','1912'] yry =[-3.1 ,-4.85 ,-3.03 ,-2.63 ,-3.32 ,-6.0 ,-5. ,-3.42 ,$ -3.81 ,-4.15 ,-3.52 ,-3.35 ,-3.05 ,-5.03 ,-3.5 ,-3.6 ,-4.1] yrx =[1440.5,1453 ,1495 ,1544 ,1587 ,1589 ,1642 ,1658 ,$ 1675 ,1695 ,1711 ,1740 ,1783 ,1817 ,1837 ,1884 ,1912] nlab=n_elements(yrlab) for i = 0 , nlab-1 do xyouts,yrx(i),yry(i),/data,yrlab(i),size=0.6,align=0.5 ; keepcol=!p.color !p.color=53 yrvei=[1450,1452,1471,1477,1480,1482,1580,1586,1593,1600,1640,1641,1660,$ 1663,1667,1673,1680,1707,1739,1800,1815,1835,1853,1854,1883,1886,1902,$ 1907,1912,1932,1956,1980,1982,1991] yrsiz=[5,6,5,5,5,5,6,5,5,6,5,6,6,5,5,5,5,5,5,5,7,5,5,5,6,5,6,5,6,5,5,5,5,6] nvei=n_elements(yrvei) ym=!y.crange(0) for i = 0 , nvei-1 do begin z=yrvei(i) y=yrsiz(i) case y of 5: plots,[z,z,z-1,z,z+1,z],[ym,ym+0.4,ym+0.4,ym+0.7,ym+0.4,ym+0.4] 6: begin plots,[z,z],[ym,ym+1.1],thick=3 polyfill,[z-2.,z,z+2.,z-2.],[ym+1.1,ym+1.6,ym+1.1,ym+1.1] end 7: begin plots,[z,z],[ym,ym+2.2],thick=5 polyfill,[z-3.,z,z+3.,z-3.],[ym+2.2,ym+2.7,ym+2.2,ym+2.2] end endcase endfor ; xyouts,1449.,ym+0.3,align=1.,'?',size=0.6 xyouts,1579.,ym+0.3,align=1.,'?',size=0.6 xyouts,1659.,ym+0.3,align=1.,'?',size=0.6 ; xyouts,1420.,ym+0.3,orient=90.,'VEI',size=0.9 plots,[1432.,1432.],[ym,ym+2.7],thick=3 plots,[1432.,1438.],[ym+2.7,ym+2.7],thick=3 plots,[1432.,!x.crange(1)],[ym+2.7,ym+2.7],linestyle=1 plots,[1432.,1438.],[ym+1.6,ym+1.6],thick=3 plots,[1432.,!x.crange(1)],[ym+1.6,ym+1.6],linestyle=1 plots,[1432.,1438.],[ym+0.7,ym+0.7],thick=3 plots,[1432.,!x.crange(1)],[ym+0.7,ym+0.7],linestyle=1 xyouts,1431.,ym+2.5,'7',size=0.6,align=1. xyouts,1431.,ym+1.4,'6',size=0.6,align=1. xyouts,1431.,ym+0.5,'5',size=0.6,align=1. !p.color=keepcol ; axis,yaxis=1,ystyle=1,yrange=!y.crange*0.116626,$ ;-0.105023,$ yminor=4,yticks=5,ytickv=findgen(6)*0.2-0.8,$ ytitle=' Estimated Northern Hemisphere temperature anomaly (!Uo!NC)' ; ;plot,x,n,ytitle='Sample size',/xstyle ; ;plot,x,rbar_td,/xstyle,ytitle='rbar' ; ;plot,x,neff,/xstyle,ytitle='Effective sample size' ; ; Now save the weighted mean timeseries, for further analysis ; densadj=meanmxd save,filename='densadj_MEAN.idlsave',x,densadj,n,neff ; ; Also plot the fixed grid correlations below ; xsk=[1400,1450,1500,1550,1600,1620,1640,1660,1680,1700,1720,1740,1760,1780,$ 1800,1820,1840,1860,1880,1901,1970,1980,1990,1994] ysk=[0.32,0.37,0.33,0.34,0.35,0.38,0.39,0.40,0.46,0.46,0.48,0.50,0.52,0.54,$ 0.55,0.55,0.55,0.55,0.55,0.57,0.57,0.51,0.39,0.29] ;plot,xsk,ysk,/xstyle,xtitle='Year of fixed grid',$ ; /ystyle,yrange=[0.,1.],ytitle='Fixed-grid correlation over 1881-1960' ;ysk=[0.51,0.65,0.60,0.62,0.61,0.61,0.60,0.62,0.68,0.67,0.69,0.71,0.73,0.74,$ ; 0.76,0.77,0.77,0.77,0.78,0.80,0.80,0.74,0.58,0.40] ;oplot,xsk,ysk,thick=5 ;oplot,!x.crange,[0.2,0.2],linestyle=1 ;oplot,!x.crange,[0.4,0.4],linestyle=1 ;oplot,!x.crange,[0.6,0.6],linestyle=1 ;oplot,!x.crange,[0.8,0.8],linestyle=1 ; plotoff,/gv ; end