MODULE model_mod (for CM1)

DART project logo

Jump to DART Documentation Main Index
version information for this file:
$Id: model_mod.html 11626 2017-05-11 17:27:50Z nancy@ucar.edu $

NAMELIST / INTERFACES / FILES / REFERENCES / ERRORS / PLANS / PRIVATE COMPONENTS / TERMS OF USE

Overview

Cloud Model 1 (CM1) version 18 (CM1r18) is compatible with the DART. CM1 is a non-hydrostatic numerical model in Cartesian 3D coordinates designed for the study of micro- to mesoscale atmospheric phenomena in idealized to semi-idealized simulations. The CM1 model was developed and is maintained by George Bryan at the National Center for Atmospheric Research (NCAR) Mesoscale and Microscale Meteorology Laboratory (MMM). The model code is freely available from the CM1 website: http://www2.mmm.ucar.edu/people/bryan/cm1 and must be downloaded and compiled outside of DART.

This model interface and scripting support were created by Luke Madaus. Thanks Luke!

Several modifications to the CM1 namelist namelist.input are required to produce model output compatible with DART. The values are described here and and example is shown below. Using CM1 output files as a prior ensemble state in DART requires each ensemble member to produce a restart file in netCDF format (&param9: restart_format=2) only containing output at the analysis time (&param9: restart_filetype=2). The only required state variable to be updated is potential temperature (theta) and requires the following namelist setting to ensure this is present in the CM1 restart files: (&param9: restart_file_theta = .true., restart_use_theta = .true.).

Additional state variables that have been tested within DART include ua, va, wa, ppi, u0, v0, u10, v10, t2, th2, tsk, q2, psfc, qv, qc, qr, qi qs, and qg. At present, observation times are evaluated relative to the date and time specified in section &param11. Observation locations are specified in meters relative to the domain origin as defined in &param2: iorigin.

&param9 
   restart_format     = 2         restart needs to be netCDF
   restart_filetype   = 2         restart must be the analysis time - ONLY
   restart_file_theta = .true.    make sure theta is in restart file
   restart_use_theta  = .true.
  /

About testing CM1 and DART

There are two sets of scripts in the shell_scripts directory. Luke contributed a set written in python, and the DART team had a set written in csh. The csh scripts have not been tested in quite some time, so use with the understanding that they will need work. Those csh scripts and some unfinished python scripts reside in a shell_scripts/unfinished directory and should be used with the understanding that they will need work.

Strategy and Instructions for using the python scripts.

There are some prerequisites:

  1. CM1 is required to use netCDF restart files.
  2. A collection of CM1 model states for initial conditions will be available.
  3. There is a separate observation sequence file for each assimilation time.
  4. The DART input.nml file has some required values as defined below.
  5. Each time CM1 is advanced, it will start from the same filename, and the restart number in that filename will be 000001 - ALWAYS. That filename will be a link to the most current model state.

Testing a cycling experiment.

The big picture: three scripts ( setup_filter.py, run_filter.py, and advance_ensemble.py ) are alternated to configure an experiment, perform an assimilation on a set of restart files, and make the ensemble forecast. Time management is controlled through command-line arguments.

It is required that you have generated the DART executables before you test. The term {centraldir} refers to a filesystem and directory that will be used to run the experiment, the working directory. {centraldir} should have a lot of capacity, as ensemble data assimilation will require lots of disk. The term {dart_dir} will refer to the location of the DART source code.

The data referenced in the directories (the initial ensemble, etc.) are provided as a compressed tar file cm1r18_3member_example_data.tar.gz You will have to download the tar file, uncompress it, and modify the scripts to use these directories instead of the example directories in the scripts. You will also have to compile your own cm1 executable.

  1. Set some variables in both shell_scripts/setup_filter.py and shell_scripts/advance_ensemble.py as described below.

  2. In the {dart_dir}/models/cm1/shell_scripts directory, run ./setup_filter.py -d YYYYmmDDHHMMSS -i where YYYYmmDDHHMMSS is the date and time of the first assimilation cycle (the -i option indicates this is the initial setup and extra work will be performed). This will create the working directory {centraldir}, link in required executables, copy in the initial conditions for each member from some predetermined location, copy in the observation sequence file for this assimilation time from some predetermined location, modify namelists, and build a queue submission script in the {centraldir}: run_filter.py.

  3. Change into {centraldir} and verify the contents of run_filter.py. Ensure the assimilation settings in input.nml are correct. Once you are satisfied, submit run_filter.py to the queue to perform an assimilation.

  4. After the assimilation job completes, check to be sure that the assimilation completed successfully, and the archived files requested in the setup_filter.py files_to_archive variable are in {centraldir}/archive/YYYYmmDDHHMMSS

  5. Change into {dart_dir}/models/cm1/shell_scripts and advance the ensemble to the next assimilation time by running ./advance_ensemble.py -d YYYYmmDDHHMMSS -l nnnn where YYYYmmDDHHMMSS is the date of the COMPLETED analysis (the start time for the model) and nnnn is the length of model integration in seconds (the forecast length). (The forecast length option is specified by 'hypen ell' - the lowercase letter L, not the number one.) advance_ensemble.py will submit jobs to the queue to advance the ensemble.

  6. After all ensemble members have successfully completed, run ./setup_filter.py -d YYYYmmDDHHMMSS where YYYYmmDDHHMMSS is the new current analysis time. Note the -i flag is NOT used here, as we do not need to (should not!) re-initialize the entire directory structure.

  7. Change into {centraldir} and submit the run_filter.py script to perform the assimilation.

  8. Go back to step 4 and repeat steps 4-7 for each assimilation cycle until the end of the experiment.

Within the setup_filter.py and advance_ensemble.py scripts, the following variables need to be set between the "BEGIN USER-DEFINED VARIABLES" and "END USER-DEFINED VARIABLES" comment blocks:

variable description
jobname A name for this experiment, will be included in the working directory path.
ens_size Number of ensemble members.
restart_filename The filename for each ensemble member's restart. Highly recommended to leave this as 'cm1out_rst_000001.nc'
window_mins The assimilation window width (in minutes) for each assimilation cycle.
copy The copy command with desired flags for this system.
link The link command with desired flags for this system.
remove The remove command with desired flags for this system.
files_to_archive A list of DART output files to archive for each assimilation cycle. Note that any inflation files generated are automatically carried over.
centraldir Directory (which will be created if setup_filter.py is run in intialization mode) where the assimilation and model advances will take place. Should be on a system with enough space to allow for several assimilation cycles of archived output.
dart_dir Path to the cm1 subdirectory of DART.
cm1_dir Path to the cm1 model executable (cm1.exe)
icdir Path to the ensemble of initial conditions. It is assumed that within this directory, each ensemble member has a subdirectory (m1, m2, m3, ...) that contains:
  • a restart file for cm1 at the desired start time and having the filename defined in restart_filename above
  • a namelist.input file compatible with the generation of that restart file.
obsdir Path to a directory containing observation sequence files to be assimilated. It is assumed that the observation sequence files are named following the convention YYYYmmDDHHMMSS_obs_seq.prior, where the date of the analysis time whose observations are contained in that file is the first part of the file name.
setup_filter.py and advance_ensemble.py assume that mpi queue submissions are required to run cm1.exe and filter. These variables control how that is handled.
queue_system The name of the queueing system
mpi_run_command The command used in a submitted script to execute an mpi task in the queue, including any required flags
queue_sub_command The command used to submit a script to the queue
job_sub_info A dictionary of all flags required to execute a job in the queue, with the key being the flag and the value being the variable. e.g. {'-P' : 'PROJECT CODE HERE', '-W' : '00:20'}, etc.

Additionally, {dart_dir}/work/input.nml should be modified with the desired assimilation settings. Some of the variables listed above will override the values in {dart_dir}/work/input.nml should be modified

[top]

NAMELIST

This namelist is read from the file input.nml. Namelists start with an ampersand '&' and terminate with a slash '/'. Character strings that contain a '/' must be enclosed in quotes to prevent them from prematurely terminating the namelist.

&model_nml 
   assimilation_period_days     = 0
   assimilation_period_seconds  = 21600
   model_perturbation_amplitude = 0.2
   cm1_template_file            = 'null'
   calendar                     = 'Gregorian'
   periodic_x                   = .true.
   periodic_y                   = .true.
   periodic_z                   = .false.
   debug                        = 0
   model_variables              = ' '
  /
Item Type Description
assimilation_period_[days,seconds]   integer This specifies the width of the assimilation window. The current model time is used as the center time of the assimilation window. All observations in the assimilation window are assimilated. BEWARE: if you put observations that occur before the beginning of the assimilation_period, DART will error out because it cannot move the model 'back in time' to process these observations.
model_perturbation_amplitude real(r8) unsupported
cm1_template_file character(len=256) filename used to read the variable sizes, location metadata, etc.
calendar character(len=256) Character string to specify the calendar in use. Usually 'Gregorian' (since that is what the observations use).
model_variables character(:,5) Strings that identify the CM1 variables, their DART quantity, the minimum & maximum possible values, and whether or not the posterior values should be written to the output file. The DART QUANTITY must be one found in the DART/obs_kind/obs_kind_mod.f90 AFTER it gets built by preprocess.
model_variables(:,1) Specifies the CM1 variable name in the netCDF file.
model_variables(:,2) Specifies the DART quantity for that variable.
model_variables(:,3) Specifies a minimum bound (if any) for that variable.
model_variables(:,4) Specifies a maximum bound (if any) for that variable.
model_variables(:,5) Specifies if the variable should be updated in the restart file. The value may be "UPDATE" or anything else.
periodic_x logical a value of .true. means the 'X' dimension is periodic.
periodic_y logical a value of .true. means the 'Y' dimension is periodic.
periodic_z logical unsupported
debug integer switch to control the amount of run-time output is produced. Higher values produce more output. 0 produces the least.


Note: the values above are the default values. A more realistic (useful?) example is shown below and closely matches the values in the default input.nml.

&model_nml 
   assimilation_period_days     = 0
   assimilation_period_seconds  = 60
   cm1_template_file            = 'cm1out_rst_000001.nc'
   calendar                     = 'Gregorian'
   periodic_x                   = .true.
   periodic_y                   = .true.
   periodic_z                   = .false.
   debug                        = 0
   model_variables = 'ua'   , 'QTY_U_WIND_COMPONENT'      , 'NULL', 'NULL', 'UPDATE',
                     'va'   , 'QTY_V_WIND_COMPONENT'      , 'NULL', 'NULL', 'UPDATE',
                     'wa'   , 'QTY_VERTICAL_VELOCITY'     , 'NULL', 'NULL', 'UPDATE',
                     'theta', 'QTY_POTENTIAL_TEMPERATURE' , 0.0000, 'NULL', 'UPDATE',
                     'ppi'  , 'QTY_PRESSURE'              , 'NULL', 'NULL', 'UPDATE',
                     'u10'  , 'QTY_10M_U_WIND_COMPONENT'  , 'NULL', 'NULL', 'UPDATE',
                     'v10'  , 'QTY_10M_V_WIND_COMPONENT'  , 'NULL', 'NULL', 'UPDATE',
                     't2'   , 'QTY_2M_TEMPERATURE'        , 0.0000, 'NULL', 'UPDATE',
                     'th2'  , 'QTY_POTENTIAL_TEMPERATURE' , 0.0000, 'NULL', 'UPDATE',
                     'tsk'  , 'QTY_SURFACE_TEMPERATURE'   , 0.0000, 'NULL', 'UPDATE',
                     'q2'   , 'QTY_SPECIFIC_HUMIDITY'     , 0.0000, 'NULL', 'UPDATE',
                     'psfc' , 'QTY_SURFACE_PRESSURE'      , 0.0000, 'NULL', 'UPDATE',
                     'qv'   , 'QTY_VAPOR_MIXING_RATIO'    , 0.0000, 'NULL', 'UPDATE',
                     'qc'   , 'QTY_CLOUD_LIQUID_WATER'    , 0.0000, 'NULL', 'UPDATE',
                     'qr'   , 'QTY_RAINWATER_MIXING_RATIO', 0.0000, 'NULL', 'UPDATE',
                     'qi'   , 'QTY_CLOUD_ICE'             , 0.0000, 'NULL', 'UPDATE',
                     'qs'   , 'QTY_SNOW_MIXING_RATIO'     , 0.0000, 'NULL', 'UPDATE',
                     'qg'   , 'QTY_GRAUPEL_MIXING_RATIO'  , 0.0000, 'NULL', 'UPDATE'

  /

[top]

OTHER MODULES USED

types_mod
time_manager_mod
location_mod (multiple choices here)
utilities_mod
POSSIBLY MANY OTHERS DEPENDING ON MODEL DETAILS
[top]

PUBLIC INTERFACES

use model_mod, only : get_model_size
 adv_1step
 get_state_meta_data
 model_interpolate
 get_model_time_step
 static_init_model
 init_time
 init_conditions
 nc_write_model_atts
 nc_write_model_vars
 pert_model_copies
 get_close_maxdist_init
 get_close_obs_init
 get_close_obs
 get_close_state_init
 get_close_state
 query_vert_localization_coord
 vert_convert
 read_model_time
 write_model_time
 end_model

A namelist interface &model_nml may be defined by the module, in which case it will be read from file input.nml. The details of the namelist are always model-specific (there are no generic namelist values).

A note about documentation style. Optional arguments are enclosed in brackets [like this].


model_size = get_model_size( )
integer :: get_model_size

Returns the length of the model state vector. Required.

model_size The length of the model state vector.


call adv_1step(x, time)
real(r8), dimension(:), intent(inout) :: x
type(time_type),        intent(in)    :: time

Does a single timestep advance of the model. The input value of the vector x is the starting condition and x must be updated to reflect the changed state after a timestep. The time argument is intent in and is used for models that need to know the date/time to compute a timestep, for instance for radiation computations. This interface is only called if the namelist parameter async is set to 0 in perfect_model_obs or filter or if the program integrate_model is to be used to advance the model state as a separate executable. If one of these options is not going to be used (the model will only be advanced as a separate model-specific executable), this can be a NULL INTERFACE. (The subroutine name must still exist, but it can contain no code and it will not be called.)

x State vector of length model_size.
time    Specifies time of the initial model state.


call get_state_meta_data (state_handle, index_in, location, [, var_type] )
type(ensemble_type), intent(in)  :: state_handle
integer,             intent(in)  :: index_in
type(location_type), intent(out) :: location
integer, optional,   intent(out) ::  var_type 

Given a handle to a state structure and an integer index into the state vector, returns the associated location. A second intent(out) optional argument returns the generic quantity of this item, e.g. QTY_TEMPERATURE, QTY_DENSITY, QTY_SALINITY, QTY_U_WIND_COMPONENT. This interface is required to be functional for all applications.

state_handle    The handle to the state structure containing the state vector about which information is requested.
index_in Index of state vector element about which information is requested.
location The location of state variable element.
var_type The generic kind of the state variable element.


call model_interpolate(state_handle, x, location, itype, obs_val, istatus)
type(ensemble_type),    intent(in)  :: state_handle
real(r8), dimension(:), intent(in)  :: x
type(location_type),    intent(in)  :: location
integer,                intent(in)  :: itype
real(r8),               intent(out) :: obs_val
integer,                intent(out) :: istatus

Given a handle containing information for a state vector, a state vector, a location, and a model state variable kind interpolates the state variable field to that location and returns the value in obs_val. The istatus variable should be returned as 0 unless there is some problem in computing the interpolation in which case a positive value should be returned. The itype variable is one of the KIND parameters defined in the obs_kind_mod.f90 file and defines which generic kind of item is being interpolated. In low-order models that have no notion of kinds of variables this argument may be ignored. For applications in which only perfect model experiments with identity observations (i.e. only the value of a particular state variable is observed), this can be a NULL INTERFACE. Otherwise it is required (which is the most common case).

state_handle    The handle to the state structure containing information about the state vector about which information is requested.
x A model state vector.
location    Location to which to interpolate.
itype Quantity of state field to be interpolated.
obs_val The interpolated value from the model.
istatus Integer value returning 0 for success. Other values can be defined for various failures.


var = get_model_time_step()
type(time_type) :: get_model_time_step

Returns the time step (forecast length) of the model; the smallest increment in time that the model is capable of advancing the state in a given implementation. The actual value may be set by the model_mod namelist (depends on the model). This interface is required for all applications.

var    Smallest time step of model.


call static_init_model()

Called to do one time initialization of the model. As examples, might define information about the model size or model timestep. In models that require pre-computed static data, for instance spherical harmonic weights, these would also be computed here. Can be a NULL INTERFACE for the simplest models.



call init_time(time)
type(time_type), intent(out) :: time

Companion interface to init_conditions. Returns a time that is somehow appropriate for starting up a long integration of the model. At present, this is only used if the perfect_model_obs namelist parameter read_input_state_from_file = .false. If this option should not be used in perfect_model_obs, calling this routine should issue a fatal error.

time    Initial model time.


call init_conditions(x)
real(r8), dimension(:), intent(out) :: x

Returns a model state vector, x, that is some sort of appropriate initial condition for starting up a long integration of the model. At present, this is only used if the perfect_model_obs namelist parameter read_input_state_from_file = .false. If this option should not be used in perfect_model_obs, calling this routine should issue a fatal error.

x    Initial conditions for state vector.


ierr = nc_write_model_atts(ncFileID)
integer, intent(in)  :: ncFileID
logical, intent(out) :: model_mod_writes_state_variables
integer              :: nc_write_model_atts

This routine writes the model-specific attributes to netCDF files that DART creates. This includes coordinate variables and any metadata, but NOT the actual model state vector. models/template/model_mod.f90 contains code that can be used for any model as-is.

The typical sequence for adding new dimensions, variables, attributes:

NF90_OPEN             ! open existing netCDF dataset               
   NF90_redef         ! put into define mode                       
   NF90_def_dim       ! define additional dimensions (if any)     
   NF90_def_var       ! define variables: from name, kind, and dims
   NF90_put_att       ! assign attribute values                    
NF90_ENDDEF           ! end definitions: leave define mode         
   NF90_put_var       ! provide values for variable                
NF90_CLOSE            ! close: save updated netCDF dataset        
ncFileID Integer file descriptor to previously-opened netCDF file.
model_mod_writes_state_variables    logical descriptor to flag whether or not the model_mod will write the actual state variables or the DART-intrinsic routines will write them.
ierr Returns a 0 for successful completion.


ierr = nc_write_model_vars(ncFileID, statevec, copyindex, timeindex)
integer                            :: nc_write_model_vars
integer,                intent(in) :: ncFileID
real(r8), dimension(:), intent(in) :: statevec
integer,                intent(in) :: copyindex
integer,                intent(in) :: timeindex

This routine may be used to write the model-specific state vector (data) to a netCDF file. Only used if model_mod_writes_state_variables = .true.

Typical sequence for adding new dimensions,variables,attributes:

NF90_OPEN             ! open existing netCDF dataset               
   NF90_redef         ! put into define mode                       
   NF90_def_dim       ! define additional dimensions (if any)      
   NF90_def_var       ! define variables: from name, kind, and dims
   NF90_put_att       ! assign attribute values                    
NF90_ENDDEF           ! end definitions: leave define mode         
   NF90_put_var       ! provide values for variable                
NF90_CLOSE            ! close: save updated netCDF dataset         
ncFileID file descriptor to previously-opened netCDF file.
statevec A model state vector.
copyindex    Integer index of copy to be written.
timeindex The timestep counter for the given state.
ierr Returns 0 for normal completion.


call pert_model_copies(state_ens_handle, ens_size, pert_amp, interf_provided)
type(ensemble_type), intent(inout) :: state_ens_handle
integer,             intent(in)    :: ens_size
real(r8),            intent(in)    :: pert_amp
logical,             intent(out)   :: interf_provided

Given an ensemble handle, the ensemble size, and a perturbation amplitude; perturb the ensemble. Used to generate initial conditions for spinning up ensembles. If the model_mod does not want to do this, instead allowing the default algorithms in filter to take effect, interf_provided =&nbps;.false. and the routine can be trivial. Otherwise, interf_provided must be returned as .true.

state_ens_handle handle containing and ensemble of state vectors to be perturbed.
ens_size the number of ensemble members to perturb.
pert_amp the amplitude of the perturbations. The interpretation is based on the model-specific implementation.
interf_provided    Returns false if model_mod cannot do this, else true.


call get_close_maxdist_init(gc, maxdist [, maxdist_list])
type(get_close_type), intent(inout) :: gc
real(r8),             intent(in)    :: maxdist
real(r8), optional,   intent(in)    :: maxdist_list(:)

In distance computations any two locations closer than the given maxdist will be considered close by the get_close_obs() routine. In general this is a PASS-THROUGH ROUTINE. It is listed on the use line for the locations_mod, and in the public list for this module, but has no subroutine declaration and no other code in this module:

use location_mod, only: get_close_maxdist_init

public :: get_close_maxdist_init

The location module has code which stores maxdist in the gc derived type. However, if the model needs to alter the value or wants to supply an alternative implementation it can intercept the call like so:

use location_mod, only: &
        lm_get_close_maxdist_init => get_close_maxdist_init

public :: get_close_maxdist_init

In this case a local get_close_maxdist_init() routine must be supplied. To call the original code in the location module use:

call lm_get_close_maxdist_init(gc, mymaxdist)

This subroutine will be called before get_close_obs_init and get_close_obs.

In most cases the PASS-THROUGH ROUTINE will be used.

gc The get_close_type which stores precomputed information about the locations to speed up searching
maxdist Anything closer than this will be considered close.
maxdist_list(:)    An array of distances - one for each specific type. This allows different types to have different localization behavior.


call get_close_obs_init(gc, num, obs)
type(get_close_type), intent(inout) :: gc
integer,              intent(in)    :: num
type(location_type),  intent(in)    :: obs(num)

In general this is a PASS-THROUGH ROUTINE. The default routine in the location module precomputes information to accelerate the distance computations done by get_close_obs(). Like the other PASS-THROUGH ROUTINES it is listed on the use line for the locations_mod, and in the public list for this module, but has no subroutine declaration and no other code in this module:

use location_mod, only: get_close_obs_init

public :: get_close_obs_init

The location module code bins the list of locations and precomputes maximum possible distances between bins. However, if the model needs to alter the values or wants to supply an alternative implementation it can intercept the call like so:

use location_mod, only: &
        lm_get_close_obs_init => get_close_obs_init
        
public :: get_close_obs_init

In this case a local get_close_obs_init() routine must be supplied. To call the original code in the location module use:

call lm_get_close_obs_init(gc, num, obs)

This subroutine will be called after get_close_maxdist_init and before get_close_obs.

In most cases the PASS-THROUGH ROUTINE will be used.

gc The get_close_type which stores precomputed information about the locations to speed up searching
num The number of items in the third argument
obs    A list of locations which will be part of the subsequent distance computations


call get_close_obs(gc, base_obs_loc, base_obs_kind, obs_loc, obs_kind, num_close, close_ind, dist, state_handle)
type(get_close_type), intent(in)  :: gc
type(location_type),  intent(in)  :: base_obs_loc
integer,              intent(in)  :: base_obs_kind
type(location_type),  intent(in)  :: obs_loc(:)
integer,              intent(in)  :: obs_kind(:)
integer,              intent(out) :: num_close
integer,              intent(out) :: close_ind(:)
real(r8),             intent(out) :: dist(:)
type(ensemble_type),  intent(in)  :: state_handle

Given a location and kind, compute the distances to all other locations in the obs list. The return values are the number of items which are within maxdist of the base, the index numbers in the original obs list, and optionally the distances. The gc contains precomputed information to speed the computations.

In general this is a PASS-THROUGH ROUTINE. It is listed on the use line for the locations_mod, and in the public list for this module, but has no subroutine declaration and no other code in this module:

use location_mod, only: get_close_obs

public :: get_close_obs

However, if the model needs to alter the values or wants to supply an alternative implementation it can intercept the call like so:

use location_mod, only: &
        lm_get_close_obs => get_close_obs
        
public :: get_close_obs

In this case a local get_close_obs() routine must be supplied. To call the original code in the location module use:

call lm_get_close_obs(gc, base_obs_loc, ...)

This subroutine will be called after get_close_maxdist_init and get_close_obs_init.

In most cases the PASS-THROUGH ROUTINE will be used, but some models need to alter the actual distances depending on the observation or state vector kind, or based on the observation or state vector location. It is reasonable in this case to leave get_close_maxdist_init() and get_close_obs_init() as pass-through routines and intercept only get_close_obs(). The local get_close_obs() can first call the location mod routine and let it return a list of values, and then inspect the list and alter or remove any entries as needed. See the CAM and WRF model_mod files for examples of this use.

gc The get_close_type which stores precomputed information about the locations to speed up searching
base_obs_loc Reference location. The distances will be computed between this location and every other location in the obs list
base_obs_kind    The kind of base_obs_loc
obs_loc(:) Compute the distance between the base_obs_loc and each of the locations in this list
obs_kind(:) The corresponding kind of each item in the obs list
num_close The number of items from the obs list which are within maxdist of the base location
close_ind(:) The list of index numbers from the obs list which are within maxdist of the base location
dist(:) If present, return the distance between each entry in the close_ind list and the base location. If not present, all items in the obs list which are closer than maxdist will be added to the list but the overhead of computing the exact distances will be skipped.
state_handle    The handle to the state structure containing information about the state vector about which information is requested.


call get_close_state_init(gc, num, obs)
type(get_close_type), intent(inout) :: gc
integer,              intent(in)    :: num
type(location_type),  intent(in)    :: obs(num)

In general this is a PASS-THROUGH ROUTINE. The default routine in the location module precomputes information to accelerate the distance computations done by get_close_state(). Like the other PASS-THROUGH ROUTINES it is listed on the use line for the locations_mod, and in the public list for this module, but has no subroutine declaration and no other code in this module:

use location_mod, only: get_close_state_init

public :: get_close_state_init

The location module code bins the list of locations and precomputes maximum possible distances between bins. However, if the model needs to alter the values or wants to supply an alternative implementation it can intercept the call like so:

use location_mod, only: &
        loc_get_close_state_init => get_close_state_init
        
public :: get_close_state_init

In this case a local get_close_state_init() routine must be supplied. To call the original code in the location module use:

call loc_get_close_state_init(gc, num, obs)

This subroutine will be called after get_close_maxdist_init and before get_close_state.

In most cases the PASS-THROUGH ROUTINE will be used.

gc The get_close_type which stores precomputed information about the locations to speed up searching
num The number of items in the third argument
obs    A list of locations which will be part of the subsequent distance computations


call get_close_state(gc, base_obs_loc, base_obs_kind, state_loc, state_kind, num_close, close_ind, dist, state_handle)
type(get_close_type), intent(in)  :: gc
type(location_type),  intent(in)  :: base_obs_loc
integer,              intent(in)  :: base_obs_kind
type(location_type),  intent(in)  :: state_loc(:)
integer,              intent(in)  :: state_kind(:)
integer,              intent(out) :: num_close
integer,              intent(out) :: close_ind(:)
real(r8),             intent(out) :: dist(:)
type(ensemble_type),  intent(in)  :: state_handle

Given a location and kind, compute the distances to all other locations in the state_loc list. The return values are the number of items which are within maxdist of the base, the index numbers in the original state_loc list, and optionally the distances. The gc contains precomputed information to speed the computations.

In general this is a PASS-THROUGH ROUTINE. It is listed on the use line for the locations_mod, and in the public list for this module, but has no subroutine declaration and no other code in this module:

use location_mod, only: get_close_state

public :: get_close_state

However, if the model needs to alter the values or wants to supply an alternative implementation it can intercept the call like so:

use location_mod, only: &
        lm_get_close_state => get_close_state
        
public :: get_close_state

In this case a local get_close_state() routine must be supplied. To call the original code in the location module use:

call loc_get_close_state(gc, base_obs_loc, ...)

This subroutine will be called after get_close_maxdist_init and get_close_state_init.

In most cases the PASS-THROUGH ROUTINE will be used, but some models need to alter the actual distances depending on the observation or state vector kind, or based on the observation or state vector location. It is reasonable in this case to leave get_close_maxdist_init() and get_close_state_init() as pass-through routines and intercept only get_close_state(). The local get_close_state() can first call the location mod routine and let it return a list of values, and then inspect the list and alter or remove any entries as needed. See the CAM and WRF model_mod files for examples of this use.

gc The get_close_type which stores precomputed information about the locations to speed up searching
base_obs_loc Reference location. The distances will be computed between this location and every other location in the obs list
base_obs_kind    The kind of base_obs_loc
state_loc(:) Compute the distance between the base_obs_loc and each of the locations in this list
state_kind(:) The corresponding kind of each item in the state_loc list
num_close The number of items from the state_loc list which are within maxdist of the base location
close_ind(:) The list of index numbers from the state_loc list which are within maxdist of the base location
dist(:) If present, return the distance between each entry in the close_ind list and the base location. If not present, all items in the state_loc list which are closer than maxdist will be added to the list but the overhead of computing the exact distances will be skipped.
state_handle    The handle to the state structure containing information about the state vector about which information is requested.


ivert = query_vert_localization_coord()

Returns the integer code for the vertical coordinate system used for localization.

ivert The integer code for the vertical coordinate system.


call vert_convert(state_handle, location, obs_kind, istatus)
type(ensemble_type), intent(in)  :: state_handle
type(location_type), intent(in)  :: location
integer,             intent(in)  :: obs_kind
integer,             intent(out) :: istatus

Converts the state to the desired vertical localization coordinate system. Some models (toy models with no 'real' observations) will not need this. Most (real) models have observations in one or more coordinate systems (pressure, height) and the model is generally represented in only one coordinate system. To be able to interpolate the model state to the observation location, or to compute the true distance between the state and the observation, it is necessary to convert everything to one coodinate system.

state_handle    The handle to that state.
location the desired location
obs_kind the quantity defining the part of state to convert.
istatus Specifies the success or failure of the vertical conversion. If istatus = 0, the conversion was a sucess. Any other value is a failure.


model_time = read_model_time(filename)
character(len=*), intent(in) :: filename
type(time_type)              :: model_time

Reads the valid time of the model state in a netCDF file. There is a default routine in assimilation_code/modules/io/dart_time_io_mod.f90 that can be used as a pass-through. That routine will read the last timestep of a 'time' variable - which is the same strategy used for reading netCDF files that have multiple timesteps in them. If your model has some other representation of time (i.e. it does not use a netCDF variable named 'time') - you will have to write this routine.

ncid handle to an open netCDF file
dart_time    Specifies the (last) time of the model state.


call write_model_time(ncid, dart_time)
integer,          intent(in) :: ncid
type(time_type),  intent(in) :: dart_time

Writes the assimilation time to a netCDF file. There is a default routine in assimilation_code/modules/io/dart_time_io_mod.f90 that can be used as a pass-through. If your model has some other representation of time (i.e. it does not use a netCDF variable named 'time') - you will have to write this routine.

ncid handle to an open netCDF file
dart_time    Specifies the time of the assimilation (the current time step).


call end_model()

Does any shutdown and clean-up needed for model. Can be a NULL INTERFACE if the model has no need to clean up storage, etc.


[top]

FILES

[top]

REFERENCES

  1. Madaus, L. and G. Hakim (2017): "Constraining Ensemble Forecasts of Discrete Convective Initiation with Surface Observations." Mon. Wea. Rev. doi: 10.1175/MWR-D-16-0395.1

[top]

ERROR CODES and CONDITIONS

KNOWN BUGS

none at this time

[top]

FUTURE PLANS

It is likely that a number of additional optional interfaces will be added to the model_mod structure. For instance, hints about how to divide the state vector into regions for parallel assimilation will need to be obtained from the model. It is planned that the interf_provided mechanism used in pert_model_copies will allow those who do not wish to support enhanced interfaces to add NULL interfaces by simply pasting in an interface block.

[top]

PRIVATE COMPONENTS

N/A

[top]

Terms of Use

DART software - Copyright UCAR. This open source software is provided by UCAR, "as is", without charge, subject to all terms of use at http://www.image.ucar.edu/DAReS/DART/DART_download

Contact: dart@ucar.edu
Revision: $Revision: 11626 $
Source: $URL: https://svn-dares-dart.cgd.ucar.edu/DART/releases/Manhattan/models/cm1/model_mod.html $
Change Date: $Date: 2017-05-11 11:27:50 -0600 (Thu, 11 May 2017) $
Change history:  try "svn log" or "svn diff"