etiennesky / gdal-frmts-netcdf-issues Goto Github PK
View Code? Open in Web Editor NEWgdal netcdf driver (not full clone)
gdal netcdf driver (not full clone)
There are several cases where the NetCDF driver's import capabilities aren't properly dealing with the NetCDF projection info -> OGC WKT format, especially in cases where NetCDF uses the one projection name to deal with multiple 'variants' - eg Mercator_1SP and Mercator_2SP, and possibly Lambert Cylindrical Equal Area having different variants.
I just committed a fix to handle Mercator 1SP vs 2SP, but should do an audit of other projections too.
...to take into account current thinking and ideas, and send to ET.
This is then basis to ask for input/feedback on mailing list.
Relevant from chat logs:
ust a comment/question: I feel that the two sections named "Continue saving NetCDF custom attribute tags by default, but add option not to save" and "Add saving of Datum reference Ellipsoid parameters in CF-1 compliant format, in addition to within GDAL's WKT" are redundant. perhaps they should be merged or re-written
so I would see this in 2 parts 1) just mention the options to save wkt, and drop the NN,SE,WE,EE and 2) upon import test for GDAL projection+ datum vs CF projection and dump the GDAL if it conflicts with the CF one
I mean parts that are in 1) should go in 2)
also re. my comments on that wiki. A) there shouldn't be an option to force assume south-regular B) about the testing of datums, that should remain an open question left to debate C) I added a section about new attributes (history, GDAL version and CF-version)
Re: C) Perhaps these should be a top-level heading? I don't think they'd be controversial, but they're a bit hard to find currently.
1- Gdal 1.8.1 not setting GeoTransform properly when reading from (new) exported NetCDFs
Just testing and found a back-compatibility issue with current GDAL 1.8.1 code:
(1) exported files using new GDAL NetCDF driver developed here (eg as done during the netcdf_cf tests)
(2) Try to do a gdalinfo on them with gdal 1.8.1, and the GeoTransform isn't detected correctly:
pds@pds-mbp-laptop:~/GIS_Projects/GDAL-NetCDF-Testing/gdal_trans_testing_all_cf1_weld/savedResults/Melb_USGS/refactored_gdal_gtiff$ /usr/local/gdal_1.8.1/bin/gdalinfo melb-small_AEA.nc
Driver: netCDF/Network Common Data Format
Files: melb-small_AEA.nc
Size is 199, 101
Coordinate System is:
PROJCS["GDA94 / Australian Albers",
GEOGCS["GDA94",
DATUM["Geocentric_Datum_of_Australia_1994",
SPHEROID["GRS 1980",6378137,298.2572221010042,
AUTHORITY["EPSG","7019"]],
TOWGS84[0,0,0,0,0,0,0],
AUTHORITY["EPSG","6283"]],
PRIMEM["Greenwich",0],
UNIT["degree",0.0174532925199433],
AUTHORITY["EPSG","4283"]],
PROJECTION["Albers_Conic_Equal_Area"],
PARAMETER["standard_parallel_1",-18],
PARAMETER["standard_parallel_2",-36],
PARAMETER["latitude_of_center",0],
PARAMETER["longitude_of_center",132],
PARAMETER["false_easting",0],
PARAMETER["false_northing",0],
UNIT["metre",1,
AUTHORITY["EPSG","9001"]],
AUTHORITY["EPSG","3577"]]
Origin = (0.000000000000000,0.000000000000000)
Pixel Size = (-0.000000000000000,0.000000000000000)
Where origin and pixel size should be:
Origin = (1232736.866607666946948,-4118493.958211950026453)
Pixel Size = (30.080459243366306,-30.080459243366494)
Reason why this is happening:
The logic in the old 1.8.1 code inside SetProjection() has a flaw: it first (correctly) reads the GeoTransform from X and Y coordinate variables, but then if the GDAL-added OGC WKT string is added, first looks for 'GeoTransform' string to override the above, and then the 4 corner attributes. And even if the last of these isn't found, it goes ahead and sets blank values.
So I guess there are 2 options:
(1) Update the new code to deal with the issue: and write the "GeoTransform" string whenever WRITE_GDAL_TAGS is enabled;
(2) Publish the back-compat issue, and update the 1.8 driver.
I know it's not ideal, but I slightly favour the former .... to reduce compatibility problems for installed userbase of people using gdal 1.8 and who will want to be able to process files written with the new driver.
We count how many examples before start executing example, we used to use --dry--run option to achieve this in Rspec 1. Now dry run is removed in Rspec 2, how do I count examples without executing examples?
Thanks.
As discussed in Google chat.
Aim for:
Send note to ET for comment when done.
me: 1) what about my suggestion (well yours originally) of refactoring createcopy() so that its parts can be used in create(), and also for optimizing file creation (header first, data last)
Sent at 11:34 PM on Thursday
Patrick: Yes, I think that's the right approach. NetCDF really is quite modular with respect to this, so if we have clearly deliniated stages (header first, then data), then code should be mostly reusable between the two, and also IWriteBlock().
Probably the trick is with Projection attributes - and I'm not sure if with the Create() method, you're told the projection to use at create stage - think that you are.
me: there's going to be a little bit or code duplication (mainly for the projected/geographic separation), but I think it's worth it.
Sent at 11:37 PM on Thursday
me: actually the projection is not known during the Create(). So we have to make sure (or hope that) users do the following 1) Create 2) SetProjection 3) IO. look at the API examlpes, that is what is done usually
what we could do is print a warning, if IO is done before SetProjection()
Sent at 11:39 PM on Thursday
me: I think that a helper class might be useful for the Create() and CreateCopy() so that the many variables (currently in Createcopy()) are not passed as parameters for each of the helper functions
Sent at 11:41 PM on Thursday
Patrick: Ok - just a note that on http://www.unidata.ucar.edu/software/netcdf/workshops/most-recent/classic_perf/FileFormat.html , it says you can reserve extra space in the header, which could be a good approach given the CF-1 projection attributes are of fairly constant size presumably
Sent at 11:42 PM on Thursday
me: that's a good point - how much space could we loose you think?
in case the header is much smaller than a reasonable space
Sent at 11:43 PM on Thursday
Patrick: I think it's pretty small - IE we're talking probably a few hundred bytes for the grid_mapping variable, vs KB to GB of raster data in most realistic cases.
me: ok, let's wtick to that for now
see gdal bug#4221
requires netcdf-4 and hdf5, should test for them (already code in the driver for file detection)
should set filetype also, perhaps need an option? e.g. FILETYPE=NC/NC2/NC4
COMPRESS=NONE/DEFLATE/PACKING and ZLEVEL=[1-9]
http://www.unidata.ucar.edu/software/netcdf/workshops/2011/nc4chunking/CompressionAdvice.html
http://www.unidata.ucar.edu/software/netcdf/docs/netcdf-c/nc_005fdef_005fvar_005fdeflate.html
http://www.unidata.ucar.edu/software/netcdf/workshops/2008/nc4chunking/Compression.html
example
ncdump -h MCD12Q1.IGBP.SA.2002.wgs84.mosaic_0p5.nc
netcdf MCD12Q1.IGBP.SA.2002.wgs84.mosaic_0p5 {
dimensions:
lon = 94 ;
lat = 138 ;
variables:
float Band1(lat, lon) ;
Band1:long_name = "GDAL Band Number 1" ;
Band1:_FillValue = -1.e+10f ;
Band1:vegtype_id = 1 ;
Band1:vegtype_name = "Evergreen Needleleaf forest" ;
Band1:grid_mapping = "crs" ;
Band1:coordinates = "lon lat" ;
...
float Band18(lat, lon) ;
Band18:long_name = "GDAL Band Number 18" ;
Band18:_FillValue = -1.e+10f ;
Band18:vegtype_id = 255 ;
Band18:vegtype_name = "Fill Value" ;
Band18:grid_mapping = "crs" ;
Band18:coordinates = "lon lat" ;
char crs ;
crs:grid_mapping_name = "latitude_longitude" ;
crs:longitude_of_prime_meridian = 0. ;
crs:semi_major_axis = 6378137. ;
crs:inverse_flattening = 298.257223563 ;
crs:spatial_ref = "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\"]]" ;
crs:GeoTransform = "-81.4999999545239 0.5000000000000004 0 12.99999999820093 0 -0.5000000000000004 " ;
double lat(lat) ;
lat:standard_name = "latitude" ;
lat:long_name = "latitude" ;
lat:units = "degrees_north" ;
double lon(lon) ;
lon:standard_name = "longitude" ;
lon:long_name = "longitude" ;
lon:units = "degrees_east" ;
netcdf-4 w/zip and also szip for input
packing: int variables with add_offset and scale_factor (i.e. ERA40)
COMPRESS=NONE/DEFLATE/PACKING and ZLEVEL=[1-9]
http://www.gdal.org/frmt_gtiff.html
http://www.unidata.ucar.edu/software/netcdf/workshops/2011/nc4chunking/CompressionAdvice.html
A declarative, efficient, and flexible JavaScript library for building user interfaces.
๐ Vue.js is a progressive, incrementally-adoptable JavaScript framework for building UI on the web.
TypeScript is a superset of JavaScript that compiles to clean JavaScript output.
An Open Source Machine Learning Framework for Everyone
The Web framework for perfectionists with deadlines.
A PHP framework for web artisans
Bring data to life with SVG, Canvas and HTML. ๐๐๐
JavaScript (JS) is a lightweight interpreted programming language with first-class functions.
Some thing interesting about web. New door for the world.
A server is a program made to process requests and deliver data to clients.
Machine learning is a way of modeling and interpreting data that allows a piece of software to respond intelligently.
Some thing interesting about visualization, use data art
Some thing interesting about game, make everyone happy.
We are working to build community through open source technology. NB: members must have two-factor auth.
Open source projects and samples from Microsoft.
Google โค๏ธ Open Source for everyone.
Alibaba Open Source for everyone
Data-Driven Documents codes.
China tencent open source team.