; ; Creates a regional timeseries and adjust it for changing sample size ; plot,[0,1] multi_plot,nrow=6 if !d.name eq 'X' then begin window,ysize=900 initx endif else begin device,xoffset=2,xsize=17 endelse ; ncid=ncdf_open('longvolc.nc') ncdf_diminq,ncid,'time',dummy,ntime ncdf_diminq,ncid,'station',dummy,nstat ncdf_varget,ncid,'year',x ncdf_varget,ncid,'density',density ncdf_attget,ncid,'density','missing_value',valmiss ncdf_varget,ncid,'fraction',weight ncdf_close,ncid ; misslist=where(density eq valmiss,nmiss) density(misslist)=!values.f_nan wt=weight(*,*) ; nms=['Polar Urals','Tornetrask','Jasper'] for i = 0 , 2 do begin cpl_barts,x,density(*,i),title='longvolc record for '+nms(i),$ xtitle='Year',/xstyle,$ zeroline=0.,yrange=[-4,4] plot,x,wt(*,i),title='effective sample size for '+nms(i),$ xtitle='Year',/xstyle endfor ; restore,filename='rbar_longvolc.idlsave' ; dens=density(*,*) wt=weight(*,*) ; n=total(dens*0.+1.,2,/nan) ; compute timeseries of number of chronols ncore=total(dens*0.+wt,2,/nan) ; compute timeseries of summed core fraction totwdens=total(dens*wt,2,/nan) totudens=total(dens,2,/nan) totwvals=total(float(finite(dens))*wt,2) totuvals=total(float(finite(dens)),2) wdens=totwdens/float(totwvals) ; this is the weighted anomaly udens=totudens/float(totuvals) ; this is the unweighted anomaly ; pause cpl_barts,x,wdens,title='longvolc record: simple weighted mean',$ xtitle='Year',/xstyle,$ zeroline=0.,yrange=[-4,4] cpl_barts,x,udens,title='longvolc record: simple mean',$ xtitle='Year',/xstyle,$ zeroline=0.,yrange=[-4,4] ; ; Now need to adjust for varying sample size ; xadj=mkeffective(udens,n,rbar=allrbar,neff=neff,/local) keepold=udens udens=xadj ; ; Now normalise w.r.t. 1901-1970 ; mknormal,keepold,x,refperiod=[1901,1970],refmean=refmean,refsd=refsd mknormal,udens,x,refperiod=[1901,1970],refmean=refmean,refsd=refsd print,refmean,refsd ; ; Now plot them ; ;filter_cru,20,tsin=keepold,tslow=tslow,/nan cpl_barts,x,keepold,title='Mean record normalised over 1901-70',$ xtitle='Year',/xstyle,$ zeroline=tslow,yrange=[-6,6] ; oplot,x,tslow,thick=2 oplot,!x.crange,[0.,0.] oplot,!x.crange,[-2.,-2.],linestyle=1 oplot,!x.crange,[2.,2.],linestyle=1 ; ;filter_cru,20,tsin=udens,tslow=tslow,/nan cpl_barts,x,udens,title='Mean record adjusted for sample size and then normalised over 1901-70',$ xtitle='Year',/xstyle,$ zeroline=0.,yrange=[-6,6] ;oplot,x,tslow,thick=2 ;moms=moment(udens(where(finite(udens))),sdev=sdev) oplot,!x.crange,[0.,0.] oplot,!x.crange,[-2.,-2.],linestyle=1 oplot,!x.crange,[2.,2.],linestyle=1 ; plot,x,n,title='Number of chronologies',yrange=[-0.5,3.5],/xstyle plot,x,neff,title='Effective number of independent chronologies',$ yrange=[-0.5,3.5],/xstyle ; ; Now save the weighted mean timeseries, for further analysis ; densadj=udens save,filename='longvolcadj.idlsave',x,densadj,n,neff,ncore ; end