Opened 7 years ago
Closed 7 years ago
#3414 closed defect (fixed)
v.import: ERROR: Invalid e-w resolution field: 0 when importing EU_LAEA point into Sinusoidal location
| Reported by: | neteler | Owned by: | |
|---|---|---|---|
| Priority: | normal | Milestone: | 7.4.0 |
| Component: | Projections/Datums | Version: | svn-releasebranch74 |
| Keywords: | v.import | Cc: | |
| CPU: | Unspecified | Platform: | Linux |
Description
When trying to import a vector point in EU LAEA projection (see attached SHP file) into a location in MODIS Sinusoidal projection, an error occurs while v.proj does the job correctly:
GRASS 7.3.svn (modis_sinusoidal):~ > v.import in=point_eu_laea.shp out=point_eu_laea2 ERROR: Invalid e-w resolution field: 0
With OGR debugging on:
GRASS 7.3.svn (modis_sinusoidal):~ > PROJ_DEBUG=ON CPL_DEBUG=ON v.import in=point_eu_laea.shp out=point_eu_laea2 GDAL: In GDALDestroy - unloading GDAL shared library. Shape: DBF Codepage = LDID/87 for point_eu_laea.shp Shape: Treating as encoding 'ISO-8859-1'. GDAL: GDALOpen(point_eu_laea.shp, this=0x255f710) succeeds as ESRI Shapefile. GDAL: In GDALDestroy - unloading GDAL shared library. GDAL: force close of point_eu_laea.shp (0x255f710) in GDALDriverManager cleanup. GDAL: GDALClose(point_eu_laea.shp, this=0x255f710) Shape: DBF Codepage = LDID/87 for point_eu_laea.shp Shape: Treating as encoding 'ISO-8859-1'. GDAL: GDALOpen(point_eu_laea.shp, this=0x9d57e0) succeeds as ESRI Shapefile. ERROR: Invalid e-w resolution field: 0 GDAL: In GDALDestroy - unloading GDAL shared library. GDAL: force close of point_eu_laea.shp (0x9d57e0) in GDALDriverManager cleanup. GDAL: GDALClose(point_eu_laea.shp, this=0x9d57e0) ERROR: Unable to create location from OGR datasource <point_eu_laea.shp>
The projection of the location:
GRASS 7.3.svn (modis_sinusoidal):~ > g.proj -p -PROJ_INFO------------------------------------------------- name : unnamed a : 6371007.181 es : 0 proj : sinu lon_0 : 0 x_0 : 0 y_0 : 0 no_defs : defined -PROJ_UNITS------------------------------------------------ unit : Meter units : Meters meters : 1
I am a bit clueless here.
Attachments (4)
Change History (24)
by , 7 years ago
| Attachment: | point_eu_laea.zip added |
|---|
follow-up: 3 comment:2 by , 7 years ago
Interestingly:
g.proj georef=point_eu_laea.shp location=test2 Essai d'ouverture à l'aide d'OGR... ...réussi. Location <test2> created You can switch to the new location by `g.mapset mapset=PERMANENT location=test2`
but
v.in.ogr point_eu_laea.shp location=test
ATTENTION: All available OGR layers will be imported into vector map
<point_eu_laea>
Location <test> created
ERREUR :Invalid e-w resolution field: 0
Both are created, but:
diff -u test/PERMANENT/ test2/PERMANENT/
diff -u test/PERMANENT/.bashrc test2/PERMANENT/.bashrc
--- test/PERMANENT/.bashrc 2017-09-06 15:54:44.389963786 +0200
+++ test2/PERMANENT/.bashrc 2017-09-06 15:54:55.545841900 +0200
@@ -1,5 +1,5 @@
test -r ~/.alias && . ~/.alias
-PS1='GRASS 7.3.svn (test):\w > '
+PS1='GRASS 7.3.svn (test2):\w > '
grass_prompt() {
LOCATION="`g.gisenv get=GISDBASE,LOCATION_NAME,MAPSET separator='/'`"
if test -d "$LOCATION/grid3/G3D_MASK" && test -f "$LOCATION/cell/MASK" ; then
diff -u test/PERMANENT/DEFAULT_WIND test2/PERMANENT/DEFAULT_WIND
--- test/PERMANENT/DEFAULT_WIND 2017-09-06 15:51:39.891978100 +0200
+++ test2/PERMANENT/DEFAULT_WIND 2017-09-06 15:52:13.227614355 +0200
@@ -10,9 +10,9 @@
n-s resol: 0
top: 1.000000000000000
bottom: 0.000000000000000
-cols3: 1
-rows3: 1
+cols3: 20
+rows3: 20
depths: 1
-e-w resol3: 1
-n-s resol3: 1
+e-w resol3: 0
+n-s resol3: 0
t-b resol: 1
Common subdirectories: test/PERMANENT/sqlite and test2/PERMANENT/sqlite
Only in test/PERMANENT/: .tmp
diff -u test/PERMANENT/WIND test2/PERMANENT/WIND
--- test/PERMANENT/WIND 2017-09-06 15:51:39.891978100 +0200
+++ test2/PERMANENT/WIND 2017-09-06 15:52:13.227614355 +0200
@@ -10,9 +10,9 @@
n-s resol: 0
top: 1.000000000000000
bottom: 0.000000000000000
-cols3: 1
-rows3: 1
+cols3: 20
+rows3: 20
depths: 1
-e-w resol3: 1
-n-s resol3: 1
+e-w resol3: 0
+n-s resol3: 0
t-b resol: 1
and in both I get:
g.proj -p ERREUR :Invalid e-w resolution field: 0
because DEFAULT_WIND/WIND contains:
proj: 99 zone: 0 north: 1205000 south: 1205000 east: 3674000 west: 3674000 cols: 20 rows: 20 e-w resol: 0 n-s resol: 0 top: 1.000000000000000 bottom: 0.000000000000000 cols3: 1 rows3: 1 depths: 1 e-w resol3: 1 n-s resol3: 1 t-b resol: 1
which seems logical as a point does not have dimensions, but the default region should still be defined with one cell and a resolution of 1....
So, yes, there seems to be a link to #2878...
follow-up: 4 comment:3 by , 7 years ago
Replying to mlennert:
[...] I get:
g.proj -p ERREUR :Invalid e-w resolution field: 0because DEFAULT_WIND/WIND contains:
proj: 99 zone: 0 north: 1205000 south: 1205000 east: 3674000 west: 3674000 [...]which seems logical as a point does not have dimensions, but the default region should still be defined with one cell and a resolution of 1....
A simple solution would be to add/subtract 0.5 if north == south and/or east == west with rows and/or cols set to 1. For latlong, this is only guaranteed to work in trunk.
So, yes, there seems to be a link to #2878...
It seems so, because g.proj uses code copied from v.in.ogr.
follow-up: 5 comment:4 by , 7 years ago
Replying to mmetz:
Replying to mlennert:
[...] I get:
g.proj -p ERREUR :Invalid e-w resolution field: 0because DEFAULT_WIND/WIND contains:
proj: 99 zone: 0 north: 1205000 south: 1205000 east: 3674000 west: 3674000 [...]which seems logical as a point does not have dimensions, but the default region should still be defined with one cell and a resolution of 1....
A simple solution would be to add/subtract 0.5 if north == south and/or east == west with rows and/or cols set to 1. For latlong, this is only guaranteed to work in trunk.
I think that the default region created when using a vector as input for georeference should always be the default n=1 s=0 e=1 w=0 res=1 region (i.e. the same that is created when creating a location without determining a default region). One could argue that the user would like to have a default region with the extents of the vector used for creating the location, but as there is no way to determine a sensible resolution, the risk is high that you get an arbitrarily large region with a resolution of 1, i.e. a very large number of pixels. Or we would have to decide on a reasonable default number of pixels and adjust the region accordingly. But I think that having an obviously non-sense region is better as it forces the user to explicitly set the region.
follow-up: 6 comment:5 by , 7 years ago
Replying to mlennert:
Replying to mmetz:
Replying to mlennert:
[...] I get:
g.proj -p ERREUR :Invalid e-w resolution field: 0because DEFAULT_WIND/WIND contains:
proj: 99 zone: 0 north: 1205000 south: 1205000 east: 3674000 west: 3674000 [...]which seems logical as a point does not have dimensions, but the default region should still be defined with one cell and a resolution of 1....
A simple solution would be to add/subtract 0.5 if north == south and/or east == west with rows and/or cols set to 1. For latlong, this is only guaranteed to work in trunk.
I think that the default region created when using a vector as input for georeference should always be the default n=1 s=0 e=1 w=0 res=1 region (i.e. the same that is created when creating a location without determining a default region). One could argue that the user would like to have a default region with the extents of the vector used for creating the location, but as there is no way to determine a sensible resolution, the risk is high that you get an arbitrarily large region with a resolution of 1, i.e. a very large number of pixels.
I aggree.
Or we would have to decide on a reasonable default number of pixels and adjust the region accordingly.
I don't see a way to determine a reasonable resolution/number of pixels from vector extents.
But I think that having an obviously non-sense region is better as it forces the user to explicitly set the region.
Usually, users need to set (or at least review) the region anyway.
That means that duplicate code in g.proj/v.in.ogr where new region settings are estimated needs to be removed.
comment:6 by , 7 years ago
Replying to mmetz:
Replying to mlennert:
Replying to mmetz:
Replying to mlennert:
[...] I get:
g.proj -p ERREUR :Invalid e-w resolution field: 0because DEFAULT_WIND/WIND contains:
proj: 99 zone: 0 north: 1205000 south: 1205000 east: 3674000 west: 3674000 [...]which seems logical as a point does not have dimensions, but the default region should still be defined with one cell and a resolution of 1....
A simple solution would be to add/subtract 0.5 if north == south and/or east == west with rows and/or cols set to 1. For latlong, this is only guaranteed to work in trunk.
I think that the default region created when using a vector as input for georeference should always be the default n=1 s=0 e=1 w=0 res=1 region (i.e. the same that is created when creating a location without determining a default region). One could argue that the user would like to have a default region with the extents of the vector used for creating the location, but as there is no way to determine a sensible resolution, the risk is high that you get an arbitrarily large region with a resolution of 1, i.e. a very large number of pixels.
I aggree.
Or we would have to decide on a reasonable default number of pixels and adjust the region accordingly.
I don't see a way to determine a reasonable resolution/number of pixels from vector extents.
But I think that having an obviously non-sense region is better as it forces the user to explicitly set the region.
Usually, users need to set (or at least review) the region anyway.
That means that duplicate code in g.proj/v.in.ogr where new region settings are estimated needs to be removed.
+1
follow-up: 8 comment:7 by , 7 years ago
I have attached the 0.5 suggestion being tested as dirty hack.
Seems to do the job while a correct solution in v.in.ogr and g.proj would obviously be preferred... (not sure how to implement it).
comment:8 by , 7 years ago
Replying to neteler:
I have attached the 0.5 suggestion being tested as dirty hack.
Seems to do the job while a correct solution in v.in.ogr and g.proj would obviously be preferred... (not sure how to implement it).
I don't know if you might have to make sure that cellhd.ns_res is big enough to avoid resolutions which might cause precision issues ? Also, in #2878 I mention issues with vertical resolution.
Personally, I think it better to go with the idea that when a location is defined based on a vector, no attempt is made to define the region based on the vector.
by , 7 years ago
| Attachment: | no_vector_default_region.diff added |
|---|
diff for not basing region on vector extent
follow-ups: 10 11 comment:9 by , 7 years ago
I just attached a diff which I think should work for not attempting to create a region based on vector extent. The region will simply be n=1 s=0 e=1 w=0 res=1.
I don't have time for much testing, so if others could, this would be good. Not sure though that this is secure enough to introduce now into the 7.2.2 release. Probably should be noted under known bugs and backported to release72 after the 7.2.2 release.
follow-up: 12 comment:10 by , 7 years ago
Replying to mlennert:
I just attached a diff which I think should work for not attempting to create a region based on vector extent. The region will simply be n=1 s=0 e=1 w=0 res=1.
For v.in.ogr, some more changes are needed, see attached patch no_vector_default_region_2.diff (the patch also solves some TODOs).
I don't have time for much testing, so if others could, this would be good. Not sure though that this is secure enough to introduce now into the 7.2.2 release. Probably should be noted under known bugs and backported to release72 after the 7.2.2 release.
+1
comment:11 by , 7 years ago
Replying to mlennert:
I don't have time for much testing, so if others could, this would be good. Not sure though that this is secure enough to introduce now into the 7.2.2 release. Probably should be noted under known bugs and backported to release72 after the 7.2.2 release.
... this means that we will have to release 7.2.3 for sure because that's a relevant bug.
follow-up: 13 comment:12 by , 7 years ago
Replying to mmetz:
For v.in.ogr, some more changes are needed, see attached patch no_vector_default_region_2.diff (the patch also solves some TODOs).
I get a compilation error:
[mneteler@oboe v.in.ogr ]$ make
gcc -O2 -march=native -std=gnu99 -fexceptions -fstack-protector -m64 -fdiagnostics-color -I/home/mneteler/software/grass73/dist.x86_64-pc-linux-gnu/include -I/home/mneteler/software/grass73/dist.x86_64-pc-linux-gnu/include -I/usr/include/gdal -I/usr/include -DPACKAGE=\""grassmods"\" -I/usr/include/gdal -I/home/mneteler/software/grass73/dist.x86_64-pc-linux-gnu/include -I/home/mneteler/software/grass73/dist.x86_64-pc-linux-gnu/include -DRELDIR=\"vector/v.in.ogr\" -o OBJ.x86_64-pc-linux-gnu/main.o -c main.c
main.c: In function ‘main’:
main.c:886:35: error: ‘have_ogr_extents’ undeclared (first use in this function); did you mean ‘have_ogr_extent’?
if (flag.no_clean->answer || !have_ogr_extents) {
^~~~~~~~~~~~~~~~
have_ogr_extent
main.c:886:35: note: each undeclared identifier is reported only once for each function it appears in
make: *** [../../include/Make/Compile.make:32: OBJ.x86_64-pc-linux-gnu/main.o] Error 1
by , 7 years ago
| Attachment: | no_vector_default_region_2.diff added |
|---|
fixing another diff for not basing region on vector extent
comment:13 by , 7 years ago
Replying to neteler:
Replying to mmetz:
For v.in.ogr, some more changes are needed, see attached patch no_vector_default_region_2.diff (the patch also solves some TODOs).
I get a compilation error:
Oops, fixed in the updated attached patch no_vector_default_region_2.diff.
About backporting: v.in.ogr in trunk is already quite different from relbr72, and with this patch, it would be more different. IMHO, the fix for OSM import and this fix are both bug fixes, i.e. once confirmed working, patched v.in.ogr in trunk should be backported as is to relbr72.
comment:14 by , 7 years ago
Thanks for the new updated no_vector_default_region_2.diff, tested here:
GRASS 7.3.svn (modis_sinusoidal):~ > v.import --v in=~/point_eu_laea.shp out=point_eu_laea2Creating temporary location for </home/mundialis/point_eu_laea.shp>... Using OGR driver 'ESRI Shapefile' Projection of dataset does not appear to match current location. GRASS LOCATION PROJ_INFO is: name: unnamed a: 6371007.181 es: 0 proj: sinu lon_0: 0 x_0: 0 y_0: 0 no_defs: defined Import dataset PROJ_INFO is: name: ETRS89_ETRS_LAEA datum: etrs89 ellps: grs80 proj: laea lat_0: 52 lon_0: 10 x_0: 4321000 y_0: 3210000 no_defs: defined In case of no significant differences in the projection definitions, use the -o flag to ignore them and use current location definition. Consider generating a new location with 'location' parameter from input data set. switch to temp location -PROJ_INFO------------------------------------------------- name : ETRS89_ETRS_LAEA datum : etrs89 ellps : grs80 proj : laea lat_0 : 52 lon_0 : 10 x_0 : 4321000 y_0 : 3210000 no_defs : defined -PROJ_UNITS------------------------------------------------ unit : Meter units : Meters meters : 1 Importing </home/mundialis/point_eu_laea.shp> ... Using OGR driver 'ESRI Shapefile' Projection of input dataset and current location appear to match Check if OGR layer <point_eu_laea> contains polygons... 100% Importing 1 features Using native format Default driver / database set to: driver: sqlite database: $GISDBASE/$LOCATION_NAME/$MAPSET/sqlite/sqlite.db Column name <cat> renamed to <cat_> Importing 1 features (OGR layer <point_eu_laea>)... 100% ----------------------------------------------------- Building topology for vector map <point_eu_laea2@PERMANENT>... Registering primitives... One primitive registered One vertex registered Building areas... 100% 0 areas built 0 isles built Attaching islands... Attaching centroids... 100% Topology was built Number of nodes: 0 Number of primitives: 1 Number of points: 1 Number of lines: 0 Number of boundaries: 0 Number of centroids: 0 Number of areas: 0 Number of isles: 0 Reprojecting <point_eu_laea2>... Input Projection Parameters: +proj=laea +lat_0=52 +lon_0=10 +x_0=4321000 +y_0=3210000 +no_defs +a=6378137 +rf=298.257222101 +towgs84=0.000,0.000,0.000 Input Unit Factor: 1 Output Projection Parameters: +proj=sinu +lon_0=0 +x_0=0 +y_0=0 +no_defs +a=6371007.181 +b=6371007.181 Output Unit Factor: 1 Using native format Reprojecting primitives ... Building topology for vector map <point_eu_laea2@Germany_LST_1km_daily_gdd_2016>... Registering primitives... One primitive registered One vertex registered Building areas... 100% 0 areas built 0 isles built Attaching islands... Attaching centroids... 100% Topology was built Number of nodes: 0 Number of primitives: 1 Number of points: 1 Number of lines: 0 Number of boundaries: 0 Number of centroids: 0 Number of areas: 0 Number of isles: 0
Now the file is smoothly imported :-)
follow-up: 17 comment:16 by , 7 years ago
Replying to neteler:
Markus M, do you plan to submit the patch?
Yes, I plan to submit the patch, but the patch affects the options spatial=, where= and the -r and -e flags which I want to test a bit more before submitting.
comment:17 by , 7 years ago
comment:20 by , 7 years ago
| Milestone: | 7.2.3 → 7.4.0 |
|---|---|
| Resolution: | → fixed |
| Status: | new → closed |
| Version: | svn-trunk → svn-releasebranch74 |
Replying to martinl:
What is status of this issue? Backport needed?
Since it was fixed in trunk 6 months ago it should be part of 7.4.0, that's fine. Closing.
Test point in EU LAEA (EPSG3035) projection - SHAPE file