; ; Reads in Mauna Loa and Barrow CO2 observations ; headst=strarr(10) footst=strarr(5) yrst=1958 yren=1995 nyr=yren-yrst+1 nmon=12 nt=nyr*nmon co2dat=fltarr(2,nt)*!values.f_nan timey=indgen(nt)/nmon + yrst timem=indgen(nt) mod nmon xtime=float(timey)+float(timem)/12. ; openr,1,'mauna.dat' ; read in Mauna Loa readf,1,headst print,headst nlin=38 nval=nlin*12 rawdat=fltarr(14,nlin) readf,1,rawdat,format='(I4,13F7.2)' co2dat(0,nt-nval:nt-1)=reform(rawdat(1:12,*),12*nlin) readf,1,footst ; read in Barrow ;readf,1,headst ;print,headst ;nlin=22 ;nval=nlin*12 ;rawdat=fltarr(14,nlin) ;readf,1,rawdat,format='(I4,13F7.2)' ;co2dat(1,nt-nval:nt-1)=reform(rawdat(1:12,*),12*nlin) ; close,1 ; ; get Barrow from other file ; openr,1,'co2.txt' headst=strarr(301) readf,1,headst print,headst(297:300) nval=272 rawdat=fltarr(nval) readf,1,rawdat,format='(19X,F11.2)' co2dat(1,nt-24-nval:nt-1-24)=rawdat close,1 ; misslist=where(co2dat eq -99.99,nmiss) if nmiss gt 0 then co2dat(misslist)=!values.f_nan ; ;plot full timeseries ; multi_plot,nrow=3,layout='large' if !d.name eq 'X' then window,ysize=820 ; plot,xtime,co2dat(1,*),/ynozero,thick=1,$ xtitle='Year',ytitle='[CO!D2!N] (ppmv)',$ title='Observed timeseries' oplot,xtime,co2dat(0,*),thick=4 legend,['Mauna Loa','Barrow'],thick=[4,1] ; ; compute and plot mean seasonal cycles (only an approximation, since we ; are ignoring the underlying trend when computing the values) ; co2seas=fltarr(2,nmon)*!values.f_nan for j = 0 , nmon-1 do begin wantlist=where(timem eq j,nwant) for i = 0 , 1 do begin vals=co2dat(i,wantlist) totvals=total(vals,/nan) numvals=total(finite(vals)) print,j,i,numvals if numvals gt 2 then co2seas(i,j)=totvals/float(numvals) endfor endfor plot,indgen(12)+1,co2seas(1,*),thick=1,/ynozero,$ yrange=[330,355],/ystyle,psym=10,$ xrange=[0,13],/xstyle,$ title='Observed mean seasonal cycle (biased by long term trend!!)',$ xtitle='Month',ytitle='[CO!D2!N] (ppmv)' oplot,indgen(12)+1,co2seas(0,*),thick=4,psym=10 legend,['Mauna Loa','Barrow'],thick=[4,1] ; ; compute and plot means, mins and maxs and seasonal amplitude ; (from mean seasonal cycle, it appears reasonable to search for maximum ; during months 0 to 8, and for a minimum during months 5 to 11) ; co2mean=fltarr(2,nyr)*!values.f_nan co2max=fltarr(2,nyr)*!values.f_nan co2min=fltarr(2,nyr)*!values.f_nan for j = 0 , nyr-1 do begin wantlist=where(timey eq j+yrst,nwant) if nwant gt 0 then begin for i = 0 , 1 do begin vals=co2dat(i,wantlist) totvals=total(vals,/nan) numvals=total(finite(vals)) print,j,i,numvals if numvals eq 12 then begin co2mean(i,j)=totvals/float(numvals) co2max(i,j)=max(vals(0:8),/nan) co2min(i,j)=min(vals(5:11),/nan) endif endfor endif endfor co2amp=co2max-co2min ; xx=indgen(nyr)+yrst plot,xx,co2amp(1,*),thick=1,$ yrange=[4,18],/ystyle,$ psym=10,$ title='Drawdown',$ xtitle='Year',ytitle='[CO!D2!N] (ppmv)' oplot,xx,co2amp(0,*),thick=4,psym=10 legend,['Mauna Loa','Barrow'],thick=[4,1] ; ; Now do first differences ; fdmean=co2mean*!values.f_nan fdmean(*,1:nyr-1)=co2mean(*,1:nyr-1)-co2mean(*,0:nyr-2) fdmax=co2max*!values.f_nan fdmax(*,1:nyr-1)=co2max(*,1:nyr-1)-co2max(*,0:nyr-2) fdmin=co2min*!values.f_nan fdmin(*,1:nyr-1)=co2min(*,1:nyr-1)-co2min(*,0:nyr-2) pause plot,xx,fdmean(1,*),thick=1,$ psym=10,$ title='First difference of annual mean values',$ xtitle='Year',ytitle='[CO!D2!N] (ppmv)' oplot,xx,fdmean(0,*),thick=4,psym=10 legend,['Mauna Loa','Barrow'],thick=[4,1] plot,xx,fdmax(1,*),thick=1,$ yrange=[-3,5],/ystyle,$ psym=10,$ title='First difference of annual maximum values',$ xtitle='Year',ytitle='[CO!D2!N] (ppmv)' oplot,xx,fdmax(0,*),thick=4,psym=10 legend,['Mauna Loa','Barrow'],thick=[4,1] plot,xx,fdmin(1,*),thick=1,$ psym=10,$ title='First difference of annual minimum values',$ xtitle='Year',ytitle='[CO!D2!N] (ppmv)' oplot,xx,fdmin(0,*),thick=4,psym=10 legend,['Mauna Loa','Barrow'],thick=[4,1] ; ; Now get regional data for correlating! ; regname='ALL' restore,'densadj_'+regname+'.idlsave' restore,'rwidadj_'+regname+'.idlsave' restore,'instradj_'+regname+'.idlsave' ;refper=[yrst,1990] refper=[1982,1990] kpall=where( (x ge refper(0)) and (x le refper(1)) , nkp) regts=fltarr(6,nkp) regts(0,*)=densadj(kpall) regts(1,*)=rwidadj(kpall) regts(2,*)=instradj(kpall) regname='NORTH' restore,'densadj_'+regname+'.idlsave' restore,'rwidadj_'+regname+'.idlsave' restore,'instradj_'+regname+'.idlsave' kpall=where( (x ge refper(0)) and (x le refper(1)) , nkp) regts(3,*)=densadj(kpall) regts(4,*)=rwidadj(kpall) regts(5,*)=instradj(kpall) ; kpco2=where( (xx ge refper(0)) and (xx le refper(1)) , nkp) co2amp=co2amp(*,kpco2) ; thalf=5. allr=fltarr(6,2) for i = 0 , 5 do begin for j = 0 , 1 do begin cmpls=where( finite(regts(i,*)) and finite(co2amp(j,*)) ) a=regts(i,cmpls) b=co2amp(j,cmpls) filter_cru,thalf,tsin=reform(a),tshigh=aa,/nan filter_cru,thalf,tsin=reform(b),tshigh=bb,/nan allr(i,j)=correlate(aa,bb) endfor endfor print,'ALL: MXD TRW INSTR NORTH: MXD TRW INSTR' print,'MAUNA' print,allr print,'BARROW' ; pause multi_plot,nrow=3,ncol=2,layout='large' !p.multi(4)=1 t=['ALL: MXD','ALL: TRW','ALL: INSTR','NORTH: MXD','NORTH: TRW','NORTH: INSTR'] for i = 0 , 5 do begin a1=regts(i,*) filter_cru,thalf,tsin=reform(a1),tshigh=aa1,/nan mknormal,aa1,refperiod=[0,nkp-1] plot,x(kpall),aa1,thick=2,$ title=t(i)+' r(vs. barrow,mauna)='+string(allr(i,1),format='(F5.2)')+$ ' '+string(allr(i,0),format='(F5.2)') a1=co2amp(0,*) filter_cru,thalf,tsin=reform(a1),tshigh=aa1,/nan mknormal,aa1,refperiod=[0,nkp-1] oplot,xx(kpco2),aa1,linestyle=1,thick=3 a1=co2amp(1,*) filter_cru,thalf,tsin=reform(a1),tshigh=aa1,/nan mknormal,aa1,refperiod=[0,nkp-1] oplot,xx(kpco2),aa1,thick=5 endfor !p.multi(4)=0 ; end