[Dart-dev] DART/branches Revision: 13185

[dart at ucar.edu]

thoar at ucar.edu
2019-06-20 21:30:09 -0600 (Thu, 20 Jun 2019)
1027
Added missing vertical coordinate information for (?,<anylevelcoordinate>,lat,lon) variables.
There are multiple vertical coordinate possibilities, levgrnd, levsoi, levtot, ...
Added code to use any variable explicitly stated as QTY_SOIL_MOISTURE as opposed to
adding QTY_SOIL_LIQUID_WATER and QTY_SOIL_ICE together.   

The gridded history variable H2OSOI is QTY_SOIL_MOISTURE and has units mm3/mm3
The restart variables H2OSOI_LIQ and H2OSOI_ICE  have units kg/m2

The satellite soil moistures typically come as m3/m3 so there would seem to
be no units conversion necessary if using H2OSOI for QTY_SOIL_MOISTURE.
Some satellite soil moistures come as percentages, these have to be divided by 100.0
to have the same units as H2OSOI. 

One of the benefits here is that the COSMOS operator is written in m3/m3 ...
It was always a pain to convert kg/m2 to m3/m3

The downside is that H2OSOI is only available at the gridcell level.
If more accurate forward operators are desired (for the column level)
then more work is needed.




Modified: DART/branches/cesm_clm/models/clm/model_mod.f90
===================================================================
--- DART/branches/cesm_clm/models/clm/model_mod.f90	2019-06-20 23:29:26 UTC (rev 13184)
+++ DART/branches/cesm_clm/models/clm/model_mod.f90	2019-06-21 03:30:09 UTC (rev 13185)
@@ -841,22 +841,17 @@
 
    elseif (rank == 3) then
  
-      ! restart file variables never have 3 dimensions
-      ! vector_history variables   may have 3 dimensions [time, lat, lon]
-      ! history file variables always have 3 dimensions [time, lat, lon]
-      !     exception is float H2OSOI(time, levgrnd, lat, lon) ... but we
-      !     float       H2OSOI(time, levsoi, lat, lon) ;
-      !     float        HR_vr(time, levsoi, lat, lon) ;
-      !     float    LITR2C_vr(time, levsoi, lat, lon) ;
+      ! restart file   variables  never have 3 dimensions
+      ! vector_history variables    may have 3 dimensions [time, lat, lon]
+      ! history file   variables always have 3 dimensions [time, lat, lon]
+      !     float      H2OSOI(time, levgrnd, lat, lon) ;
+      !     float    TSOI_ICE(time, levgrnd, lat, lon) ;
       !     float PCT_GLC_MEC(time, glc_nec, lat, lon) ;
       !     float  PCT_LANDUNIT(time, ltype, lat, lon) ;
       !     float  PCT_NAT_PFT(time, natpft, lat, lon) ;
       !     float        TLAKE(time, levlak, lat, lon) ;
-      !     float    TSOI_ICE(time, levgrnd, lat, lon) ;
       !     float       VEGWP(time, nvegwcs, lat, lon) ;
 
-      !     have access to restart file h2osoi_[liq,ice]
-
       if ((debug > 8) .and. do_output()) then
          write(*,*)
          write(*,*)'variable ',trim(varname)
@@ -878,9 +873,22 @@
                  source, revision, revdate, text2=string2)
       endif
 
+      !>@todo  extend fill_levels to support all the vertical levels
+
+      SELECT CASE ( trim(dimnames(3)) )
+         CASE ("levgrnd")
+            levtot(1:dimlens(3)) = LEVGRND;
+         CASE DEFAULT
+            write(string1,*)'(3d) unsupported vertical dimension name "'//trim(dimnames(3))//'"'
+            write(string2,*)' while trying to create metadata for "'//trim(varname)//'"'
+            call error_handler(E_ERR, routine, string1, &
+                       source, revision, revdate, text2=string2)
+      END SELECT
+
       do k = 1, dimlens(3)
          do j = 1, dimlens(2)
             do i = 1, dimlens(1)
+               levels(  indx) = levtot(k)
                lonixy(  indx) = i
                latjxy(  indx) = j
                landarea(indx) = AREA2D(i,j) * LANDFRAC2D(i,j)
@@ -1767,6 +1775,20 @@
 
    case ( QTY_SOIL_MOISTURE )
 
+      ! only clm (history) variable 'H2OSOI' is actually soil moisture.
+      ! If this is part of the state - get it, if not, construct it from
+      ! QTY_SOIL_LIQUID_WATER and QTY_SOIL_ICE ... beware units.
+      ! In the history file, the units are mm3/mm3
+      ! In the restart file, the units are kg/m2 ... thanks ...
+
+      call get_grid_vertval(state_handle, ens_size, location, QTY_SOIL_MOISTURE, interp_val_liq, istatus)
+      if (any(istatus == 0)) then
+         where(istatus == 0) expected_obs = interp_val_liq
+         return
+      endif
+
+      istatus = 0
+
       ! TJH FIXME : make sure this is consistent with the COSMOS operator
       ! This is terrible ... the COSMOS operator wants m3/m3 ... CLM is kg/m2
       call get_grid_vertval(state_handle, ens_size, location, QTY_SOIL_LIQUID_WATER, interp_val_liq, istatus_liq)
@@ -2147,6 +2169,8 @@
 below         = 0.0_r8
 area_below    = 0.0_r8
 area_above    = 0.0_r8
+value_above   = MISSING_R8
+value_below   = MISSING_R8
 
 ELEMENTS : do indexi = index1, indexN
 
@@ -2199,7 +2223,7 @@
    write(string1, *)'Variable '//trim(varname)//' has peculiar interpolation problems.'
    write(string2, *)'uneven number of values "above" and "below"'
    write(string3, *)'counter_above == ',counter_above,' /= ',counter_below,' == counter_below'
-   call error_handler(E_MSG,routine, string1, &
+   call error_handler(E_ERR,routine, string1, &
                   text2=string2,text3=string3)
    deallocate(counter_above, counter_below, above, below, area_above, area_below)
    istatus = 22
@@ -2216,11 +2240,13 @@
       area_above(imem, :) = area_above(imem, :) / total_area(imem)
       value_above(imem) = sum(above(imem, :) * area_above(imem, :))
    else
-      write(string1, *)'Variable '//trim(varname)//' had no viable data above'
-      write(string2, *)'at gridcell lon/lat/level = (',gridloni,',',gridlatj,',',levelabove,')'
-      call write_location(0,location,charstring=string3)