INTERFACE / PUBLIC COMPONENTS / NAMELIST / FILES / REFERENCES / ERRORS / BUGS / PLANS / PRIVATE COMPONENTS

MODULE obs_def_mod

Contact: Jeff Anderson
Reviewers:  
Revision: $Revision: 1.4 $
Release Name: $Name: HEAD $
Change Date: $Date: 2005/06/11 14:24:44 $
Change history: see CVS log

OVERVIEW

Provides an abstraction of the definition of an observation. An observation sequence at a higher level is composed of observation definitions associated with observed values. An observation definition can become quite complicated by using a variety of lower order modules that define metadata associated with an observation. This is the only module in DART that requires preprocessing to exclude unnecessary code. Code for reading, writing, getting the expected value of, and interactive definition of (yes, this sentence is not parallel) obs_def_type's for different kinds of observations need to be turned off when they are not in use since they may involve model variables that are not available for all models. So, obs_def_mod.F90 is stored in the repository and preprocessed to form obs_def_mod.f90 which is compiled by mkmf. The syntax #IFDEF observation_type starts a section of optional code and it is concluded by #ENDIF. See the preprocess program for more details of how this is done. Preprocess reads the obs_def_nml to see what types of observations are being assimilated or evaluated to determine what sections of code to keep.





OTHER MODULES USED

types_mod
utilities_mod
location_mod (depends on model choice)
time_manager_mod
assim_model_mod
Other obs_def_kind modules as required



PUBLIC INTERFACE

use assim_tools_mod, only : obs_def_type
 init_obs_def
 get_obs_def_location
 get_obs_def_kind
 get_obs_def_time
 get_obs_def_error_variance
 set_obs_def_location
 set_obs_def_kind
 set_obs_def_time
 set_obs_def_error_variance
 interactive_obs_def
 write_obs_def
 read_obs_def
 set_radar_obs_def
 destroy_obs_def
 copy_obs_def
 assignment(=)




PUBLIC COMPONENTS





type obs_def_type
   private
   type(location_type)    :: location
   type(obs_kind_type)    :: kind
   type(time_type)        :: time
   real(r8)               :: error_variance
   integer                :: key
end type obs_def_type

Description

Models all that is known about an observation except for actual values. Includes a location, kind, time and error variance.

Component Description
location Location of the observation.
kind The kind of the observation.
time Time of the observation.
error_variance Error variance of the observation.
key Unique identifier for observations of a particular kind.


call init_obs_def(obs_def,location,kind,time,error_variance)
 type(obs_def_type), intent(out) :: obs_def
 type(location_type), intent(in) :: location
 type(obs_kind_type), intent(in) :: kind
 type(time_type), intent(in)     :: time
 real(r8), intent(in)            :: error_variance
 

Description

Creates an obs_def type with location, kind, time and error_variance specified.

obs_def    The obs_def that is created
location    Location for this obs_def
kind    Observation kind for obs_def
time    Time for obs_def
error_variance    Error variance of this observation


var = get_obs_def_location(obs_def)
 type(location_type), intent(out) :: get_obs_def_location
 type(obs_def_type), intent(in)   :: obs_def
 

Description

Returns the location from an obs_def.

get_obs_def_location    Returns location from an obs_def
obs_def    An obs_def


var = get_obs_def_kind(obs_def)
 type(obs_kind_type)            :: get_obs_def_kind
 type(obs_def_type), intent(in) :: obs_def
 

Description

Returns obs_kind from an obs_def.

get_obs_def_kind    Returns the obs_kind from an obs_def
obs_def    An obs_def


var = get_obs_def_time(obs_def)
 type(time_type)                :: get_obs_def_time
 type(obs_def_type), intent(in) :: obs_def
 

Description

Returns time from an obs_def.

get_obs_def_time    Returns time from an obs_def
obs_def    An obs_def


var = get_obs_def_error_variance(obs_def)
 real(r8)                       :: get_obs_def_error_variance
 type(obs_def_type), intent(in) :: obs_def
 

Description

Returns error variance from an obs_def.

get_obs_def_error_variance    Error variance from an obs_def
obs_def    An obs_def


call set_obs_def_location(obs_def, location)
 type(obs_def_type), intent(inout) :: obs_def
 type(location_type), intent(in)   :: location
 

Description

Set the location in an obs_def.

obs_def    An obs_def
location    A location


call set_obs_def_kind(obs_def, kind)
 type(obs_def_type), intent(inout) :: obs_def
 type(obs_kind_type), intent(in)   :: kind
 

Description

Set the kind of observation in an observation definition.

obs_def    An obs_def
kind    An obs_kind


call set_obs_def_time(obs_def, time)
 type(obs_def_type), intent(inout) :: obs_def
 type(time_type), intent(in)       :: time
 

Description

Sets time for an obs_def.

obs_def    An obs_def
time    Time to set


call set_obs_def_error_variance(obs_def, error_variance)
 type(obs_def_type), intent(inout) :: obs_def
 real(r8), intent(in)              :: error_variance
 

Description

Set error variance for an obs_def.

obs_def    An obs_def
error_variance    Error variance


call interactive_obs_def(obs_def)
 type(obs_def_type), intent(inout) :: obs_def
 

Description

Creates an obs_def via input from standard in.

obs_def    An obs_def


call write_obs_def(ifile,obs_def [,fform])
 integer, intent(in)                    :: ifile
 type(obs_def_type), intent(in)         :: obs_def
 character(len=*), optional, intent(in) :: fform
 

Description

Writes an obs_def to file open on channel ifile. Uses format specified in fform or FORMATTED if fform is not present.

ifile    File unit open to output file
obs_def    Observation definition to be written
fform    File format specifier: FORMATTED or UNFORMATTED; default FORMATTED


call read_obs_def(ifile,obs_def [,fform])
 integer, intent(in)                    :: ifile
 type(obs_def_type), intent(inout)      :: obs_def
 character(len=*), optional, intent(in) :: fform
 

Description

Reads an obs_def from file open on channel ifile. Uses format specified in fform or FORMATTED if fform is not present.

ifile    File unit open to output file
obs_def    Observation definition to be read
fform    File format specifier: FORMATTED or UNFORMATTED; default FORMATTED


call set_radar_obs_def(rad_loc,rgate,raz,elev_rad,ae,dir,var,obs_def)
 type(location_type), intent(in)    :: rad_loc
 real(r8), intent(in)               :: rgate
 real(r8), intent(in)               :: raz
 real(r8), intent(in)               :: elev_rad
 real(r8), intent(in)               :: ae
 real(r8), dimension(3), intent(in) :: dir
 real(r8), intent(in)               :: var
 type(obs_def_type), intent(inout)  :: obs_def
 

Description

Creates a radar observation obs_def: currently not fully supported.

rad_loc    Location of radar
rgate    ???
raz    Radar azimuth?
elev_rad    Radar elevation?
ae    ???
dir    Direction of beam??
var    Error variance??
obs_def    Obs_def being created


call destroy_obs_def(obs_def)
 type(obs_def_type), intent(inout) :: obs_def
 

Description

Releases all storage associated with an obs_def and its subcomponents.

obs_def    An obs_def to be released.


call copy_obs_def(obs_def1,obs_def2)
 type(obs_def_type), intent(out) :: obs_def1
 type(obs_def_type), intent(in)  :: obs_def2
 

Description

Copies obs_def2 to obs_def1, overloaded with assignment (=).

obs_def1    obs_def to be copied into
obs_def2    obs_def to be copied from




NAMELIST

We adhere to the F90 standard of starting a namelist with an ampersand '&' and terminating with a slash '/'.

 namelist / obs_def_nml /  &
 assimilate_these_obs, evaluate_these_obs
 

Discussion

Control what observation types are to be assimilated and what types are to be evaluated (forward operators are computed for diagnostics but observations don't impact the state). The observation kinds are selected by using the string values found in the data structure obs_kind_info. Any number of these can be selected to be assimilated or evaluated but not both. This namelist is also used by the preprocess program to select those parts of the code in obs_def_mod.F90 that should be included in the obs_def_mod.f90 that is generated.

This namelist is read in a file called input.nml

Contents Type Description
assimilate_these_obs character(len=129), dimension(max_obs_kinds) Types of observations to be assimilated. Default: 'null'
evaluate_these_obs character(len=129), dimension(max_obs_kinds) Types of observations to be evaluated only. Default: 'null'




FILES




REFERENCES


ERROR CODES and CONDITIONS

RoutineMessageComment
read_obs_def Expected location header "obdef" in input file The format of the input file is not consistent.



KNOWN BUGS




FUTURE PLANS