DART icon indicating copy to clipboard operation
DART copied to clipboard

WRF model_mod - Why not allow last level?

Open lkugler opened this issue 2 years ago • 7 comments

Problem

My WRF simulation has 50 levels, (51 staggered) but debugging the obs_def_rttov_mod.f90 showed that it sees 49 levels and 48 layers. (see image here)

Reason

rttov_mod.f90 determines the number of levels in the COUNTLEVELS loop (https://github.com/NCAR/DART/blob/9c7d829a507311625edcaf0c145a580811fcd513/observations/forward_operators/obs_def_rttov_mod.f90#L3526) and this loop stops at level 49 (wrong) since interpolate returns istatus 99 This is because interpolate checks the bounds and does NOT allow to retrieve the value at the last level (https://github.com/NCAR/DART/blob/9c7d829a507311625edcaf0c145a580811fcd513/models/wrf/model_mod.f90#L2042)

Solution

So I propose to change the code so that it also allows to access the last level (https://github.com/NCAR/DART/commit/25efb6c8ac49e5322300e0a5dd2244dc61e157bb)

This change will be part of the https://github.com/lkugler/DART/tree/rttov13 push request

lkugler avatar Jun 29 '22 09:06 lkugler

when you say "last" level is that the top or bottom of the atmosphere? is it different in some way from the others? not interpolating at the boundary of the model space may have a good reason for being there. asking craig s or glen about this might be a good idea.

nancycollins avatar Jun 29 '22 12:06 nancycollins

Hi, I mean the top-of-atmosphere level, specifically wrf%dom(id)%bt and wrf%dom(id)%bts. https://github.com/NCAR/DART/blob/9c7d829a507311625edcaf0c145a580811fcd513/models/wrf/model_mod.f90#L6704

There is a remark in the code that an observation at integer index value = bt is not allowed. But why? A possible reason I could imagine is that increments on a level with boundary conditions is problematic but whether I want to update this level or not is still my own job to decide, right?

Allowable REAL-VALUED Index Values ==> INTEGER Index Values
  ! In vertical (z) direction
  !                  M_grid ==> [1 bt]         ==> [1 bt)
  !                  W_grid ==> [1.5 bt+0.5]   ==> [1 bts)

lkugler avatar Jun 29 '22 12:06 lkugler

you'll have to look at the code that calls this routine to see if ind == bt is now allowed, would any code try to index ind+1 as it's starting to interpolate. sometimes the code is written to find the lower and upper level bounds and then figure out where the location is between 0 and 1. boundary cases are tedious to test.

nancycollins avatar Jun 29 '22 15:06 nancycollins

The rttov code is using

do until fail:
   model_interpotate(QTY_PRESSURE) 

https://github.com/NCAR/DART/blob/9c7d829a507311625edcaf0c145a580811fcd513/observations/forward_operators/obs_def_rttov_mod.f90#L3528

to get the number of levels.

QTY_PRESSURE in WRF is MU which is a 2D variable, which is then used to calculate the 3D pressure. There is a call to model_pressure_t_distrib for the k levels and Interpolation for the P field at level k+1. This seems overkill for counting the number of levels.

https://github.com/NCAR/DART/blob/9c7d829a507311625edcaf0c145a580811fcd513/models/wrf/model_mod.f90#L2052-L2069

  1. general solution - rttov calls a "get_num_levels" function for a variable.
  2. Hack solution, for model_interpolate() for the level check call a variable on the bts (k+1 level stagger).

hkershaw-brown avatar Jun 29 '22 16:06 hkershaw-brown

sorry I clicked close by accident half way through my comment.

@mjs2369 Marlee if you are following this discussion, this is similar to what you brought up yesterday for TIEGCM, if there are levs and ilevs, how many levels do we have?

hkershaw-brown avatar Jun 29 '22 16:06 hkershaw-brown

to implement a general 'get_num_levels()' routine, you'd have to add a new required entry point to the model_mod interface. the model_mods for all models would have to dummy it out (use some_generic_mod : get_num_levels) or add support for it. there's also the question of models with vertical stagger, and ocean models which may have different numbers of levels at different locations based on the ocean bottom. gets tricky fast. that's why we've been avoiding it up to now. maybe rttov could use the state structure info to query for a real 3d variable before calling interpolate to count levels. if it doesn't already, it could have a count_num_levels() routine that it calls once and then subsequent code could query the saved value.

nancycollins avatar Jun 29 '22 16:06 nancycollins

the other thing that might work is to have a default state variable quantity in the rttov code (temperature, pressure, whatever) for probing for number of model levels. then have a namelist item for a different quantity that a model can set to override it. if changed, they should select one that's cheap to return (e.g. maybe not pressure if it's expensive to compute, or avoid temperature if there is a conversion to sensible T, etc), isn't staggered, etc.

nancycollins avatar Jun 29 '22 22:06 nancycollins