MODULE model_mod (GITM)

DART project logo

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



The Global Ionosphere Thermosphere Model (GITM) is a 3-dimensional spherical code that models the Earth's thermosphere and ionosphere system using a stretched grid in latitude and altitude.

The GITM interface for Data Assimilation Research Testbed (DART) is under development. If you wish to use GITM, you are urged to contact us. The original scripts were configured to run on the University of Michigan machine NYX using the Portable Batch System (PBS). We have attempted to extend the scripts to work with both PBS and LSF and are only partway through the process at the time of the Lanai release.

DART does not come with the GITM code. You need to get that on your own. The normal procedure of building GITM creates some resource files that are subsequently needed by DART - just to compile. These include:

  1. GITM/src/ModSize.f90,
  2. GITM/src/ModTime.f90, and (for example)
  3. GITM/src/ModEarth.f90 .

GITM uses binary files for their restart mechanisms, so no metadata is available to confirm the number and order of fields in the file. Care must be used to make sure the namelist-controlled set of variables to be included in the DART state vector is consistent with the restart files. Each variable must also correspond to a DART "KIND"; required for the DART interpolate routines. For example, this is nowhere close to being correct:

   gitm_state_variables = 'Temperature',      'QTY_TEMPERATURE',
                          'eTemperature',     'QTY_TEMPERATURE_ELECTRON',
                          'ITemperature',     'QTY_TEMPERATURE_ION',
                          'iO_3P_NDensityS',  'QTY_DENSITY_NEUTRAL_O3P',
                          'iO2_NDensityS',    'QTY_DENSITY_NEUTRAL_O2',
                          'iN2_NDensityS',    'QTY_DENSITY_NEUTRAL_N2',
                               ...                    ...

These variables are then adjusted to be consistent with observations and stuffed back into the same netCDF restart files. Since DART is an ensemble algorithm, there are multiple restart files for a single restart time: one for each ensemble member. Creating the initial ensemble of states is an area of active research.

DART reads grid information for GITM from several sources. The file specifies the number of latitudes/longitudes per block, and the number of blocks comes from the GITM2/src/ModSize.f90 module. Internal to the DART code, the following variables exist:

real(r8) :: LON(:) longitude array [0, 360)
real(r8) :: LAT(:) latitude array (-90,90)
real(r8) :: ALT(:) altitude array (0,~inf)
integer :: NgridLon the length of the longitude array
integer :: NgridLat the length of the latitude array
integer :: NgridAlt the length of the altitude array


was compiled with the gfortran 4.2.3 compilers and run on a Mac.

The DART components were built with the following mkmf.template settings:

      FC = gfortran
      LD = gfortran
      NETCDF = /Users/thoar/GNU
      INCS = -I${NETCDF}/include
      LIBS = -L${NETCDF}/lib -lnetcdf -lcurl -lhdf5_hl -lhdf5 -lz  -lm
      FFLAGS = -O0 -fbounds-check -frecord-marker=4 -ffpe-trap=invalid $(INCS)

Converting between DART files and GITM restart files

The binary GITM files contain no metadata, so care is needed when converting between DART state variables and GITM files.

There are two programs - both require the list of GITM variables to use in the DART state vector: the model_nml:gitm_state_variables variable in the input.nml file.

gitm_to_dart.f90 converts a set of GITM restart files (there is one restart file per block) bxxxx.rst into a DART-compatible file normally called dart_ics . We usually wind up linking to this static filename.
dart_to_gitm.f90 inserts the DART output into existing GITM restart files. There are two different types of DART output files, so there is a namelist option to specify if the DART file has two time records or just one. If there are two, the first one is the 'advance_to' time, followed by the 'valid_time' of the ensuing state. If there is just one, it is the 'valid_time' of the ensuing state. dart_to_gitm determines the GITM restart file name from the input.nml model_nml:gitm_restart_dirname. If the DART file contains an 'advance_to' time, dart_to_gitm creates a DART_GITM_time_control.txt file which can be used to control the length of the GITM integration.

Simple test.

The simplest way to test the converter is to compile GITM and run a single model state forward using work/ To build GITM ... download GITM and unpack the code into DART/models/gitm/GITM2 and follow these instructions:

cd models/gitm/GITM2
./ -install -compiler=ifortmpif90 -earth
cd ../work
./ 1 1 0 150.0 170.0 1.0


We adhere to the F90 standard of starting a namelist with an ampersand '&' and terminating with a slash '/' for all our namelist input. Consider yourself forewarned that character strings that contain a '/' must be enclosed in quotes to prevent them from prematurely terminating the namelist.

# The list of variables to put into the state vector is here:
# The definitions for the DART kinds are in DART/observations/forward_operators/obs_def*f90
# The order doesn't matter to DART. It may to you.

   gitm_restart_dirname         = 'advance_temp_e1/UA/restartOUT',
   assimilation_period_days     = 0,
   assimilation_period_seconds  = 1800,
   model_perturbation_amplitude = 0.2,
   output_state_vector          = .false.,
   calendar                     = 'Gregorian',
   debug                        = 0,
   gitm_state_variables  = 'Temperature',            'QTY_TEMPERATURE',
                           'eTemperature',           'QTY_TEMPERATURE_ELECTRON',
                           'ITemperature',           'QTY_TEMPERATURE_ION',
                           'iO_3P_NDensityS',        'QTY_DENSITY_NEUTRAL_O3P',
                           'iO2_NDensityS',          'QTY_DENSITY_NEUTRAL_O2',
                           'iN2_NDensityS',          'QTY_DENSITY_NEUTRAL_N2',
                           'iN_4S_NDensityS',        'QTY_DENSITY_NEUTRAL_N4S',
                           'iNO_NDensityS',          'QTY_DENSITY_NEUTRAL_NO',
                           'iN_2D_NDensityS',        'QTY_DENSITY_NEUTRAL_N2D',
                           'iN_2P_NDensityS',        'QTY_DENSITY_NEUTRAL_N2P',
                           'iH_NDensityS',           'QTY_DENSITY_NEUTRAL_H',
                           'iHe_NDensityS',          'QTY_DENSITY_NEUTRAL_HE',
                           'iCO2_NDensityS',         'QTY_DENSITY_NEUTRAL_CO2',
                           'iO_1D_NDensityS',        'QTY_DENSITY_NEUTRAL_O1D',
                           'iO_4SP_IDensityS',       'QTY_DENSITY_ION_O4SP',
                           'iO2P_IDensityS',         'QTY_DENSITY_ION_O2P',
                           'iN2P_IDensityS',         'QTY_DENSITY_ION_N2P',
                           'iNP_IDensityS',          'QTY_DENSITY_ION_NP',
                           'iNOP_IDensityS',         'QTY_DENSITY_ION_NOP',
                           'iO_2DP_IDensityS',       'QTY_DENSITY_ION_O2DP',
                           'iO_2PP_IDensityS',       'QTY_DENSITY_ION_O2PP',
                           'iHP_IDensityS',          'QTY_DENSITY_ION_HP',
                           'iHeP_IDensityS',         'QTY_DENSITY_ION_HEP',
                           'ie_IDensityS',           'QTY_DENSITY_ION_E',
                           'U_Velocity_component',   'QTY_VELOCITY_U',
                           'V_Velocity_component',   'QTY_VELOCITY_V',
                           'W_Velocity_component',   'QTY_VELOCITY_W',
                           'U_IVelocity_component',  'QTY_VELOCITY_U_ION',
                           'V_IVelocity_component',  'QTY_VELOCITY_V_ION',
                           'W_IVelocity_component',  'QTY_VELOCITY_W_ION',
                           'iO_3P_VerticalVelocity', 'QTY_VELOCITY_VERTICAL_O3P',
                           'iO2_VerticalVelocity',   'QTY_VELOCITY_VERTICAL_O2',
                           'iN2_VerticalVelocity',   'QTY_VELOCITY_VERTICAL_N2',
                           'iN_4S_VerticalVelocity', 'QTY_VELOCITY_VERTICAL_N4S',
                           'iNO_VerticalVelocity',   'QTY_VELOCITY_VERTICAL_NO',
                           'f107',                   'QTY_1D_PARAMETER',
                           'Rho',                    'QTY_DENSITY',

This namelist is read in a file called input.nml. This namelist provides control over the assimilation period for the model. All observations within (+/-) half of the assimilation period are assimilated. The assimilation period is the minimum amount of time the model can be advanced, and checks are performed to ensure that the assimilation window is a multiple of the model dynamical timestep.

Contents Type Description
gitm_restart_dirname   character(len=256)   The name of the directory containing the GITM restart files and runtime control information.
assimilation_period_days integer The number of days to advance the model for each assimilation.
assimilation_period_seconds integer In addition to assimilation_period_days, the number of seconds to advance the model for each assimilation.
model_perturbation_amplitude real(r8) Reserved for future use.
output_state_vector logical The switch to determine the form of the state vector in the output netCDF files. If .true. the state vector will be output exactly as DART uses it ... one long array. If .false., the state vector is parsed into prognostic variables and output that way -- much easier to use with 'ncview', for example.
calendar character(len=32) Character string specifying the calendar being used by GITM.
debug integer The switch to specify the run-time verbosity.
0 is as quiet as it gets.
> 1 provides more run-time messages.
> 5 provides ALL run-time messages.
gitm_state_variables character(len=NF90_MAX_NAME)::  
The table that relates the GITM variables to use to build the DART state vector, and the corresponding DART kinds for those variables.






Only a select number of interfaces used are discussed here. Each module has its own discussion of their routines.

Required Interface Routines

use model_mod, only : get_model_size

Unique Interface Routines

use model_mod, only :
use location_mod, only : get_close_obs

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

Required Interface Routines

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

adv_1step is not used for the gitm model. Advancing the model is done through the advance_model script. This is a NULL_INTERFACE, provided only for compatibility with the DART requirements.

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

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 

get_state_meta_data returns metadata about a given element of the DART representation of the model state vector. Since the DART model state vector is a 1D array and the native model grid is multidimensional, get_state_meta_data returns information about the native model state vector representation. Things like the location, or the type of the variable (for instance: temperature, u wind component, ...). 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 of state vector element about which information is requested.
location Returns the 3D location of the indexed state variable. The location_ type comes from DART/location/threed_sphere/location_mod.f90. Note that the lat/lon are specified in degrees by the user but are converted to radians internally.
var_type Returns the type of the indexed state variable as an optional argument. The type is one of the list of supported observation types, found in the block of code starting ! Integer definitions for DART TYPES in DART/assimilation_code/modules/observations/obs_kind_mod.f90

The list of supported variables in DART/assimilation_code/modules/observations/obs_kind_mod.f90 is created by preprocess.

call model_interpolate(x, location, itype, obs_val, istatus)
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 model state, model_interpolate returns the value of the desired observation type (which could be a state variable) that would be observed at the desired location. The interpolation method is either completely specified by the model, or uses some standard 2D or 3D scalar interpolation routines. Put another way, model_interpolate will apply the forward operator H to the model state to create an observation at the desired location.

If the interpolation is valid, istatus = 0. In the case where the observation operator is not defined at the given location (e.g. the observation is below the lowest model level, above the top level, or 'dry'), interp_val is returned as 0.0 and istatus = 1.

x A model state vector.
location    Location to which to interpolate.
itype Integer indexing which type of observation is desired.
obs_val The interpolated value from the model.
istatus Integer flag indicating the success of the interpolation.
success == 0, failure == anything else

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

get_model_time_step returns the forecast length to be used as the "model base time step" in the filter. This is the minimum amount of time the model can be advanced by filter. This is also the assimilation window. All observations within (+/-) one half of the forecast length are used for the assimilation. In the GITM case, this is set from the namelist values for input.nml&model_nml:assimilation_period_days, assimilation_period_seconds.

var    Smallest time step of model.

call static_init_model()

static_init_model is called for runtime initialization of the model. The namelists are read to determine runtime configuration of the model, the grid coordinates, etc. There are no input arguments and no return values. The routine sets module-local private attributes that can then be queried by the public interface routines.

See the GITM documentation for all namelists in gitm_in . Be aware that DART reads the GITM &grid_nml namelist to get the filenames for the horizontal and vertical grid information as well as the topography information.

The namelists (all mandatory) are:
gitm_in&domain_nml, and

call end_model()

end_model is used to clean up storage for the model, etc. when the model is no longer needed. There are no arguments and no return values. The grid variables are deallocated.

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

init_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 GITM model. The only time this routine would get called is if the input.nml&perfect_model_obs_nml:start_from_restart is .false., which is not supported in the GITM model.

time    the starting time for the model if no initial conditions are to be supplied. This is hardwired to 0.0

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

init_conditions returns default initial conditions for model; generally used for spinning up initial model states. For the GITM model it is just a stub because the initial state is always provided by the input files.

x    Initial conditions for state vector. This is hardwired to 0.0

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

nc_write_model_atts writes model-specific attributes to an opened netCDF file: In the GITM case, this includes information like the coordinate variables (the grid arrays: ULON, ULAT, TLON, TLAT, ZG, ZC, KMT, KMU), information from some of the namelists, and either the 1D state vector or the prognostic variables (SALT,TEMP,UVEL,VVEL,PSURF). All the required information (except for the netCDF file identifier) is obtained from the scope of the model_mod module. Both the input.nml and gitm_in files are preserved in the netCDF file as variables inputnml and gitm_in, respectively.

ncFileID    Integer file descriptor to previously-opened netCDF file.
ierr Returns a 0 for successful completion.

nc_write_model_atts is responsible for the model-specific attributes in the following DART-output netCDF files:,, and

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

nc_write_model_vars 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. Whether the state vector is parsed into prognostic variables (SALT, TEMP, UVEL, VVEL, PSURF) or simply written as a 1D array is controlled by input.nml&model_mod_nml:output_state_vector. If output_state_vector = .true. the state vector is written as a 1D array (the simplest case, but hard to explore with the diagnostics). If output_state_vector = .false. the state vector is parsed into prognostic variables before being written.

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_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, pert_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. Since the DART state vector for the GITM model contains both 'wet' and 'dry' cells, it is imperative to provide an interface to perturb just the wet cells (interf_provided == .true.).

The magnitude of the perturbation is wholly determined by input.nml&model_mod_nml:model_perturbation_amplitude and utterly, completely fails.

A more robust perturbation mechanism is needed. Until then, avoid using this routine by using your own ensemble of initial conditions. This is determined by setting input.nml&filter_nml:start_from_restart = .false.

state State vector to be perturbed.
pert_state The perturbed state vector.
interf_provided    Because of the 'wet/dry' issue discussed above, this is always .true., indicating a model-specific perturbation is available.

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. See 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 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), dimension(:), intent(in ) :: obs
integer,             dimension(:), intent(in ) :: obs_kind
integer,                           intent(out) :: num_close
integer,             dimension(:), intent(out) :: close_ind
real(r8), optional,  dimension(:), intent(out) :: dist

Given a DART location (referred to as "base") and a set of locations, and a definition of 'close' - return a subset of locations that are 'close', as well as their distances to the DART location and their indices. This routine intentionally masks a routine of the same name in location_mod because we want to be able to discriminate against selecting 'dry land' locations.

Given a single location and a list of other locations, returns the indices of all the locations close to the single one along with the number of these and the distances for the close ones. The list of locations passed in via the obs argument must be identical to the list of obs passed into the most recent call to get_close_obs_init(). If the list of locations of interest changes, get_close_obs_destroy() must be called and then the two initialization routines must be called before using get_close_obs() again.

For vertical distance computations, the general philosophy is to convert all vertical coordinates to a common coordinate. This coordinate type is defined in the namelist with the variable "vert_localization_coord".

gc Structure to allow efficient identification of locations 'close' to a given location.
base_obs_loc Single given location.
base_obs_kind Kind of the single location.
obs List of candidate locations.
obs_kind Kind associated with candidate locations.
num_close Number of locations close to the given location.
close_ind Indices of those locations that are close.
dist Distance between given location and the close ones identified in close_ind.

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

ens_mean_for_model normally saves a copy of the ensemble mean to module-local storage. This is a NULL_INTERFACE for the GITM model. At present there is no application which requires module-local storage of the ensemble mean. No storage is allocated.

ens_mean State vector containing the ensemble mean.


Unique Interface Routines

call get_gridsize( num_x, num_y, num_z )
integer, intent(out) :: num_x, num_y, num_z

get_gridsize returns the dimensions of the compute domain. The horizontal gridsize is determined from

num_x The number of longitudinal gridpoints.
num_y The number of latitudinal gridpoints.
num_z The number of vertical gridpoints.

call restart_file_to_statevector(dirname, state_vector, model_time)
character(len=*),       intent(in)    :: dirname
real(r8), dimension(:), intent(inout) :: state_vector
type(time_type),        intent(out)   :: model_time

restart_file_to_statevector Reads a GITM netCDF format restart file and packs the desired variables into a DART state vector. The desired variables are specified in the input_nml model_nml namelist.

dirname The name of the directory containing the GITM restart files. The restart files have names like b0001.rst. There are as many restart files as blocks.
state_vector the 1D array containing the concatenated GITM variables.
model_time the time of the model state as read from <dirname>/header.rst.

call statevector_to_restart_file(state_vector, dirname, statedate)
real(r8), dimension(:), intent(in) :: state_vector
character(len=*),       intent(in) :: dirname
type(time_type),        intent(in) :: statedate

statevector_to_restart_file updates the variables in the GITM restart files with values from the DART vector state_vector. The time in the GITM files must match the statedate.

state_vector the 1D array containing the DART state vector.
dirname The name of the directory containing the GITM restart files. The restart files have names like b0001.rst. There are as many restart files as blocks.
statedate the 'valid_time' of the DART state vector.

call get_gitm_restart_dirname( dirname )
character(len=*), intent(out) :: dirname

get_gitm_restart_dirname returns the name of the directory containing the GITM restart files.

filename    The name of the GITM restart directory as specified by the input.nml:model_nml:gitm_restart_dirname varible.

time = get_state_time( dirname )
character(len=*), intent(in)  :: dirname
type(time_type),  intent(out) :: time

get_state_time extracts the time of the model state as contained in the <gitm_restart_dirname>/header.rst file.



filename purpose
input.nml to read the model_mod namelist
Several GITM source modules: ModConstants, ModSizeGitm, ModEarth ... provides grid dimensions, model state, and 'valid_time' of the model state
header.rst, bNNNN.rst provides the 'valid_time' of the model state and the model state, respectively the time-history of the "true" model state from an OSSE the time-history of the model state before assimilation  the time-history of the model state after assimilation
dart_log.out [default name] the run-time diagnostic output
dart_log.nml [default name] the record of all the namelists actually USED - contains the default values





restart_file_to_statevector cannot open file "xxxx" for reading The GITM restart file "xxxx" does not exist.
restart_file_to_statevector 'WARNING!!! year 0 not supported; setting to year 1 year 0 ... is not supported in a Gregorian calendar. Our intent here is to do data assimilation, normally 'real' observations have 'real' dates.
statevector_to_restart_file current time /= model time. FATAL error. The DART time does not match the time of the GITM restart file. This message is preceeded by several lines indicating the expected times of both DART and GITM.


The scripts to advance the model and run an experiment are in a state of flux in response to changing hardware.



Provide a better mechanism for generating a set of perturbed initial conditions - pert_model_state()





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

Contact: DART core group
Revision: $Revision: 11626 $
Source: $URL: $
Change Date: $Date: 2017-05-11 11:27:50 -0600 (Thu, 11 May 2017) $
Change history:  try "svn log" or "svn diff"