Re: [thredds] GeoTIFF response from THREDDS WCS seems displaced by 1/2 grid cell

  • To: John Caron <caron@xxxxxxxxxxxxxxxx>
  • Subject: Re: [thredds] GeoTIFF response from THREDDS WCS seems displaced by 1/2 grid cell
  • From: Tom Kunicki <tkunicki@xxxxxxxx>
  • Date: Fri, 15 Mar 2013 15:01:32 -0500
Hi John and Rich,

From what Rich has described I would argue this is an issue with THREDDS WCS 
GeoTIFF output.

If the grid data represents an area (GTRasterTypeGeoKey=PixelIsArea),  as 
opposed to a point (GTRasterTypeGeoKey=PixelIsPoint) the origin (value stored 
with ModelTiepointTag) needs to be offset a half pixel.  That is, the origin 
(aka "tie point") is the outer corner of the (0,0) pixel and not the center  
(per http://www.remotesensing.org/geotiff/spec/geotiff2.5.html).

From what Rich has shown below the PixelIsArea flag is used but the tie point 
value is equal to the cell center at (0,0).   The tie point needs to be the 
outer corner and not the center of the cell at (0,0).

An alternate resolution would be to use the PixelIsPoint flag as this would 
allow one to set the tiepoint using the cell center.  But PixelIsPoint probably 
doesn't truly represent the pixel value.

Tom Kunicki
Center for Integrated Data Analytics
U.S. Geological Survey
8505 Research Way
Middleton, WI  53562

On Mar 14, 2013, at 5:42 PM, John Caron <caron@xxxxxxxxxxxxxxxx> wrote:


> Hi Rich:
> 
> it appears geotiff can do either midpoint or vertex,
> 
> see"PixelIsArea" or "PixelIsPoint":
> 
>  http://www.remotesensing.org/geotiff/spec/geotiff2.5.html
> 
> it sems likely that whatever underlies the "bathy/gom03_v1_0" dataset might 
> be getting this wrong. Or our geotiff writer is ignoring the distinction.
> 
> might be good to look at the underlying file.
> 
> John
> 
> On 3/14/2013 4:24 PM, Signell, Richard wrote:
>> TDS folks,
>> 
>> I think the GeoTIFF response from THREDDS WCS is displaced by 1/2 grid
>> cell, caused by center values being written into what should be corner
>> values in the GeoTIFF file.
>> 
>> Here's what I did:
>> 1. Request a small chunk of bathymetry (with ~ 90m grid spacing) as NetCDF3:
>> 
>> curl -o gom3_tiny.nc
>> 'http://geoport.whoi.edu/thredds/wcs/bathy/gom03_v1_0?SERVICE=wcs&VERSION=1.0.0&REQUEST=GetCoverage&COVERAGE=topo&CRS=OGC:CRS84&BBOX=-71.5,39.5,-71.0,40.0&FORMAT=NetCDF3'
>> 
>> 2. Request the same chunk as GeoTIFFfloat:
>> curl -o gom3_tiny.nc
>> 'http://geoport.whoi.edu/thredds/wcs/bathy/gom03_v1_0?SERVICE=wcs&VERSION=1.0.0&REQUEST=GetCoverage&COVERAGE=topo&CRS=OGC:CRS84&BBOX=-71.5,39.5,-71.0,40.0&FORMAT=NetCDF3'
>> 
>> 3. Drag and drop both datasets into ArcGIS and contour them at 100 m 
>> intervals.
>> 
>> 4. Zoom in and notice that the contours do not lie on top of each
>> other.   The GeoTIFF contours are shifted 1/2 grid cell (~45 m) to the
>> north and east.
>> 
>> I think the reason is that coordinates values of lon and lat in NetCDF
>> represent the centers of the grid cell, but coordinates in GeoTIFF are
>> measured from the *corner* of the grid cell, not the grid cell center.
>> 
>> Yet if we look at the NetCDF file, the minimum lon:
>> 
>> $ ncks -v lon -d lon,0  gom3_tiny.nc | tail -2
>> lon[0]=-71.4995833333 degrees_east
>> 
>> and the maximum lat:
>> 
>> $ ncks -v lat -d lat,599 gom3_tiny.nc | tail -2
>> lat[599]=39.9995833333 degrees_north
>> 
>> are what end up getting assigned to the lower bound for lon and the
>> upper bound for lat in the GeoTIFF:
>> 
>> c:\Users\rsignell\Downloads>gdalinfo gom3_tiny.tif
>> Warning 1: TIFFFetchNormalTag:ASCII value for tag "Software" does not
>> end in null byte
>> Warning 1: TIFFFetchNormalTag:Incompatible type for "GDALNoDataValue";
>> tag ignored
>> Driver: GTiff/GeoTIFF
>> Files: gom3_tiny.tif
>> Size is 600, 600
>> Coordinate System is:
>> GEOGCS["WGS 84",
>>     DATUM["WGS_1984",
>>         SPHEROID["WGS 84",6378137,298.257223563,
>>             AUTHORITY["EPSG","7030"]],
>>         AUTHORITY["EPSG","6326"]],
>>     PRIMEM["Greenwich",0],
>>     UNIT["degree",0.0174532925199433],
>>     AUTHORITY["EPSG","4326"]]
>> Origin = (-71.499580383300781,39.999583333336069)
>> Pixel Size = (0.000833333333333,-0.000833333333332)
>> Metadata:
>>   TIFFTAG_SOFTWARE=nc2geotiff
>>   AREA_OR_POINT=Area
>> Image Structure Metadata:
>>   INTERLEAVE=BAND
>> Corner Coordinates:
>> Upper Left  ( -71.4995804,  39.9995833) ( 71d29'58.49"W, 39d59'58.50"N)
>> Lower Left  ( -71.4995804,  39.4995833) ( 71d29'58.49"W, 39d29'58.50"N)
>> Upper Right ( -70.9995804,  39.9995833) ( 70d59'58.49"W, 39d59'58.50"N)
>> Lower Right ( -70.9995804,  39.4995833) ( 70d59'58.49"W, 39d29'58.50"N)
>> Center      ( -71.2495804,  39.7495833) ( 71d14'58.49"W, 39d44'58.50"N)
>> Band 1 Block=600x1 Type=Float32, ColorInterp=Gray
>> 
>> 
>> I guess the good news is that not many people use the WCS response.
>> And the other good news is that it will take 1 minute to fix.
>> 
>> I'm attaching a screenshot of my ArcGIS session.
>> 
>> -Rich
>> 
>> 
>> 
>> 
>> _______________________________________________
>> thredds mailing list
>> thredds@xxxxxxxxxxxxxxxx
>> For list information or to unsubscribe,  visit: 
>> http://www.unidata.ucar.edu/mailing_lists/
>> 
> 
> _______________________________________________
> thredds mailing list
> thredds@xxxxxxxxxxxxxxxx
> For list information or to unsubscribe,  visit: 
> http://www.unidata.ucar.edu/mailing_lists/ 

  • 2013 messages navigation, sorted by:
    1. Thread
    2. Subject
    3. Author
    4. Date
    5. ↑ Table Of Contents
  • Search the thredds archives: