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

MODULE obs_def_1d_state_mod

Contact: Jeff Anderson
Reviewers:  
Revision: $Revision: 1.1 $
Release Name: $Name: $
Change Date: $Date: 2005/10/17 19:47:35 $
Change history: see CVS log

OVERVIEW

Beginning with the I-release of DART, a more flexible. powerful (and complicated) mechanism for incorporating new types of observations is part of DART. The speical obs_def module being described here is used by the program preprocess.f90 to insert appropriate code sections into a DEFAULT_obs_def module and a DEFAULT_obs_kind module to generate the obs_def_mod.f90 and obs_kind_mod.f90 that are used by filer programs.

This is an extended format Fortran90 module that provides the definition for observation types designed for use with idealized low-order models that use the 1D location module and can be thought of as having a state vector that is equally spaced on a 1D cyclic domain. The two observation types are: a straight linear interpolation to a point on a [0,1] domain called a RAW_STATE_VARIABLE; an area-weighted 'integral' of the state variable over some part of the cylic 1D domain. The second type can be used to do idealized studies related to remote sensing observations that are best thought of as weighted integrals of some quantity over a finite volume.

The second observation type is referred to as a RAW_STATE_1D_INTEGRAL. It has an associated half_width and localization type (see cov_cutoff module) and a number of points at which to compute the associated integral by quadrature. The location of the observation defines the center of mass of the integral. The integral is centered around the location and extends outward on each side to 2*half_width. The weight associated with the integral is defined by the weight of the localization function (for instance Gaspari Cohn) selected by the localization in the cov_cutoff module. The number of points are used to equally divide the range for computing the integral by quadrature.

This module was partly motivated to provide a prototype for those developing more complicated specialized observation definition modules. Special observation modules like this contain Fortran90 code PLUS additional specially formatted commented code that is used to guide the preprocessor in constructing the obs_def_mod.f90 and obs_kind_mod.f90. These specially formatted comments are most conveniently placed at the beginning of the module and comprise six sections, each beginning and ending with a special F90 comment line that must be included VERBATIM. The six sections and there specific instance for the 1d_raw_state_mod are:

1. A list of all observation kinds defined by this module and their associated generic kinds (see obs_kind module). The header line is followed by lines that have the observation kind name (an all caps Fortran 90 identifier) and their associated generic kind identifier from the obs_kind module.

! BEGIN DART PREPROCESS KIND LIST
! RAW_STATE_VARIABLE, KIND_RAW_STATE_VARIABLE
! RAW_STATE_1D_INTEGRAL, KIND_1D_INTEGRAL
! END DART PREPROCESS KIND LIST

2. A list of all the use statements that the completed obs_def_mod.f90 must have in order to use the public interfaces provided by this special obs_def module.

! BEGIN DART PREPROCESS USE OF SPECIAL OBS_DEF MODULE
! use obs_def_1d_state_mod, only : write_1d_integral, read_1d_integral, ! interactive_1d_integral, get_expected_1d_integral
! END DART PREPROCESS USE OF SPECIAL OBS_DEF MODULE

3. Case statement entries for each observation kind defined by this special obs_def module stating how to compute the forward observation operator. Note that there must be a case statement entry for each type obs observation defined even if no special action (just put in a continue) is required.

! BEGIN DART PREPROCESS GET_EXPECTED_OBS_FROM_DEF
! case(RAW_STATE_VARIABLE)
! call interpolate(state, location, 1, obs_val, istatus)
! case(RAW_STATE_1D_INTEGRAL)
! call get_expected_1d_integral(state, location, obs_def%key, obs_val, istatus)
! END DART PREPROCESS GET_EXPECTED_OBS_FROM_DEF

4. Case statement entries for each observation kind defined by this special obs_def module stating how to read any extra required information from an obs sequence file. Note that there must be a case statement entry for each type obs observation defined even if no special action (just put in a continue) is required.

! BEGIN DART PREPROCESS READ_OBS_DEF
! case(RAW_STATE_VARIABLE)
! continue
! case(RAW_STATE_1D_INTEGRAL)
! call read_1d_integral(obs_def%key, ifile, fileformat)
! END DART PREPROCESS READ_OBS_DEF

5. Case statement entries for each observation kind defined by this special obs_def module stating how to write any extra required information to an obs sequence file. Note that there must be a case statement entry for each type obs observation defined even if no special action (just put in a continue) is required.
! BEGIN DART PREPROCESS WRITE_OBS_DEF
! case(RAW_STATE_VARIABLE)
! continue
! case(RAW_STATE_1D_INTEGRAL)
! call write_1d_integral(obs_def%key, ifile, fileformat)
! END DART PREPROCESS WRITE_OBS_DEF

6. Case statement entries for each observation kind defined by this special obs_def module stating how to interactively create any extra required information. Note that there must be a case statement entry for each type obs observation defined even if no special action (just put in a continue) is required.

! BEGIN DART PREPROCESS INTERACTIVE_OBS_DEF
! case(RAW_STATE_VARIABLE)
! continue
! case(RAW_STATE_1D_INTEGRAL)
! call interactive_1d_integral(obs_def%key)
! END DART PREPROCESS INTERACTIVE_OBS_DEF





OTHER MODULES USED

types_mod
utilities_mod
location_mod (1d_location_mod_only)
time_manager_mod
assim_model_mod
cov_cutoff_mod



PUBLIC INTERFACE

use obs_def_mod, only : write_1d_integral
 read_1d_integral
 interactive_1d_integral
 get_expected_1d_integral




PUBLIC COMPONENTS



call write_1d_integral(key,ifile,fileformat)
 integer, intent(in)           :: key
 integer, intent(in)           :: ifile
 character(len=32), intent(in) :: fileformat
 

Description

Writes out the extra information for observation with unique identifier key for a 1d_integral observation kind. This information includes a header which identifies the total numbe of 1d_integral observations in the whole sequence (this is only written on the first call) followed by the half-width, localization type and number of quadrature points for this observation,

key    Unique integer key associated with the observation being processed.
ifile    Unit number on which observation sequence file is open
fileformat    String noting whether file is opened for formatted or unformatted IO.


call read_1d_integral(key,ifile,fileformat)
 integer, intent(out)           :: key
 integer, intent(in)           :: ifile
 character(len=32), intent(in) :: fileformat
 

Description

Reads the extra information for observation with unique identifier key for a 1d_integral observation kind. This information includes a header which identifies the total numbe of 1d_integral observations in the whole sequence (this is only read on the first call) followed by the half-width, localization type and number of quadrature points for this observation,

key    Unique integer key associated with the observation being processed.
ifile    Unit number on which observation sequence file is open
fileformat    String noting whether file is opened for formatted or unformatted IO.


call interactive_1d_integral(key)
 integer, intent(out) :: key
 

Description

Uses input from standard in to define the characteristics of a 1D integral observation. The key that is returned is uniquely associated with the definition that has been created and can be used by this module to keep track of the associated parameters (half_width, localization kind, number of quadrature points) for this key.

key    Unique identifier associated with the created observation definition in the obs sequence.


call get_expected_1d_integral(state,location,key,val,istatus)
 real(r8), intent(in)            :: state
 type(location_type), intent(in) :: location
 integer, intent(in)             :: key
 real(r8), intent(out)           :: val
 integer, intent(out)            :: istatus
 

Description

Computes the forward observation operator for a 1d integral observation.

state    Model state vector (or extended state vector).
location    Location of this observation.
key    Unique integer key associated with this observation.
val    Returned value of forward observation operator.
istatus    Returned as 0 if forward operator was successfully computed, else retur ned non-zero.




FILES




REFERENCES


ERROR CODES and CONDITIONS

RoutineMessageComment
interactive_1d_integral Not enough space for a 1d_integral_obs. Can only have max_1d_integral_obs (currently _____). There is only room for a fixed number of 1d integral observations. The max number is defined by max_1d_integral_obs. Set this to a larger value if more are needed.



KNOWN BUGS




FUTURE PLANS