; ; Plots distributions and cutoffs of trends and means of different length ; segments. ; ; Get pre-analysed data ; restore,filename='trendmean_distrib.idlsave' ; ntry,tlen,nts,alltit,tren95,mean95,binmean,bintren,freqmean,freqtren,$ ; nbin,binsc ; ; Get observed 20th century data ; restore,filename='obs_ts.idlsave' ; hccits,timey,nts,nyr obsts=hccits obsyr=timey obsnyr=nyr ; ; Reference the observed data to the same 'pre-industrial' mean ; ; (a) adjust series of those whose palaeo series were not referenced to 61-90 x=reform(obsts(*,1)) mkanomaly,x,obsyr,refperiod=[1902,1980] obsts(*,1)=x(*) ; (b) adjust to pre-20th century mean refmean=[-0.5518,-0.1709,-0.3272,-0.6815,-0.7504,-0.0401,$ 9.0556,12.9851,5.1216] for its = 0 , nts-1 do obsts(*,its)=obsts(*,its)-refmean(its) ; ; Get perturbed run ensemble mean data ; restore,filename='gsax_ts.idlsave' ; hccits,timey,nts,nyr ; ; Reference model data to the same 'pre-industrial' (control-run) mean ; refmean=[293.073,287.785,289.273,277.088,277.939,277.607,$ 282.143,285.193,279.093]-273.15 for its = 0 , nts-1 do hccits(*,its)=hccits(*,its)-refmean(its) ; ; Define end-points of perturbed trends/means ; ep=[1997,1998] ; in fact ep(>1) are starting points! nep=n_elements(ep) ; ; Compute trends/means over various length segments, ending at end-points ; ghgtren=fltarr(ntry,nts,nep)*!values.f_nan ghgmean=fltarr(ntry,nts,nep)*!values.f_nan obstren=fltarr(ntry,nts)*!values.f_nan obsmean=fltarr(ntry,nts)*!values.f_nan for its = 0 , nts-1 do begin ; do observed trends ending in 1997 iep=0 iend=where(obsyr eq ep(iep)) iend=iend(0) for itry = 0 , ntry-1 do begin t1=tlen(itry) ist=iend-t1+1 if ist gt 0 then begin x=obsyr(ist:iend) y=obsts(ist:iend,its) acoeff=linfit(x,y) obstren(itry,its)=acoeff(1)*float(t1) obsmean(itry,its)=total(y)/float(t1) endif endfor endfor for its = 0 , nts-1 do begin ; do simulated trends ending in 1997 iep=0 iend=where(timey eq ep(iep)) iend=iend(0) for itry = 0 , ntry-1 do begin t1=tlen(itry) ist=iend-t1+1 if ist gt 0 then begin x=timey(ist:iend) y=hccits(ist:iend,its) acoeff=linfit(x,y) ghgtren(itry,its,iep)=acoeff(1)*float(t1) ghgmean(itry,its,iep)=total(y)/float(t1) endif endfor endfor for its = 0 , nts-1 do begin ; do simulated trends starting in 1998 for iep = 1 , nep-1 do begin ist=where(timey eq ep(iep)) ist=ist(0) for itry = 0 , ntry-1 do begin t1=tlen(itry) iend=ist+t1-1 x=timey(ist:iend) y=hccits(ist:iend,its) acoeff=linfit(x,y) ghgtren(itry,its,iep)=acoeff(1)*float(t1) ghgmean(itry,its,iep)=total(y)/float(t1) endfor endfor endfor ; ; Prepare window/paper for plotting ; loadct,39 def_1color,20,color='red' def_1color,21,color='blue' def_1color,22,color='purple' multi_plot,nrow=3 if !d.name eq 'X' then begin window, ysize=800 !p.font=-1 endif else begin !p.font=0 device,/helvetica,/bold,font_size=12 endelse ; ; Select which distributions to show ; iuse=[2,6,12] nuse=n_elements(iuse) ; ; Compute distribution totals for scaling by ; fmtot=total(freqmean,1)/100. fttot=total(freqtren,1)/100. ; ; Repeat for each time series ; for its = 2 , 5 , 3 do begin ; ; plot,binmean*binsc(its),freqmean(*,iuse(0),its,0),/nodata,$ ; /xstyle,xtitle='Mean anomaly (!Uo!NC)',$ ; /ystyle,ytitle='Frequency of occurrence!C!C(%)',yrange=[0.,50.] ; oplot,[0.,0.],!y.crange,linestyle=1 ; ; for j = 0 , nuse-1 do begin ; oplot,binmean*binsc(its),freqmean(*,iuse(j),its,0)/fmtot(iuse(j),its,0),$ ; thick=4 ; oplot,binmean*binsc(its),freqmean(*,iuse(j),its,1)/fmtot(iuse(j),its,1) ; endfor ; ; xyouts,/normal,0.5,0.99,align=0.5,alltit(its) ; ; plot,bintren*binsc(its),freqtren(*,iuse(0),its,0),/nodata,$ ; /xstyle,xtitle='Trend (!Uo!NC per segment length)',$ ; /ystyle,ytitle='Frequency of occurrence!C!C(%)',yrange=[0.,50.] ; oplot,[0.,0.],!y.crange,linestyle=1 ; ; for j = 0 , nuse-1 do begin ; oplot,bintren*binsc(its),freqtren(*,iuse(j),its,0)/fttot(iuse(j),its,0),$ ; thick=4 ; oplot,bintren*binsc(its),freqtren(*,iuse(j),its,1)/fttot(iuse(j),its,1) ; endfor ; if its le 2 then yr=[0.,1.1] else yr=[-0.3,1.9] plot,tlen,mean95(*,its,0),thick=4,psym=-def_sym(10),$ /xstyle,xtitle='Segment length (years)',$ ytitle='95th centile of segment means!C!C(!Uo!NC)',$ /ystyle,yrange=yr oplot,tlen,mean95(*,its,1),psym=-def_sym(15) oplot,!x.crange,[0.,0.],linestyle=1 ; oplot,tlen,ghgmean(*,its,0),color=20,thick=4 oplot,tlen,ghgmean(*,its,1),color=22,thick=4 ; oplot,tlen,ghgmean(*,its,2),color=20,thick=5 ; oplot,tlen,obsmean(*,its),color=21,thick=4 ; ; if its le 2 then yr=[-0.1,1.6] else yr=[-0.5,3.] ; plot,tlen,30.*tren95(*,its,0)/tlen,thick=4,psym=-def_sym(10),$ ; /xstyle,xrange=[10,100],xtitle='Segment length (years)',$ ; ytitle='95th centile of segment trends!C!C(!Uo!NC per 30 yr)',$ ; /ystyle,yrange=yr ; oplot,tlen,30.*tren95(*,its,1)/tlen,psym=-def_sym(15) ; oplot,!x.crange,[0.,0.],linestyle=1 ; ; oplot,tlen,30.*ghgtren(*,its,0)/tlen,color=20,thick=4 ; oplot,tlen,30.*ghgtren(*,its,1)/tlen,color=22,thick=4 ; oplot,tlen,30.*ghgtren(*,its,2)/tlen,color=20,thick=5 ; ; oplot,tlen,30.*obstren(*,its)/tlen,color=21,thick=4 ; xyouts,/normal,0.5,0.33+0.33*!p.multi(0),align=0.5,alltit(its) ; endfor ; end