3

I'm trying to take a set of NetCDFs and turn them into georeferenced datasets using GDAL.

However, the only Geospatial data is a set of 2-D lat/lon arrays. The projection type is not lat/lon. THe arrays give the lat/lon of every individual cell.

I can't use a geotransform to reference the data because of the irregular gridding.

How can I create a dataset (such as a GeoTiff) that properly references the data?

Ideally I would then warp to some even grid.

I tried using an even spread of Ground Control Points, but the resulting GeoTiff did not display correctly... Despite giving GCPs up to 90 deg lat, ArcMap claims the GTiff has an extent ending at 70 deg latitude...

I was able to solve the issue using gdalwarp -geoloc as suggested. I created 3 VRTs (1 with cell values, 1 with latitude, 1 with longitude), and wrote the lat/lon vrts in a geoarrays in the data vrt metadata. Then using gdalwarp -geoloc worked like a charm.

The issue with extent was unrelated and due simply to the large difference between the original coordinate system and lat/lon. It was solved by explicitely stating the extent instead of having GDAL try to guess what it should be. (i.e. I added -te -180 0 180 90 to the gdalwarp line)

PolyGeo
  • 65,136
  • 29
  • 109
  • 338
JFerg
  • 347
  • 1
  • 2
  • 9
  • Could someone, @AndreJ or JFerg or ..., please explain this: "the only Geospatial data is a set of 2-D lat/lon arrays. The projection type is not lat/lon. THe arrays give the lat/lon of every individual cell."? I cannot grasp the situation at all. – Martin F Mar 03 '15 at 00:51
  • There are some netcdf and hdf files out there that do not follow a ceratin CRS, but consist of several subdatasets, where the geolocation in lat/lon of a cell is stored in two subdatasets X and Y, while the actual cell information is stored in the other subdatasets. This is useful for satellite gathered information, like MODIS SWATH data. Gdalinfo is able to report this geolocation information stored in the files metadata. – AndreJ Mar 03 '15 at 06:08

1 Answers1

1

gdalwarp -geoloc allows you to use the complete 2d-array of latlon as georeference. With that, you can use any target CRS to reproject your data to a commonly used projection.

See my answer to this question for an example:

How to match a raster NetCDF data with a vector layer in QGIS?

Community
  • 1
AndreJ
  • 76,698
  • 5
  • 86
  • 162
  • When I attempt to run the code in the question you linked, I get the following error: ERROR 1: Unable to compute a GEOLOC_ARRAY based transformation between pixel/line and georeferenced coordinates. I get this using the files from the question too. Any idea what it means? Thanks for your response. – JFerg Mar 02 '15 at 20:35
  • Have you selected the subdatasets correctly? Try gdalinfo to get the correct subdataset name. – AndreJ Mar 03 '15 at 06:09
  • I created a set of vrts that were interconnected (1 for data, 1 for lat, 1 for lon) and used gdalwarp -geoloc on the data vrt. This worked weel enough but the issue with the cutting off of the data above 70 deg lat remains. I think something may be wrong with the initial file. Thank you for your help - geoloc is very useful! – JFerg Mar 03 '15 at 16:06
  • For the data loss, this might help: http://gis.stackexchange.com/questions/134239/gdal-warping-raster-to-new-domain-resolution-projection-causes-data-loss – AndreJ Mar 03 '15 at 16:44
  • 1
    @JFerg - Could you briefly explain how you created the "set of vrts that were interconnected"? Can't figure out how you finally managed to use gdalwarp -geoloc after getting the error you mentioned in the comment above – s6hebern Mar 05 '18 at 07:54