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/
