clamp: ameriflux/15.plot.ncl@4be95183fbcd

     1 ;************************************************************

     2 ; find year_ob >= 4

     3 ; add year_ob_i and year_ob_f

     4 ; change long_name

     5 ;************************************************************

     6 load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"

     7 load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"

     8 load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"

     9 ;************************************************************

    10 procedure set_line(lines:string,nline:integer,newlines:string)

    11 begin

    12 ; add line to ascci/html file

    14   nnewlines = dimsizes(newlines)

    15   if(nline+nnewlines-1.ge.dimsizes(lines))

    16     print("set_line: bad index, not setting anything.")

    17     return

    18   end if

    19   lines(nline:nline+nnewlines-1) = newlines

    20 ;  print ("lines = " + lines(nline:nline+nnewlines-1))

    21   nline = nline + nnewlines

    22   return

    23 end

    24 ;*************************************************************

    25 begin

    27   plot_type     = "ps"

    28   plot_type_new = "png"

    30 ; for 6 fields, 12-monthly

    31   nmon      = 12

    32   nfield    = 6

    34 ;************************************************

    35 ; read model data

    36 ;************************************************

    37   model = "cn"

    38 ; model = "casa"

    40   model_name = "i01.10" + model

    42   dirm  = "/fis/cgd/cseg/people/jeff/surface_data/"

    43   film  = "lnd_T42.nc"

    44   fm    = addfile(dirm+film,"r")

    46   xm    = fm->lon

    47   ym    = fm->lat

    48 ;------------------------------------------------

    49   nlat = dimsizes(ym)

    50   nlon = dimsizes(xm)

    52   data_mod0 = new ((/nfield,nmon,nlat,nlon/),float)

    54 ; change to unit of observed (u mol/m2/s)

    55 ; Model_units [=] gC/m2/s

    56 ; 12. = molecular weight of C

    57 ; u mol = 1e-6 mol

    58   factor = 1e6 /12.

    60   ENERGY ="new"

    62 ;************************************************

    63 ; read data: observed

    64 ;************************************************

    66  station = (/"ARM_Oklahoma" \

    67             ,"ARM_Oklahoma_burn" \

    68             ,"ARM_Oklahoma_control" \

    69             ,"Atqasuk" \

    70             ,"Audubon" \

    71             ,"AustinCary" \

    72             ,"Bartlett" \

    73             ,"Bondville" \

    74             ,"Brookings" \

    75             ,"Donaldson" \

    76             ,"Duke_Forest_Hardwoods" \

    77             ,"Duke_Forest_Open_Field" \

    78             ,"Duke_Forest_Pine" \

    79             ,"Fermi_Ag" \

    80             ,"Fermi_Prairie" \

    81             ,"Flagstaff_Managed" \

    82             ,"Flagstaff_Unmanaged" \

    83             ,"Flagstaff_Wildfire" \

    84             ,"FortPeck" \

    85             ,"FreemanRanch_mesquite" \

    86             ,"Goodwin_Creek" \

    87             ,"HarvardForest" \

    88             ,"HarvardForestHemlock" \

    89             ,"HowlandForestMain" \

    90             ,"HowlandForestWest" \

    91             ,"Ivotuk" \

    92             ,"KendallGrasslands" \

    93             ,"KennedySpaceCenterPine" \

    94             ,"KennedySpaceCenterScrub" \

    95             ,"LittleProspect" \

    96             ,"LostCreek" \

    97             ,"Mead-irrigated" \

    98             ,"Mead-irrigated-rotation" \

    99             ,"Mead-rainfed" \

   100             ,"Metolius_2nd_YoungPonderosaPine" \

   101             ,"MetoliusEyerly" \

   102             ,"MetoliusIntermediatePine" \

   103             ,"MetoliusOldPonderosaPine" \

   104             ,"MissouriOzark" \

   105             ,"Mize" \

   106             ,"MorganMonroe" \

   107             ,"NiwotRidge" \

   108             ,"NorthCarolina_cc" \

   109             ,"NorthCarolina_lp" \

   110             ,"ParkFalls" \

   111             ,"Rayonier" \

   112             ,"SantaRita" \

   113             ,"SkyOaks_Old" \

   114             ,"SkyOaks_PostFire" \

   115             ,"SkyOaks_Young" \

   116             ,"SylvaniaWilderness" \

   117             ,"Toledo" \

   118             ,"Tonzi" \

   119             ,"UCI_1850" \

   120             ,"UCI_1930" \

   121             ,"UCI_1964" \

   122             ,"UCI_1964wet" \

   123             ,"UCI_1981" \

   124             ,"UCI_1989" \

   125             ,"UCI_1998" \

   126             ,"UMBS" \

   127             ,"Vaira" \

   128             ,"WalkerBranch" \

   129             ,"WillowCreek" \

   130             ,"WindRiver" \

   131             ,"Wisconsin_ihw" \

   132             ,"Wisconsin_irp" \

   133             ,"Wisconsin_mrp" \

   134             ,"Wisconsin_myjp" \

   135             ,"Wisconsin_pb" \

   136             ,"Wisconsin_rpcc" \

   137             ,"Wisconsin_yhw" \

   138             ,"Wisconsin_yjp" \

   139             ,"Wisconsin_yrp" \

   140             /)

   142  year_ob = (/"2003-2006" \

   143             ,"2005-2006" \

   144             ,"2005-2006" \

   145             ,"1999-2006" \

   146             ,"2002-2006" \

   147             ,"2000-2005" \

   148             ,"2004-2005" \

   149             ,"1996-2006" \

   150             ,"2004-2006" \

   151             ,"1999-2004" \

   152             ,"2003-2005" \

   153             ,"2001-2005" \

   154             ,"2001-2005" \

   155             ,"2005-2006" \

   156             ,"2004-2006" \

   157             ,"2005-2006" \

   158             ,"2005-2006" \

   159             ,"2005-2006" \

   160             ,"2000-2006" \

   161             ,"2004-2006" \

   162             ,"2002-2006" \

   163             ,"1991-2004" \

   164             ,"2004-2004" \

   165             ,"1996-2004" \

   166             ,"1999-2004" \

   167             ,"2003-2006" \

   168             ,"2004-2006" \

   169             ,"2002-2002" \

   170             ,"2000-2006" \

   171             ,"2002-2005" \

   172             ,"2001-2005" \

   173             ,"2001-2005" \

   174             ,"2001-2005" \

   175             ,"2001-2005" \

   176             ,"2004-2005" \

   177             ,"2004-2005" \

   178             ,"2003-2005" \

   179             ,"1996-2000" \

   180             ,"2004-2006" \

   181             ,"1998-2004" \

   182             ,"1999-2005" \

   183             ,"1999-2003" \

   184             ,"2005-2006" \

   185             ,"2005-2006" \

   186             ,"1996-2003" \

   187             ,"1998-1998" \

   188             ,"2004-2006" \

   189             ,"1997-2006" \

   190             ,"2004-2006" \

   191             ,"1997-2006" \

   192             ,"2002-2006" \

   193             ,"2004-2005" \

   194             ,"2001-2006" \

   195             ,"2004-2005" \

   196             ,"2001-2005" \

   197             ,"2001-2005" \

   198             ,"2002-2004" \

   199             ,"2001-2005" \

   200             ,"2001-2005" \

   201             ,"2002-2005" \

   202             ,"1999-2003" \

   203             ,"2001-2006" \

   204             ,"1995-1999" \

   205             ,"1999-2005" \

   206             ,"1999-2004" \

   207             ,"2003-2003" \

   208             ,"2003-2003" \

   209             ,"2002-2005" \

   210             ,"2004-2004" \

   211             ,"2002-2002" \

   212             ,"2005-2005" \

   213             ,"2002-2002" \

   214             ,"2004-2005" \

   215             ,"2002-2002" \

   216             /)

   218  year_ob_i = (/2003 \

   219             ,2005 \

   220             ,2005 \

   221             ,1999 \

   222             ,2002 \

   223             ,2000 \

   224             ,2004 \

   225             ,1996 \

   226             ,2004 \

   227             ,1999 \

   228             ,2003 \

   229             ,2001 \

   230             ,2001 \

   231             ,2005 \

   232             ,2004 \

   233             ,2005 \

   234             ,2005 \

   235             ,2005 \

   236             ,2000 \

   237             ,2004 \

   238             ,2002 \

   239             ,1991 \

   240             ,2004 \

   241             ,1996 \

   242             ,1999 \

   243             ,2003 \

   244             ,2004 \

   245             ,2002 \

   246             ,2000 \

   247             ,2002 \

   248             ,2001 \

   249             ,2001 \

   250             ,2001 \

   251             ,2001 \

   252             ,2004 \

   253             ,2004 \

   254             ,2003 \

   255             ,1996 \

   256             ,2004 \

   257             ,1998 \

   258             ,1999 \

   259             ,1999 \

   260             ,2005 \

   261             ,2005 \

   262             ,1996 \

   263             ,1998 \

   264             ,2004 \

   265             ,1997 \

   266             ,2004 \

   267             ,1997 \

   268             ,2002 \

   269             ,2004 \

   270             ,2001 \

   271             ,2004 \

   272             ,2001 \

   273             ,2001 \

   274             ,2002 \

   275             ,2001 \

   276             ,2001 \

   277             ,2002 \

   278             ,1999 \

   279             ,2001 \

   280             ,1995 \

   281             ,1999 \

   282             ,1999 \

   283             ,2003 \

   284             ,2003 \

   285             ,2002 \

   286             ,2004 \

   287             ,2002 \

   288             ,2005 \

   289             ,2002 \

   290             ,2004 \

   291             ,2002 \

   292             /)

   294  year_ob_f = (/2006 \

   295             ,2006 \

   296             ,2006 \

   297             ,2006 \

   298             ,2006 \

   299             ,2005 \

   300             ,2005 \

   301             ,2006 \

   302             ,2006 \

   303             ,2004 \

   304             ,2005 \

   305             ,2005 \

   306             ,2005 \

   307             ,2006 \

   308             ,2006 \

   309             ,2006 \

   310             ,2006 \

   311             ,2006 \

   312             ,2006 \

   313             ,2006 \

   314             ,2006 \

   315             ,2004 \

   316             ,2004 \

   317             ,2004 \

   318             ,2004 \

   319             ,2006 \

   320             ,2006 \

   321             ,2002 \

   322             ,2006 \

   323             ,2005 \

   324             ,2005 \

   325             ,2005 \

   326             ,2005 \

   327             ,2005 \

   328             ,2005 \

   329             ,2005 \

   330             ,2005 \

   331             ,2000 \

   332             ,2006 \

   333             ,2004 \

   334             ,2005 \

   335             ,2003 \

   336             ,2006 \

   337             ,2006 \

   338             ,2003 \

   339             ,1998 \

   340             ,2006 \

   341             ,2006 \

   342             ,2006 \

   343             ,2006 \

   344             ,2006 \

   345             ,2005 \

   346             ,2006 \

   347             ,2005 \

   348             ,2005 \

   349             ,2005 \

   350             ,2004 \

   351             ,2005 \

   352             ,2005 \

   353             ,2005 \

   354             ,2003 \

   355             ,2006 \

   356             ,1999 \

   357             ,2005 \

   358             ,2004 \

   359             ,2003 \

   360             ,2003 \

   361             ,2005 \

   362             ,2004 \

   363             ,2002 \

   364             ,2005 \

   365             ,2002 \

   366             ,2005 \

   367             ,2002 \

   368             /)

   370  year_ob_long = ind((year_ob_f - year_ob_i) .ge. 4 .and. year_ob_i .le. 2001)

   371  print (year_ob_long)

   372  print (year_ob_i(year_ob_long))

   375  field   = (/"NEE Flux" \

   376             ,"Net Radiation" \

   377             ,"Latent Heat" \

   378             ,"Sensible Heat" \

   379             ,"GPP Flux" \

   380             ,"Respiration" \

   381             /)

   383  field_unit = (/"u mol/m2/s" \

   384                ,"W/m2" \

   385                ,"W/m2" \

   386                ,"W/m2" \

   387                ,"u mol/m2/s" \

   388                ,"u mol/m2/s" \

   389                /)

   391 ;nstation  = dimsizes(station)

   392  nstation  = dimsizes(year_ob_long)

   394  data_mod  = new ((/nstation, nfield, nmon/),float)

   395  data_ob   = new ((/nstation, nfield, nmon/),float)

   396  lat_ob    = new ((/nstation/),float)

   397  lon_ob    = new ((/nstation/),float)

   398  unit_ob   = new ((/nfield/),string)

   400  diro_root  = "/fis/cgd/cseg/people/jeff/clamp/ameriflux/data/"

   401  dirm_root  = "/fis/cgd/cseg/people/jeff/clamp_data/model/"

   403  do n = 0,nstation-1

   405 ;-------------------------------------------------

   406 ;   get ob data

   408 ;   nn = n

   409     nn = year_ob_long(n)

   411     diro = diro_root + station(nn)+"/"

   412     filo = year_ob_i(nn)+"-"+year_ob_f(nn)+"_L4_m.nc"

   413     fo   = addfile (diro+filo,"r")

   415     print (filo)

   417     lon_ob(n) = fo->lon

   418     lat_ob(n) = fo->lat

   420     data      = fo->NEE_or_fMDS

   421     data_ob(n,0,:) = dim_avg(data(month|:,year|:))

   422     delete (data)

   424     data      = fo->Rg_f

   425     data_ob(n,1,:) = dim_avg(data(month|:,year|:))

   426     delete (data)

   428     data      = fo->LE_f

   429     data_ob(n,2,:) = dim_avg(data(month|:,year|:))

   430     delete (data)

   432     data      = fo->H_f

   433     data_ob(n,3,:) = dim_avg(data(month|:,year|:))

   434     delete (data)

   436     data      = fo->GPP_or_MDS

   437     data_ob(n,4,:) = dim_avg(data(month|:,year|:))

   438     delete (data)

   440     data      = fo->Reco_or

   441     data_ob(n,5,:) = dim_avg(data(month|:,year|:))

   442     delete (data)

   444     delete (fo)

   445 ;---------------------------------------------------

   446 ;   get model data

   448 ;   film = model_name+"_"+ year_ob(n)+"_MONS_climo.nc"

   449     film = model_name+"_"+"1990-2004"+"_MONS_climo.nc"

   450     fm   = addfile (dirm_root+film,"r")

   452 ;   print (film)

   454 if (ENERGY .eq. "old") then

   456   data = fm->NEE

   457   data_mod0(0,:,:,:) = data(:,:,:) * factor

   458   delete (data)

   460 ; data  = fm->LATENT

   461   data1 = fm->FCEV

   462   data2 = fm->FCTR

   463   data3 = fm->FGEV

   464   data_mod0(2,:,:,:) = data1(:,:,:)+data2(:,:,:)+data3(:,:,:)

   465   delete (data1)

   466   delete (data2)

   467   delete (data3)

   469 ; data = fm->SENSIBLE

   470   data  = fm->FSH

   471   data_mod0(3,:,:,:) = data(:,:,:)

   472   delete (data)

   474 ; data  = fm->NETRAD

   475   data1 = fm->FSA

   476   data2 = fm->FIRA

   477   data_mod0(1,:,:,:) = data1(:,:,:)-data2(:,:,:)-data_mod0(2,:,:,:)-data_mod0(3,:,:,:)

   478   delete (data1)

   479   delete (data2)

   481 else

   483   data = fm->NEE

   484   data_mod0(0,:,:,:) = data(:,:,:) * factor

   485   delete (data)

   487   data = fm->NETRAD

   488   data_mod0(1,:,:,:) = data(:,:,:)

   489   delete (data)

   491   data = fm->LATENT

   492   data_mod0(2,:,:,:) = data(:,:,:)

   493   delete (data)

   495 ; data = fm->SENSIBLE

   496   data = fm->FSH

   497   data_mod0(3,:,:,:) = data(:,:,:)

   498   delete (data)

   499 end if

   501   data = fm->GPP

   502   data_mod0(4,:,:,:) = data(:,:,:) * factor

   503   delete (data)

   505   if (model .eq. "cn") then

   506      data = fm->ER

   507   else

   508      data1 = fm->AR

   509      data2 = fm->HR

   510      data  = data1 + data2

   512      delete (data1)

   513      delete (data2)

   514   end if

   516   data_mod0(5,:,:,:) = data(:,:,:) * factor

   517   delete (data)

   519   delete (fm)

   521 ;************************************************************

   522 ; interpolate model data into observed station

   523 ; note: model is 0-360E, 90S-90N

   524 ;************************************************************

   526 ; to be able to handle observation at (-89.98,-24.80)

   527   ym(0) = -90.

   529   yy = linint2_points_Wrap(xm,ym,data_mod0,True,lon_ob(n),lat_ob(n),0)

   531 ; print (yy)

   533   data_mod(n,:,:) = yy(:,:,0)

   535   delete (yy)

   537  end do

   538 ;************************************************************

   539 ; compute correlation coef and M score

   540 ;************************************************************

   542  score_max = 5.

   544  ccr     = new ((/nstation, nfield/),float)

   545  M_score = new ((/nstation, nfield/),float)

   547  do n=0,nstation-1

   548  do m=0,nfield-1

   549     ccr(n,m) = esccr(data_ob(n,m,:),data_mod(n,m,:),0)

   550     bias = sum(abs(data_mod(n,m,:)-data_ob(n,m,:))/(abs(data_mod(n,m,:))+abs(data_ob(n,m,:))))

   551     M_score(n,m) = (1. -(bias/nmon)) * score_max

   552  end do

   553  end do

   555  M_co2 = avg(M_score(:,0))

   556  M_rad = avg(M_score(:,1))

   557  M_lh  = avg(M_score(:,2))

   558  M_sh  = avg(M_score(:,3))

   559  M_gpp = avg(M_score(:,4))

   560  M_er  = avg(M_score(:,5))

   561  M_all = M_co2+ M_rad +M_lh + M_sh + M_gpp + M_er

   563  M_energy_co2 = sprintf("%.2f", M_co2)

   564  M_energy_rad = sprintf("%.2f", M_rad)

   565  M_energy_lh  = sprintf("%.2f", M_lh )

   566  M_energy_sh  = sprintf("%.2f", M_sh )

   567  M_energy_gpp = sprintf("%.2f", M_gpp)

   568  M_energy_er  = sprintf("%.2f", M_er )

   569  M_energy_all = sprintf("%.2f", M_all)

   571 ;*******************************************************************

   572 ; for station line plot

   573 ;*******************************************************************

   575 ; for x-axis in xyplot

   576   mon = ispan(1,12,1)

   577   mon@long_name = "month"

   579   res                   = True               ; plot mods desired

   580   res@xyLineThicknesses = (/2.0,2.0/)        ; make 2nd lines thicker

   581   res@xyLineColors      = (/"blue","red"/)   ; line color (ob,model)

   582 ;-------------------------------------------------------------------------

   583 ; Add a boxed legend using the more simple method

   585   res@pmLegendDisplayMode    = "Always"

   586 ; res@pmLegendWidthF         = 0.1

   587   res@pmLegendWidthF         = 0.08

   588   res@pmLegendHeightF        = 0.06

   589 ; res@pmLegendOrthogonalPosF = -1.17

   590 ; res@pmLegendOrthogonalPosF = -1.00  ;(downward)

   591   res@pmLegendOrthogonalPosF = -0.30  ;(downward)

   593 ; res@pmLegendParallelPosF   =  0.18

   594   res@pmLegendParallelPosF   =  0.23  ;(rightward)

   596 ; res@lgPerimOn             = False

   597   res@lgLabelFontHeightF     = 0.015

   598   res@xyExplicitLegendLabels = (/"observed",model_name/)

   599 ;-------------------------------------------------------------------

   600 ; for panel plot

   601   res@gsnFrame     = False                   ; Do not draw plot

   602   res@gsnDraw      = False                   ; Do not advance frame

   604   pres                            = True     ; panel plot mods desired

   605   pres@gsnPanelYWhiteSpacePercent = 5        ; increase white space around

   606                                              ; indiv. plots in panel

   607   pres@gsnMaximize                = True     ; fill the page

   608 ;-------------------------------------------------------------------

   610   plot_data   = new((/2,12/),float)

   611   plot_data!0 = "case"

   612   plot_data!1 = "month"

   614   do n = 0,nstation-1

   616 ;    nn = n

   617      nn = year_ob_long(n)

   619 ;----------------------------

   620 ; for observed

   622     plot_name = station(nn)+"_ob"

   623     title = station(nn)+"("+sprintf("%5.2f",lat_ob(n))+","+sprintf("%5.2f",lon_ob(n))+")"

   624     res@tiMainString = title

   626     wks = gsn_open_wks (plot_type,plot_name)

   627     plot=new(nfield,graphic)                        ; create graphic array

   629     do i = 0,nfield-1

   630        plot_data(0,:) = (/data_ob (n,i,:)/)

   631        plot_data@long_name = field(i)+" ("+field_unit(i)+")"

   632        plot(i)=gsn_csm_xy(wks,mon,plot_data(0,:),res)    ; create plot

   633     end do

   635     gsn_panel(wks,plot,(/3,2/),pres)                 ; create panel plot

   637     system("convert "+plot_name+"."+plot_type+" "+plot_name+"."+plot_type_new+";"+ \

   638            "rm "+plot_name+"."+plot_type)

   640     clear (wks)

   641     delete (plot)

   642 ;----------------------------

   643 ; for model_vs_ob

   645     plot_name = station(nn)+"_model_vs_ob"

   646     title = station(nn)+"("+sprintf("%5.2f",lat_ob(n))+","+sprintf("%5.2f",lon_ob(n))+")"

   647     res@tiMainString = title

   649     wks = gsn_open_wks (plot_type,plot_name)

   650     plot=new(nfield,graphic)                         ; create graphic array

   652     do i = 0,nfield-1

   653        plot_data(0,:) = (/data_ob (n,i,:)/)

   654        plot_data(1,:) = (/data_mod(n,i,:)/)

   655        plot_data@long_name = field(i)+" ("+field_unit(i)+")"

   656        plot(i)=gsn_csm_xy(wks,mon,plot_data,res)     ; create plot

   657     end do

   659     gsn_panel(wks,plot,(/3,2/),pres)                 ; create panel plot

   661     system("convert "+plot_name+"."+plot_type+" "+plot_name+"."+plot_type_new+";"+ \

   662            "rm "+plot_name+"."+plot_type)

   664     clear (wks)

   665     delete (plot)

   666  end do

   668 ;*******************************************************************

   669 ; html table of site: observed

   670 ;*******************************************************************

   671   output_html = "line_ob.html"

   673   header = (/"<HTML>" \

   674             ,"<HEAD>" \

   675             ,"<TITLE>CLAMP metrics</TITLE>" \

   676             ,"</HEAD>" \

   677             ,"<H1>Energy at Site: Observation</H1>" \

   678             /)

   679   footer = "</HTML>"

   681   table_header = (/ \

   682         "<table border=1 cellspacing=0 cellpadding=3 width=60%>" \

   683        ,"<tr>" \

   684        ,"   <th bgcolor=DDDDDD >Site Name</th>" \

   685        ,"   <th bgcolor=DDDDDD >Latitude</th>" \

   686        ,"   <th bgcolor=DDDDDD >Longitude</th>" \

   687        ,"   <th bgcolor=DDDDDD >Observed</th>" \

   688        ,"</tr>" \

   689        /)

   690   table_footer = "</table>"

   691   row_header = "<tr>"

   692   row_footer = "</tr>"

   694   lines = new(50000,string)

   695   nline = 0

   697   set_line(lines,nline,header)

   698   set_line(lines,nline,table_header)

   699 ;-----------------------------------------------

   700 ; row of table

   702   do n = 0,nstation-1

   704 ;    nn = n

   705      nn = year_ob_long(n)

   707      set_line(lines,nline,row_header)

   709      txt0 = station(nn)

   710      txt1 = sprintf("%5.2f", lat_ob(n))

   711      txt2 = sprintf("%5.2f", lon_ob(n))

   712      txt3 = year_ob(nn)

   714      set_line(lines,nline,"<th><a href="+txt0+"_ob.png>"+txt0+"</a></th>")

   715      set_line(lines,nline,"<th>"+txt1+"</th>")

   716      set_line(lines,nline,"<th>"+txt2+"</th>")

   717      set_line(lines,nline,"<th>"+txt3+"</th>")

   719      set_line(lines,nline,row_footer)

   720   end do

   721 ;-----------------------------------------------

   722   set_line(lines,nline,table_footer)

   723   set_line(lines,nline,footer)

   725 ; Now write to an HTML file.

   726   idx = ind(.not.ismissing(lines))

   727   if(.not.any(ismissing(idx))) then

   728     asciiwrite(output_html,lines(idx))

   729   else

   730    print ("error?")

   731   end if

   732   delete (idx)

   734 ;*******************************************************************

   735 ; score and line table : model vs observed

   736 ;*******************************************************************

   737   output_html = "score+line_vs_ob.html"

   739   header = (/"<HTML>" \

   740             ,"<HEAD>" \

   741             ,"<TITLE>CLAMP metrics</TITLE>" \

   742             ,"</HEAD>" \

   743             ,"<H1>Energy at Site: Model "+model_name+"</H1>" \

   744             /)

   745   footer = "</HTML>"

   747   delete (table_header)

   748   table_header = (/ \

   749         "<table border=1 cellspacing=0 cellpadding=3 width=100%>" \

   750        ,"<tr>" \

   751        ,"   <th bgcolor=DDDDDD >Site Name</th>" \

   752        ,"   <th bgcolor=DDDDDD >Latitude</th>" \

   753        ,"   <th bgcolor=DDDDDD >Longitude</th>" \

   754        ,"   <th bgcolor=DDDDDD >Observed</th>" \

   755        ,"   <th bgcolor=DDDDDD >CO2 Flux</th>" \

   756        ,"   <th bgcolor=DDDDDD >Net Radiation</th>" \

   757        ,"   <th bgcolor=DDDDDD >Latent Heat</th>" \

   758        ,"   <th bgcolor=DDDDDD >Sensible Heat</th>" \

   759        ,"   <th bgcolor=DDDDDD >GPP Glux</th>" \

   760        ,"   <th bgcolor=DDDDDD >Respiration</th>" \

   761        ,"   <th bgcolor=DDDDDD >Average</th>" \

   762        ,"</tr>" \

   763        /)

   764   table_footer = "</table>"

   765   row_header = "<tr>"

   766   row_footer = "</tr>"

   768   lines = new(50000,string)

   769   nline = 0

   771   set_line(lines,nline,header)

   772   set_line(lines,nline,table_header)

   773 ;-----------------------------------------------

   774 ; row of table

   776   do n = 0,nstation-1

   778 ;    nn = n

   779      nn = year_ob_long(n)

   781      set_line(lines,nline,row_header)

   783      txt0  = station(nn)

   784      txt1  = sprintf("%5.2f", lat_ob(n))

   785      txt2  = sprintf("%5.2f", lon_ob(n))

   786      txt3  = year_ob(nn)

   787      txt4  = sprintf("%5.2f", M_score(n,0))

   788      txt5  = sprintf("%5.2f", M_score(n,1))

   789      txt6  = sprintf("%5.2f", M_score(n,2))

   790      txt7  = sprintf("%5.2f", M_score(n,3))

   791      txt8  = sprintf("%5.2f", M_score(n,4))

   792      txt9  = sprintf("%5.2f", M_score(n,5))

   793      txt10 = sprintf("%5.2f", avg(M_score(n,:)))

   795      set_line(lines,nline,"<th><a href="+txt0+"_model_vs_ob.png>"+txt0+"</a></th>")

   796      set_line(lines,nline,"<th>"+txt1+"</th>")

   797      set_line(lines,nline,"<th>"+txt2+"</th>")

   798      set_line(lines,nline,"<th>"+txt3+"</th>")

   799      set_line(lines,nline,"<th>"+txt4+"</th>")

   800      set_line(lines,nline,"<th>"+txt5+"</th>")

   801      set_line(lines,nline,"<th>"+txt6+"</th>")

   802      set_line(lines,nline,"<th>"+txt7+"</th>")

   803      set_line(lines,nline,"<th>"+txt8+"</th>")

   804      set_line(lines,nline,"<th>"+txt9+"</th>")

   805      set_line(lines,nline,"<th>"+txt10+"</th>")

   807      set_line(lines,nline,row_footer)

   808   end do

   810 ; last row, summary

   811   set_line(lines,nline,row_header)

   813   txt0  = "All_"+sprintf("%.0f", nstation)

   814   txt1  = "-"

   815   txt2  = "-"

   816   txt3  = "-"

   817   txt4  = M_energy_co2

   818   txt5  = M_energy_rad

   819   txt6  = M_energy_lh

   820   txt7  = M_energy_sh

   821   txt8  = M_energy_gpp

   822   txt9  = M_energy_er

   823   txt10 = M_energy_all

   825   set_line(lines,nline,"<th>"+txt0+"</th>")

   826   set_line(lines,nline,"<th>"+txt1+"</th>")

   827   set_line(lines,nline,"<th>"+txt2+"</th>")

   828   set_line(lines,nline,"<th>"+txt3+"</th>")

   829   set_line(lines,nline,"<th>"+txt4+"</th>")

   830   set_line(lines,nline,"<th>"+txt5+"</th>")

   831   set_line(lines,nline,"<th>"+txt6+"</th>")

   832   set_line(lines,nline,"<th>"+txt7+"</th>")

   833   set_line(lines,nline,"<th>"+txt8+"</th>")

   834   set_line(lines,nline,"<th>"+txt9+"</th>")

   835   set_line(lines,nline,"<th>"+txt10+"</th>")

   837   set_line(lines,nline,row_footer)

   838 ;-----------------------------------------------

   839   set_line(lines,nline,table_footer)

   840   set_line(lines,nline,footer)

   842 ; Now write to an HTML file.

   843   idx = ind(.not.ismissing(lines))

   844   if(.not.any(ismissing(idx))) then

   845     asciiwrite(output_html,lines(idx))

   846   else

   847    print ("error?")

   848   end if

   849   delete (idx)

   851 ;***************************************************************************

   852 ; output plots

   853 ;***************************************************************************

   854   output_dir = model_name+"/ameriflux"

   856   system("mv *.png *.html " + output_dir)

   857 ;***************************************************************************

   858 end

author	Forrest Hoffman <forrest@climatemodeling.org>
	Thu, 26 Mar 2009 14:02:21 -0400
changeset 1	4be95183fbcd
permissions	-rw-r--r--