MMAB C++ Class Library Descriptions
Robert W. Grumbine
17 April 2000
MMAB Technical Note 186
Description of the Classes
points.h - classes point3, ijpt, fijpt
Data members of point3<T>:
These classes declare useful 'point' classes. A point3 is an ordered
triple of things of type T, with member values i, j, k.
ijpt and fijpt are restricted to be integer and floating point point3's.
A structure, latpt, is defined in this include. This is a structure
rather than a class member, so operations cannot be performed on latpts
(the reason being that the latpts name members differently).
One can set an fijpt = latpt, and operate on the fijpts.
T i, j, k;
Operations in the class are:
Declaration in your program:
Arithmetic with a single element of type T (say a float) (=, +, -, *, /
) Operation is applied to each element of the point3. x = y
Arithmetic with another point3 (=, +, -, +=, -=, *=, element by element)
x += y;
logical operations (element-by element. ==, !=, >=, <=, >, <
). >, < are true only if every element is > or <.
type casting (to int, float, double. The value cast is the L2 norm
of the 3 vector.)
point3<type> x, x(5.0), x(1., 2.), x(1., 2., 3.);
where 'type' is the type of variable you want x to contain. May
be any declared type, including
int, float, double, unsigned char.
vector.h: classes vector, metricvector, time_series
metricvector : public vector
time_series : public vector
Each of these classes is templated.
The vector class is a base class for vectors considered as mathematical,
1 dimensional, objects. This contrasts with the 'vector' which may be defined
in a C++ installation, which is a somewhat different creature.
The metricvector is a pair of vectors, one which carries values, and
one which defines a coordinate associated with the values. This might
be, for example, a vector of the depths in a sounding and the associated
The time_series class is experimental. As the name suggests, it
is directed towards time series operations.
vector : Data members
No public data members.
You may find and change elements of the vector by
value = x[int]; where x is your vector, and value takes the
value of the vector at the [int] point.
x[int] = value; which will put 'value' in to the vector.
vector : Operations
construction, destruction: vector<T> x, x(number of elements)
resizing: x.resize(new number of elements) - this is a destructive resize
Finding number of elements in a vector: nx = x.xpoints();
referencing elements x[int] = y[int]
Are two vectors equal: x == y; (returns bool)
Vector to scalar x = T;
Vector to vector x = y;
void binout(FILE *fout) - unformatted binary output (C-style)
void binin(FILE *fin) - unformatted binary input (C-style)
void printer(FILE *fout) - formatted output (C-style in index, value pairs)
+=, -=, *=, /=
Note that the division has no protection against division by zero
The right hand side may be a vector, in which case the operation is applied
element by element
The rhs may be a scalar (Thing of type T) in which case the same value
is applied to every element of the vector
Vector to scalar functions
x.norm(L-norm value, integer) -- Compute the L norm. If invoked as
x.norm() (no argument provided) L2 is computed
x.rms() -- compute the rms of the vector
x.maximum(), x.minimum(), x.average() -- as you expect
x.average(maskval ) -- Compute the average of the values in the vector
for those points which are equal to the mask value (use the usual care
if you're working on floating point numbers)
x.complete(maskval ) -- Compute the completness (fraction of the vector
which does not have a mask value) of a vector
Vector level functions
x.fill(maskval) -- Fill in a vector by interpolating (linearly) across
masked points. If an end point (first, last) is masked, then the
first (last) interior point is continued to the end of the vector
x.rescale(scaling) -- divide vector by scaling (redundant with /=)
x.normalize() -- Normalize (L2 sense) the invoking vector
vector<T> x, y(number of elements);
metricvector : public vector -> this means that anything you can do
to/with a vector, you can do to/with a metricvector
No public data. The metric may be read or set by set_metric
metricvector: Additional Operations:
Set the metric of the metricvector x to be y.
vector<T> y1 = x.get_metric(vector<T> y2);
both puts the metric of x into the vector y2, and returns it to y1.
May also be invoked as x.get_metric(vector<T> y2); with the return value
Interpolate from the metricvector y onto the metricvector x -- c.f. there
are two different depth sets of interest and you have data on one (the
y metricvector) and want it on the other (the x)
metricvector<T> x, y(500);
time_series : public vector -> this means that anything you can
do to/with a vector, you can do to/with a time_series
time_series : Additional Operations:
y.autocovary(int lag) -- compute autocovariance at lag
y.autocovary(T maskval, int lag) -- compute the autocovariance for 'lag'
but ignore pairs for which one or the other has the mask value. Note
that the autocovary and crossvary functions are overloaded.
y.crossvary(time_series<T> x, int lag) -- compute the cross-covariance
between time series x and y. lag has the default
value of 1.
y.crossvary(time_series<T> x, T maskval, int lag) -- compute the cross-covariance
between time series x and y, skipping points where either series is masked.
lag has the default value of 1.
y.fft(vector<T> real_part, vector<T> imaginary_part) -- compute the
fft of a data set, returning the real and imaginary parts as vectors.
Pads to 2^N data points and uses Numerical Recipes algorithm for
y.ifft(vector<T> real_part, vector<T> imaginary_part) -- compute
the inverse fft given the real and imaginary parts, and place values into
time_series<T> x, y(500);
templated classes for working with 2d arrays.
Each class here inherits
from the class above. Anything a grid2_base can do, so can a grid2,
a metricgrid, or any of the many ncepgrids. Metricgrids add even
more operations. The ncepgrids, in general, don't add any new operations.
ncepgrids specialize the data segment of the class. That is, the
global_ice grid doesn't provide any new operations. It is a latitude-longitude
grid and can do whatever it is that such grids can do. The new thing
is the data -- the grid resolution, orientation, etc. -- are specified.
grid2 (grid2_base with mathematics added)
metric (grid2's which have a mapping between
i,j and latitude-longitude, can also be written out in grib) The metric
grid is a pure virtual class with subclasses:
psgrid (polar stereographic)
mercator (as name, experimental)
cofs (native cofs grid)
chal_native (native Chalikov cubical earth grid)
ncepgrids (metric grids whose parameters describe
specific metricgrids which are in use at ncep. They inherit from
the appropriate metricgrid subclass)
Note here a side benefit. I just said that
the global_ice grid was a latitude-longitude grid. Today, this is
true. Tomorrow, perhaps it won't be. But, regardless,
you won't have to rewrite your program at all (you would have to
recompile, as things stand at the moment. In the future even this
may not be necessary.). What you write is global_ice<float> x;
to let x represent a global_ice grid filled with floating point values.
To read it in (for example), you write x.binin(input file). You don't
have to change this in your program, as your logic is unchanged.
What goes on behind the scenes to carry out your request, you don't care
about as long as it works. And the behind the scenes is the part
that changes to the definition of global_ice affect. If next week,
the resolution changes, again, you don't change your program, just recompile.
(If you write your program making explicit use of the resolution, shame
on you. What you should use is x.dlat, x.dlon, the grid spacings
of the grid. Better yet, use the locate function since grid types
other than latitude-longitude may not have these values.)
No public data.
Number of grid points in x direction is returned by x.xpoints()
Number of grid points in y direction is returned by x.ypoints()
Element at an ijpoint is given by x[ijpt] and can be used to access or
change values at a point
grid2_base<T> x(nx); (ny has default value of 30 for no particular reason)
grid2_base<T> x(nx, ny); (Build a grid2_base which has nx elements in
the x direction and ny in the y direction. x direction is first index,
y direction is second, Fortran ordering used.)
x.resize(int, int) - changes the size of a grid2_base by deleting the current
grid2_base and reconstructing a blank one
x[fijpt] - return the value of x nearest to the floating point location
x[ijpt] - return the value of x at the given location
x[int ] - return the value of x at the given element number (deprecated,
use the corresponding ijpt)
iscyclicx - returns true if the grid is cyclic in x (false by default,
but descendant grids may be cyclic)
iscyclicy - returns true if the grid is cyclic in y (false by default,
but descendant grids may be cyclic)
Setting or exporting
x.set(T value); - set the entire grid to the given value
x.set(grid2_base<T> y) - copy the values in y to the grid x (deprecated,
use x = y)
x.strip(T *y) - y is a pointer to an array of things type T. The
values of x get copied over to y (i.e, the data are stripped out of the
object). (deprecated. a function
to do whatever it is should either exist in the class itself, or should
be written in to a library)
IO on grid2_base :
void binout(FILE *fout);
void binin(FILE *fin); -- deprecated. Use
x.read(FILE *fin) instead.
void ftnin(FILE *fin); -- do a 'fortran' binary read. Assumes ieee
int read(grid2_base<T> y, FILE *fin); -- Read from file 'fin' in to
grid2_base y , returns number of elements read - (deprecated,
use y.read(FILE *fin);
int read(FILE *fin); -- read x in from file 'fin', returns number of elements
int read(char *fname); -- open file 'fname' and read the values in to x.
Returns number of elements read
Logical grid/point operations
int anyof(T type, int range, ijpt ¢er); - returns the number
of points in a 'range' of grid points around the 'center' whose value is
bool in(fijpt &x); -- returns a boolean denoting
whether 'x' is (or rounds to being) 'in' the grid2_base.
bool in(ijpt &i); -- returns a boolean denoting
whether 'i' is 'in' the grid2_base.
bool mask(T maskval, fijpt &loc); -- is the location
'loc' a 'masked' point?
bool mask(T maskval, ijpt ¢er); -- does the
location 'center' have a mask value?
T nonmask(grid2_base<T>&mask, T maskval, T flagval, fijpt &y);
-- returns a value T corresponding to the flagval if all the points are
masked, else it returns the value for the linear interpolator excluding
= grid2_base<T> = T -- set grid to a value
= grid2_base<T> = grid2_base<T> -- set grid to
== grid2_base<T> x == grid2_base<T> y; -- return true
if every value in grid x is the same as the corresponding location in y
put_vector(vector<T> &x, int j); put vector x in to row j of the
put_transect(ijpt &ll, ijpt &ur, vector<T> &x); put vector
x into grid2_base from point ll to point ur (need not correspond to a full
get_transect(ijpt &ll, ijpt &ur, vector<T> &x); get vector
x from grid2_base from point ll to point ur (need not correspond to a full
Logical grid operations
grid2_base<T> & shift(ijpt delta ) -- shift (without wraparound)
the values in the grid2_base by delta.i points in the x, and delta.j points
in the y direction.
grid2_base<T> & subset(int x1, int y1, int x2, int y2) -- deprecated
(use subset(ijpt ll, ijpt ur)) -- return a grid2_base<T> which has only
the values from x1,y1 to x2,y2 of the original grid
grid2_base<T> & subset(ijpt &ll, ijpt &ur) -- return a grid2_base
which has values of the original grid from the lower left (ll) specified
point to the upper right (ur).
grid2_base<T> & magnify(int N); -- return a grid2_base which is
magnified by an integer factor N times of the original grid2_base.
The values of the original are simply repeated in to NxN boxes.
grid2_base<T> x, y(5), z(250, 728);
grid2 : public grid2_base
Note that below, even though we're now in the grid2 class, some of the
arguments are given as being grid2_base. This is because a grid2
_is_ a grid2_base. (Think of grid2_base as being 'mammals' and grid2
as being 'carnivore'. 'Carnivore' can fill any spot that a 'mammal'
Construction: as for grid2_base
Unitary Mathematical operators:
+=, -=, *=, /= with either grid2 or scalars. Operations are
applied element by element. *= is not a matrix multiplying a matrix.
Binary mathematical operators
+, -, *, / with either grid2 or scalars. Use is deprecated
(take two steps, z = x; z += y; rather than z = x + y. The reason
for deprecation is memory and execution efficiency. The 'efficiency'
involved is that in z = x + y, the + operation constructs a grid2 that
gets copied in to z, but is never deleted. If you use this kind of
statement a lot, you can easily overrun system memory limits. This
usage is also slower than using z = x; z += y;)
= grid2_base<T> (be sure that the T of the grid2_base is something that
mathematical operations are defined for. Obviously 'int' is fine.
But 'buoy_report' may not be, even though it may make sense to have a 2d
grid of 'buoy_report's)
= T (assign a scalar to every element of the array)
Grid level functions:
T gridmax(), T gridmin(), T average (), T rms() -- as defined for math.
T average(T maskval) -- average excluding points that are masked
laplace(grid2<T> &x) -- puts the laplacean (cartesian grid assumed
in this class) of the invoking grid2 into x
grid2<T> & laplace() -- returns the laplacean (cartesian grid assumed
in this class) of the invoking grid2
grid2<T> & gradsq() -- returns the square of the gradient (cartesian
grid assumed in this class)
void gradsq(grid2<T> &x) -- puts the square of the gradient (cartesian
grid assumed in this class) of the invoking grid2 into x
void scale() -- rescales the invoking grid into the range 0:127, potentially
useful for graphics
void grib_scale(int &prec, int &nbit, grid2<float> &round)
-- perform rescaling of the invoking grid, multiplying by 10^prec, rounding,
and putting the result into the 'round' grid2. nbit, the number of
bits required to represent the data range, is also computed. NOTE: This
is probably a function which should be protected (from the user).
The user is interested in writing out grib files. This is part of
the sausage-making process involved in doing so. Since 'metricgrids'
have a function to write grib files, users should use that, with the metricgrid
process (gribit) taking care of the grib scaling.
Logical grid operations:
void crit (T cutoff) -- set every element which is not above the cutoff
void crit (T cutoff, grid2<T> &x) -- copy the elements of the invoking
grid to x for every point which is above the cutoff and set to zero otherwise
grid2<T> & reduce(grid2_base<T> &x) -- perform an averaging
reduction of the grid x to the size of the invoking grid.
grid2<T> & reduce(grid2_base<T> &x, grid2_base<T> &mask,
T maskval) -- perform an averaging reduction of the grid x to the size
of the invoking grid, _excluding_ points which are masked (if all points
are masked, zero is returned).
IO on grid2's:
void printer(FILE *fout) -- assumes a conversion to float exists for all
elements of the grid, writes out i, j, value for every element of the grid
void reader(FILE *fin) -- reads in the result of a 'printer' execution
void xpm(char *fname, int resolution, grid2<T> &land, palette<unsigned
char> gg) -- Deprecated, use next xpm and colorproc instead.
void xpm(char *fname, int resolution, palette<unsigned char> gg) --
construct an .xpm file named 'fname' by dividing the values in the invoking
grid by 'resolution' and coloring according to 'palette'. It may
be useful to run colorproc with a function of your design to do the color
void colorproc(grid2<T> &x, int arg1, int arg2, T (*fn)(T, T, int,
int)) -- will apply the function fn to every element of the invoking grid
and place the result in to x, passing arg1 and arg2 to your function.
The last argument is a function which returns a T. Declare as (for
example): unsigned char (*fn)(unsigned char, unsigned char, int, int) ;
fn = &fun1; -- where fun1 is the function you wrote with these arguments.
You'll then invoke the xpm function with a 'resolution' of 1 (you've already
scaled the values to the right range in your function.)
grid2<T> x, x(5), x(70, 900);
The metricgrid class itself is a pure virtual class. That is,
there are certain operations which are common to all grids which have metrics,
but they can't be written out because the operation depends on the actual
metricgrid : Data Elements:
llgrid (latitude-longitude grids)
rotllgrid (rotated latitude-longitude grids, experimental)
psgrid (polar stereographic grids)
mercator (Mercator grids, experimental)
cofsnative -- cofs grid, detailed under ncepgrids
chal_native -- native cubical earth grid
metricgrid : Operations:
grib_pds pds -- all metricgrids have the data for a grib PDS, in a variable
grib_gds gds -- all metricgrids have the data for a grib GDS, in a variable
(Note: This usage will probably be changed in the future, to having metricgrids
inherit from grib_pds and grib_gds classes.)
others which vary by grid type
llgrid: public metricgrid
metricgrid(grid2<T> &) -- promote a grid2 to a metricgrid
Required operators -- You can construct your own grids as descendants of
the metricgrid class, but you have to write these operations in to your
class. They are required for the grid to grid interpolation to work.
latpt locate(ijpt &x) -- return the latitude-longitude location of
the integer location x
latpt locate(fijpt &x) -- return the latitude-longitude location of
the floating point x
fijpt locate(latpt &loc) -- return the location in your grid (as a
floating point location, i.e., i = 2.53, j = 3.45) of the given latitude-longitude
bool llin (latpt &loc) -- Is the given location in the lat-long space
of the grid?
T llget (latpt &loc) -- return the value of the grid point nearest
the given lat-long location
void set_gds() -- Set up the GRIB GDS (default values used, but will generally
be overridded by the descendant)
double integrate() -- integrate the grid over physical space. Currently
this is a dummy function which returns zero except for those grids which
have written their own integration routine (lat-long and polar stereographic
have them). In the future, this will be a requirement for all descendants.
int gribit(int parmno, int depth, int fcst_lead, char *grib, int lgrib,
int mxbit) -- return the grib-encoded binary in 'grib' (you need not have
already allocated space for it), for a given physical parameter (parmno),
at forecast lead time (fcst_lead), subject to a limited number of bits
precision (mxbit), and placing the length of the grib message in to lgrib.
The return 'int' is the error code from W3FI72.
void fromall(metricgrid<T> &input, metricgrid<T> &mask, T
landvalkval, T nonval) -- interpolates from the 'input' grid to the invoking
grid, subject to there being a data mask (grid 'mask', with values 'maskval'
being the mask), and putting in the value 'nonval' if it is impossible
to find a valid value for interpolation on the original grid. [N.B.
This function and its support were the genesis of much of what is present
in the class library.]
float dlat, dlon, firstlat, firstlon -- the grid
spacing in latitude, longitude, and the location in latitude and longitude
of the first grid point (0,0). Southern latitudes are negative.
Longitude is measured positive to the east, i.e., 1 is 1 degree east.
Deltas may be negative (meaning that as the index increases,
you move south, or west).
Note: These parameters should, in the future, be
made private so as to avoid accidental overwriting of their values.
(Though with accessors so that you can find out their values!)
llgrid : Operations: -- required locate functions, and:
llgrid(int nx, int ny, float dlat , float dlon, float firstlat, float firstlon);
-- construct a lat-long grid with these parameters
void set_gds() -- overrides the default metricgrid operation and specialized
the gds for a latitude-longitude grid
double integrate() -- integral over the grid for a spherical earth
double cellarea(ijpt &x) -- area of the grid cell at point x
llgrid<T> x, y(360, 180, 1.0, 1.0, -89.5, 0.5);
(The second declares a 1x1 degree global grid centered on the half-integer
points and starting from the south pole. This happens to be the Levitus
psgrid : public metricgrid
No additional public data
psgrid : Operations: -- required locate functions, and:
psgrid(int nx, int ny,float polei, float polej, float slat, float slon,
float sign, float dx, float dy) -- construct a polar stereo grid with these
parameters. See the declaration section for more.
void subset(psgrid<T> &x, ijpt &ll, ijpt &ur) -- copy original
grid's data out from the lower-left point (ll) to the upper right point
(ur) and put it in to the polar-stereographic grid x. X will be resized
to fit this subset.
void set_gds() -- overrides the default metricgrid operation and specialized
the gds for the current polar stereographic grid
double integrate() -- integral over the grid for a spherical earth
double cellarea(ijpt &x) -- area of the grid cell at point x
psgrid<T> x, y(385, 465, 190.0, 230.0, 60.0, -10.0, +1.0,
The latter specifies the northern hemisphere sea ice analysis grid.
The grid is 385x465, with its pole at 190, 230. The standard latitude
of the projection is 60 N, and the standard longitude is 10 degrees east
of 90 W (I agree this is awkward). +1.0 denotes that this grid is
in the northern hemisphere (-1.0 for the southern). The grid resolution
is 25.4e3 meters in x and y. The standard latitude and longitude
here are taken from the NCEP Handbook (standard longitude is to be 80 W).
mercator : public metricgrid
rotllgrid : public llgrid
float rotation -- rotation angle in degrees.
Operations: only required locates are present
x is rotated by 30 degrees. (Currently assumes dimensions and
values from Chalikov global ocean model.)
derived from llgrids:
chalgrid - Global Chalikov ocean model grid
chalreg - Regional Chalikov ocean model grid
mrf1deg - 1 degree grid in global atmospheric model standard form (start
at north pole, 0 E, increase by 1 degree East, then South
global_otis - Global OTIS analysis grid
reg_otis - Regional OTIS analysis grid
global_wave - Global Wave model grid (2.5x2.5, operational at this writing)
NEW_global_wave - NWW3 Global wave model grid (1.25x1.0 degrees, will be
ECG_wave - East Coast and Gulf of Mexico wave model
ECG_hur - Hurricane version of East Coast and Gulf of Mexico wave model
nh_ocean_weather - Northern hemisphere ocean weather grid (fog, icing)
nh_hazard - Northern hemisphere hazard map grid (20 S to 90 N)
global_ice - Global ice analysis grid
global_sst - Reynolds analysis grid
levitus_atlas - Levitus atlas grid
global_ll - provides a default (1 degree) global lat-long grid
cfsreg - Regular grid representing (at least) the ROFS domain
nopp - Regular grid representing the NOPP/CMDP domain
derived from psgrids:
northgrid - northern hemisphere sea ice analysis grid (25.4 km)
southgrid - southern hemisphere sea ice analysis grid (25.4 km)
northhigh - northern hemisphere high resolution sea ice grid (12.7 km)
southhigh - southern hemisphere high resolution sea ice grid (12.7 km)
northmodel - northern hemisphere sea ice model grid (50.8 km)
southmodel - southern hemisphere sea ice model grid (50.8 km)
nsidcnorth - NASA/National Snow and Ice Data Center grid for northern hemisphere,
nsidcsouth - NASA/National Snow and Ice Data Center grid for southern hemisphere,
north_nesdis_eighth - 1/8th bedient (47.6 km) NESDIS mastermap grid
cfsgrid : public grid2
Cofs native grids:
-- specifies a 181 by 101 grid2. No metric, but it is the right
size for cofs native grids
cfsnative : public metricgrid
-- specifies the native cfs grid, with metric and correct size
cfsnative: Data Elements
cfsgrid<float> alat, alon, depth; -- specifies the latitude, longitude,
and ocean depth on the cofs native grid. Must be read in from outside
cfsgrid<float> t; --dummy. Currently does nothing.
cfsreg<float> fi, fj; -- stores the location (floating point) in the
native grid that corresponds to points on the regular grid. Must
be read in from outside files
vector<float> z, zz, dz, dzz; -- Depth grid coordinate information.
In same file as alat, alon, depth.
Required locate functions
int land(ijpt &x) -- Is the point x a land point (depth <= 10 meters)?
Grid space conversions
void from(llgrid<T> &input, llgrid<T> &mask, T maskval, T
nonval, cfsgrid<T> & );
void from(psgrid<T> &input, psgrid<T> &mask, T maskval, T
nonval, cfsgrid<T> &);
void from(mercator<T> &input, mercator<T> &mask, T maskval,
T nonval, cfsgrid<T> & );
void from(cfslevel<T> &input, cfsgrid<T> &mask, T maskval,
T, llgrid<T> &);
Chalikov Native Grid
chal_native : public metricgrid:
chal_native : Data elements:
chal_native : Operations
grid2<float> alat, alon, depth -- grids of the i,j to latitude, longitude,
and depth mapping
Required locate functions
class grib_pds declares data and operations for the PDS portion of a
grib message. This is used by
metricgrids for the purpose of writing out grib messages.
grib_pts : Data
pds : array of the 25 elements which describe a PDS message
grib_pds : Operations
get_precision() -- returns the precision to which data is recorded
show_date() -- prints the date from the PDS
to the stdout
fix_year() -- checks and year and century for consistency with WMO standards
fix_century(int ) -- given a year, checks and year and century for consistency
with WMO standards
set_time(int, int, int, int, int = 0) -- set the year, month, day, hour,
and optionally, minute. Minute may be omitted, in which case it is
set to zero.
set_precision(float ) -- set the data precision, pass in the decimal
precision desired, e.g. 0.01
set_table (default = 3) -- set the grib table to be used. 3 is the
NCEP standard table, and is the default (value given for no argument).
128 is the NCEP MMAB grib table.
set_center (default = 7) -- declare generating center. Default is
set_process(int ) -- set the process generating number
set_gridid(int ) -- set the grid ID (255 is 'undefined')
set_parmno(int ) -- set the parameter number of your physical parameter
set_layer(int ) -- set the vertical layer type (experimental)
set_depth(int ) -- set the depth of your layer (experimental)
set_fcst_lead(int ) -- set the forecast lead
grib_gds : Data
gdslen, gds; // Length of the gds and the gds
All of these are declared by default, and must be overwritten in the
construction of the data grids. Experimental. This is an undesirable
way to manage the type.
grid3 : public grid2
grid3 is an impoverished class. At the moment is it almost strictly
a placeholder of sorts. It's purpose is to declare a class suitable
for working on 3 dimensional grids such as the global ocean temperature
field. Ultimate implementation should be quite different than the
present as this version makes no real inheritance from the grid2 class.
grid3 : Data elements
int nz; // Number of levels in the vertical, data stored as nz levels of
grid3 : Operations
[ijpt] -- return the value at grid point i, j, k
[int] -- return the i'th element
set(T ) -- set every element of the grid to T
binout(FILE *) -- binary output of the 3d field.
binin(FILE *) -- binary input of a 3d field
printer(FILE *) -- ASCII output of a 3d field: i, j, k, value;
vector<T> get_sounding(ijpt ) -- get a sounding from coordinates i,
j and return as a vector<T>
void put_sounding(ijpt, vector<T> &) -- take a sounding from the
input vector and put it in to the grid3 at location i,j.
grid2<T> get_layer(int) -- get layer N (numbers from 0 to nz-1) from
the 3d field and return as a grid2
put_layer(int, grid2<T> ) -- take a grid2 and enter it in to the grid3
at layer N
get_transect(ijpt, ijpt, grid2<T> ) -- construct a transect on native
levels between the two points given as arguments.
grid3<T> & operator=(grid3<T> &) -- equate two grid3's
void operator+=(grid3<T> &) -- add the argument grid3 to the invoking
void operator-=(grid3<T> &) -- subtract the argument grid3 from
the invoking grid3
void operator/=(T) -- divide every element of the grid by T
void operator*=(T) -- multiply every element of the grid by T
buoy is a class developed for working with buoy reports, particularly
from the CPC buoy data files. Note that as of 1 March 1997, the CPC
coding changed, as did the permissable classes of reporting platform.
This date is a private part of the class, and functions like isbuoy will
permit users to determine whether a report is a buoy report (as opposed
to ship, e.g.) without remembering this date.
buoy_report : Public Data
buoy_report : Operations
float slp, wind_dir, wind_speed, air_temp, dewpt_depression;
float latitude, longitude
time_t obs_secs; // Time in seconds since 1 January 1970 0:00:00
of the observation. This is limited by the size of the type 'time_t'
to approximately 60 years.
int year, month, day
Adjusting class defaults:
set_space_range(float ) -- Reset the default spatial separation inside
which things are considered 'close'. Default is 50 km. Units
set_time_range(time_t ) -- Temporal separation, in seconds, which is considered
'close'. Default is 3600 seconds (1 hour).
Working with time
set_time -- ensure that the time information is completely initialized
tm get_time() -- get the observation time 'tm' structure. See /usr/include/time.h.
set_secs(time_t &x) -- set the observation time in seconds since 1
January 1970 0:00:00
time_t get_secs() -- get the observation time in seconds since 1 January
buoy_report& read(FILE *) -- read a CPC buoyreport from the file
buoy_report& iabpread(FILE *) -- read an IABP buoy report to a buoy_report
write(FILE *) -- write a buoy_report out in CPC format
Tests on report type
bool same_id(buoy_report &x) -- true if the invoking and the
x buoy have the same identifier
bool isbuoy() -- true if the report is from a fixed or drifting buoy (as
opposed to a ship or CMAN)
bool isdrifter() -- true if the report is from a drifting buoy.
bool isfixed() -- true if report is from a fixed buoy
bool isship() -- true if the report is from a ship
bool is_moving_ship() -- true if the report is from a moving ship (Archaic:
As of 1 March 1997 the CPC files ceased distinguishing between moving and
fixed ships, as well as changing the encoding for permitted types)
bool is_other_fixed() -- true if the report is from fixed platform other
than a fixed ship or buoy.
Tests on whether a buoy_report is 'near' to something
near(tm &x, time_t toler) -- true if the buoy report is within toler
seconds of time x
near(buoy_report &, const char *, time_t toler) -- true if near
in time and has the same ID
poleward(float ) -- true if latitude of report is closer to the
pole than the argument
synoptic(float) -- true if report is within float hours of 00 UTC.
A good additional function would be one which let the user specify the
near(buoy_report &x, time_t toler) -- wrapper for near(tm &, time_t),
true if the invoking buoy report is within toler seconds of the argument
near(latpt &x, float toler) -- true if the report is within toler km
of the given (x) latitude-longitude point.
near(buoy_report &x, float toler) -- note overloading, float here,
time_t above. This one tests whether the buoy is within toler km
of the argument x.
near(buoy_report &x, time_t toler, float distance) -- true if the buoy
x is within toler seconds of the invoking buoy and is within distance km.
class avbuoy : public buoy_report
This class makes it easier to average a number of buoy_reports.
As currently written, it only truly averages latitude/longitude information.
avbuoy : Public data
None (averaging is done privately)
avbuoy : Operations
avbuoy& operator=(buoy_report &x) -- give the avbuoy the same
data as the buoy_report x
avbuoy& operator=(avbuoy &x) -- equate two avbuoy's
avbuoy& operator+=(buoy_report &x) -- add a buoy_report to
avbuoy& average() -- average the avbuoy
avbuoy& read(FILE *) -- read in an avbuoy (obsolete)
void & write(FILE *) -- write out an avbuoy (in buoy_report
It is intended that this class be sufficient for handling color graphic
processes. The basic objects are RGB colors as ordered triples of
point3<T> type. T is unsigned char for 8 bit graphics.
palette : Public data
ncol, cbase -- ncol is number of colors, and cbase is the base character
for .xpm file output. Cbase should be 65.
point3<T> *pal -- a pointer to point3<T>'s (the colors).
Better would be a vector of point3's.
palette : Operations
palette() -- default color palette, no colors, cbase = 65
palette(int N) -- N colors, cbase = 65, palette initialized
palette(int , int) -- set the number of colors and the character base
resize(int ) -- destructive resizing of the color palette
color_bar(int I, int J) -- print out blocks of pixels with each of
the colors in the table, size of block is I by J
print() -- print out index and R, G, B values for each color
set_color(point3<T> x) -- set entire palette to x.
lighten(float ) -- lighten (increase magnitude of) all the colors
in the table by the specified percentage
darken(float ) -- darken as above
invert() -- invert the color values (new = 255 - old)
set_color(int I, point3<T> x) -- set color I in the palette
to be x
set_color(int I, int R, int G, int B) -- set color I in the palette
to have RGB values R, G, B
The MMAB C++ library should be usable as any ordinary library, with
only a couple of points of particular care needed. First, and most
significant, is that the MMAB library declared a vector class of its
own. It is necessary that it and the C++ standard library vector
class not be used in the same program. It is also necessary to be
careful of defining "vector.h" rather than <vector.h>.
Second, in order to use the cfsnative and chal_native classes, it will
be necessary to get the ASCII data files for them and write them back out
in your machine's native binary format. More detail is given in the
README file in interp.doc.tar.
Third, it may be necessary to edit the files for platform dependances.
This is mostly a reflection of the fact that some of the omblib is in Fortran,
and calling Fortran from C++ is platform dependant. Defined
platforms are SGI (Origin workstations tested), LINUX (Intel hardware and
g++), and IBM (IBM SP). If this is needed, the error message
you'll get on attempts to link are things like arcdis_ not found (there
is a Fortran function called arcdis, and it gets silently renamed by the
compiler according to local conventions).
Straightforward, but necessary, is that the sources in omblib.tar need
to be compiled and placed in a library to which your programs can and will
Standard type in C++ (not present on all systems,
in particular our IBM SP). Has values of true and false.
Functions in C++ may be specified with default values.
Where they exist, the corresponding argument may be omitted. The
vector norm function has the default value of 2, so x.norm() and
are the same. The defaults must be at the end of argument lists if
there are non-default arguments present.
Various operations in the classes are noted as being
deprecated. This is a reflection of the fact that the class library
was developed over a signficant period of time, and some things that seemed
like good ideas at the time, no longer do. If you're writing new
code, don't use these. Alternatives are noted in the documentation.
If you have old code that uses these, on your next build you should replace
these calls. Deprecated features will be removed from the library
in the future. They are retained solely for backwards compatibility.
The grid classes make extensive use of this feature
of object oriented programming. The basic idea is that an ncepgrid
(for example) is some special type of grid -- it can do all the normal
kinds of 'grid' things _as well as_ doing a few additional things which
_only_ an ncepgrid can do.
We often encounter situations in which we want to
do the same basic thing, but there may well be many different ways of doing
it. In Fortran, this sets up an undesirable situation. Either
we write a separate function -- with a different name! -- for each and
every method, or we write a single mammoth routine which handles each and
every case. Neither is desirable. The first requires the user
to track down which of a dozen different routines he wants, and use a different
name every time (the different names thereby obscuring what exactly is
going on as it won't be so obvious from just the name). The second
requires the user, even in a case where he really does want the simplest
version with the fewest options, to pass the full set of arguments (many
of which may make no sense to someone who wants a default behavior).
C++ manages this situation with overloaded functions.
The function name is left unchanged. The compiler selects
which version of the function you want by reading the argument list.
This permits the naming to remain clear, as in x.crossvary(y), which
computes the cross-covariance between time_series x and y (default set
up, no options), and x.crossvary(y, maskval, lag), which computes
the cross-covariance between x and y again, but now skips masked points
and works at a user-selected time lag. The program invocation remains
clear as to what is happening -- cross-covariances are being computed --
and the compiler, rather than the user, selects the right routine to carry
out the desired action.
Through the documentation, I've noted that something
is 'public'. Part of the design in object orientation is to avoid
some types of errors by making things 'private'. The private things
can be done inside the class, but cannot be done directly by the user.
This prevents the user from making certain kinds of mistakes. For
example, although users may read and modify the contents of arrays and
can _find_ the dimensions of the array, they may not change on their
own the dimensions of the array. The information on array size
is packaged with the array contents into the object. The user may
well want to operate on all the elements (or the first half in x and last
half in y, or some other portion) of the array, but the operations themselves
don't depend on being element 23, rather on being half way through the
array. nx needs to be knowable, but its value doesn't need to be
changed by the user.
If the array size itself were public, one has the
problem (even if nobody else has made this mistake, I have) of passing
the array size to a routine and accidentally changing it inside the routine
(i.e, nx is the array size. You mean to write limit = nx, but actually
write nx = limit). The fact that nx is not public means that you couldn't
make this mistake even if you tried. If you want to change
the size of the array, you _can_ do so, by asking the object to resize
Where classes are listed as being 'public someotherclass',
this means that all the operations that you could do with an object of
the 'someotherclass', you can do explicitly with the objects of the descendant
class. For example, grid2 inherits publicly from grid2_base, so all
the grid2_base operations can be done by a user to a grid2 object.
All the classes currently in the MMAB class library use public inheritance.
There are reasons why one may want to use other types of inheritance, but
they haven't yet occurred for classes we're working with.
Pure virtual class:
This is used to declare a class that you'll never actually have an
object from, but which specifies things that all its descendants have to
be able to do. The metricgrid class is a pure virtual class.
All the descendants, latitude-longitude grids, polar stereographic grids,
the ROFS native grid, etc., have to write their own functions to provide
mapping between grid coordinates and latitude-longitude coordinates.
We simply don't know ahead of time how this will be done for new grid classes,
but any metric grid _must_ do this. The benefit is that _if_ there
is such a mapping, we can do other things and write (and debug!) the algorithm
for it in the metricgrid class. One example is grid to grid interpolation.
The _logic_ is the same, regardless of the two grids involved. That
is, do a bilinear interpolation from the data source grid onto the destination
grid given the location of the destination point in the i,j space of the
data grid. Bilinear interpolation is a known and solved problem,
so we write the routine once, debug it. Then regardless of
what new grid we try to work on, once we write the functions for mapping
between physical (lat-long) and grid (i,j) space, we can interpolate from
grid to grid.
Abstraction is a major reason for working in C++. When you have
a function which performs the same logical operations regardless of the
data types it is working on, you want to write (and debug!) only one such
function. In Fortran, for example, if you want to compute the average
value of an array (say) you need different routines for an array of integers,
an array of floating point numbers, and an array of complex numbers.
Logically, this makes no sense, as averaging is defined without regard
to data type. In C++ you would use a template function:
template <class T>
T average(grid<T> x)
and inside the function, simply write the logic for how to average
a grid x which has values of type
T. And you don't, in writing the function, need to consider what
T is. This example is taken from the grid_math class.