Even Sphere
It is frequently useful to generate a series of synthetic observations located at roughly evenly-distributed locations [1] on a sphere.
There are three methods described here.
A Matlab script and the standard DART observation generation utilities.
A csh script with all the parts of 1. (Not available for all models).
A stand-alone Fortran program.
The Fortran program does not generate the nice plots that the Matlab process does, but it may be faster and easier to automate for generating a large number of obs.
Matlab Scripts Plus Standard DART Observation Executables
This involves multiple steps:
determine how many locations are needed
determine the vertical levels needed
run the MATLAB function even_sphere.m to create the text file containing the input for
create_obs_sequence
run program create_obs_sequence to create an observation sequence (usually
set_def.out
, although it is possible to createobs_seq.out
files directly if you don’t really care about the observation values).if desired, run program create_fixed_network_seq to create a longer observation sequence file.
run program perfect_model_obs to harvest the synthetic observations from a chosen model.
This directory contains a MATLAB function (even_sphere.m
)
that generates input for the
program create_obs_sequence .
After executing create_obs_sequence
, the resulting observation sequence file
will have a template for ‘RADIOSONDE_TEMPERATURE’,’RADIOSONDE_U_WIND_COMPONENT’,and
‘RADIOSONDE_V_WIND_COMPONENT’ observations at specified pressure levels and roughly
evenly-distributed locations across the entire globe. Optionally, bogus observation
values may also inserted; which may be useful in certain circumstances.
even_sphere.m
even_sphere.m
has many optional arguments to tailor its behavior.
It has exactly 1 required argument - the number of horizontal locations desired.
It will create a text file
even_create_input
to be used as input tocreate_obs_sequence
The choice of pressure levels is described here.
The default number of pressure levels is 21. Argument
nlevels
specifies how many levels to use from the beginning of the levels list.The default observation error variances for each observation type are level-dependent and are consistent with
DART/observations/obs_converters/obs_error/ncep_obs_err_mod.f90
The default is to create ‘empty’ observation sequences - i.e. they have no actual observation values and are suitable to be used with
perfect_model_obs
The default date of the observations is 2017-12-25 00:00:00
A plot of the locations will be created. The number of gridlines is configurable but defaults to 288 in longitude and 192 in latitude.
All the defaults can be changed by specifying ‘variable-value’ pairs of options, as described below. Examples of some options are also available via the normal MATLAB help facility. (Documenting all of them in the help makes the help page too long.)
Note that the number of observations will be the number of locations *
the number of vertical levels * the number of variables (i.e. 3)
even_sphere.m
also takes observation error variances
and includes them in the observation sequences.
Optional Argument Variable-Value pairs
The optional variable-value pairs can appear in any order.
optional variable |
example value |
Description |
---|---|---|
‘nlevels’ |
5 |
number of pressure levels to use. May be less than the length of the ‘levels’ array, but cannot be more. |
levels |
[1000 500 300 200 100] |
pressure levels desired. see Levels section for discussion. |
T_error_var |
[1.44 0.64 0.81 1.44 0.64] |
level-specific Temperature error variances. see Levels section for discussion. |
W_error_var |
[1.96 4.41 9.00 7.29 4.41] |
level-specific error variances for both U, V wind components. see Levels section for discussion. |
‘YMD’ |
‘2017-12-25’ |
Date required for |
fill_obs |
false |
‘true’ inserts a bogus observation value of 1.0
and a bogus QC value of 0.’false’ does not insert
bogus values and essentially creates an empty
obs sequence file (typically |
‘nlon’ |
288 |
number of longitude grid lines in plot |
‘nlat’ |
192 |
number of latitude grid lines in plot |
Examples
30 horizontal locations at 6 pressure levels:
nprofiles = 30;
levels = [1000 850 500 300 200 100];
T_error_var = [1.44 0.64 0.64 0.81 1.44 0.64];
W_error_var = [1.96 2.25 4.41 9.00 7.29 4.41];
even_sphere(nprofiles, 'levels', levels, ...
'T_error_var', T_error_var, 'W_error_var', W_error_var)
30 horizontal locations at 3 pressure levels. Note that the nlevels argument specifies that only the first 3 pressure levels are used even though there are 6 potential pressure levels. Similarly, only the matching error variances are used.
nprofiles = 30;
nlevels = 3 ;
levels = [1000 850 500 300 200 100];
T_error_var = [1.44 0.64 0.64 0.81 1.44 0.64];
W_error_var = [1.96 2.25 4.41 9.00 7.29 4.41];
even_sphere(nprofiles, 'nlevels', nlevels, 'levels', levels, ...
'T_error_var', T_error_var, 'W_error_var', W_error_var)
Levels
Attention
If you need realistic error variances attached to your observations, be careful to align your levels and variances.
The default levels that this program generates are the mandatory pressure levels defined in the AMS glossary. The corresponding error variances are from ncep_obs_err_mod. See Obs Error. Levels at the top can be excluded by setting nprofiles < 21 (size(levels)).
levels = [1000 925 850 700 500 400 300 250 200 150 100 70 50 30 20 10 7 5 3 2 1];
T_error_var = [1.44 1.00 0.64 0.64 0.64 0.64 0.81 1.44 1.44 1.00 0.64 0.64 0.81 1.00 1.69 2.25 2.25 2.25 2.25 2.25 2.25];
W_error_var = [1.96 2.25 2.25 2.56 4.41 6.76 9.00 10.24 7.29 5.76 4.41 4.41 4.41 4.41 4.41 4.41 4.41 4.41 4.41 4.41 4.41];
Here’s an example of replacing the AMS levels with a set that has more levels near the surface and none above 150 hPa. Note that the error variances should change to be consistent with the levels.
levels = [1000 950 900 850 800 750 700 650 600 550 500 400 300 200 150];
T_error_var = [1.44 1.21 0.81 0.64 0.64 0.64 0.64 0.64 0.64 0.64 0.64 0.64 0.81 1.44 1.00];
W_error_var = [1.96 2.25 2.25 2.25 2.56 2.56 2.56 3.24 3.61 4.00 4.41 6.76 9.00 7.29 5.76];
Running Matlab in Batch Mode
If you would prefer to run even_sphere.m
in batch mode (i.e. from within a shell script),
here is an example syntax that worked for me. The script ran in the same directory
as even_sphere.m
. There are many ways to construct the input, naturally - but you don’t have
to explicitly edit even_sphere.m
this way.
#!/bin/csh
\rm -rf matlab_input.m
cat >> matlab_input.m << EndOfInput
nprofiles = 30;
levels = [1000 850 500 300 200 100];
T_error_var = [1.44 0.64 0.64 0.81 1.44 0.64];
W_error_var = [1.96 2.25 4.41 9.00 7.29 4.41];
even_sphere(nprofiles, 'levels', levels, ...
'T_error_var', T_error_var, 'W_error_var', W_error_var)
fname = sprintf('even_sphere_%d_profiles',nprofiles);
orient landscape
print(fname,'-dpdf')
EndOfInput
matlab -nosplash -nodesktop -r "try; cd $PWD; matlab_input; catch; end; exit";
Automation Scripts
Here there are also scripts (run_fixed_network_*.csh
) which use the
output from create_obs_sequence
and the
program create_fixed_network_seq
to generate a series of observation sequence files.
run_fixed_network_seq.csh
Calls create_fixed_network_seq
to create a separate file for each time period.
By default, it makes 2 files/day, 12 hours apart, single time per file.
The intervals and dates can be changed by editing the script.
It assumes that create_fixed_network
has any model-specific files it needs in this directory.
It requires a set_def.out
file (usually created by create_obs_sequence
).
run_fixed_network_daily.csh
Calls create_fixed_network_seq
to create a separate file for each time period.
By default it makes 1 file/day, single time (noon) per file.
The dates and time can be changed by editing the script.
It assumes that create_fixed_network
has any model-specific files it needs in this directory.
It requires a set_def.out
file (usually created by create_obs_sequence
).
The process, end to end:
MATLAB:
Set the number of profiles, the levels, etc. and run even_sphere.m
in
MATLAB. It creates the necessary text file even_create_input
for the next step.
It will also make a plot - which you can save.
DART:
Then you have a choice about building and running the create_obs_sequence
and create_fixed_network_seq
programs:
building them in the
models/template/work
directoryusing the ones which were built in
models/your_model/work
directory byquickbuild.sh
.
Choice A uses programs which have no model specific file dependencies, but may involve more separate steps than B.
A
Build the programs in
template/work
Link (or copy) these files to the directory in which you want to create obs_seq files.
./even_create_sequence
./run_fixed_network_{seq or daily}.csh
models/template/work/create_fixed_network_seq
models/template/work/create_obs_sequence
models/template/work/input.nml
In your obs_seq directory, run create_obs_sequence, which creates a
set_def.out
file.
./create_obs_sequence < even_create_input > /dev/null
Edit and run your choice of
run_fixed_network_*.csh
for the desired dates. These call create_fixed_network_seq, which creates anobs_seq.in
file for each specified date.
B
This choice may involve fewer steps, if there is a model specific script
which combines the steps in A).
See the cam-fv example (models/cam-fv/shell_scripts/synth_obs_locs_to_seqs.csh).
If there is not a script like that for your model,
you can follow the steps in A),
substituting your model name for the “template” in the pathnames.
NOTE: you may need to link any additional input files which your model requires
into the directory where you will run the programs.
These typically contain grid information and are found in your_model/work
.
For example, cam-fv needs a caminput.nc
and cam_phis.nc
.
Fortran program for generating obs directly
cd into the work directory and run quickbuild.sh
.
This builds the create_even_sphere
executable. Edit the input.nml
to set the number of obs to generate and the date in the namelist. Run
the program and the output file will be generated.
DETAILS of generating points evenly-distributed on a sphere
This is the algorithm that’s being used [1]:
N := the number of profiles you want
dlong := pi*(3-sqrt(5)) /* ~2.39996323 */
dy := 2.0/N
phi := 0
y := 1 - dy/2
for k := 0 .. N-1
r := sqrt(1-y*y)
node[k] := (cos(phi)*r, sin(phi)*r, y)
y := y - dy
phi := phi + dlong
For the geometric and visually minded:
Picture a unit sphere in cartesion space (x,y,z).
Choose a value -1 < y < 1, which defines an x-z plane. That plane intersects with the unit sphere to form a circle whose center is on the y axis. (The circle radius is small near y = +/-1 and is 1 at y=0.)
Choose an angle (“phi”) and draw a ray from the center of the circle to a point on the circle using this angle relative to the x positive direction. Where the ray intersects the circle (and sphere) is one of the evenly distributed points on the sphere which we want.
Its x and z coordinates can then be combined with the already defined y coordinate to define the cartesian location of the point.
The choice of the y and angle for each point is where the magic enters the algorithm. They are derived from the Fibonacci or Golden Spiral formula (derived elsewhere).