[GRASS-user] Issue with netcdf import into GRASS 7

All:

Using the help found here: http://grasswiki.osgeo.org/wiki/NetCDF I have been able to import a netcdf file into GRASS 7, but the data is not scaled correctly.

My file is found here:
https://drive.google.com/file/d/0B-Dmtz0DFEiubG5aYTk4TDZERXc/view?usp=sharing

I created a GRASS LOCATION based on this:

char Lambert;
  :grid_mapping_name = "lambert_conformal_conic";
  :latitude_of_projection_origin = 37.99999237060547; // double
  :longitude_of_central_meridian = -106.0; // double
  :standard_parallel = 34.0, 42.0; // double
  :earth_radius = 6371229.0; // double
  :_CoordinateTransformType = "Projection";
  :_CoordinateAxisTypes = "GeoX GeoY";

which I obtained using panoply from here: [http://www.giss.nasa.gov/tools/panoply/](http://www.giss.nasa.gov/tools/panoply/)
Looking at the file netcdf header information, it's clear the dx=dy=1000m; but, what imports is x=dy=1m resolution even though I set the x & y resolution to 1000m
This is my import command:
r.in.gdal input=NETCDF:"[geo_em.d02.nc](http://geo_em.d02.nc)":GREENFRAC output=greenfrac -o -e --overwrite

which produces..
Warning 1: No UNIDATA NC_GLOBAL:Conventions attribute
WARNING: Over-riding projection check
 100%
Raster map <greenfrac.1> created.
 100%
Raster map <greenfrac.2> created.
 100%
Raster map <greenfrac.3> created.
 100%
Raster map <greenfrac.4> created.
 100%
Raster map <greenfrac.5> created.
 100%
Raster map <greenfrac.6> created.
 100%
Raster map <greenfrac.7> created.
 100%
Raster map <greenfrac.8> created.
 100%
Raster map <greenfrac.9> created.
 100%
Raster map <greenfrac.10> created.
 100%
Raster map <greenfrac.11> created.
 100%
Raster map <greenfrac.12> created.
Region for the current mapset updated
r.in.gdal complete.

I am virtually certain the lower left corner is placed correctly, because I can zoom to the lower-left location and it falls about where it should – however the spatial extent of the data is not what it should be, but when zoomed-in looks correct.

Any thoughts?

Thank you,
Tom

On Tue, Jun 9, 2015 at 4:47 AM, Thomas Adams <tea3rd@gmail.com> wrote:

All:

Using the help found here: http://grasswiki.osgeo.org/wiki/NetCDF I have
been able to import a netcdf file into GRASS 7, but the data is not scaled
correctly.

...

This is my import command:

r.in.gdal input=NETCDF:"geo_em.d02.nc":GREENFRAC output=greenfrac -o -e
--overwrite

which produces..

Warning 1: No UNIDATA NC_GLOBAL:Conventions attribute

...

I am virtually certain the lower left corner is placed correctly, because I
can zoom to the lower-left location and it falls about where it should --
however the spatial extent of the data is not what it should be, but when
zoomed-in looks correct.

Tom, I believe that the gdal list is the right place here since the
netCDF file is not read properly.
Once solved, r.in.gdal will run out of the box.

Perhaps they could give you hints how to convert the file with
gdalwarp into e.g. a GeoTIFF which you can then import without
problems.

Markus

Markus,

Thank you for the suggestion.

Best,
Tom

···

On Tue, Jun 9, 2015 at 3:45 AM, Markus Neteler <neteler@osgeo.org> wrote:

On Tue, Jun 9, 2015 at 4:47 AM, Thomas Adams <tea3rd@gmail.com> wrote:

All:

Using the help found here: http://grasswiki.osgeo.org/wiki/NetCDF I have
been able to import a netcdf file into GRASS 7, but the data is not scaled
correctly.

This is my import command:

r.in.gdal input=NETCDF:“geo_em.d02.nc”:GREENFRAC output=greenfrac -o -e
–overwrite

which produces…

Warning 1: No UNIDATA NC_GLOBAL:Conventions attribute

I am virtually certain the lower left corner is placed correctly, because I
can zoom to the lower-left location and it falls about where it should –
however the spatial extent of the data is not what it should be, but when
zoomed-in looks correct.

Tom, I believe that the gdal list is the right place here since the
netCDF file is not read properly.
Once solved, r.in.gdal will run out of the box.

Perhaps they could give you hints how to convert the file with
gdalwarp into e.g. a GeoTIFF which you can then import without
problems.

Markus

Hi Tom:

These geo_em* netcdf files do not have any “proper” CRS information in the headers. So you indeed get a simple X-Y matrix with cellsize 1 when you try to import directly using GDAL.

One of the ways to work around this (explained in NCAR’s user manual, chapter 4, about p. 59) is to export to an ARCInfo ASCII file, then manually edit the ASCII file header rows and change the xllcorner, yllcorner, and cellsize to the correct values.

So first:
gdal_translate -of AAIGrid NETCDF:“geo_em.d02.nc”:GREENFRAC greenfrac.asc
Now do:
ncdump geo_em.d02.nc | grep corner
:corner_lats = 36.19099f, 39.06561f, 39.08329f, 36.20801f, 36.19086f, 39.06547f, 39.08319f, 36.20791f, 36.18649f, 39.07011f, 39.0878f, 36.20351f, 36.18636f, 39.06997f, 39.0877f, 36.2034f ;
:corner_lons = -108.7585f, -108.8684f, -104.0023f, -104.0788f, -108.7641f, -108.8742f, -103.9965f, -104.0732f, -108.7584f, -108.8686f, -104.0021f, -104.0789f, -108.7639f, -108.8744f, -103.9963f, -104.0733f ;

The above “corner_lats” and “corner_lons” are the correct Lon/Lat for the corners and center for each of the four corner pixels in the domain 2. The lower left corner of the lower left pixel is at -108.7585,36.19099. You enter these into the AAI header rows, and set the cellsize to 1000, then import that altered AAIGrid into GRASS in a WGS84 location .

One caveat: this method assumes that the whole geo_em.d03 is projected in a Lon/Lat WGS84 coordinate system, which is wrong, as you know. The way I know of to get an accurate raster is to parse out the XLAT and XLONG variables, together with the GREENFRAC variable you are interested in as a CSV list of lon,lat,greenfac values, import the values as vector points, and do an interpolation.

Another possibility comes to mind, but I’ve never tried it: You should be able to convert the above xllcorner and yllcorner values to your Lambert Conformal Conic projection in advance, and use the converted values in the AAI header rows. Then import the AAI file directly into the correct LCC location. You would do this by:

  1. create a point vector of a single point (v.in.ascii) from the Long/Lat values above in an WGS84 location
  2. project that point to your LCC location (v.proj)
  3. get the LCC X-Y coordinates (v.out.ascii)
  4. Edit the header rows with the LCC xllcorner and yllcorner values
  5. Import into an LCC location (with r.in.gdal -o as you did)

I’d be interested to know if this works :slight_smile:

Regards,
Micha

···

On 06/09/2015 05:47 AM, Thomas Adams wrote:

All:

Using the help found here: http://grasswiki.osgeo.org/wiki/NetCDF I have been able to import a netcdf file into GRASS 7, but the data is not scaled correctly.

My file is found here:
https://drive.google.com/file/d/0B-Dmtz0DFEiubG5aYTk4TDZERXc/view?usp=sharing

I created a GRASS LOCATION based on this:

char Lambert;
  :grid_mapping_name = "lambert_conformal_conic";
  :latitude_of_projection_origin = 37.99999237060547; // double
  :longitude_of_central_meridian = -106.0; // double
  :standard_parallel = 34.0, 42.0; // double
  :earth_radius = 6371229.0; // double
  :_CoordinateTransformType = "Projection";
  :_CoordinateAxisTypes = "GeoX GeoY";

which I obtained using panoply from here: [http://www.giss.nasa.gov/tools/panoply/](http://www.giss.nasa.gov/tools/panoply/)
Looking at the file netcdf header information, it's clear the dx=dy=1000m; but, what imports is x=dy=1m resolution even though I set the x & y resolution to 1000m
This is my import command:
r.in.gdal input=NETCDF:"[geo_em.d02.nc](http://geo_em.d02.nc)":GREENFRAC output=greenfrac -o -e --overwrite

which produces..
Warning 1: No UNIDATA NC_GLOBAL:Conventions attribute
WARNING: Over-riding projection check
 100%
Raster map <greenfrac.1> created.
 100%
Raster map <greenfrac.2> created.
 100%
Raster map <greenfrac.3> created.
 100%
Raster map <greenfrac.4> created.
 100%
Raster map <greenfrac.5> created.
 100%
Raster map <greenfrac.6> created.
 100%
Raster map <greenfrac.7> created.
 100%
Raster map <greenfrac.8> created.
 100%
Raster map <greenfrac.9> created.
 100%
Raster map <greenfrac.10> created.
 100%
Raster map <greenfrac.11> created.
 100%
Raster map <greenfrac.12> created.
Region for the current mapset updated
r.in.gdal complete.

I am virtually certain the lower left corner is placed correctly, because I can zoom to the lower-left location and it falls about where it should – however the spatial extent of the data is not what it should be, but when zoomed-in looks correct.

Any thoughts?

Thank you,
Tom

This mail was received via Mail-SeCure System.

_______________________________________________
grass-user mailing list
[grass-user@lists.osgeo.org](mailto:grass-user@lists.osgeo.org)
[http://lists.osgeo.org/mailman/listinfo/grass-user](http://lists.osgeo.org/mailman/listinfo/grass-user)
This mail was received via Mail-SeCure System.

Micha,

Thank you very much for your detailed help; I’ll give these approaches a try and I’ll let you know what I get. This is very helpful; thank you! I hope all is well. BTW, thank you for the book contribution!!

Best,

Tom

···

On Tue, Jun 9, 2015 at 4:22 PM, Micha Silver <micha@arava.co.il> wrote:

Hi Tom:

These geo_em* netcdf files do not have any “proper” CRS information in the headers. So you indeed get a simple X-Y matrix with cellsize 1 when you try to import directly using GDAL.

One of the ways to work around this (explained in NCAR’s user manual, chapter 4, about p. 59) is to export to an ARCInfo ASCII file, then manually edit the ASCII file header rows and change the xllcorner, yllcorner, and cellsize to the correct values.

So first:
gdal_translate -of AAIGrid NETCDF:“geo_em.d02.nc”:GREENFRAC greenfrac.asc
Now do:
ncdump geo_em.d02.nc | grep corner
:corner_lats = 36.19099f, 39.06561f, 39.08329f, 36.20801f, 36.19086f, 39.06547f, 39.08319f, 36.20791f, 36.18649f, 39.07011f, 39.0878f, 36.20351f, 36.18636f, 39.06997f, 39.0877f, 36.2034f ;
:corner_lons = -108.7585f, -108.8684f, -104.0023f, -104.0788f, -108.7641f, -108.8742f, -103.9965f, -104.0732f, -108.7584f, -108.8686f, -104.0021f, -104.0789f, -108.7639f, -108.8744f, -103.9963f, -104.0733f ;

The above “corner_lats” and “corner_lons” are the correct Lon/Lat for the corners and center for each of the four corner pixels in the domain 2. The lower left corner of the lower left pixel is at -108.7585,36.19099. You enter these into the AAI header rows, and set the cellsize to 1000, then import that altered AAIGrid into GRASS in a WGS84 location .

One caveat: this method assumes that the whole geo_em.d03 is projected in a Lon/Lat WGS84 coordinate system, which is wrong, as you know. The way I know of to get an accurate raster is to parse out the XLAT and XLONG variables, together with the GREENFRAC variable you are interested in as a CSV list of lon,lat,greenfac values, import the values as vector points, and do an interpolation.

Another possibility comes to mind, but I’ve never tried it: You should be able to convert the above xllcorner and yllcorner values to your Lambert Conformal Conic projection in advance, and use the converted values in the AAI header rows. Then import the AAI file directly into the correct LCC location. You would do this by:

  1. create a point vector of a single point (v.in.ascii) from the Long/Lat values above in an WGS84 location
  2. project that point to your LCC location (v.proj)
  3. get the LCC X-Y coordinates (v.out.ascii)
  4. Edit the header rows with the LCC xllcorner and yllcorner values
  5. Import into an LCC location (with r.in.gdal -o as you did)

I’d be interested to know if this works :slight_smile:

Regards,
Micha

On 06/09/2015 05:47 AM, Thomas Adams wrote:

All:

Using the help found here: http://grasswiki.osgeo.org/wiki/NetCDF I have been able to import a netcdf file into GRASS 7, but the data is not scaled correctly.

My file is found here:
https://drive.google.com/file/d/0B-Dmtz0DFEiubG5aYTk4TDZERXc/view?usp=sharing

I created a GRASS LOCATION based on this:

char Lambert;
  :grid_mapping_name = "lambert_conformal_conic";
  :latitude_of_projection_origin = 37.99999237060547; // double
  :longitude_of_central_meridian = -106.0; // double
  :standard_parallel = 34.0, 42.0; // double
  :earth_radius = 6371229.0; // double
  :_CoordinateTransformType = "Projection";
  :_CoordinateAxisTypes = "GeoX GeoY";

which I obtained using panoply from here: [http://www.giss.nasa.gov/tools/panoply/](http://www.giss.nasa.gov/tools/panoply/)
Looking at the file netcdf header information, it's clear the dx=dy=1000m; but, what imports is x=dy=1m resolution even though I set the x & y resolution to 1000m
This is my import command:
r.in.gdal input=NETCDF:"[geo_em.d02.nc](http://geo_em.d02.nc)":GREENFRAC output=greenfrac -o -e --overwrite

which produces..
Warning 1: No UNIDATA NC_GLOBAL:Conventions attribute
WARNING: Over-riding projection check
 100%
Raster map <greenfrac.1> created.
 100%
Raster map <greenfrac.2> created.
 100%
Raster map <greenfrac.3> created.
 100%
Raster map <greenfrac.4> created.
 100%
Raster map <greenfrac.5> created.
 100%
Raster map <greenfrac.6> created.
 100%
Raster map <greenfrac.7> created.
 100%
Raster map <greenfrac.8> created.
 100%
Raster map <greenfrac.9> created.
 100%
Raster map <greenfrac.10> created.
 100%
Raster map <greenfrac.11> created.
 100%
Raster map <greenfrac.12> created.
Region for the current mapset updated
r.in.gdal complete.

I am virtually certain the lower left corner is placed correctly, because I can zoom to the lower-left location and it falls about where it should – however the spatial extent of the data is not what it should be, but when zoomed-in looks correct.

Any thoughts?

Thank you,
Tom

This mail was received via Mail-SeCure System.

_______________________________________________
grass-user mailing list
[grass-user@lists.osgeo.org](mailto:grass-user@lists.osgeo.org)
[http://lists.osgeo.org/mailman/listinfo/grass-user](http://lists.osgeo.org/mailman/listinfo/grass-user)
This mail was received via Mail-SeCure System.

Thomas E Adams, III
718 McBurney Drive
Lebanon, OH 45036

1 (513) 739-9512 (cell)

Micha,

It looks like your first suggestion worked. I did, however, have to add -b 1 to the end of your line:

gdal_translate -of AAIGrid NETCDF:“geo_em.d02.nc”:GREENFRAC greenfrac.asc -b 1

because I was getting the error:

Warning 1: No UNIDATA NC_GLOBAL:Conventions attribute
Input file size is 420, 320
ERROR 6: AAIG driver doesn’t support 12 bands. Must be 1 band.

There are 12 bands (1 for each month) in the file ‘geo_em.d02.nc’, so I had to specify which band number. Then I just needed to change the dx=dy from 1.0 to 1000. in the greenfrac.asc file. The imported file seems to be placed correctly with the right resolution.

Thank you so much for your help!

Cheers!
Tom

···

On Tue, Jun 9, 2015 at 5:22 PM, Micha Silver <micha@arava.co.il> wrote:

Hi Tom:

These geo_em* netcdf files do not have any “proper” CRS information in the headers. So you indeed get a simple X-Y matrix with cellsize 1 when you try to import directly using GDAL.

One of the ways to work around this (explained in NCAR’s user manual, chapter 4, about p. 59) is to export to an ARCInfo ASCII file, then manually edit the ASCII file header rows and change the xllcorner, yllcorner, and cellsize to the correct values.

So first:
gdal_translate -of AAIGrid NETCDF:“geo_em.d02.nc”:GREENFRAC greenfrac.asc
Now do:
ncdump geo_em.d02.nc | grep corner
:corner_lats = 36.19099f, 39.06561f, 39.08329f, 36.20801f, 36.19086f, 39.06547f, 39.08319f, 36.20791f, 36.18649f, 39.07011f, 39.0878f, 36.20351f, 36.18636f, 39.06997f, 39.0877f, 36.2034f ;
:corner_lons = -108.7585f, -108.8684f, -104.0023f, -104.0788f, -108.7641f, -108.8742f, -103.9965f, -104.0732f, -108.7584f, -108.8686f, -104.0021f, -104.0789f, -108.7639f, -108.8744f, -103.9963f, -104.0733f ;

The above “corner_lats” and “corner_lons” are the correct Lon/Lat for the corners and center for each of the four corner pixels in the domain 2. The lower left corner of the lower left pixel is at -108.7585,36.19099. You enter these into the AAI header rows, and set the cellsize to 1000, then import that altered AAIGrid into GRASS in a WGS84 location .

One caveat: this method assumes that the whole geo_em.d03 is projected in a Lon/Lat WGS84 coordinate system, which is wrong, as you know. The way I know of to get an accurate raster is to parse out the XLAT and XLONG variables, together with the GREENFRAC variable you are interested in as a CSV list of lon,lat,greenfac values, import the values as vector points, and do an interpolation.

Another possibility comes to mind, but I’ve never tried it: You should be able to convert the above xllcorner and yllcorner values to your Lambert Conformal Conic projection in advance, and use the converted values in the AAI header rows. Then import the AAI file directly into the correct LCC location. You would do this by:

  1. create a point vector of a single point (v.in.ascii) from the Long/Lat values above in an WGS84 location
  2. project that point to your LCC location (v.proj)
  3. get the LCC X-Y coordinates (v.out.ascii)
  4. Edit the header rows with the LCC xllcorner and yllcorner values
  5. Import into an LCC location (with r.in.gdal -o as you did)

I’d be interested to know if this works :slight_smile:

Regards,
Micha

On 06/09/2015 05:47 AM, Thomas Adams wrote:

All:

Using the help found here: http://grasswiki.osgeo.org/wiki/NetCDF I have been able to import a netcdf file into GRASS 7, but the data is not scaled correctly.

My file is found here:
https://drive.google.com/file/d/0B-Dmtz0DFEiubG5aYTk4TDZERXc/view?usp=sharing

I created a GRASS LOCATION based on this:

char Lambert;
  :grid_mapping_name = "lambert_conformal_conic";
  :latitude_of_projection_origin = 37.99999237060547; // double
  :longitude_of_central_meridian = -106.0; // double
  :standard_parallel = 34.0, 42.0; // double
  :earth_radius = 6371229.0; // double
  :_CoordinateTransformType = "Projection";
  :_CoordinateAxisTypes = "GeoX GeoY";

which I obtained using panoply from here: [http://www.giss.nasa.gov/tools/panoply/](http://www.giss.nasa.gov/tools/panoply/)
Looking at the file netcdf header information, it's clear the dx=dy=1000m; but, what imports is x=dy=1m resolution even though I set the x & y resolution to 1000m
This is my import command:
r.in.gdal input=NETCDF:"[geo_em.d02.nc](http://geo_em.d02.nc)":GREENFRAC output=greenfrac -o -e --overwrite

which produces..
Warning 1: No UNIDATA NC_GLOBAL:Conventions attribute
WARNING: Over-riding projection check
 100%
Raster map <greenfrac.1> created.
 100%
Raster map <greenfrac.2> created.
 100%
Raster map <greenfrac.3> created.
 100%
Raster map <greenfrac.4> created.
 100%
Raster map <greenfrac.5> created.
 100%
Raster map <greenfrac.6> created.
 100%
Raster map <greenfrac.7> created.
 100%
Raster map <greenfrac.8> created.
 100%
Raster map <greenfrac.9> created.
 100%
Raster map <greenfrac.10> created.
 100%
Raster map <greenfrac.11> created.
 100%
Raster map <greenfrac.12> created.
Region for the current mapset updated
r.in.gdal complete.

I am virtually certain the lower left corner is placed correctly, because I can zoom to the lower-left location and it falls about where it should – however the spatial extent of the data is not what it should be, but when zoomed-in looks correct.

Any thoughts?

Thank you,
Tom

This mail was received via Mail-SeCure System.

_______________________________________________
grass-user mailing list
[grass-user@lists.osgeo.org](mailto:grass-user@lists.osgeo.org)
[http://lists.osgeo.org/mailman/listinfo/grass-user](http://lists.osgeo.org/mailman/listinfo/grass-user)
This mail was received via Mail-SeCure System.

Thomas E Adams, III
718 McBurney Drive
Lebanon, OH 45036

1 (513) 739-9512 (cell)