Contact: | Hui Liu, Tim Hoar, Jeff Anderson |
Revision: | $Revision: 2801 $ |
Source: | $URL: http://subversion.ucar.edu/DAReS/DART/trunk/obs_def/obs_def_radar_mod.html $ |
Change Date: | $Date: 2007-04-04 22:11:48 -0600 (Wed, 04 Apr 2007) $ |
Change history: | try "svn log" or "svn diff" |
types_mod utilities_mod location_mod (most likely threed_sphere) assim_model_mod obs_kind_mod
use obs_def_radar_mod, only : | write_rad_vel |
read_rad_vel | |
set_rad_vel | |
interactive_rad_vel | |
get_expected_rad_vel | |
get_expected_rad_ref |
Optional namelist interface &obs_def_radar_mod_nml may be read from file input.nml.
integer, intent(in) :: velkey integer, intent(in) :: ifile character(len=*), optional, intent(in) :: fform
Writes the obs definition part required for radial velocity.
valkey | Unique identifier associated to this radial velocity observation |
ifile | File unit open to output file |
fform | File format specifier: FORMATTED or UNFORMATTED; default FORMATTED |
integer, intent(in) :: velkey integer, intent(in) :: ifile character(len=*), optional, intent(in) :: fform
Read companion of write_rad_vel.
velkey | Unique identifier associated to this radial velocity observation |
ifile | File unit open to output file |
fform | File format specifier: FORMATTED or UNFORMATTED; default FORMATTED |
integer, intent(in) :: velkey type(location_type), intent(in) :: location real(r8), intent(in) :: orientation
Set the radar location and the orientation of the radar beam at the observation location in the private arrays of the module.
velkey | Unique identifier associated to this radial velocity observation |
location | Radar location |
orientation | Orientation of the radar beam at the observation location |
integer, intent(in) :: velkey
Initializes the specialized part of a radial velocity observation from information provided interactively by the user.
velkey | Unique identifier associated to this radial velocity observation |
real(r8), intent(in) :: state_vector(:) type(location_type), intent(in) :: location integer, intent(in) :: velkey real(r8), intent(out) :: vr integer, intent(out) :: istatus
Calculates the radial velocity (m/s).
state_vector | A one dimensional representation of the model state vector |
location | Location of this obs |
velkey | Unique identifier associated to this radial velocity observation |
vr | The returned radial velocity value |
istatus | Returned integer describing problems with applying forward operator |
real(r8), intent(in) :: state_vector(:) type(location_type), intent(in) :: location real(r8), intent(out) :: ref integer, intent(out) :: istatus
Calculates the radar reflectivity (dBZ), where Z (mm6 m-3).
state_vector | A one dimensional representation of the model state vector |
location | Location of this obs |
vr | The returned radar reflectivity value |
istatus | Returned integer describing problems with applying forward operator |
get_reflectivity | |
write_orientation | |
read_orientation |
real(r8), intent(in) :: qr real(r8), intent(in) :: qg real(r8), intent(in) :: qs real(r8), intent(in) :: rho real(r8), intent(in) :: temp real(r8), intent(out) :: ref
Calculates the radar reflectivity Z (mm6 m-3).
qr | Rain water content (kg kg-1) |
qg | Graupel/hail content (kg kg-1) |
qs | Snow content (kg kg-1) |
rho | Air density (kg m-3) |
temp | Air temperature (K) |
ref | The returned radar reflectivity value |
integer, intent(in) :: ifile real(r8), intent(in) :: orientation(3) character(len=*), optional, intent(in) :: fform
Writes the orientation required for radial velocity.
ifile | File unit open to output file |
orientation | Three components of the velocity orientation |
fform | File format specifier: FORMATTED or UNFORMATTED; default FORMATTED |
real(r8), dimension(3) :: read_orientation integer, intent(in) :: ifile character(len=*), optional, intent(in) :: fform
Reads orientation from file that was written by write_orientation.
orientation | Returns 3 real values for the orientation. This is specific to the observation kind 'radial velocity'. |
ifile | File unit open to output file |
fform | File format specifier: FORMATTED or UNFORMATTED; default FORMATTED |
We adhere to the F90 standard of starting a namelist with an ampersand '&' and terminating with a slash '/'.
namelist / obs_def_radar_mod_nml / & convert_to_dbz, dbz_threshold, & apply_ref_limit_to_obs, reflectivity_limit_obs, lowest_reflectivity_obs, & apply_ref_limit_to_state, reflectivity_limit_state, lowest_reflectivity_state
This namelist is read in a file called input.nml
Contents | Type | Description |
---|---|---|
convert_to_dbz | logical | Convert final reflectivity to dB. Default .TRUE. Note that for now this is forced to .TRUE. in initialize_module because it is not clear how to assimilate Z. |
dbz_threshold | real(r8) |
Applied if convert_to_dbz is .TRUE.
Reflectivity below dbz_threshold is set to dbz_threshold, which is chosen
such to convert to the smallest possible dBZ.
The purpose is mainly computational: We want to avoid the occurence of log(0).
Default 0.001 .
Some useful values: Z = 3.163 -> 5.0 dBZ Z = 2.512 -> 4.0 dBZ Z = 1.259 -> 1.0 dBZ Z = 1.0233 -> 0.1 dBZ Z = 0.1 -> -10.0 dBZ Z = 0.01 -> -20.0 dBZ Z = 0.001 -> -30.0 dBZ |
apply_ref_limit_to_obs | logical | If .TRUE., replace reflectivities below "reflectivity_limit" with "lowest_reflectivity". If FALSE, ignore "reflectivity_limit" and "lowest_reflectivity". Default .FALSE. |
reflectivity_limit_obs | real(r8) | Reflectivity below this value is set to "lowest_reflectivity". Default 0.0 |
lowest_reflectivity_obs | real(r8) | If reflectivity < reflectivity_limit, reflectivity is set to this value. The default value of 'missing' is useful for perfect model experiments. Suggested options: lowest_reflectivity_obs = missing_r8 or lowest_reflectivity_obs = any real value. Default MISSING_R8. |
apply_ref_limit_to_state | logical | Similar to "apply_ref_limit_to_obs" but applied to the state (forward operator). Default .FALSE. |
reflectivity_limit_state | real(r8) | Similar to "reflectivity_limit_obs" for state. Default 0.0 |
lowest_reflectivity_state | real(r8) | Similar to "lowest_reflectivity_obs" for state. Default MISSING_R8. |
Routine | Message | Comment |
---|---|---|
read_rad_vel | Expected location header "platform" in input file | The format of the input file is not consistent. |
read_rad_vel set_rad_vel |
'key (',key,') exceed max_rad_vel_obs (',max_rad_vel_obs,')' | The number of radial velocity observations exceeds the array size set in the module. The module has to be edited and recompiled. |
read_orientation | Expected orientation header "dir3d" in input file | The format of the input file is not consistent. |