beta/06.biome_ob.ncl
author Forrest Hoffman <forrest@climatemodeling.org>
Thu, 26 Mar 2009 14:02:21 -0400
changeset 1 4be95183fbcd
permissions -rw-r--r--
Modifications to scoring and graphics production for the final version of code for the C-LAMP paper in GCB.
forrest@0
     1
;********************************************************
forrest@0
     2
; required command line input parameters:
forrest@0
     3
;  ncl 'model_name="10cn" model_grid="T42" dirm="/.../ film="..."' 01.npp.ncl
forrest@0
     4
;
forrest@0
     5
; histogram normalized by rain and compute correleration
forrest@0
     6
;**************************************************************
forrest@0
     7
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl.test"
forrest@0
     8
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl.test"
forrest@0
     9
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
forrest@0
    10
;**************************************************************
forrest@0
    11
procedure set_line(lines:string,nline:integer,newlines:string) 
forrest@0
    12
begin
forrest@0
    13
; add line to ascci/html file
forrest@0
    14
    
forrest@0
    15
  nnewlines = dimsizes(newlines)
forrest@0
    16
  if(nline+nnewlines-1.ge.dimsizes(lines))
forrest@0
    17
    print("set_line: bad index, not setting anything.") 
forrest@0
    18
    return
forrest@0
    19
  end if 
forrest@0
    20
  lines(nline:nline+nnewlines-1) = newlines
forrest@0
    21
;  print ("lines = " + lines(nline:nline+nnewlines-1))
forrest@0
    22
  nline = nline + nnewlines
forrest@0
    23
  return 
forrest@0
    24
end
forrest@0
    25
;**************************************************************
forrest@0
    26
; Main code.
forrest@0
    27
begin
forrest@0
    28
 
forrest@0
    29
 plot_type     = "ps"
forrest@0
    30
 plot_type_new = "png"
forrest@0
    31
forrest@0
    32
;************************************************
forrest@0
    33
; read data: model       
forrest@0
    34
;************************************************
forrest@0
    35
 co2_i = 283.1878
forrest@0
    36
 co2_f = 364.1252
forrest@0
    37
forrest@0
    38
 model_grid = "T42"
forrest@0
    39
forrest@0
    40
 model_name_i = "i01.07cn"
forrest@0
    41
 model_name_f = "i01.10cn"
forrest@0
    42
forrest@0
    43
;model_name_i = "i01.07casa"
forrest@0
    44
;model_name_f = "i01.10casa"
forrest@0
    45
forrest@0
    46
 model_name = model_name_f
forrest@0
    47
forrest@0
    48
 dirm = "/fis/cgd/cseg/people/jeff/clamp_data/model/"
forrest@0
    49
 film_i = model_name_i + "_1990-2004_ANN_climo.nc"
forrest@0
    50
 film_f = model_name_f + "_1990-2004_ANN_climo.nc"
forrest@0
    51
forrest@0
    52
 fm_i   = addfile (dirm+film_i,"r")
forrest@0
    53
 fm_f   = addfile (dirm+film_f,"r")
forrest@0
    54
  
forrest@0
    55
 xm     = fm_f->lon  
forrest@0
    56
 ym     = fm_f->lat
forrest@0
    57
forrest@0
    58
 npp_i  = fm_i->NPP
forrest@0
    59
 npp_f  = fm_f->NPP
forrest@0
    60
forrest@0
    61
;Units for these variables are:
forrest@0
    62
;npp_i: g C/m^2/s
forrest@0
    63
forrest@0
    64
 nsec_per_year = 60*60*24*365
forrest@0
    65
  
forrest@0
    66
 npp_i = npp_i *  nsec_per_year
forrest@0
    67
 npp_f = npp_f *  nsec_per_year 
forrest@0
    68
   
forrest@0
    69
;===================================================
forrest@0
    70
; read data: observed -station
forrest@0
    71
;===================================================
forrest@0
    72
forrest@0
    73
 station = (/"DukeFACE" \
forrest@0
    74
            ,"AspenFACE" \
forrest@0
    75
            ,"ORNL-FACE" \
forrest@0
    76
            ,"POP-EUROFACE" \
forrest@0
    77
            /)
forrest@0
    78
forrest@0
    79
 lat_ob = (/ 35.58,  45.40,  35.54, 42.22/)
forrest@0
    80
 lon_ob = (/-79.05, -89.37, -84.20, 11.48/)
forrest@0
    81
 lon_ob = where(lon_ob.lt.0.,lon_ob+360.,lon_ob)
forrest@0
    82
;print (lon_ob)
forrest@0
    83
forrest@0
    84
 n_sta  = dimsizes(station)
forrest@0
    85
 beta_4_ob = new((/n_sta/),float)
forrest@0
    86
 beta_4_ob = 0.60
forrest@0
    87
forrest@0
    88
;===================================================
forrest@0
    89
; get model data at station 
forrest@0
    90
;===================================================
forrest@0
    91
forrest@0
    92
 npp_i_4  =linint2_points(xm,ym,npp_i,True,lon_ob,lat_ob,0)
forrest@0
    93
forrest@0
    94
 npp_f_4  =linint2_points(xm,ym,npp_f,True,lon_ob,lat_ob,0)
forrest@0
    95
forrest@0
    96
;print (npp_i_4)
forrest@0
    97
;print (npp_f_4)
forrest@0
    98
forrest@0
    99
;============================
forrest@0
   100
;compute beta_4
forrest@0
   101
;============================
forrest@0
   102
forrest@0
   103
 beta_4 = new((/n_sta/),float)
forrest@0
   104
forrest@0
   105
 beta_4 = ((npp_f_4/npp_i_4) - 1.)/log(co2_f/co2_i)
forrest@0
   106
forrest@0
   107
 beta_4_avg = avg(beta_4)
forrest@0
   108
forrest@0
   109
;print (beta_4)
forrest@0
   110
;print (beta_4_avg)
forrest@0
   111
forrest@0
   112
;M_beta = abs((beta_4_avg/beta_4_ob) - 1.)* 3.
forrest@0
   113
forrest@0
   114
 bias = sum(abs(beta_4-beta_4_ob)/(abs(beta_4)+abs(beta_4_ob))) 
forrest@0
   115
 M_beta  = (1. - (bias/n_sta))*3.
forrest@0
   116
 
forrest@0
   117
 print (M_beta)
forrest@0
   118
forrest@0
   119
;=========================
forrest@0
   120
; for html table - station
forrest@0
   121
;=========================
forrest@0
   122
forrest@0
   123
  output_html = "table_station.html"
forrest@0
   124
forrest@0
   125
; column (not including header column)
forrest@0
   126
forrest@0
   127
  col_head = (/"Latitude","Longitude","CO2_i","CO2_f","NPP_i","NPP_f","Beta_model","Beta_ob"/)
forrest@0
   128
forrest@0
   129
  ncol = dimsizes(col_head)
forrest@0
   130
forrest@0
   131
; row (not including header row)
forrest@0
   132
  row_head = (/"DukeFACE" \
forrest@0
   133
            ,"AspenFACE" \
forrest@0
   134
            ,"ORNL-FACE" \
forrest@0
   135
            ,"POP-EUROFACE" \
forrest@0
   136
            ,"All Station" \                
forrest@0
   137
            /)  
forrest@0
   138
  nrow = dimsizes(row_head)                  
forrest@0
   139
forrest@0
   140
; arrays to be passed to table. 
forrest@0
   141
  text4 = new ((/nrow, ncol/),string )
forrest@0
   142
forrest@0
   143
 do i=0,nrow-2
forrest@0
   144
  text4(i,0) = sprintf("%.1f",lat_ob(i))
forrest@0
   145
  text4(i,1) = sprintf("%.1f",lon_ob(i))
forrest@0
   146
  text4(i,2) = sprintf("%.1f",co2_i)
forrest@0
   147
  text4(i,3) = sprintf("%.1f",co2_f)
forrest@0
   148
  text4(i,4) = sprintf("%.1f",npp_i_4(0,i))
forrest@0
   149
  text4(i,5) = sprintf("%.1f",npp_f_4(0,i))
forrest@0
   150
  text4(i,6) = sprintf("%.2f",beta_4(i))
forrest@0
   151
  text4(i,7) = "-"
forrest@0
   152
 end do
forrest@0
   153
  text4(nrow-1,0) = "-"
forrest@0
   154
  text4(nrow-1,1) = "-"
forrest@0
   155
  text4(nrow-1,2) = "-"
forrest@0
   156
  text4(nrow-1,3) = "-"
forrest@0
   157
  text4(nrow-1,4) = "-"
forrest@0
   158
  text4(nrow-1,5) = "-"
forrest@0
   159
  text4(nrow-1,6) = sprintf("%.2f",beta_4_avg)
forrest@0
   160
  text4(nrow-1,7) = sprintf("%.2f",avg(beta_4_ob))
forrest@0
   161
forrest@0
   162
;-----------
forrest@0
   163
; html table
forrest@0
   164
;-----------
forrest@0
   165
forrest@0
   166
  header_text = "<H1>Beta Factor: Model "+model_name+"</H1>" 
forrest@0
   167
forrest@0
   168
  header = (/"<HTML>" \
forrest@0
   169
            ,"<HEAD>" \
forrest@0
   170
            ,"<TITLE>CLAMP metrics</TITLE>" \
forrest@0
   171
            ,"</HEAD>" \
forrest@0
   172
            ,header_text \
forrest@0
   173
            /) 
forrest@0
   174
  footer = "</HTML>"
forrest@0
   175
forrest@0
   176
  table_header = (/ \
forrest@0
   177
        "<table border=1 cellspacing=0 cellpadding=3 width=80%>" \
forrest@0
   178
       ,"<tr>" \
forrest@0
   179
       ,"   <th bgcolor=DDDDDD >Station</th>" \
forrest@0
   180
       ,"   <th bgcolor=DDDDDD >"+col_head(0)+"</th>" \
forrest@0
   181
       ,"   <th bgcolor=DDDDDD >"+col_head(1)+"</th>" \
forrest@0
   182
       ,"   <th bgcolor=DDDDDD >"+col_head(2)+"</th>" \
forrest@0
   183
       ,"   <th bgcolor=DDDDDD >"+col_head(3)+"</th>" \
forrest@0
   184
       ,"   <th bgcolor=DDDDDD >"+col_head(4)+"</th>" \
forrest@0
   185
       ,"   <th bgcolor=DDDDDD >"+col_head(5)+"</th>" \
forrest@0
   186
       ,"   <th bgcolor=DDDDDD >"+col_head(6)+"</th>" \
forrest@0
   187
       ,"   <th bgcolor=DDDDDD >"+col_head(7)+"</th>" \
forrest@0
   188
       ,"</tr>" \
forrest@0
   189
       /)
forrest@0
   190
  table_footer = "</table>"
forrest@0
   191
  row_header = "<tr>"
forrest@0
   192
  row_footer = "</tr>"
forrest@0
   193
forrest@0
   194
  lines = new(50000,string)
forrest@0
   195
  nline = 0
forrest@0
   196
forrest@0
   197
  set_line(lines,nline,header)
forrest@0
   198
  set_line(lines,nline,table_header)
forrest@0
   199
;-----------------------------------------------
forrest@0
   200
;row of table
forrest@0
   201
forrest@0
   202
  do n = 0,nrow-1
forrest@0
   203
     set_line(lines,nline,row_header)
forrest@0
   204
forrest@0
   205
     txt1  = row_head(n)
forrest@0
   206
     txt2  = text4(n,0)
forrest@0
   207
     txt3  = text4(n,1)
forrest@0
   208
     txt4  = text4(n,2)
forrest@0
   209
     txt5  = text4(n,3)
forrest@0
   210
     txt6  = text4(n,4)
forrest@0
   211
     txt7  = text4(n,5)
forrest@0
   212
     txt8  = text4(n,6)
forrest@0
   213
     txt9  = text4(n,7)
forrest@0
   214
forrest@0
   215
     set_line(lines,nline,"<th>"+txt1+"</th>")
forrest@0
   216
     set_line(lines,nline,"<th>"+txt2+"</th>")
forrest@0
   217
     set_line(lines,nline,"<th>"+txt3+"</th>")
forrest@0
   218
     set_line(lines,nline,"<th>"+txt4+"</th>")
forrest@0
   219
     set_line(lines,nline,"<th>"+txt5+"</th>")
forrest@0
   220
     set_line(lines,nline,"<th>"+txt6+"</th>")
forrest@0
   221
     set_line(lines,nline,"<th>"+txt7+"</th>")
forrest@0
   222
     set_line(lines,nline,"<th>"+txt8+"</th>")
forrest@0
   223
     set_line(lines,nline,"<th>"+txt9+"</th>")
forrest@0
   224
forrest@0
   225
     set_line(lines,nline,row_footer)
forrest@0
   226
  end do
forrest@0
   227
;-----------------------------------------------
forrest@0
   228
  set_line(lines,nline,table_footer)
forrest@0
   229
  set_line(lines,nline,footer) 
forrest@0
   230
forrest@0
   231
; Now write to an HTML file.
forrest@0
   232
  idx = ind(.not.ismissing(lines))
forrest@0
   233
  if(.not.any(ismissing(idx))) then
forrest@0
   234
    asciiwrite(output_html,lines(idx))
forrest@0
   235
  else
forrest@0
   236
   print ("error?")
forrest@0
   237
  end if
forrest@0
   238
forrest@0
   239
  delete (col_head)
forrest@0
   240
  delete (row_head)
forrest@0
   241
  delete (text4)
forrest@0
   242
  delete (table_header)
forrest@0
   243
  delete (idx)
forrest@0
   244
      
forrest@0
   245
;************************************************
forrest@0
   246
; read data: observed-2
forrest@0
   247
;************************************************
forrest@0
   248
forrest@0
   249
  ob_name = "MODIS MOD 15A2 2000-2005"
forrest@0
   250
forrest@0
   251
  diro  = "/fis/cgd/cseg/people/jeff/clamp_data/lai/ob/"
forrest@0
   252
  filo  = "land_class_"+model_grid+".nc"
forrest@0
   253
forrest@0
   254
  fo = addfile(diro+filo,"r")
forrest@0
   255
 
forrest@0
   256
  classob = tofloat(fo->LAND_CLASS)
forrest@0
   257
forrest@0
   258
; observed data has 20 land-type classes
forrest@0
   259
  nclass = 20
forrest@0
   260
                
forrest@0
   261
;*******************************************************************
forrest@0
   262
; Calculate "nice" bins for binning the data in equally spaced ranges
forrest@0
   263
;********************************************************************
forrest@0
   264
forrest@0
   265
  nclassn     = nclass + 1
forrest@0
   266
  range       = fspan(0,nclassn-1,nclassn)
forrest@0
   267
; print (range)
forrest@0
   268
forrest@0
   269
; Use this range information to grab all the values in a
forrest@0
   270
; particular range, and then take an average.
forrest@0
   271
forrest@0
   272
  nr           = dimsizes(range)
forrest@0
   273
  nx           = nr-1
forrest@0
   274
  xvalues      = new((/2,nx/),float)
forrest@0
   275
  xvalues(0,:) = range(0:nr-2) + (range(1:)-range(0:nr-2))/2.
forrest@0
   276
  dx           = xvalues(0,1) - xvalues(0,0)       ; range width
forrest@0
   277
  dx4          = dx/4                              ; 1/4 of the range
forrest@0
   278
  xvalues(1,:) = xvalues(0,:) - dx/5.
forrest@0
   279
; get data
forrest@0
   280
forrest@0
   281
  base_1D  = ndtooned(classob)
forrest@0
   282
  data1_1D = ndtooned(npp_i)
forrest@0
   283
  data2_1D = ndtooned(npp_f)
forrest@0
   284
forrest@0
   285
; output
forrest@0
   286
forrest@0
   287
  yvalues = new((/2,nx/),float)
forrest@0
   288
  count   = new((/2,nx/),float)
forrest@0
   289
forrest@0
   290
  do nd=0,1
forrest@0
   291
forrest@0
   292
;   See if we are doing data1 (nd=0) or data2 (nd=1).
forrest@0
   293
forrest@0
   294
    base = base_1D
forrest@0
   295
forrest@0
   296
    if(nd.eq.0) then
forrest@0
   297
      data = data1_1D
forrest@0
   298
    else
forrest@0
   299
      data = data2_1D
forrest@0
   300
    end if
forrest@0
   301
forrest@0
   302
; Loop through each range, using base.
forrest@0
   303
forrest@0
   304
    do i=0,nr-2
forrest@0
   305
      if (i.ne.(nr-2)) then
forrest@0
   306
;        print("")
forrest@0
   307
;        print("In range ["+range(i)+","+range(i+1)+")")
forrest@0
   308
         idx = ind((base.ge.range(i)).and.(base.lt.range(i+1)))
forrest@0
   309
      else
forrest@0
   310
;        print("")
forrest@0
   311
;        print("In range ["+range(i)+",)")
forrest@0
   312
         idx = ind(base.ge.range(i))
forrest@0
   313
      end if
forrest@0
   314
forrest@0
   315
;     Calculate average 
forrest@0
   316
forrest@0
   317
      if(.not.any(ismissing(idx))) then
forrest@0
   318
        yvalues(nd,i) = avg(data(idx))
forrest@0
   319
        count(nd,i)   = dimsizes(idx)
forrest@0
   320
      else
forrest@0
   321
        yvalues(nd,i) = yvalues@_FillValue
forrest@0
   322
        count(nd,i)   = 0
forrest@0
   323
      end if
forrest@0
   324
forrest@0
   325
;#############################################################
forrest@0
   326
; set the following 4 classes to _FillValue:
forrest@0
   327
; Water Bodies(0), Urban and Build-Up(13),
forrest@0
   328
; Permenant Snow and Ice(15), Unclassified(17)
forrest@0
   329
forrest@0
   330
      if (i.eq.0 .or. i.eq.13 .or. i.eq.15 .or. i.eq.17) then
forrest@0
   331
         yvalues(nd,i) = yvalues@_FillValue
forrest@0
   332
         count(nd,i)   = 0
forrest@0
   333
      end if
forrest@0
   334
;############################################################# 
forrest@0
   335
forrest@0
   336
;     print(nd + ": " + count(nd,i) + " points, avg = " + yvalues(nd,i))
forrest@0
   337
forrest@0
   338
; Clean up for next time in loop.
forrest@0
   339
forrest@0
   340
      delete(idx)
forrest@0
   341
    end do
forrest@0
   342
forrest@0
   343
    delete(data)
forrest@0
   344
  end do
forrest@0
   345
forrest@0
   346
;============================
forrest@0
   347
;compute beta
forrest@0
   348
;============================
forrest@0
   349
forrest@0
   350
 u       = yvalues(0,:)
forrest@0
   351
 v       = yvalues(1,:)
forrest@0
   352
 u_count = count(0,:)
forrest@0
   353
 v_count = count(1,:)
forrest@0
   354
forrest@0
   355
 good = ind(.not.ismissing(u) .and. .not.ismissing(v))
forrest@0
   356
forrest@0
   357
 uu       = u(good)
forrest@0
   358
 vv       = v(good)
forrest@0
   359
 uu_count = u_count(good)
forrest@0
   360
 vv_count = v_count(good) 
forrest@0
   361
forrest@0
   362
 n_biome = dimsizes(uu)
forrest@0
   363
 beta_biome = new((/n_biome/),float)
forrest@0
   364
forrest@0
   365
 beta_biome = ((vv/uu) - 1.)/log(co2_f/co2_i)
forrest@0
   366
forrest@0
   367
;beta_biome_avg = avg(beta_biome)
forrest@0
   368
 beta_biome_avg = (sum(vv*vv_count)/sum(uu*uu_count) - 1.)/log(co2_f/co2_i)
forrest@0
   369
  
forrest@0
   370
 print (beta_biome_avg)
forrest@0
   371
forrest@0
   372
;===========================
forrest@0
   373
; for html table - biome
forrest@0
   374
;===========================
forrest@0
   375
forrest@0
   376
  output_html = "table_biome.html"
forrest@0
   377
forrest@0
   378
; column (not including header column)
forrest@0
   379
forrest@0
   380
  col_head = (/"CO2_i","CO2_f","NPP_i","NPP_f","Beta_model"/)
forrest@0
   381
forrest@0
   382
  ncol = dimsizes(col_head)
forrest@0
   383
forrest@0
   384
; row (not including header row)
forrest@0
   385
; 3 classes removed: Water Bodies, Urban and Build-Up, Unclassified
forrest@0
   386
; function "good" removed the last 2 classes
forrest@0
   387
; text4(i,2) = sprintf("%.2f",uu(i+1)) remove the first class
forrest@0
   388
forrest@0
   389
  row_head  = (/"Evergreen Needleleaf Forests" \
forrest@0
   390
               ,"Evergreen Broadleaf Forests" \
forrest@0
   391
               ,"Deciduous Needleleaf Forest" \
forrest@0
   392
               ,"Deciduous Broadleaf Forests" \
forrest@0
   393
               ,"Mixed Forests" \                      
forrest@0
   394
               ,"Closed Bushlands" \                   
forrest@0
   395
               ,"Open Bushlands" \                     
forrest@0
   396
               ,"Woody Savannas (S. Hem.)" \           
forrest@0
   397
               ,"Savannas (S. Hem.)" \                 
forrest@0
   398
               ,"Grasslands" \                         
forrest@0
   399
               ,"Permanent Wetlands" \                 
forrest@0
   400
               ,"Croplands" \                                           
forrest@0
   401
               ,"Cropland/Natural Vegetation Mosaic" \             
forrest@0
   402
               ,"Barren or Sparsely Vegetated" \                             
forrest@0
   403
               ,"Woody Savannas (N. Hem.)" \           
forrest@0
   404
               ,"Savannas (N. Hem.)" \
forrest@0
   405
               ,"All Biome" \                
forrest@0
   406
               /)  
forrest@0
   407
  nrow = dimsizes(row_head)                  
forrest@0
   408
forrest@0
   409
; arrays to be passed to table. 
forrest@0
   410
  text4 = new ((/nrow, ncol/),string )
forrest@0
   411
 
forrest@0
   412
 do i=0,nrow-2
forrest@0
   413
  text4(i,0) = sprintf("%.1f",co2_i)
forrest@0
   414
  text4(i,1) = sprintf("%.1f",co2_f)
forrest@0
   415
  text4(i,2) = sprintf("%.1f",uu(i))
forrest@0
   416
  text4(i,3) = sprintf("%.1f",vv(i))
forrest@0
   417
  text4(i,4) = sprintf("%.2f",beta_biome(i))
forrest@0
   418
 end do
forrest@0
   419
  text4(nrow-1,0) = "-"
forrest@0
   420
  text4(nrow-1,1) = "-"
forrest@0
   421
  text4(nrow-1,2) = "-"
forrest@0
   422
  text4(nrow-1,3) = "-"
forrest@0
   423
  text4(nrow-1,4) = sprintf("%.2f",beta_biome_avg)
forrest@0
   424
forrest@0
   425
;**************************************************
forrest@0
   426
; html table
forrest@0
   427
;**************************************************
forrest@0
   428
forrest@0
   429
  header_text = "<H1>Beta Factor: Model "+model_name+"</H1>" 
forrest@0
   430
forrest@0
   431
  header = (/"<HTML>" \
forrest@0
   432
            ,"<HEAD>" \
forrest@0
   433
            ,"<TITLE>CLAMP metrics</TITLE>" \
forrest@0
   434
            ,"</HEAD>" \
forrest@0
   435
            ,header_text \
forrest@0
   436
            /) 
forrest@0
   437
  footer = "</HTML>"
forrest@0
   438
forrest@0
   439
  table_header = (/ \
forrest@0
   440
        "<table border=1 cellspacing=0 cellpadding=3 width=80%>" \
forrest@0
   441
       ,"<tr>" \
forrest@0
   442
       ,"   <th bgcolor=DDDDDD >Biome Class</th>" \
forrest@0
   443
       ,"   <th bgcolor=DDDDDD >"+col_head(0)+"</th>" \
forrest@0
   444
       ,"   <th bgcolor=DDDDDD >"+col_head(1)+"</th>" \
forrest@0
   445
       ,"   <th bgcolor=DDDDDD >"+col_head(2)+"</th>" \
forrest@0
   446
       ,"   <th bgcolor=DDDDDD >"+col_head(3)+"</th>" \
forrest@0
   447
       ,"   <th bgcolor=DDDDDD >"+col_head(4)+"</th>" \
forrest@0
   448
       ,"</tr>" \
forrest@0
   449
       /)
forrest@0
   450
  table_footer = "</table>"
forrest@0
   451
  row_header = "<tr>"
forrest@0
   452
  row_footer = "</tr>"
forrest@0
   453
forrest@0
   454
  lines = new(50000,string)
forrest@0
   455
  nline = 0
forrest@0
   456
forrest@0
   457
  set_line(lines,nline,header)
forrest@0
   458
  set_line(lines,nline,table_header)
forrest@0
   459
;-----------------------------------------------
forrest@0
   460
;row of table
forrest@0
   461
forrest@0
   462
  do n = 0,nrow-1
forrest@0
   463
     set_line(lines,nline,row_header)
forrest@0
   464
forrest@0
   465
     txt1  = row_head(n)
forrest@0
   466
     txt2  = text4(n,0)
forrest@0
   467
     txt3  = text4(n,1)
forrest@0
   468
     txt4  = text4(n,2)
forrest@0
   469
     txt5  = text4(n,3)
forrest@0
   470
     txt6  = text4(n,4)
forrest@0
   471
forrest@0
   472
     set_line(lines,nline,"<th>"+txt1+"</th>")
forrest@0
   473
     set_line(lines,nline,"<th>"+txt2+"</th>")
forrest@0
   474
     set_line(lines,nline,"<th>"+txt3+"</th>")
forrest@0
   475
     set_line(lines,nline,"<th>"+txt4+"</th>")
forrest@0
   476
     set_line(lines,nline,"<th>"+txt5+"</th>")
forrest@0
   477
     set_line(lines,nline,"<th>"+txt6+"</th>")
forrest@0
   478
forrest@0
   479
     set_line(lines,nline,row_footer)
forrest@0
   480
  end do
forrest@0
   481
;-----------------------------------------------
forrest@0
   482
  set_line(lines,nline,table_footer)
forrest@0
   483
  set_line(lines,nline,footer) 
forrest@0
   484
forrest@0
   485
; Now write to an HTML file.
forrest@0
   486
  idx = ind(.not.ismissing(lines))
forrest@0
   487
  if(.not.any(ismissing(idx))) then
forrest@0
   488
    asciiwrite(output_html,lines(idx))
forrest@0
   489
  else
forrest@0
   490
   print ("error?")
forrest@0
   491
  end if
forrest@0
   492
forrest@0
   493
end
forrest@0
   494