Re: CF-1.4 conventions for lat/lon datum.

Ethan and John,

We are are using WGS84 in our software internally and would request GeoTIFF 
referenced with WGS84 if using WCS.  I would love to implement a complete 
GeoTIFF IOSP implementation supporting all GeoTIFF supported projections and 
CRS but thats not going to happen.  There are too many options and with the CF 
CRS definition requirements I would need to dig up the parameters for each CRS. 
 I really wish that CF would have allowed to use of the 'crs_id' attribute in 
lieu of defining each parameters explicitly...

The NAD83 vs. WGS84 issue really isn't an issue.  I was just confused as to the 
meaning  'longitude_of_prime_meridian'.  Is the prime meridian the nominal 
prime meridian (line on the ground in Greenwich) or the IERS reference?  The CF 
example has 0.0, the ProjWKT for WGS84 lists 0.0.  It's just confusing... 

As far as handling of these params in NetCDF-Java.  I just wanted to make sure 
that if we wanted to do a CRS transform it wasn't also being done in 
NetCDF-Java.

Tom

On Jun 24, 2010, at 11:51 PM, Ethan Davis wrote:

> Hi all,
> 
> First off, when you say you can "explicitly request [WGS84 and NAD83]
> from the WCS server", does that mean you are getting the data on a
> lat/lon grid based on the WGS84 or NAD83 datum?
> 
> Second, is the data you are using of high enough resolution that the
> differences between these datums will be noticeable?
> 
> The more I read about datums and such the more I realize I don't
> understand. From a quick google search, I found a few pages that are
> both enlightening and baffling. From one of them:
> 
>> http://sci.tech-archive.net/Archive/sci.geo.satellite-nav/2006-09/msg00307.html
> 
> it appears that WGS84 has had a number of "realizations" over the years.
> One of them changed the location of the origin (center of the earth) by
> a few meters and another changed the location of the PM by 100 meters.
> (I'm really not sure I understand how or why the origin was changed.)
> Bottom line, the reference system of the original WGS84 was very close
> to that of NAD83 but the differences between current versions are more
> noticeable.
> 
> Here's what the "Remarks" in the EPSG datum table says about the WGS84:
> 
>    EPSG's WGS 84 datum has been the then current realisation. No
>    distinction is made between the original WGS 84 frame, WGS 84
>    (G730), WGS 84 (G873) and WGS 84 (G1150). Since 1997, WGS 84
>    has been maintained within 10cm of the then current ITRF.
> 
> As for describing them in CF, I'd say that if you know which realization
> of the datum actually applies to your data then you should use the
> appropriate numbers. On the other hand, none of us involved in writing
> this part of CF really know all the ins-and-outs of datums and such and
> I now wonder if there aren't some holes in the CF spec. In particular,
> there is nothing in the CF spec that would allow the data reader to
> figure out that the datums used on two different datasets place the
> origin of the reference system 2 meters apart. Though from the above
> quote it look like the EPSG doesn't either (and for that matter I'm not
> sure WKT and such have places to contain that information either).
> 
> As to what defining these CF parameters actually does, basically it
> captures the information in a way that hopefully allows client
> applications to use it in appropriate ways.
> 
> What the netCDF-Java library does with this information John can answer
> better than I. However, to give my guess, if the data is on lat/lon
> coordinates, I'm not sure if it even looks at the datum information
> underlying those coordinates when deciding if the data are on the same grid.
> 
> As always, whenever I spend time thinking or reading about datum, CRS,
> and such I end up a bit more uncertain about how it all works than when
> I started.
> 
> Ethan
> 
> 
> On 6/15/2010 5:36 PM, John Caron wrote:
>> Hi Dave, Tom:
>> 
>> Sorry Im having trouble getting to this and other emails. Im leaving for
>> vacation on thursday. I will be on email, so i may be able to respond
>> more coherently later. anyway, im also cc'ing ethan who actually knows
>> more about these things than I do. But I also guess you guys know more
>> than we both do, so probably just do your best guess. I would want to
>> get a hold of some sample data to test things out.
>> 
>> glad you got the shwag. I think Ben Dominico from Unidata is at the OGC
>> this week, you should say hi to him. Hes not technical but hes a nice
>> guy anyway ;^)
>> 
>> On 6/14/2010 6:48 AM, David Blodgett wrote:
>>> My understanding is that NAD83 will have a slightly different PM than
>>> WGS84, but the difference is difficult to define. I'm still foggy on
>>> exactly how NAD83 changes... (it does, staying static over the "mean?"
>>> position of the north american continental plate). WGS84 is static
>>> relative to the center of mass of the earth and a number of star
>>> observations.
>>> 
>>> The important piece is that differences are quite small, and there is
>>> no great method of converting between the two (since NAD83 changes
>>> periodically you would need a bunch more metadata and a whole bunch of
>>> tectonic data to do the shift accurately). 
>>> 
>>> ahh,,, they don't have the same geoid exactly... they are whole
>>> micrometers different! heh... something about needing the accuracy for
>>> sattelite orbital calcs. (yay wikipedia!)
>>> 
>>> So, take home, I think we can ignore the distinction between NAD83 and
>>> WGS84 (for calculations on the GDP). I also think we shouldn't support
>>> data in NAD27 or any of the boutique state based reference systems.
>>> 
>>> Dave Blodgett
>>> Center for Integrated Data Analytics (CIDA)
>>> USGS WI Water Science Center
>>> 8505 Research Way Middleton WI 53562
>>> 608-821-3899 | 608-628-5855 (cell)
>>> http://cida.usgs.gov
>>> 
>>> On Jun 14, 2010, at 12:27 AM, Tom Kunicki wrote:
>>> 
>>>> 
>>>> Hey John,
>>>> 
>>>> I'm trying to tidy up the GeoTIFF IOSP.  For now we're going to only
>>>> support WGS84 or NAD83 (as we can explicitly request these from the
>>>> WCS server) and I want to make sure I declare these datums
>>>> correctly..  I am looking at the example for WGS84 (4326) in the
>>>> CF-1.4 spec [1].  When looking at the the explanation of the the
>>>> "longitude_of_prime_meridian" it is stated as "0.0" in the example.  
>>>> The definition of the attribute states "Specifies the longitude,
>>>> with respect to Greenwich, of the prime meridian associated with the
>>>> geodetic datum." [2]   It was my understanding that WGS84 used the
>>>> IERS reference meridian, which is ~5.?? arcsec (~100m at 51N lat)
>>>> east of the prime meridian [3].   If this is true, shouldn't the
>>>> "longitude_of_prime_meridian" be non-zero, something around 0.001475?
>>>> Just want to make sure I understand these parameters correctly.
>>>> 
>>>> Out of curiosity, what does defining these CF CRS parameters do?  Is
>>>> it really even necessary?
>>>> 
>>>> Dave, do NAD83 and WGS84 share the same prime meridian?  Between the
>>>> they appear to share the same semi-major-axis with a minor deviation
>>>> in inverse-flattening.
>>>> 
>>>> btw-  Thanks for the shwag.  I will be sporting my Unidata T-Shirt at
>>>> the OGC meetings this week.  Let me know if you'll be there by
>>>> chance.  I noticed a CF-netCDF SWG meeting on Thursday, I will pop in
>>>> to see what's shaking out.
>>>> 
>>>> Below is a small example of the JAI GeoTIFF IOSP.readData.
>>>> Performance looks good but a potential 2-3x memory overhead on the
>>>> read using the currently logic below (defensive copies in the JAI code)
>>>> 
>>>> [1] http://cf-pcmdi.llnl.gov/documents/cf-conventions/1.4/ch05s06.html
>>>> [2] http://cf-pcmdi.llnl.gov/documents/cf-conventions/1.4/apf.html
>>>> [3] http://en.wikipedia.org/wiki/Prime_Meridian#IERS_Reference_Meridian
>>>> 
>>>> Tom Kunicki
>>>> Center for Integrated Data Analytics
>>>> U.S. Geological Survey
>>>> 8505 Research Way
>>>> Middleton, WI  53562
>>>> 
>>>> [...]
>>>>    @Override
>>>>    public Array readData(Variable variable, Section section)
>>>>            throws IOException, InvalidRangeException {
>>>> 
>>>>        int image = Integer.MIN_VALUE;
>>>>        int band = Integer.MIN_VALUE;
>>>> 
>>>>        for(Attribute attribute : variable.getAttributes()) {
>>>>            if (ATTRIBUTE_IMAGE.equals(attribute.getName())) {
>>>>                image = attribute.getNumericValue().intValue();
>>>>            } else if (ATTRIBUTE_BAND.equals(attribute.getName())) {
>>>>                band = attribute.getNumericValue().intValue();
>>>>            }
>>>>        }
>>>> 
>>>>        if (image < 0 || band < 0) {
>>>>            throw new IllegalArgumentException("ERROR:  Can't extract
>>>> image "
>>>>                    + "or band number from variable.");
>>>>        }
>>>> 
>>>>        ImageReader imageReader = null;
>>>>        ImageInputStream imageInputStream = null;
>>>>        try {
>>>> 
>>>>            Range yRange = section.getRange(0);
>>>>            Range xRange = section.getRange(1);
>>>> 
>>>>            imageReader = createImageReader();
>>>> 
>>>>            if (section.getStride(1) > 1 &&
>>>>                variable.getDataType().isFloatingPoint() &&
>>>> 
>>>> IMAGEIO_JAI_READER.equals(imageReader.getClass().getName())) {
>>>>                throw new IllegalArgumentException(
>>>>                        "ERROR:  Data corruption bug in Java JAI
>>>> ImageIO TIFF "
>>>>                        + " ImageReader with floating point data and
>>>> X strides > "
>>>>                        + " 1.  To utilize this functionality please
>>>> obatain another"
>>>>                        + "JAI TIFF implementation such as
>>>> ImageIO-Ext and place"
>>>>                        + "in the application classpath");
>>>>            }
>>>> 
>>>>            imageInputStream = ImageIO.createImageInputStream(new
>>>> File(location));
>>>>            imageReader.setInput(imageInputStream, true, true);
>>>> 
>>>>            ImageReadParam imageReadParam =
>>>> imageReader.getDefaultReadParam();
>>>>            imageReadParam.setSourceRegion(new Rectangle(
>>>>                    xRange.first(),
>>>>                    yRange.first(),
>>>>                    xRange.last() - xRange.first() + 1,
>>>>                    yRange.last() - yRange.first() + 1));
>>>>            imageReadParam.setSourceSubsampling(
>>>>                    xRange.stride(),
>>>>                    yRange.stride(),
>>>>                    0,
>>>>                    0);
>>>>            imageReadParam.setSourceBands(new int[] { band } );
>>>> 
>>>>            BufferedImage bufferedImage = imageReader.read(image,
>>>> imageReadParam);
>>>> 
>>>>            Raster raster = bufferedImage.getData();
>>>> 
>>>>            Object arrayAsObject = createJavaArray(
>>>>                    variable.getDataType(), (int)section.computeSize());
>>>> 
>>>>            Array array = Array.factory(
>>>>                    variable.getDataType(),
>>>>                    section.getShape(),
>>>>                    arrayAsObject);
>>>> 
>>>>            raster.getDataElements(
>>>>                    0, 0,
>>>>                    bufferedImage.getWidth(), bufferedImage.getHeight(),
>>>>                    arrayAsObject);
>>>> 
>>>>            return array;
>>>> 
>>>>        } finally {
>>>>            if (imageReader != null) {
>>>>                imageReader.dispose();
>>>>            }
>>>>            if (imageInputStream != null) {
>>>>                try { imageInputStream.close(); } catch (IOException
>>>> e) { }
>>>>            }
>>>>        }
>>>>    }
>>>> [...]
>>>> 
>>>> 
>>> 
>>