Re: More Songs about Coordinate Systems and Buildings
John Sheldon (
jps@gfdl.gov)
Fri, 1 Aug 1997 19:39:29 -0400 (EDT)
Hi again John-
> I presume that the recent lull in the firestorm over coordinate systems
> is due to our patient rumination over the issues, and not the exhaustion
> of the fuel or the firefighters...:^}
No, you were right - it's exhaustion :-). But also the fact that my
other eight projects currently "in progress" weren't making much
progress :-)
> Its been very useful to follow
> other discussions, and I've been influenced by most of them, whether
> I've agreed or disagreed, commented on them or not...
That it has, and you're to be commended for all your dedication...and
for riding the bucking bronco longer and better than most of us =|:-}
> Anyway, I have come to some new impasses in my own thinking, so I've
> decided to expose my current "works in progress", despite the fact that
> it's not complete or satisfactory. I have new web pages:
...
> However it does start to touch on the issues of data display, which we
> should admit that we want to solve.
>
> >From a recent post from Stephen.Walker:
>
> > Here we need to be careful to distinguish between the coordinates which
> > locate the data in the physical world, and the coordinates which locate
> > a given line element or pixel on a plot. They are not always the same.
> > As far as I understand it, most of the discussion about coordinates
> > has been about coordinates which locate the data in the real physical world,
> > not plot coordinates (which may or may not be related to real-world coordinates
> > in a simple way).
>
> Obviously we aren't going to provide coordinates into "display space"
> (pixels on a display or whatever). However we do want to provide
> coordinates in "real world" space clearly enough that applications can
> easily map to display space. So my solution of mapping to Rn is an
> attempt to do this.
Yes. The projection onto the display device is a totally separate
issue...an *old* one, already addressed by all plotting/visualization
packages.
> As John Sheldon insists, we still haven't made it yet to the "real
> world", and so I've proposed a "geodetic" convention (still very
> incomplete), much indebted to Steve Emmerson's posts, and tried to say
> what the semantics of that coordinate system are. Since we are mostly an
> earth science community, I havent proposed other coordinate systems like
> Steve does, that would be useful in solid modeling applications. But if
> some think that's important...
If we do it right, it seems like there ought to be a natural niche for
other coordinate systems...
> ...
>
> Proposed Conventions for Coordinate Systems in Netcdf
>
> draft 7/29/97 by John Caron
...
> In order to facilitate implementation, these proposals will mark some of the
> more advanced features as "level 2". Two conventions are thus
> proposed: "level 1", and "level 2" that includes all of level 1.
>
> 1. General Coordinate Systems
>
> 1.1 Definitions. A coordinate system is best thought of as a function which
> maps a point from an index domain (really just the set of array indices for
> the variable) to a location in some user-defined coordinate space. We assume
> that the maker of the netcdf file is recommending that the user of the file
> view the data in the specified coordinate system.
Absolutely...
> ...
> 1.2 Specification. A coordinate system is specified in a netcdf file by an
> attribute whose name begins with the keyword "coordinates" and whose value
> is a blank or comma delimited list of coordinate functions. For example:
>
> :coordinates = "xpos ypos time wavelength";
I've already discussed (and won't belabor here) my preference for Gary
Granger's approach of defining global "coordinate maps", which
variables may reference. In my rendition, however, the prefix of the
coordinate map will indicate the general coordinate system (cartesian,
cylindrical, spherical, etc.) Further down, I will try to show how
this might work with your motivating examples. But, in general, I
think we are on a surprisingly similar track! :-))
> ...
> 1.4 Coordinate Values.
>
> (Level 1) A coordinate function that is a scalar netcdf variable is
> considered a point along its coordinate axis.
>
> (Level 2) A coordinate function that is multiply valued is specified by a
> list of scalar netcdf variables enclosed in parenthesis:
>
> :coordinates = "lat lon (lev_upper lev_lower lev_midpoint lev_label)";
Hmmm...like some of the other comments made about this issue, I don't
think that "boundary" information really constitutes "coordinate"
status. I prefer GDT's approach, although I'd modify it so that
different variables could use different boudaries along the same basic
axis. (eg, a 3-hr vs. a 6-hr window for data collection.)
> ...
> 1.6 Coordinate Variables. A variable with the same name as a dimension, with
> that dimension in its domain, is the coordinate variable for that dimension.
> It is recommended that coordinate variables be one-dimensional. In this case
> the value of the variable must be strictly increasing or decreasing with
> respect to the dimension index (the function is then said to be monotonic).
>
> (Level 2) Multi-dimensional coordinate variables, while allowed, may be
> misleading in their association with a single dimension, and are considered
> less desirable. Multi-dimensional coordinate variables are not in general
> monotonic, but they must satisfy the uniqueness property of any coordinate
> systems of which they are part.
Well, here's our old, familiar disagreement. Your admonitions aside, I
would really like to avoid multi-dimensional coordinate variables. I
think there are other ways to convey the necessary information that
will not run counter to the long-standing convention of 1-D coordinate
variables. For the sake of prospective users and application
developers, let's retain the simplicity of the current assumptions.
> Coordinate variables define an "implicit" coordinate system for any variable
> that uses those dimensions. For such a variable, find all dimensions with
> coordinate variables: the implicit coordinate system is the one composed of
> that list of coordinate variables. This can be turned off by adding a global
> or variable attribute with name "coordinates.implicit" and value equal to a
> blank string or the string "none". It is recommended that when multiple
> coordinate systems are intended, that all be explicitly defined.
I'm not sure what purpose your "coordinates.implicit" serves. Could
you give a CDL example? Wouldn't a global "coordinates" attribute
imply a default? Or, by "implicit" do you mean that the coordinate
system should be inferred somehow from the qualities of the coordinate
variables without actually having to list them (a la COARDS I)?
> ...
> ----------------------------------------------------------------------------
>
> Examples
>
> Shorthand "solutions" of the motivating examples :
I'll combine parts of your mail and offer an alternative set of
solutions here. The two approaches are so close that I think little
additional explanation is required. Basically, my approach differs
in the following ways:
a) define global "coordinate maps", a la Gary Granger
b) the "coordinate map" name should have a prefix indicating the
basic coordinate system (cartesian, cylindrical, spherical,
geodetic, ...)
c) the specific variables referred to in the "coordinate map"
should be assigned to the generic axes of the basic coordinate
system
cartesian : "X", "Y", and "Z"
cylindrical : "rho", "theta", and "z"
spherical : "theta", "phi", and "r"
geodetic : "latitude", "longitude", and "altitude"
In general, my solutions will look longer, but keep in mind that a
defined coordinate map can be referenced (more concisely) by many
variables, and they are somewhat more descriptive.
> ----------------------------------------------------------------------------
>
> 1. classic coordinate variables:
> var(lat, lon)
> lat(lat)
> lon(lon)
>
> 2. scattered points; assign them a location:
> var(npoints)
> lat(npoints)
> lon(npoints)
Your solution:
> :coordinates = "lat lon";
My solution:
:coordinates = "geodetic_1";
:geodetic_1="latitude=lat, longitude=lon";
Note: partial system definition is OK; it may be that no
variable needs the missing axis.
>
> 3. projective geometry:
> var(lat, lon)
> lat(lat,lon)
> lon(lat,lon)
Your solution:
> :coordinates.latlon = "lat lon";
My solution:
2-D coordinate variables not allowed.
> however, I think we could be clearer about this. var(lat, lon) isn't
> really correct; its really a function of (x,y) on a projection surface.
> So better is:
>
> var(x, y)
> x(x)
> y(y)
> lat(x,y)
> lon(x,y)
Your solution:
> :coordinates.xy = "x y";
My solution:
:coordinates="cartesian_1 geodetic_1";
:cartesian_1="X=x, Y=y";
:geodetic_1="latitude=lat, longitude=lon";
> 4. hybrid coordinates:
> var(lon, lat, lev)
> lat(lat)
> lon(lon)
> lev(lon, lat, lev)
3-D coordinate variables not allowed.
>
> again, we're being vague; really we are describing two alternative lev
> coordinates: hybrid(lev) and pressure(lon, lat, lev) so we can rewrite
> as:
>
> var(lon, lat, hybrid)
> lon(lon)
> lat(lat)
> hybrid(hybrid)
> pressure(lon, lat, hybrid)
Your solution:
> :coordinates.geodetic.hybrid = "lon lat hybrid";
> :coordinates.geodetic.pressure = "lon lat pressure";
My solution:
:coordinates="geodetic_1 geodetic_2";
:geodetic_1"=latitude=lat, longitude=lon, altitude=hybrid";
:geodetic_2"=latitude=lat, longitude=lon, altitude=pressure";
>
> 5. edge coordinates, level is some kind of altitude coordinate:
>
> var(level)
> lev_bottom(level)
> lev_top(level)
Your solution:
> :coordinates = "(lev_bottom lev_top)";
> float lev_bottom( level);
> lev_bottom:component = "lower_bound";
> float lev_top( level);
> lev_top:component = "upper_bound";
My solution:
These are not really "coordinates", but "bounds" of a "cell"
around the actual data point.
> 6. non- georeferencing coordinates. Up to now, all the examples have been
> georeferencing coordinates. A different coordinate system that we use
> in radiative transfer codes is
>
> var(lev, wavelength)
Your solution:
> :coordinates = "lev wavelength";
My solution:
This is not really a spatio-temporal coodinates problem, but we'll
stretch it a bit here:
:coordinates = "cartesian_1";
:cartesian_1="X=wavelength, Z=lev";
> Steve Emmerson's wire example for a spiral wire in a cylindrical
> coordinate system (CCS), which is a spatial, but not georeferencing
> coordinate system:
>
> temp(s) // temperature along spiral
> rho(s) // distance from CCS center axis
> theta(s) // CCS azimuth
> z(s) // CCS height
Your solution:
> :coordinates = "rho theta z";
My solution:
:coordinates = "cylindrical_1";
:cylindrical_1 = "rho=rho, theta=theta, z=z";
>
> 7. vector valued variables, for example:
> vector(lev, 3)
> velocity(lat,lon,component)
> component(3) = "u", "v", "w"
Your solution:
> :coordinates = "lev";
> :coordinates = "lat lon components";
My solution:
This is not a "coordinates" issue. "Associating" separate scalar
quantities into a vector quantity is a useful goal, but I haven't
had time to give it enough thought.
>
> 8. correlations, using a dimension more than once:
> precip(time, npoints)
> precip.correlation( npoints, npoints)
> lat(npoints)
> lon(npoints)
Your solution:
> - will need some notation not yet formally proposed, eg:
> :coordinates = "lat(npoints,) lon(npoints,) lat(,npoints)
My solution:
Not applicable..."correlation" does not possess "coordinates" the
way we think of them. I suppose you could try to document how
the correlation was calculated (ie, the lat-lon location
of the points being correlated), but I think that is probably
something to be addressed using variable attributes (perhaps
GDT's "associate").
> 9. a balloon or airplane trajectory:
>
> temperature(time)
> ch4(time)
> lat(time)
> lon(time)
> elevation(time)
Your solution:
> :coordinates.geodetic = "lon lat elevation";
My solution:
:coordinates="geodetic_1";
:geodetic_1"=latitude=lat, longitude=lon, altitude=elevation";
>
> 10. moving coordinate system. For example, satellite images tracking a
> tropical hurricane keep the coordinate system centered on the hurricane
> as it moves:
> pressure(time, x, y)
> lat(time, x, y)
> lon(time, x, y)
Your solution:
> :coordinates.geodetic = "lon lat pressure";
My solution:
What is the "data" variable here? It appears to be "pressure",
as if you are monitoring sea-level pressure. It that case,
"pressure" should probably not be in your list of coordinates,
and my solution would be:
:coordinates="geodetic_1";
:geodetic_1"=latitude=lat, longitude=lon";
>
> 11. multiple time coordinates:
>
> var(time)
> year(time)
> day_of_year(time)
> second_of_day(time)
Your solution:
> :coordinates = "year, day_of_year, second_of_day";
? This makes it appear that "var" is occupying a 3-dimensional
space. Or are you assuming that that fact that the units on
"year", "day_of_year", and "second_of_day" will imply that
there is really only 1 axis?
My solution:
short time(time); // simply index values, unless you
// have a traditional value, too
time:component="year, day_of_year, second_of_day"; // a la GDT
> and the famous NUWG case:
> var(time)
> generate_time(time)
> valid_time(time)
Your solution:
> :coordinates = "generate_time, valid_time";
My solution:
time:time="generate_time, valid_time";
or
time:associate="generate_time, valid_time";
>
> 12. sparse data variables. Adapted from Harvey Davies' example:
> soil.temperature(time, land_point)
> land_index(lat, lon) // if over land = land_point index, else -1
> lat(lat)
> lon(lon)
> time(time)
>
> or another way of getting a similar result:
> soil.temperature(time, land_point)
> latidx(land_point) // latitude index of ith land point
> lonidx(land_point) // longitude index of ith land point
> lat(lat)
> lon(lon)
> time(time)
Your solution:
> :coordinates = "lat*latidx, lon*lonidx, time"
I am a bit apprehensive of adding a "language" with a certain
nomenclature that has to be parsed out.
My solution:
Go back to the first suggestion above...
I am not hard-over on this, but I like the idea of having a
"mask" readily available for visual referral:
float soil.temperature(time, land_point)
soil.temperature:mask="land_point=land_index";
short land_index(lat,lon); // could also be 0's and 1's
:coordinates="geodetic_1";
:geodetic_1"=latitude=lat, longitude=lon";
So, at land_point=i, "mask" tells you to look at the "land_index"
array to find the point with a value of "i" (or, if array is
simple 0's and 1's, find the i-th point that is "1"), then get
the latitude and longitude from "lat" and "lon".
>
> 13. reduced grid From CSM: the number and spacing of the longitude elements
> decreases toward the pole:
>
> Solution 1: use maximum number of longitude points, and put "missing
> data" in the extra elements (note this violates the rule "no missing
> data in coordinates", but it does seem natural):
>
> var(nlat, max_lon)
> lat(nlat)
> lon(nlat,max_lon)
Your solution:
> In solution one, we would have to allow missing values in coordinates,
> which seems ok.
My solution:
I'm stumped.
>
> Solution 2: To save storage, you could store just the real points, but
> you lose the 2-D "connectedness" of the coordinate system:
>
> var(npoints)
> lat(npoints)
> lon(npoints)
Your solution:
> :coordinates = "lon lat";
My solution:
:coordinates = ""geodetic_1";
:geodetic_1"=latitude=lat, longitude=lon";
-------------------------
That's it for now.
Now, about those other 8 projects....
John
