forrest@0: ;************************************************************
forrest@0: ; sort by latitude in decending order (N->S)
forrest@0: ; add year_ob_i and year_ob_f
forrest@0: ; change long_name
forrest@0: ;************************************************************
forrest@0: load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
forrest@0: load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
forrest@0: load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
forrest@0: ;************************************************************
forrest@0: procedure set_line(lines:string,nline:integer,newlines:string)
forrest@0: begin
forrest@0: ; add line to ascci/html file
forrest@0:
forrest@0: nnewlines = dimsizes(newlines)
forrest@0: if(nline+nnewlines-1.ge.dimsizes(lines))
forrest@0: print("set_line: bad index, not setting anything.")
forrest@0: return
forrest@0: end if
forrest@0: lines(nline:nline+nnewlines-1) = newlines
forrest@0: ; print ("lines = " + lines(nline:nline+nnewlines-1))
forrest@0: nline = nline + nnewlines
forrest@0: return
forrest@0: end
forrest@0: ;*************************************************************
forrest@0: begin
forrest@0:
forrest@0: plot_type = "ps"
forrest@0: plot_type_new = "png"
forrest@0:
forrest@0: ; for 6 fields, 12-monthly
forrest@0: nmon = 12
forrest@0: nfield = 6
forrest@0:
forrest@0: ;************************************************
forrest@0: ; read model data
forrest@0: ;************************************************
forrest@0: model = "cn"
forrest@0: ; model = "casa"
forrest@0:
forrest@0: model_name = "i01.10" + model
forrest@0:
forrest@0: dirm = "/fis/cgd/cseg/people/jeff/surface_data/"
forrest@0: film = "lnd_T42.nc"
forrest@0: fm = addfile(dirm+film,"r")
forrest@0:
forrest@0: xm = fm->lon
forrest@0: ym = fm->lat
forrest@0: ;------------------------------------------------
forrest@0: nlat = dimsizes(ym)
forrest@0: nlon = dimsizes(xm)
forrest@0:
forrest@0: data_mod0 = new ((/nfield,nmon,nlat,nlon/),float)
forrest@0:
forrest@0: ; change to unit of observed (u mol/m2/s)
forrest@0: ; Model_units [=] gC/m2/s
forrest@0: ; 12. = molecular weight of C
forrest@0: ; u mol = 1e-6 mol
forrest@0: factor = 1e6 /12.
forrest@0:
forrest@0: ENERGY ="new"
forrest@0:
forrest@0: ;************************************************
forrest@0: ; read data: observed
forrest@0: ;************************************************
forrest@0:
forrest@0: station = (/"ARM_Oklahoma" \
forrest@0: ,"ARM_Oklahoma_burn" \
forrest@0: ,"ARM_Oklahoma_control" \
forrest@0: ,"Atqasuk" \
forrest@0: ,"Audubon" \
forrest@0: ,"AustinCary" \
forrest@0: ,"Bartlett" \
forrest@0: ,"Bondville" \
forrest@0: ,"Brookings" \
forrest@0: ,"Donaldson" \
forrest@0: ,"Duke_Forest_Hardwoods" \
forrest@0: ,"Duke_Forest_Open_Field" \
forrest@0: ,"Duke_Forest_Pine" \
forrest@0: ,"Fermi_Ag" \
forrest@0: ,"Fermi_Prairie" \
forrest@0: ,"Flagstaff_Managed" \
forrest@0: ,"Flagstaff_Unmanaged" \
forrest@0: ,"Flagstaff_Wildfire" \
forrest@0: ,"FortPeck" \
forrest@0: ,"FreemanRanch_mesquite" \
forrest@0: ,"Goodwin_Creek" \
forrest@0: ,"HarvardForest" \
forrest@0: ,"HarvardForestHemlock" \
forrest@0: ,"HowlandForestMain" \
forrest@0: ,"HowlandForestWest" \
forrest@0: ,"Ivotuk" \
forrest@0: ,"KendallGrasslands" \
forrest@0: ,"KennedySpaceCenterPine" \
forrest@0: ,"KennedySpaceCenterScrub" \
forrest@0: ,"LittleProspect" \
forrest@0: ,"LostCreek" \
forrest@0: ,"Mead-irrigated" \
forrest@0: ,"Mead-irrigated-rotation" \
forrest@0: ,"Mead-rainfed" \
forrest@0: ,"Metolius_2nd_YoungPonderosaPine" \
forrest@0: ,"MetoliusEyerly" \
forrest@0: ,"MetoliusIntermediatePine" \
forrest@0: ,"MetoliusOldPonderosaPine" \
forrest@0: ,"MissouriOzark" \
forrest@0: ,"Mize" \
forrest@0: ,"MorganMonroe" \
forrest@0: ,"NiwotRidge" \
forrest@0: ,"NorthCarolina_cc" \
forrest@0: ,"NorthCarolina_lp" \
forrest@0: ,"ParkFalls" \
forrest@0: ,"Rayonier" \
forrest@0: ,"SantaRita" \
forrest@0: ,"SkyOaks_Old" \
forrest@0: ,"SkyOaks_PostFire" \
forrest@0: ,"SkyOaks_Young" \
forrest@0: ,"SylvaniaWilderness" \
forrest@0: ,"Toledo" \
forrest@0: ,"Tonzi" \
forrest@0: ,"UCI_1850" \
forrest@0: ,"UCI_1930" \
forrest@0: ,"UCI_1964" \
forrest@0: ,"UCI_1964wet" \
forrest@0: ,"UCI_1981" \
forrest@0: ,"UCI_1989" \
forrest@0: ,"UCI_1998" \
forrest@0: ,"UMBS" \
forrest@0: ,"Vaira" \
forrest@0: ,"WalkerBranch" \
forrest@0: ,"WillowCreek" \
forrest@0: ,"WindRiver" \
forrest@0: ,"Wisconsin_ihw" \
forrest@0: ,"Wisconsin_irp" \
forrest@0: ,"Wisconsin_mrp" \
forrest@0: ,"Wisconsin_myjp" \
forrest@0: ,"Wisconsin_pb" \
forrest@0: ,"Wisconsin_rpcc" \
forrest@0: ,"Wisconsin_yhw" \
forrest@0: ,"Wisconsin_yjp" \
forrest@0: ,"Wisconsin_yrp" \
forrest@0: /)
forrest@0:
forrest@0: year_ob = (/"2003-2006" \
forrest@0: ,"2005-2006" \
forrest@0: ,"2005-2006" \
forrest@0: ,"1999-2006" \
forrest@0: ,"2002-2006" \
forrest@0: ,"2000-2005" \
forrest@0: ,"2004-2005" \
forrest@0: ,"1996-2006" \
forrest@0: ,"2004-2006" \
forrest@0: ,"1999-2004" \
forrest@0: ,"2003-2005" \
forrest@0: ,"2001-2005" \
forrest@0: ,"2001-2005" \
forrest@0: ,"2005-2006" \
forrest@0: ,"2004-2006" \
forrest@0: ,"2005-2006" \
forrest@0: ,"2005-2006" \
forrest@0: ,"2005-2006" \
forrest@0: ,"2000-2006" \
forrest@0: ,"2004-2006" \
forrest@0: ,"2002-2006" \
forrest@0: ,"1991-2004" \
forrest@0: ,"2004-2004" \
forrest@0: ,"1996-2004" \
forrest@0: ,"1999-2004" \
forrest@0: ,"2003-2006" \
forrest@0: ,"2004-2006" \
forrest@0: ,"2002-2002" \
forrest@0: ,"2000-2006" \
forrest@0: ,"2002-2005" \
forrest@0: ,"2001-2005" \
forrest@0: ,"2001-2005" \
forrest@0: ,"2001-2005" \
forrest@0: ,"2001-2005" \
forrest@0: ,"2004-2005" \
forrest@0: ,"2004-2005" \
forrest@0: ,"2003-2005" \
forrest@0: ,"1996-2000" \
forrest@0: ,"2004-2006" \
forrest@0: ,"1998-2004" \
forrest@0: ,"1999-2005" \
forrest@0: ,"1999-2003" \
forrest@0: ,"2005-2006" \
forrest@0: ,"2005-2006" \
forrest@0: ,"1996-2003" \
forrest@0: ,"1998-1998" \
forrest@0: ,"2004-2006" \
forrest@0: ,"1997-2006" \
forrest@0: ,"2004-2006" \
forrest@0: ,"1997-2006" \
forrest@0: ,"2002-2006" \
forrest@0: ,"2004-2005" \
forrest@0: ,"2001-2006" \
forrest@0: ,"2004-2005" \
forrest@0: ,"2001-2005" \
forrest@0: ,"2001-2005" \
forrest@0: ,"2002-2004" \
forrest@0: ,"2001-2005" \
forrest@0: ,"2001-2005" \
forrest@0: ,"2002-2005" \
forrest@0: ,"1999-2003" \
forrest@0: ,"2001-2006" \
forrest@0: ,"1995-1999" \
forrest@0: ,"1999-2005" \
forrest@0: ,"1999-2004" \
forrest@0: ,"2003-2003" \
forrest@0: ,"2003-2003" \
forrest@0: ,"2002-2005" \
forrest@0: ,"2004-2004" \
forrest@0: ,"2002-2002" \
forrest@0: ,"2005-2005" \
forrest@0: ,"2002-2002" \
forrest@0: ,"2004-2005" \
forrest@0: ,"2002-2002" \
forrest@0: /)
forrest@0:
forrest@0: year_ob_i = (/2003 \
forrest@0: ,2005 \
forrest@0: ,2005 \
forrest@0: ,1999 \
forrest@0: ,2002 \
forrest@0: ,2000 \
forrest@0: ,2004 \
forrest@0: ,1996 \
forrest@0: ,2004 \
forrest@0: ,1999 \
forrest@0: ,2003 \
forrest@0: ,2001 \
forrest@0: ,2001 \
forrest@0: ,2005 \
forrest@0: ,2004 \
forrest@0: ,2005 \
forrest@0: ,2005 \
forrest@0: ,2005 \
forrest@0: ,2000 \
forrest@0: ,2004 \
forrest@0: ,2002 \
forrest@0: ,1991 \
forrest@0: ,2004 \
forrest@0: ,1996 \
forrest@0: ,1999 \
forrest@0: ,2003 \
forrest@0: ,2004 \
forrest@0: ,2002 \
forrest@0: ,2000 \
forrest@0: ,2002 \
forrest@0: ,2001 \
forrest@0: ,2001 \
forrest@0: ,2001 \
forrest@0: ,2001 \
forrest@0: ,2004 \
forrest@0: ,2004 \
forrest@0: ,2003 \
forrest@0: ,1996 \
forrest@0: ,2004 \
forrest@0: ,1998 \
forrest@0: ,1999 \
forrest@0: ,1999 \
forrest@0: ,2005 \
forrest@0: ,2005 \
forrest@0: ,1996 \
forrest@0: ,1998 \
forrest@0: ,2004 \
forrest@0: ,1997 \
forrest@0: ,2004 \
forrest@0: ,1997 \
forrest@0: ,2002 \
forrest@0: ,2004 \
forrest@0: ,2001 \
forrest@0: ,2004 \
forrest@0: ,2001 \
forrest@0: ,2001 \
forrest@0: ,2002 \
forrest@0: ,2001 \
forrest@0: ,2001 \
forrest@0: ,2002 \
forrest@0: ,1999 \
forrest@0: ,2001 \
forrest@0: ,1995 \
forrest@0: ,1999 \
forrest@0: ,1999 \
forrest@0: ,2003 \
forrest@0: ,2003 \
forrest@0: ,2002 \
forrest@0: ,2004 \
forrest@0: ,2002 \
forrest@0: ,2005 \
forrest@0: ,2002 \
forrest@0: ,2004 \
forrest@0: ,2002 \
forrest@0: /)
forrest@0:
forrest@0: year_ob_f = (/2006 \
forrest@0: ,2006 \
forrest@0: ,2006 \
forrest@0: ,2006 \
forrest@0: ,2006 \
forrest@0: ,2005 \
forrest@0: ,2005 \
forrest@0: ,2006 \
forrest@0: ,2006 \
forrest@0: ,2004 \
forrest@0: ,2005 \
forrest@0: ,2005 \
forrest@0: ,2005 \
forrest@0: ,2006 \
forrest@0: ,2006 \
forrest@0: ,2006 \
forrest@0: ,2006 \
forrest@0: ,2006 \
forrest@0: ,2006 \
forrest@0: ,2006 \
forrest@0: ,2006 \
forrest@0: ,2004 \
forrest@0: ,2004 \
forrest@0: ,2004 \
forrest@0: ,2004 \
forrest@0: ,2006 \
forrest@0: ,2006 \
forrest@0: ,2002 \
forrest@0: ,2006 \
forrest@0: ,2005 \
forrest@0: ,2005 \
forrest@0: ,2005 \
forrest@0: ,2005 \
forrest@0: ,2005 \
forrest@0: ,2005 \
forrest@0: ,2005 \
forrest@0: ,2005 \
forrest@0: ,2000 \
forrest@0: ,2006 \
forrest@0: ,2004 \
forrest@0: ,2005 \
forrest@0: ,2003 \
forrest@0: ,2006 \
forrest@0: ,2006 \
forrest@0: ,2003 \
forrest@0: ,1998 \
forrest@0: ,2006 \
forrest@0: ,2006 \
forrest@0: ,2006 \
forrest@0: ,2006 \
forrest@0: ,2006 \
forrest@0: ,2005 \
forrest@0: ,2006 \
forrest@0: ,2005 \
forrest@0: ,2005 \
forrest@0: ,2005 \
forrest@0: ,2004 \
forrest@0: ,2005 \
forrest@0: ,2005 \
forrest@0: ,2005 \
forrest@0: ,2003 \
forrest@0: ,2006 \
forrest@0: ,1999 \
forrest@0: ,2005 \
forrest@0: ,2004 \
forrest@0: ,2003 \
forrest@0: ,2003 \
forrest@0: ,2005 \
forrest@0: ,2004 \
forrest@0: ,2002 \
forrest@0: ,2005 \
forrest@0: ,2002 \
forrest@0: ,2005 \
forrest@0: ,2002 \
forrest@0: /)
forrest@0:
forrest@0: field = (/"NEE Flux" \
forrest@0: ,"Net Radiation" \
forrest@0: ,"Latent Heat" \
forrest@0: ,"Sensible Heat" \
forrest@0: ,"GPP Flux" \
forrest@0: ,"Respiration" \
forrest@0: /)
forrest@0:
forrest@0: field_unit = (/"u mol/m2/s" \
forrest@0: ,"W/m2" \
forrest@0: ,"W/m2" \
forrest@0: ,"W/m2" \
forrest@0: ,"u mol/m2/s" \
forrest@0: ,"u mol/m2/s" \
forrest@0: /)
forrest@0:
forrest@0: nstation = dimsizes(station)
forrest@0:
forrest@0: data_mod = new ((/nstation, nfield, nmon/),float)
forrest@0: data_ob = new ((/nstation, nfield, nmon/),float)
forrest@0: lat_ob = new ((/nstation/),float)
forrest@0: lon_ob = new ((/nstation/),float)
forrest@0: unit_ob = new ((/nfield/),string)
forrest@0:
forrest@0: diro_root = "/fis/cgd/cseg/people/jeff/clamp/ameriflux/data/"
forrest@0: dirm_root = "/fis/cgd/cseg/people/jeff/clamp_data/model/"
forrest@0:
forrest@0: do n = 0,nstation-1
forrest@0:
forrest@0: ;-------------------------------------------------
forrest@0: ; get ob data
forrest@0:
forrest@0: diro = diro_root + station(n)+"/"
forrest@0: filo = year_ob_i(n)+"-"+year_ob_f(n)+"_L4_m.nc"
forrest@0: fo = addfile (diro+filo,"r")
forrest@0:
forrest@0: print (filo)
forrest@0:
forrest@0: lon_ob(n) = fo->lon
forrest@0: lat_ob(n) = fo->lat
forrest@0:
forrest@0: data = fo->NEE_or_fMDS
forrest@0: data_ob(n,0,:) = dim_avg(data(month|:,year|:))
forrest@0: delete (data)
forrest@0:
forrest@0: data = fo->Rg_f
forrest@0: data_ob(n,1,:) = dim_avg(data(month|:,year|:))
forrest@0: delete (data)
forrest@0:
forrest@0: data = fo->LE_f
forrest@0: data_ob(n,2,:) = dim_avg(data(month|:,year|:))
forrest@0: delete (data)
forrest@0:
forrest@0: data = fo->H_f
forrest@0: data_ob(n,3,:) = dim_avg(data(month|:,year|:))
forrest@0: delete (data)
forrest@0:
forrest@0: data = fo->GPP_or_MDS
forrest@0: data_ob(n,4,:) = dim_avg(data(month|:,year|:))
forrest@0: delete (data)
forrest@0:
forrest@0: data = fo->Reco_or
forrest@0: data_ob(n,5,:) = dim_avg(data(month|:,year|:))
forrest@0: delete (data)
forrest@0:
forrest@0: delete (fo)
forrest@0: ;---------------------------------------------------
forrest@0: ; get model data
forrest@0:
forrest@0: ; film = model_name+"_"+ year_ob(n)+"_MONS_climo.nc"
forrest@0: film = model_name+"_"+"1990-2004"+"_MONS_climo.nc"
forrest@0: fm = addfile (dirm_root+film,"r")
forrest@0:
forrest@0: ; print (film)
forrest@0:
forrest@0: if (ENERGY .eq. "old") then
forrest@0:
forrest@0: data = fm->NEE
forrest@0: data_mod0(0,:,:,:) = data(:,:,:) * factor
forrest@0: delete (data)
forrest@0:
forrest@0: ; data = fm->LATENT
forrest@0: data1 = fm->FCEV
forrest@0: data2 = fm->FCTR
forrest@0: data3 = fm->FGEV
forrest@0: data_mod0(2,:,:,:) = data1(:,:,:)+data2(:,:,:)+data3(:,:,:)
forrest@0: delete (data1)
forrest@0: delete (data2)
forrest@0: delete (data3)
forrest@0:
forrest@0: ; data = fm->SENSIBLE
forrest@0: data = fm->FSH
forrest@0: data_mod0(3,:,:,:) = data(:,:,:)
forrest@0: delete (data)
forrest@0:
forrest@0: ; data = fm->NETRAD
forrest@0: data1 = fm->FSA
forrest@0: data2 = fm->FIRA
forrest@0: data_mod0(1,:,:,:) = data1(:,:,:)-data2(:,:,:)-data_mod0(2,:,:,:)-data_mod0(3,:,:,:)
forrest@0: delete (data1)
forrest@0: delete (data2)
forrest@0:
forrest@0: else
forrest@0:
forrest@0: data = fm->NEE
forrest@0: data_mod0(0,:,:,:) = data(:,:,:) * factor
forrest@0: delete (data)
forrest@0:
forrest@0: data = fm->NETRAD
forrest@0: data_mod0(1,:,:,:) = data(:,:,:)
forrest@0: delete (data)
forrest@0:
forrest@0: data = fm->LATENT
forrest@0: data_mod0(2,:,:,:) = data(:,:,:)
forrest@0: delete (data)
forrest@0:
forrest@0: ; data = fm->SENSIBLE
forrest@0: data = fm->FSH
forrest@0: data_mod0(3,:,:,:) = data(:,:,:)
forrest@0: delete (data)
forrest@0: end if
forrest@0:
forrest@0: data = fm->GPP
forrest@0: data_mod0(4,:,:,:) = data(:,:,:) * factor
forrest@0: delete (data)
forrest@0:
forrest@0: if (model .eq. "cn") then
forrest@0: data = fm->ER
forrest@0: else
forrest@0: data1 = fm->AR
forrest@0: data2 = fm->HR
forrest@0: data = data1 + data2
forrest@0:
forrest@0: delete (data1)
forrest@0: delete (data2)
forrest@0: end if
forrest@0:
forrest@0: data_mod0(5,:,:,:) = data(:,:,:) * factor
forrest@0: delete (data)
forrest@0:
forrest@0: delete (fm)
forrest@0:
forrest@0: ;************************************************************
forrest@0: ; interpolate model data into observed station
forrest@0: ; note: model is 0-360E, 90S-90N
forrest@0: ;************************************************************
forrest@0:
forrest@0: ; to be able to handle observation at (-89.98,-24.80)
forrest@0: ym(0) = -90.
forrest@0:
forrest@0: yy = linint2_points_Wrap(xm,ym,data_mod0,True,lon_ob(n),lat_ob(n),0)
forrest@0:
forrest@0: ; print (yy)
forrest@0:
forrest@0: data_mod(n,:,:) = yy(:,:,0)
forrest@0:
forrest@0: delete (yy)
forrest@0:
forrest@0: end do
forrest@0:
forrest@0: ;************************************************************
forrest@0: ; compute correlation coef and M score
forrest@0: ;************************************************************
forrest@0:
forrest@0: score_max = 5.
forrest@0:
forrest@0: ccr = new ((/nstation, nfield/),float)
forrest@0: M_score = new ((/nstation, nfield/),float)
forrest@0:
forrest@0: do n=0,nstation-1
forrest@0: do m=0,nfield-1
forrest@0: ccr(n,m) = esccr(data_ob(n,m,:),data_mod(n,m,:),0)
forrest@0: bias = sum(abs(data_mod(n,m,:)-data_ob(n,m,:))/(abs(data_mod(n,m,:))+abs(data_ob(n,m,:))))
forrest@0: M_score(n,m) = (1. -(bias/nmon)) * score_max
forrest@0: end do
forrest@0: end do
forrest@0:
forrest@0: M_co2 = avg(M_score(:,0))
forrest@0: M_rad = avg(M_score(:,1))
forrest@0: M_lh = avg(M_score(:,2))
forrest@0: M_sh = avg(M_score(:,3))
forrest@0: M_gpp = avg(M_score(:,4))
forrest@0: M_er = avg(M_score(:,5))
forrest@0: M_all = M_co2+ M_rad +M_lh + M_sh + M_gpp + M_er
forrest@0:
forrest@0: M_energy_co2 = sprintf("%.2f", M_co2)
forrest@0: M_energy_rad = sprintf("%.2f", M_rad)
forrest@0: M_energy_lh = sprintf("%.2f", M_lh )
forrest@0: M_energy_sh = sprintf("%.2f", M_sh )
forrest@0: M_energy_gpp = sprintf("%.2f", M_gpp)
forrest@0: M_energy_er = sprintf("%.2f", M_er )
forrest@0: M_energy_all = sprintf("%.2f", M_all)
forrest@0:
forrest@0: ;*******************************************************************
forrest@0: ; for station line plot
forrest@0: ;*******************************************************************
forrest@0:
forrest@0: ; for x-axis in xyplot
forrest@0: mon = ispan(1,12,1)
forrest@0: mon@long_name = "month"
forrest@0:
forrest@0: res = True ; plot mods desired
forrest@0: res@xyLineThicknesses = (/2.0,2.0/) ; make 2nd lines thicker
forrest@0: res@xyLineColors = (/"blue","red"/) ; line color (ob,model)
forrest@0: ;-------------------------------------------------------------------------
forrest@0: ; Add a boxed legend using the more simple method
forrest@0:
forrest@0: res@pmLegendDisplayMode = "Always"
forrest@0: ; res@pmLegendWidthF = 0.1
forrest@0: res@pmLegendWidthF = 0.08
forrest@0: res@pmLegendHeightF = 0.06
forrest@0: ; res@pmLegendOrthogonalPosF = -1.17
forrest@0: ; res@pmLegendOrthogonalPosF = -1.00 ;(downward)
forrest@0: res@pmLegendOrthogonalPosF = -0.30 ;(downward)
forrest@0:
forrest@0: ; res@pmLegendParallelPosF = 0.18
forrest@0: res@pmLegendParallelPosF = 0.23 ;(rightward)
forrest@0:
forrest@0: ; res@lgPerimOn = False
forrest@0: res@lgLabelFontHeightF = 0.015
forrest@0: res@xyExplicitLegendLabels = (/"observed",model_name/)
forrest@0: ;-------------------------------------------------------------------
forrest@0: ; for panel plot
forrest@0: res@gsnFrame = False ; Do not draw plot
forrest@0: res@gsnDraw = False ; Do not advance frame
forrest@0:
forrest@0: pres = True ; panel plot mods desired
forrest@0: pres@gsnPanelYWhiteSpacePercent = 5 ; increase white space around
forrest@0: ; indiv. plots in panel
forrest@0: pres@gsnMaximize = True ; fill the page
forrest@0: ;-------------------------------------------------------------------
forrest@0:
forrest@0: plot_data = new((/2,12/),float)
forrest@0: plot_data!0 = "case"
forrest@0: plot_data!1 = "month"
forrest@0:
forrest@0: do n = 0,nstation-1
forrest@0: ;----------------------------
forrest@0: ; for observed
forrest@0:
forrest@0: plot_name = station(n)+"_ob"
forrest@0: title = station(n)+"("+sprintf("%5.2f",lat_ob(n))+","+sprintf("%5.2f",lon_ob(n))+")"
forrest@0: res@tiMainString = title
forrest@0:
forrest@0: wks = gsn_open_wks (plot_type,plot_name)
forrest@0: plot=new(nfield,graphic) ; create graphic array
forrest@0:
forrest@0: do i = 0,nfield-1
forrest@0: plot_data(0,:) = (/data_ob (n,i,:)/)
forrest@0: plot_data@long_name = field(i)+" ("+field_unit(i)+")"
forrest@0: plot(i)=gsn_csm_xy(wks,mon,plot_data(0,:),res) ; create plot
forrest@0: end do
forrest@0:
forrest@0: gsn_panel(wks,plot,(/3,2/),pres) ; create panel plot
forrest@0:
forrest@0: system("convert "+plot_name+"."+plot_type+" "+plot_name+"."+plot_type_new+";"+ \
forrest@0: "rm "+plot_name+"."+plot_type)
forrest@0:
forrest@0: clear (wks)
forrest@0: delete (plot)
forrest@0: ;----------------------------
forrest@0: ; for model_vs_ob
forrest@0:
forrest@0: plot_name = station(n)+"_model_vs_ob"
forrest@0: title = station(n)+"("+sprintf("%5.2f",lat_ob(n))+","+sprintf("%5.2f",lon_ob(n))+")"
forrest@0: res@tiMainString = title
forrest@0:
forrest@0: wks = gsn_open_wks (plot_type,plot_name)
forrest@0: plot=new(nfield,graphic) ; create graphic array
forrest@0:
forrest@0: do i = 0,nfield-1
forrest@0: plot_data(0,:) = (/data_ob (n,i,:)/)
forrest@0: plot_data(1,:) = (/data_mod(n,i,:)/)
forrest@0: plot_data@long_name = field(i)+" ("+field_unit(i)+")"
forrest@0: plot(i)=gsn_csm_xy(wks,mon,plot_data,res) ; create plot
forrest@0: end do
forrest@0:
forrest@0: gsn_panel(wks,plot,(/3,2/),pres) ; create panel plot
forrest@0:
forrest@0: system("convert "+plot_name+"."+plot_type+" "+plot_name+"."+plot_type_new+";"+ \
forrest@0: "rm "+plot_name+"."+plot_type)
forrest@0:
forrest@0: clear (wks)
forrest@0: delete (plot)
forrest@0: end do
forrest@0:
forrest@0: ;###################################################################
forrest@0: ; for the following tables,
forrest@0: ; sort by latitude in decending order (N->S)
forrest@0: ; sort by lat in decending order (N->S)
forrest@0:
forrest@0: isort = dim_pqsort(lat_ob,-1)
forrest@0:
forrest@0: station_sort = station(isort)
forrest@0: lat_ob_sort = lat_ob(isort)
forrest@0: lon_ob_sort = lon_ob(isort)
forrest@0: year_ob_sort = year_ob(isort)
forrest@0: M_score_sort = M_score(isort,:)
forrest@0:
forrest@0: print(isort)
forrest@0: print(lat_ob_sort)
forrest@0: ;###################################################################
forrest@0: ;*******************************************************************
forrest@0: ; html table of site: observed
forrest@0: ;*******************************************************************
forrest@0: output_html = "line_ob.html"
forrest@0:
forrest@0: header = (/"" \
forrest@0: ,"
" \
forrest@0: ,"CLAMP metrics" \
forrest@0: ,"" \
forrest@0: ,"Energy at Site: Observation
" \
forrest@0: /)
forrest@0: footer = ""
forrest@0:
forrest@0: table_header = (/ \
forrest@0: "" \
forrest@0: ,"" \
forrest@0: ," Site Name | " \
forrest@0: ," Latitude | " \
forrest@0: ," Longitude | " \
forrest@0: ," Observed | " \
forrest@0: ,"
" \
forrest@0: /)
forrest@0: table_footer = "
"
forrest@0: row_header = ""
forrest@0: row_footer = "
"
forrest@0:
forrest@0: lines = new(50000,string)
forrest@0: nline = 0
forrest@0:
forrest@0: set_line(lines,nline,header)
forrest@0: set_line(lines,nline,table_header)
forrest@0: ;-----------------------------------------------
forrest@0: ; row of table
forrest@0:
forrest@0: do n = 0,nstation-1
forrest@0: set_line(lines,nline,row_header)
forrest@0:
forrest@0: txt0 = station_sort(n)
forrest@0: txt1 = sprintf("%5.2f", lat_ob_sort(n))
forrest@0: txt2 = sprintf("%5.2f", lon_ob_sort(n))
forrest@0: txt3 = year_ob_sort(n)
forrest@0:
forrest@0: set_line(lines,nline,""+txt0+" | ")
forrest@0: set_line(lines,nline,""+txt1+" | ")
forrest@0: set_line(lines,nline,""+txt2+" | ")
forrest@0: set_line(lines,nline,""+txt3+" | ")
forrest@0:
forrest@0: set_line(lines,nline,row_footer)
forrest@0: end do
forrest@0: ;-----------------------------------------------
forrest@0: set_line(lines,nline,table_footer)
forrest@0: set_line(lines,nline,footer)
forrest@0:
forrest@0: ; Now write to an HTML file.
forrest@0: idx = ind(.not.ismissing(lines))
forrest@0: if(.not.any(ismissing(idx))) then
forrest@0: asciiwrite(output_html,lines(idx))
forrest@0: else
forrest@0: print ("error?")
forrest@0: end if
forrest@0: delete (idx)
forrest@0:
forrest@0: ;*******************************************************************
forrest@0: ; score and line table : model vs observed
forrest@0: ;*******************************************************************
forrest@0: output_html = "score+line_vs_ob.html"
forrest@0:
forrest@0: header = (/"" \
forrest@0: ,"" \
forrest@0: ,"CLAMP metrics" \
forrest@0: ,"" \
forrest@0: ,"Energy at Site: Model "+model_name+"
" \
forrest@0: /)
forrest@0: footer = ""
forrest@0:
forrest@0: delete (table_header)
forrest@0: table_header = (/ \
forrest@0: "" \
forrest@0: ,"" \
forrest@0: ," Site Name | " \
forrest@0: ," Latitude | " \
forrest@0: ," Longitude | " \
forrest@0: ," Observed | " \
forrest@0: ," CO2 Flux | " \
forrest@0: ," Net Radiation | " \
forrest@0: ," Latent Heat | " \
forrest@0: ," Sensible Heat | " \
forrest@0: ," GPP Glux | " \
forrest@0: ," Respiration | " \
forrest@0: ," Average | " \
forrest@0: ,"
" \
forrest@0: /)
forrest@0: table_footer = "
"
forrest@0: row_header = ""
forrest@0: row_footer = "
"
forrest@0:
forrest@0: lines = new(50000,string)
forrest@0: nline = 0
forrest@0:
forrest@0: set_line(lines,nline,header)
forrest@0: set_line(lines,nline,table_header)
forrest@0: ;-----------------------------------------------
forrest@0: ; row of table
forrest@0:
forrest@0: do n = 0,nstation-1
forrest@0: set_line(lines,nline,row_header)
forrest@0:
forrest@0: txt0 = station_sort(n)
forrest@0: txt1 = sprintf("%5.2f", lat_ob_sort(n))
forrest@0: txt2 = sprintf("%5.2f", lon_ob_sort(n))
forrest@0: txt3 = year_ob_sort(n)
forrest@0: txt4 = sprintf("%5.2f", M_score_sort(n,0))
forrest@0: txt5 = sprintf("%5.2f", M_score_sort(n,1))
forrest@0: txt6 = sprintf("%5.2f", M_score_sort(n,2))
forrest@0: txt7 = sprintf("%5.2f", M_score_sort(n,3))
forrest@0: txt8 = sprintf("%5.2f", M_score_sort(n,4))
forrest@0: txt9 = sprintf("%5.2f", M_score_sort(n,5))
forrest@0: txt10 = sprintf("%5.2f", avg(M_score_sort(n,:)))
forrest@0:
forrest@0: set_line(lines,nline,""+txt0+" | ")
forrest@0: set_line(lines,nline,""+txt1+" | ")
forrest@0: set_line(lines,nline,""+txt2+" | ")
forrest@0: set_line(lines,nline,""+txt3+" | ")
forrest@0: set_line(lines,nline,""+txt4+" | ")
forrest@0: set_line(lines,nline,""+txt5+" | ")
forrest@0: set_line(lines,nline,""+txt6+" | ")
forrest@0: set_line(lines,nline,""+txt7+" | ")
forrest@0: set_line(lines,nline,""+txt8+" | ")
forrest@0: set_line(lines,nline,""+txt9+" | ")
forrest@0: set_line(lines,nline,""+txt10+" | ")
forrest@0:
forrest@0: set_line(lines,nline,row_footer)
forrest@0: end do
forrest@0:
forrest@0: ; last row, summary
forrest@0: set_line(lines,nline,row_header)
forrest@0:
forrest@0: txt0 = "All_"+sprintf("%.0f", nstation)
forrest@0: txt1 = "-"
forrest@0: txt2 = "-"
forrest@0: txt3 = "-"
forrest@0: txt4 = M_energy_co2
forrest@0: txt5 = M_energy_rad
forrest@0: txt6 = M_energy_lh
forrest@0: txt7 = M_energy_sh
forrest@0: txt8 = M_energy_gpp
forrest@0: txt9 = M_energy_er
forrest@0: txt10 = M_energy_all
forrest@0:
forrest@0: set_line(lines,nline,""+txt0+" | ")
forrest@0: set_line(lines,nline,""+txt1+" | ")
forrest@0: set_line(lines,nline,""+txt2+" | ")
forrest@0: set_line(lines,nline,""+txt3+" | ")
forrest@0: set_line(lines,nline,""+txt4+" | ")
forrest@0: set_line(lines,nline,""+txt5+" | ")
forrest@0: set_line(lines,nline,""+txt6+" | ")
forrest@0: set_line(lines,nline,""+txt7+" | ")
forrest@0: set_line(lines,nline,""+txt8+" | ")
forrest@0: set_line(lines,nline,""+txt9+" | ")
forrest@0: set_line(lines,nline,""+txt10+" | ")
forrest@0:
forrest@0: set_line(lines,nline,row_footer)
forrest@0: ;-----------------------------------------------
forrest@0: set_line(lines,nline,table_footer)
forrest@0: set_line(lines,nline,footer)
forrest@0:
forrest@0: ; Now write to an HTML file.
forrest@0: idx = ind(.not.ismissing(lines))
forrest@0: if(.not.any(ismissing(idx))) then
forrest@0: asciiwrite(output_html,lines(idx))
forrest@0: else
forrest@0: print ("error?")
forrest@0: end if
forrest@0: delete (idx)
forrest@0:
forrest@0: ;***************************************************************************
forrest@0: ; output plots
forrest@0: ;***************************************************************************
forrest@0: output_dir = model_name+"/ameriflux"
forrest@0:
forrest@0: system("mv *.png *.html " + output_dir)
forrest@0: ;***************************************************************************
forrest@0: end