MODULE model_mod (CAM)

DART project logo

Jump to DART Documentation Main Index
version information for this file:
$Id: model_mod.html 6382 2013-08-07 20:45:16Z nancy $

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

Overview

The DART system supports data assimilation into CAM -- The Community Atmosphere Model. It is being actively used by graduate students, post-graduates, and scientists at universities and research labs. Others are using analyses for their time period and resolution of interest, produced here at NCAR using DART-CAM. In addition to the standard DART features described elsewhere, current capabilities include the ability to:

In addition to the standard DART package there is a collection of initial condition files at the large file website http://www.image.ucar.edu/pub/DART/CAM/ that are helpful for interfacing CAM with DART.

Sample sets of observations, which can be used with DART-CAM assimilations, can be found at http://www.image.ucar.edu/pub/DART/Obs_sets/ of which the NCEP BUFR observations are the most widely used.

Experience on a variety of machines has shown that it is a very good idea to make sure your run-time environment has the following:

limit stacksize unlimited
limit datasize unlimited

This page contains the documentation for the DART interface module for the CAM model, using the either the Eulerian and Finite Volume dynamical core, soon the HOMME core on the cubed sphere, and possibly the Semi-Lagrangian core. It is designed to work for versions of CAM from 2.0.1 through the currently released version(s) by specifying the version number in the model_mod namelist. This implementation of DART-CAM uses the CAM initial files for transferring the model state to/from the filter/CAM. The DART interface model was designed to facilitate using any CAM version with minimal changes to the model_mod code and accompanying scripts.

The DART interfaces to CAM and the rest of the CESM components have been integrated with the CESM run scripts. Unlike previous versions of DART-CAM, CESM runs using its normal scripts, then stops and calls a DART script which runs a single assimilation step, then returns to the CESM run script to continue the model advances. See the CESM interface documentation for more information on running DART with CESM.

The DART state vector should include any variable in the CAM initial files which is a fundamental quantity and required to reconstruct any other variable that is carried along in the model files. In practice the state vector sometimes contains derived quantities to enable DART to compute forward operators (expected observation values) efficiently. The derived quantities are often overwritten when the model runs the next timestep, so the work DART does to update them is waster work.

Expected observation values on pressure, height or model levels can be requested from model_interpolate. Surface observations can not yet be interpolated, due to the difference between the model surface and the earth's surface where the observations are made. Model_interpolate can be queried for any (non-surface) variable in the state vector (which are variables native to CAM) plus pressure on height levels. The default state vector is PS, T, U, V, Q, CLDLIQ, CLDICE. Variables which are not in the initial file can be added, but minor modifications to model_mod.f90 (and perhaps CAM) may be necessary.

The 16 public interfaces in model_mod are standardized for all DART compliant models. These interfaces allow DART to get the model state and metadata describing this state, find state variables that are close to a given location, and do spatial interpolation for a variety of variables required by observational operators.

Help with several aspects of setting up assimilation experiments using CAM and parallel architecture computers is available:

[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
   output_state_vector = .false.,
   model_version       = '3.0',
   model_config_file   = 'caminput.nc',
   cam_phis            = 'cam_phis.nc',
   state_num_0d        = 0,
   state_num_1d        = 0,
   state_num_2d        = 1,
   state_num_3d        = 4,
   state_names_0d      = '',
   state_names_1d      = '',
   state_names_2d      = 'PS'
   state_names_3d      = 'T', 'U', 'V', 'Q',
   which_vert_1d       = -2 ,
   which_vert_2d       = -1 , 
   which_vert_3d       = 1  ,
   pert_names          = '',
   pert_sd             = -888888.0,
   pert_base_vals      = -888888.0,
   highest_obs_pressure_mb   = 150.0,
   highest_state_pressure_mb = 150.0,
   max_obs_lat_degree  = 90.0,
   Time_step_seconds   = 21600,
   Time_step_days      = 0,
   impact_only_same_kind = '',
   print_details       = .false.
   /


The specification of lists of names and numbers for the various dimensions enables the very flexible definition of the state vector. It can be done via the namelist, instead of recompiling DART for each different set. One hurdle that still remains is that distinct filter initial condition files are necessary for each distinct set of fields which compose the state vector.

The dimension of these lists is currently hardwired to size 100. If more fields need to be assimilated (e.g. many chemical species), look for the integer parameter MAX_STATE_NAMES in the source code and change it to a long enough value and recompile DART. Longer term we intend to investigate using 2 different namelists inside model_mod; one for setting the length of the lists and another to actually read in the data which fills the lists.

The values for which_vert_#d is described in DART/location/threed_sphere/location_mod.html.

The names of the fields to put into the state vector come from the CAM initial file field names.

Item Type Description
output_state_vector logical Controls the output to netCDF files. If .true., output the raw dart state vector. If .false., output the prognostic flavor (gridded data) for easier plotting (recommended).
model_version character(len=128) The number of the CAM version being used, i.e. '3.0.7'. (no letters allowed, so rename cam3_0_p1)
model_config_file character(len=128) CAM initial file used to provide configuration information.
state_num_#d,
#=0,1,2,3
integer Numbers of fields of various dimensions to put into the state vector.
state_names_#d,
#=0,1,2,3
character(len=8), dimension(100) Names of fields of various dimensions to put into the state vector.
which_vert_#d,
#=1,2,3
integer, dimension(100) Vertical location types of fields in state_names_#d. See the 3D sphere location documentation for the mapping of integer values to vertical location types.
pert_names character(len=8), dimension(100) To make filter generate an ensemble from a single model state by randomly perturbing it, list the field(s) to be perturbed here. To make trans_pv_sv_pert0 reset a single field to a constant value, list that field here. Trans_pv_sv_pert0 would be run as many times as there are ensemble members (see pert_base_vals), in order to provide spread to that variable.
pert_sd real(r8), dimension(100) If positive, it's the standard deviation of the perturbation for each field in the pert_names list (filter). If negative, then pert_names can contain only one entry, and that field will be set to a different constant value for each ensemble member (trans_pv_sv_pert0). Those values come from pert_base_vals. Defaults to a MISSING real value and unused unless pert_names is set.
pert_base_vals real(r8), dimension(100) If pert_sd is negative, this is the list of values to use for each ensemble member when perturbing the single field named in pert_names. Otherwise, it's the list of values to which the field(s) listed in pert_names will be reset if filter is told to create an ensemble from a single state vector. Defaults to a MISSING real value and unused unless pert_names is set and pert_base_vals is /= this missing value (-888888.0d0).
max_obs_lat_degree real(r8) Observations closer to the poles than this latitude will be ignored.
highest_obs_pressure_mb real(r8) Observations higher than this pressure are ignored.
highest_state_pressure_mb real(r8) Influence of all obs on model points higher than this is reduced.
Time_step_seconds real(r8) Minimum forecast duration (the part < 1 day)
Time_step_days real(r8) Minimum forecast duration (the part > 24*3600 sec)
impact_only_same_kind character(len=32) Name of one observation kind which can only affect state variables of the same kind.
print_details logical If true, print out detailed information about the sizes, shapes, offsets, etc of items in the CAM state vector. If false, just print out the names of the fields as they are read into the state vector.


[top]

OTHER MODULES USED

obs_kind_mod
random_nr_mod
random_seq_mod
threed_sphere/location_mod
time_manager_mod
types_mod
utilities_mod
mpi_utilities
[top]

PUBLIC INTERFACES

FILTER INTERFACE

Here they are listed alphabetically, except for the last 4 interfaces, which are only required for low-order models where advancing the model can be done by a call to a subroutine. The last 4 interfaces only appear as stubs in the CAM module. Following this list they are listed in the order in which they appear in model_mod.f90, which is separated into sections corresponding to the main functions of model_mod.f90: static_init_model, I/O, model_interpolate, vector<->field translations, get_close_obs, and utility routines.

use model_mod, only : ens_mean_for_model
 get_close_maxdist_init
 get_close_obs_init
 get_close_obs
 get_model_time_step
 get_model_size
 get_state_meta_data
 model_interpolate
 nc_write_model_atts
 nc_write_model_vars
 pert_model_state
 static_init_model
 adv_1step
 end_model
 init_time
 init_conditions

The interface pert_model_state is provided for CAM, and it allows each field of the state vector to be randomly perturbed using a separate standard deviation. Even this may not be flexible enough to handle variables such as specific humidity, which can vary by orders of magnitude from the surface to the top of the model.

OTHER PUBLIC INTERFACES

The list of interfaces above is part of the list of public routines in this module. The rest of the public list are used by programs other than filter; dart_to_cam, cam_to_dart, etc.

use model_mod, only : prog_var_to_vector
 read_cam_init
 vector_to_prog_var
 write_cam_init

Namelist interface &model_nml is read from file input.nml.



call static_init_model( )

Used for runtime initialization of the model. This is the first call made to the model by any DART compliant assimilation routine. It reads the model_mod namelist parameters, sets the calendar type (the GREGORIAN calendar is used with the CAM model), and determines the dart vector length, among other things. This subroutine requires that caminput.nc (or the name in namelist variable model_config_file) be present in the working directory to retrieve model information (grid dimensions and spacing including the vertical hybrid coordinate coefficients, time step, and Gaussian weights, etc).

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

Returns metadata about a given element, indexed by index_in, in the model state vector. The location defines where the state variable is located while the type of the variable (for instance temperature, or u wind component) is returned by var_type. The integer values used to indicate different variable types in var_type are themselves defined as public interfaces to model_mod if required.

index_in Index into the long state vector.
location Returns location of indexed state variable. The location should use a location_mod that is appropriate for the model domain. For realistic atmospheric models, for instance, a three-dimensional spherical location module that can represent height in a variety of ways is provided.
var_type Returns the type of the indexed state variable as an optional argument.

call ens_mean_for_model(ens_mean)
 real(r8), dimension(:), intent(in)  :: ens_mean
 

saves a local copy of the ensemble means which can be used for computing vertical heights, for example.

ens_mean    Ensemble mean state vector

model_size = get_model_size( )
integer                               ::  get_model_size 

Returns the length of the model state vector as an integer. This includes all nested domains.

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

Returns the forecast length to be used as the "model base time step" in the filter. The choice of initial files instead of restart files restricts the assimilation from using "time steps" which are too small. This is because CAM uses a leapfrog time scheme, but the initial files store only a single timestep. The initial forecast step is a simple forward step, so if repeated short forecasts (< ~9 model time steps) are made, the model becomes unstable. In the long run, a more general extended interface may be required that specifies the models range of time stepping possibilities.

ierr = nc_write_model_atts(ncFileId)
integer,                              ::  nc_write_model_atts 
integer,                  intent(in)  ::  ncFileId 

Function to write model specific attributes to a netCDF file. At present, DART is using the NetCDF format to output diagnostic information. This is not a requirement, and models could choose to provide output in other formats. This function writes the metadata associated with the model to a NetCDF file opened to a file identified by ncFileID.

ncFileId     Integer file descriptor opened to NetCDF file.
ierr Returned error code.

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 

Writes a copy of the state variables to a NetCDF file. Multiple copies of the state for a given time are supported, allowing, for instance, a single file to include multiple ensemble estimates of the state.

ncFileID     Integer file descriptor to opened NetCDF file.
statevec State vector.
copyindex Integer index to which copy is to be written.
timeindex Integer index of which time in the file is being written.
ierr Returned error code. success == 0, failure == -1

call read_cam_init(file_name,var)
character(len, intent(in)     :: file_name
type(model_type), intent(out) :: var

Reads state vector fields from a CAM initial file. Fields are specified in model_mod.nml.

file_name    CAM initial file name.
var    Structure to hold the state vector fields read from file_name.

call write_cam_init(file_name,var)
 character, intent(in)                     :: file_name
 type(model_type), allocatable, intent(in) :: var
 

Write fields that have been updated by assimilation to the CAM 'initial' file.

file_name    Name of CAM initial file to which var will be written.
var    Structure containing all the fields of the state vector.

call model_interpolate(x, location, obs_type, interp_val, istatus)
real(r8), dimension(:),   intent(in)  ::  x 
type(location_type),      intent(in)  ::  location 
integer,                  intent(in)  ::  obs_type 
real(r8),                 intent(out) ::  interp_val 
integer,                  intent(out) ::  istatus 

Given model state, returns the value of observation type interpolated to a given location by a method of the model's choosing. Currently observation types: KIND_U_WIND_COMPONENT, KIND_V_WIND_COMPONENT, KIND_SURFACE_PRESSURE, KIND_TEMPERATURE, KIND_SPECIFIC_HUMIDITY, KIND_PRESSURE are supported, but others can be added. KIND_PRESSURE does not have a corresponding field on CAM initial files, but is routinely calculated in CAM using its own subroutines. Interpolation of this field has been incorporated in order to facilitate assimilations of observations which require it, such as GPS radio occultation. If the interpolation is valid, istatus = 0. In the case where the observational operator is not defined at the given location (e.g. the observation is below the lowest model level or above the top level), interp_val is returned as 0.0 and istatus = 1. CAM is highly damped in the upper levels of the model, which has required the exclusion of otherwise valid observations above a certain level, which can be specified in the model_mod namelist variable highest_obs_pressure_mb. Such cases return istatus = 2, and also do the interpolation and return the value, which is NOT used by filter. Eventually such quality control may be moved to another module, but for now is performed in subroutine get_val_pressure.

x     Model state vector.
location Location to which to interpolate.
obs_type Integer indexing which type of observation is to be interpolated.
interp_val Value interpolated to location.
istatus Integer flag indicating the success of the interpolation.

call prog_var_to_vector(var,x)
 type(model_type), allocatable, intent(in)        :: var
 real(r8), allocatable, dimension(:), intent(out) :: x
 

Insert CAM fields (N-D) into DART state vector array (1-D)

var    Structure containing all the fields of the state vector.
x    DART state vector.

call vector_to_prog_var(x,var)
 real(r8), allocatable, dimension(:), intent(in) :: x
 type(model_type), allocatable, intent(out)      :: var
 

Extract CAM fields (N-D) from state vector (1-D).

x    DART state vector.
var    Structure containing all the fields of the state vector.

call get_close_obs(gc, base_obs_loc, base_obs_kind, obs, obs_kind, num_close, close_ind [, dist])
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(:)
integer,              intent(in)  :: obs_kind(:)
integer,              intent(out) :: num_close
integer,              intent(out) :: close_ind(:)
real(r8), optional,   intent(out) :: dist(:)

First calls the location module code and then updates the distances to accomodate damping CAM's highest levels. Pass-through to the 3-D sphere locations module. See get_close_obs() for the documentation of this subroutine.

CAM uses a hybrid vertical coordinate, which requires the surface pressure beneath a point in order to determine the point's vertical coordinate. When the model state vector is divided up among several regions during parallelization, surface pressure points from other regions become unavailable. The ensemble mean state is available to provide complete columns of data for these calculations.

Due to the damping at high levels (see model_interpolate) there is also code which reduces the influence of observations on model points above some altitude. Currently namelist variables highest_obs_pressure_mb and highest_state_pressure_mb control this. All observations on height and level are handled automatically, using the ensemble mean for calculations. For CAM with 26 vertical levels, the influence declines to 0 above model level 6 for highest_obs_pressure_mb = 150 ( ~ model level 12).

call get_close_maxdist_init(gc, maxdist)
type(get_close_type), intent(inout) :: gc
real(r8),             intent(in)    :: maxdist

Pass-through to the 3-D sphere locations module routine of the same name. See location_mod:get_close_maxdist_init() for the documentation of this subroutine.

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)

Pass-through to the 3-D sphere locations module. See get_close_obs_init() for the documentation of this subroutine.

call pert_model_state(state, pert_state, interf_provided)
real(r8), dimension(:),   intent(in)    ::  state 
real(r8), dimension(:),   intent(out)   ::  pert_state 
logical,                  intent(out)   ::  interf_provided 

Given a model state, produces a perturbed model state. This is used to generate ensemble initial conditions perturbed around some control trajectory state when one is preparing to spin-up ensembles. A DART compliant model can choose not to provide an implementation of this algorithm and use the default mechanism in DART by simply returning .false. as a returned value for the interf_provided argument. In this case, DART perturbs the state to generate ensemble members by adding a random sample from a N(0.0, 0.002) distribution independently to each state variable. Models should override this if some structure is required for perturbations or if the magnitude of perturbations in DART is wrong.

This implementation allows each field in the state vector to be randomly perturbed with a separate standard deviation. The fields to be perturbed, and the associated standard deviations are specified in the model_mod namelist using state_names_pert and state_names_sd. The entries in state_names_pert should be in the same order as those in state_names_#d (# = 0,1,2,3 in that order). As in the default, the perturbations sd*N[0,1] are added onto the basic state field, so sd should not be a percentage of the basic state field, but an actual physical size.

This subroutine is also used by trans_pv_sv_pert0 for the off-line perturbation of fields the model considers to be parameters. In this use a new random sequence seed must be provided via an input file ("ens_member") for each ensemble member, since trans_pv_sv_pert0 is executed once for each member and can't keep a series of different seeds intact.

state State vector to be perturbed.
pert_state Perturbed state vector is returned.
interf_provided Return false to have DART perturb state .

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

This operation is not defined for the CAM model. This interface is only required if `synchronous' model state advance is supported (the model is called directly as a Fortran90 subroutine from the assimilation programs). This is generally not the preferred method for large models and a stub for this interface is provided for the CAM model.

x State vector of length model_size.
time Gives time of the initial model state. Needed for models that have real time state requirements, for instance the computation of radiational parameters. Note that DART provides a time_manager_mod module that is used to support time computations throughout the facility.

call end_model( )

Called when use of a model is completed to clean up storage, etc. A stub is provided for the CAM model.

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

Returns the time at which the model will start if no input initial conditions are to be used. This is frequently used to spin-up models from rest, but is not meaningfully supported for the CAM model.

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

Returns default initial conditions for model; generally used for spinning up initial model states. For the CAM model it's just a stub because initial state is always to be provided from input files.

x Model state vector.

[top]

PRIVATE COMPONENTS

Here they are listed alphabetically. Following this list they are listed in the order in which they appear in model_mod.f90, which is separated into sections corresponding to the main functions of model_mod.f90: static_init_model, I/O, model_interpolate, vector<->field translations, get_close_obs, and utility routines.

use model_mod, only : convert_vert
 coord_val
 coord_index
 dcz2
 end_grid_1d_instance
 end_model_instance
 find_name
 get_interp_prof
 get_val
 get_val_level
 get_val_height
 get_val_pressure
 gph2gmh
 gravity
 index_from_grid
 init_grid_1d_instance
 init_model_instance
 map_kinds
 model_heights
 read_cam_coord
 read_cam_horiz
 read_cam_init_size
 read_topog_size
 order_state_fields
 plevs_cam
 set_ps_arrays
 trans_coord
 write_cam_coord

call read_cam_init_size(file_name, num_lons, num_lats, num_levs)
 character(len, intent(in) :: file_name
 integer, intent(out)      :: num_lons
 integer, intent(out)      :: num_lats
 integer, intent(out)      :: num_levs
 

Gets the number of lons, lats and levels from a netcdf CAM initial file

file_name    CAM initial file.
num_lons    Number of longitudes in this CAM resolution.
num_lats    Number of latitudes.
num_levs    Number of data levels.

call trans_coord(ncfileid)
 integer, intent(in) :: ncfileid
 

Rearranges sizes of the coordinates of each variable as found on caminput.nc file (stored in variable f_dim_#d) into the order model_mod wants to see them for the state (s_dim_#d); lev, lon, lat. # = 1,2,3 for the numbers of dimensions of the variables.

ncfileid    The file ID number of the caminput.nc file.

Uses the model_version from input.nml:model_nml to figure out the order of lev, lat, and lon on the caminput.nc file. Designed to be flexible enough to handle staggered and unstaggered U and V, 2D fields which may be (lev,lat) instead of (lon,lat), and unstructured grids which may come from HOMME or the CLM.

call read_topog_size(ncfileid, num_lons, num_lats)
 integer, intent(in)   :: ncfileid
 integer, intent(out)  :: num_nlons
 integer, intent(out)  :: num_nlats
 

Reads the grid size from the cam_phis.nc file, which contains the surface elevation data.

ncfileid    The file ID number of the cam_phis.nc file.
num_nlons    Number of model longitudes.
num_nlats    Number of model latitudes

call read_cam_horiz(ncfileid, var, dim1, dim2, cfield)
 integer, intent(in)                          :: ncfileid
 real(r8), dimension(dim1, dim2), intent(out) :: var
 integer, intent(in)                          :: dim1
 integer, intent(in)                          :: dim2
 character, len=8,  intent(in)                :: cfield
 

Read a 2-D field from a NetCDF file, so far just surface height (phis) from cam_phis.nc.

ncfileid    The file ID number of the caminput.nc file.
var    The surface height values returned to the calling routine.
dim1    The size of the first dimension of var.
dim2    Similarly for 2nd dimension.
cfield    Name of the field (on the file) to be read.

call nc_read_model_atts(att,att_vals,nflds)
 character, intent(in)                :: att
 character, dimension(:), intent(out) :: att_vals
 integer, intent(in)                  :: nflds
 

Reads the value of an attribute for each of the fields in cflds.

att    The name of an attribute of the fields on a CAM initial file.
att_vals    The values which that attribute has for each of the nflds fields needed for the state vector.
nflds    The number of fields in the state vector.

call read_cam_coord(var,idim,cfield)
 type(grid_1d_type), intent(out) :: var
 integer, intent(in)             :: idim
 character, intent(in)           :: cfield
 

Reads a coordinate array and metadata from a CAM initial file.

var    A coordinate array from a CAM initial file.
idim    Length of var.
cfield    Name of var, as found on the CAM initial file.

call init_grid_1d_instance(var, length, num_atts)
 type(grid_1d_type), intent(out) :: var
 integer,            intent(in ) :: length
 integer,            intent(in ) :: num_atts
 

Allocate space for the variable of defined type grid_1d_type. Put vector length and the number of attributes in structure. These are coordinate variables read from caminput.nc, and have lots of characteristics which are helpful to keep together in one structure. The rest of var is filled in read_cam_coord.

var    The coordinate variable to read from caminput.nc.
length    The length of var.
num_atts    The number of attributes of this coordinate, from caminput.nc.
   type grid_1d_type
      private
      character (len=8)            :: label          ! e.g. 'lat     '
      integer                      :: dim_id         ! NetCDF dimension ID from the file
      integer                      :: length         ! number of elements in the coordinate array
      real(r8)                     :: resolution     ! spacing between elements, or 0. for irreg.
      real(r8), pointer            :: vals(:)        ! coordinate values
      integer                      :: num_atts       ! number of NetCDF attributes
      character (len=128), pointer :: atts_names(:)  ! names of those attributes
      character (len=128), pointer :: atts_vals(:)   ! values of those attributes
   end type grid_1d_type
 
call end_grid_1d_instance(var)
 type(grid_1d_type), intent(in) :: var
 

Deallocate the array components of the variable of defined type grid_1d_type.

var    Coordinate variable and metadata.

call order_state_fields(cflds,nflds)
 character, dimension(:), intent(out) :: cflds
 integer, intent(in)                  :: nflds
 

Fills cflds with state_names for use in I/O of caminput.nc Also assigns field TYPE_ variables for use by get_state_meta_data, and other routines.

cflds    Master list of CAM fields to be incorporated in the DART state vector.
nflds    Number of CAM fields in the state vector.

call map_kinds()
 

Makes an array (dart_to_cam_kinds) of 'locations within the state vector' of all the available obs kinds that come from obs_kind_mod. Also maps the local model_mod TYPE_s onto the DART KIND_s by the same mechanism.

The obs kind that's needed will be the index into this array, the corresponding value will be the position of that field (not individual variable) within the state vector according to state_name_Xd. This subroutine is called from static_init_model, and arrays dart_to_cam_kinds and cam_to_dart_kinds are global, so they will not have to be recomputed for every obs.

call write_cam_coord_def(ncFileID, c_name, coord, dim_id, c_id)
  character (len=8),  intent(in)  :: c_name
  integer,            intent(in)  :: ncFileID
  integer,            intent(in)  :: dim_id
  type(grid_1d_type), intent(in)  :: coord
  integer,            intent(out) :: c_id
 

Define a coordinate variable on the P[oste]rior_Diag.nc file and write some of it's attributes there.

c_name    Coordinate name.
ncFileID    Integer ID number of the NetCDF file.
dim_id    The dimension ID of the coordinate on the NetCDF file.
coord    Coordinate values and metadata.
c_id    Coordinate variable ID number on the NetCDF file.

call get_val_level(val,x,lon_index,lat_index,level, obs_kind,istatus)
 real(r8), intent(out)                           :: val
 real(r8), allocatable, dimension(:), intent(in) :: x
 integer, intent(in)                             :: lon_index
 integer, intent(in)                             :: lat_index
 integer, intent(in)                             :: level
 integer, intent(in)                             :: obs_kind
 integer, intent(out)                            :: istatus
 

Gets the value on the model level for variable obs_kind at the lon_index, lat_index horizontal grid point. Some quality control of observations is done here, but may be moved later. For now, observations above the lower of highest_obs_level and the models highest level, and below the lowest model level are excluded. So are PS and Q observations.

val    The value of the obs_kind variable interpolated from state vector x.
x    DART state vector.
lon_index    Index of longitude of this val.
lat_index    Index of latitude of this val.
level    The level for which this obs_kind must be retrieved.
obs_kind    The kind of variable being interpolated, i.e. T.
istatus    Flag containing the status of the interpolation. 0=all good, 1=can't interpolate, 2=can interpolate but don't use the value.

call get_val_pressure(val,x,lon_index,lat_index,pressure, obs_kind,istatus)
 real(r8), intent(out)                           :: val
 real(r8), allocatable, dimension(:), intent(in) :: x
 integer, intent(in)                             :: lon_index
 integer, intent(in)                             :: lat_index
 real(r8), intent(in)                            :: pressure
 integer, intent(in)                             :: obs_kind
 integer0, intent(out)                           :: istatus
 

Gets the vertically interpolated value on pressure for variable obs_kind at the lon_index, lat_index horizontal grid point. Some quality control of observations is done here, but may be moved later. For now, observations above the lower of highest_obs_pressure_mb and the model's highest level, and below the lowest model level are excluded.

val    The value of the obs_kind variable interpolated from state vector x.
x    DART state vector.
lon_index    Index of longitude of this val.
lat_index    Index of latitude of this val.
pressure    The pressure to which this obs_kind must be interpolated.
obs_kind    The kind of variable being interpolated, i.e. T.
istatus    Flag containing the status of the interpolation. 0=all good, 1=can't interpolate, 2=can interpolate but don't use the value.

call get_val_height(val,x,lon_index,lat_index,height, obs_kind,istatus)
 real(r8), intent(out)                           :: val
 real(r8), allocatable, dimension(:), intent(in) :: x
 integer, intent(in)                             :: lon_index
 integer, intent(in)                             :: lat_index
 real(r8), intent(in)                            :: height
 integer, intent(in)                             :: obs_kind
 integer, intent(out)                            :: istatus
 

Gets the vertically interpolated value on height for variable obs_kind at the lon_index, lat_index horizontal grid point. Some quality control of observations is done here, but may be moved later. For now, observations above the lower of highest_obs_height_m and the model's highest level, and below the lowest model level are excluded. So are PS and Q observations.

val    The value of the obs_kind variable interpolated from state vector x.
x    DART state vector.
lon_index    Index of longitude of this val.
lat_index    Index of latitude of this val.
height    The height to which this obs_kind must be interpolated.
obs_kind    The kind of variable being interpolated, i.e. T.
istatus    Flag containing the status of the interpolation. 0=all good, 1=can't interpolate, 2=can interpolate but don't use the value.

call get_val(val,x,lon_index,lat_index,level,obs_kind,istatus)
 real(r8), intent(out)                           :: val
 real(r8), allocatable, dimension(:), intent(in) :: x
 integer, intent(in)                             :: lon_index
 integer, intent(in)                             :: lat_index
 integer, intent(in)                             :: level
 integer, intent(in)                             :: obs_kind
 integer, intent(out)                            :: istatus
 

Extracts the value of a field at a specified location from the DART state vector.

val    The value of the obs_kind variable at a grid point, from state vector x.
x    DART state vector.
lon_index    Index of longitude of this val.
lat_index    Index of latitude of this val.
level    The level of the obs_kind variable desired.
obs_kind    The kind of variable being interpolated, i.e. T.
istatus    Flag containing the status of the interpolation. 0=all good, 1=can't in terpolate, 2=can interpo but don't use the value.

call convert_vert(old_array, old_which, new_array, new_which, dart_kind)
  real(r8), dimension(3), intent(in)    :: old_array
  integer,                intent(in)    :: old_which
  real(r8), dimension(3), intent(inout) :: new_array
  integer,                intent(out)   :: new_which
  integer,                intent(in)    :: dart_kind
 

Uses model information and subroutines to convert the vertical location of an ob (prior, model state variable, or actual ob) into the standard vertical coordinate (pressure).

old_array    The location array of the incoming ob.
old_which    The vertical coordinate type of the incoming ob.
new_array    The location array of the converted ob.
new_which    The vertical coordinate type of the converted ob.
dart_kind    The location array of the incoming ob.

Called by model_mod:get_close_obs. Accommodates staggered grid US and VS of the FV core.

index = index_from_grid(lev_ind, lon_ind, lat_ind, ifld)
  integer, intent(in) :: lev_ind
  integer, intent(in) :: lon_ind
  integer, intent(in) :: lat_ind
  integer, intent(in) :: ifld
 

Function to generate the state vector index corresponding to the grid location and variable given.

lev_ind    Level of the desired variable.
lon_ind    Longitude of the desired variable.
lat_ind    Latitude of the desired variable.
ifld    TYPE_ (PS or T or ...) index of the variable whose index is needed.

Called by many routines.

name = find_name(nam, list)
  character (len=*),              intent(in) :: nam
  character (len=*), dimension(:),intent(in) :: list
 

Function to return the index of a character string as found within a list of character strings, typically variable or dimension names.

nam    The string to be found within list.
list    The list which hopefully contains nam.

call coord_val(dim_name, indx, lon_val, lat_val, lev_val)
  character (len=8), intent(in)    :: dim_name
  integer,           intent(in)    :: indx
  real(r8),          intent(inout) :: lon_val
  real(r8),          intent(inout) :: lat_val
  real(r8),          intent(inout) :: lev_val
 

Given the name of the coordinate to be searched and the index into that array, returns the coordinate value in either lon_val, lat_val, or lev_val.

dim_name    The name of the dimension being searched for the value at an index.
indx    The index of the variable whose location is needed.
lon_val    (Possibly) the location of the grid point in the longitudinal direction.
lat_val    (Possibly) the location of the grid point in the latitudinal direction.
lev_val    (Possibly) the location of the grid point in the vertical direction.

All 3 _val arguments are present so that this routine can return the value in the coordinate that the calling routine wants it to be, and search/placement doesn't have to be done there. Used by get_state_meta_data and model_interpolate.

call coord_index(dim_name, val, indx , other_indx)
  character (len=8), intent(in)  :: dim_name
  real(r8),          intent(in)  :: val
  integer,           intent(out) :: indx
  integer, optional, intent(out) :: other_indx
 

Given the name of the coordinate to be searched and the coordinate value, Returns the index of the closest coordinate value. Optionally returns the next closest index too, which may be < or > the closest.

dim_name    The name of the dimension being searched for the value at an index.
val    The coordinate value whose nearest neighbors are needed.
indx    The index whose value is closest to val.
other_indx    The index whose value is 2nd closest to val.

Used by get_state_meta_data. Uses coordinate metadata to search differently, depending on whether coordinate is regularly spaced or not.

call set_ps_arrays(vec)
  real(r8), intent(in) :: vec(:)
 

Subroutine to put the ensemble average PS into a globally defined array. Also provides PS on the FV staggered grids, if needed.

vec    The ensemble average state vector.

call plevs_cam(ncol,ncold,ps,pmid)
 integer, intent(in)                 :: ncol
 integer, intent(in)                 :: ncold
 real(r8), dimension(:), intent(in)  :: ps
 real(r8), dimension(:), intent(out) :: pmid
 

Define the pressures of the CAM data levels (layer midpoints) from the coordinate definitions and the surface pressure.

ncol    Number of columns of pressures to calculate (DART-CAM uses 1)
ncold    Dimension of ps.
ps    Surface pressure at this latitude and longitude (Pa).
pmid    Pressures at the CAM "midpoint" levels.

call model_heights( x,lon_index,lat_index,model_h,istatus)
 real(r8), dimension(:), intent(in)              :: x
 integer, intent(in)                             :: lon_index
 integer, intent(in)                             :: lat_index
 real(r8), intent(out)                           :: model_h
 integer, intent(out)                            :: istatus
 

This routine calculates geometrical height (m) at mid-layers of the CAM model

x    DART state vector.
lon_index    Index of longitude of this val.
lat_index    Index of latitude of this val.
model_h    geometrical height at midlayer (m).
istatus    Status of the calculation returned to the calling routine.

call get_interp_prof(prof, vec, num_levs, lon_index, lat_index, stagr_lon, stagr_lat, kind_cam, vstatus)
  real(r8), intent(out) :: prof(num_levs)
  real(r8), intent(in)  :: vec(:)
  integer,  intent(in)  :: num_levs
  integer,  intent(in)  :: lon_index
  integer,  intent(in)  :: lat_index
  logical,  intent(in)  :: stagr_lon
  integer,  intent(in)  :: stagr_lat
  integer,  intent(in)  :: kind_cam
  integer,  intent(out) :: vstatus
 

Interpolate the (4) nearest profiles to the needed location for the desired variable.

prof    The profile which results from the interpolation.
vec    The (ensemble average) state vector.
num_levs    Vertical size of the profile.
lon_index    The longitude index of the point requiring a profile.
lat_index    The latitude index of the point requiring a profile.
stagr_lon    Flag to say whether kind_cam must be interpolated in the longitudinal direction.
stagr_lat    Flag to say whether kind_cam must be interpolated in the latitudinal direction.
kind_cam    The variable index whose profile is needed.
vstatus    Reports the status of the calculations back to the calling routine.

So far only called by routine model_heights, for getting T and Q profiles. Handles staggered variables in either or both directions. Pole points should be handled in the calling routine by passing the correct stagr_xx, so that this program can count on having values for all the lat and lon indices referenced.

call dcz2(ps,phis0,tv,hprb,hyba,hybb,kmax,idim,imax,pmln, & hypdln,hyalph,pterm,z2)
 real(r8), dimension(idim), intent(in)           :: ps
 real(r8), dimension(idim), intent(in)           :: phis0
 real(r8), dimension(idim,kmax), intent(in)      :: tv
 real(r8), intent(in)                            :: hprb
 real(r8), dimension(2,kmax+1), intent(in)       :: hyba
 real(r8), dimension(2,kmax+1), intent(in)       :: hybb
 integer, intent(in)                             :: kmax
 integer, intent(in)                             :: idim
 integer, intent(in)                             :: imax
 real(r8), dimension(idim,kmax+1), intent(in)    :: pmln
 real(r8), dimension(idim,kmax), intent(in)      :: hypdln
 real(r8), dimension(idim,kmax), intent(in)      :: hyalph
 real(r8), dimension(idim,kmax), intent(in)      :: pterm
 real(r8), dimension(idim,kmax), intent(out)     :: z2
 

To compute geopotential height for a CCM2 hybrid coordinate vertical slice. Since the vertical integration matrix is a function of latitude and longitude, it is not explicitly computed as for sigma coordinates. The integration algorithm is derived from Boville's mods in the ibm file hybrid 1mods (6/17/88). All vertical slice arrays are oriented top to bottom as in CCM2. This field is on full model levels (aka "midpoints") not half levels. See code for more history.

ps    Surface pressure (pascals).
phis0    Surface geoptential.
tv    Virtual temperature, top to bottom.
hprb    Hybrid base pressure (pascals).
hyba    Hybrid coord coeffs for base pressure.
hybb    Hybrid coord coeffs for surf pressure (in same format as hyba).
kmax    Number of vertical levels.
idim    Longitude dimension.
imax    Num of longitude points to compute.
pmln    Vertical slice scratch space used to hold logs of midpoint pressures.
hpdln    Vertical slice scratch space used to hold log p layer thickness.
hyalph    Vertical slice scratch space used to hold distance from interface to level during vertical integration.
pterm    Vertical scratch space.
z2    Geopotential height, top to bottom.

gph2gmh = gph2gmh(h, lat)
 real(r8), intent(out)                           :: gph2gmh
 real(r8), intent(in)                            :: h
 real(r8), intent(in)                            :: lat
 

Convert a geopotential altitude to mean sea level altitude.

gph2gmh    MSL altitude, in km.
h    geopotential altitude (in km).
lat    latitude in degrees.

call gravity (xlat,alt,galt)
 real(r8), intent(in)                            :: xlat
 real(r8), intent(in)                            :: alt
 real(r8), intent(out)                           :: galt
 

This subroutine computes the Earth's gravity at any altitude and latitude. The model assumes the Earth is an oblate spheriod rotating at a the Earth's spin rate. The model was taken from "Geophysical Geodesy, Kurt Lambeck, 1988". Compute acceleration due to the Earth's gravity at any latitude/altitude author Bill Schreiner 5/95

xlat    latitude in radians.
alt    altitude above the reference ellipsiod, km.
galt    gravity at the given lat and alt, cm/sec.

call init_model_instance(var)
 type(model_type), allocatable, intent(out) :: var
 

Initializes an instance of a cam model state variable; all the fields specified for the state vector.

var    Structure which contains all the fields of various dimensions which are in the DART state vector.
 type model_type
     private
     real(r8), pointer :: vars_0d(:)           ! scalars
     real(r8), pointer :: vars_1d(:, :)        ! vectors
     real(r8), pointer :: vars_2d(:, :, :)     ! 2-D fields
     real(r8), pointer :: vars_3d(:, :, :, :)  ! 3-D fields
 end type model_type
 
call end_model_instance(var)
 type(model_type), allocatable, intent(inout) :: var
 

Ends an instance of a cam model_type state variable

var    Structure which contains all the fields included in the state vector.

[top]

FILES

[top]

REFERENCES

[top]

ERROR CODES and CONDITIONS

RoutineMessageComment
read_cam_init_size
read_cam_init
write_cam_init
nc_read_model_atts
nc_write_model_atts
nc_write_model_vars
read_cam_horiz
read_cam_coord
Various NetCDF-f90 interface error messages From one of the NetCDF calls in the named routine
prog_var_to_vector
vector_to_prog_var
get_state_meta_data
get_val
model_get_close_states
scalar and 1-D vector components of state vector are not coded into this routine Only 2D and 3D fields can be part of state vector so far.
prog_var_to_vector
vector_to_prog_var
indx # and model_size # must be equal indx was tallied during insertion of fields into state vector
model_get_close_states which_vert = # not handled in model_get_close_states See which_vert description in location/threed_sphere/location_mod.html
nc_write_model_atts Time dimension ID # must match Unlimited Dimension ID # NetCDF file writing error
order_state_fields nfld = #, nflds = # must be equal Mismatch(es) of state_names_#d and state_num_#d in model_nml

[top]

KNOWN BUGS

See description of dealing with highly damped top layers in CAM in model_interpolate and model_get_close_states.

GPS observations require that highest_obs_pressure_mb be <= 100 mb (hPa).

[top]

FUTURE PLANS

[top]

Terms of Use

DART software - Copyright 2004 - 2013 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: Kevin Raeder
Revision: $Revision: 6382 $
Source: $URL: https://svn-dares-dart.cgd.ucar.edu/DART/releases/Lanai/models/cam/model_mod.html $
Change Date: $Date: 2013-08-07 14:45:16 -0600 (Wed, 07 Aug 2013) $
Change history:  try "svn log" or "svn diff"