Create location fails using Advanced Methods in Grass 7.8.3

[bal-agates]

Creating a new location from the wxPython GRASS 7.8.3 Startup fails when using "Advanced methods" and specifying "projection code = ll" (lat/lon) and "datum code = nad83". I successfully created a location by specifying a CRS code rather than using "Advanced methods".

Popup Error Dialog:

Error ERROR 1: PROJ: proj_create: Error -13 (major axis or radius = 0 or not given) ERROR: Can't parse PROJ.4-style parameter string

Terminal message:

$ grass78 Starting GRASS GIS... Traceback (most recent call last): File "/opt/local/share/grass78/gui/wxpython/location_wizard/wizard.py", line 1079, in OnPageChanging dlg = SelectTransformDialog( File "/opt/local/share/grass78/gui/wxpython/location_wizard/dialogs.py", line 617, in __init__ transforms = '---\n\n0\nDo not apply any datum transformations\n\n' + transforms TypeError: can only concatenate str (not "bytes") to str

Grass 7.8.3 built using MacPorts [grass7 @7.8.3_1+postgresql12+proj6+python38] macOS 10.13.6

g.version -rge version=7.8.3 date=2020 revision=exported build_date=2020-06-21 build_platform=x86_64-apple-darwin17.7.0 build_off_t_size=8 libgis_revision=2020-06-22T00:02:36+00:00 libgis_date=2020-06-22T00:02:36+00:00 proj=6.3.2 gdal=3.1.0 geos=3.8.1 sqlite=3.32.3

On the surface this would appear to be a Python3 string vs bytes issue. The concatenation string seems a bit strange to me. The root cause of the failure is likely something else but it seems fragile to append what looks like a comment to the transform string. This would depend on PROJ4 ignoring everything after the first newline. The solution to this problem might be as simple as type casting the string as bytes.

[bal-agates]

I added a debug print on my system before dialogs.py:617 and transforms = b'' confirming one problem of trying to append string+bytes. Not sure what transforms should be at this point. If it should not be null then there are more problems than just with concat. I didn't see any recent changes in dialogs.py that seems related.

[petrasovaa]

This bytes/str issue was addressed by #702. The rest is still valid, this is debug output I get:

GRASS 7.9.dev (nc_spm_08_grass7):~/dev/grass/grass/gui/wxpython > GUI D1/4: gcmd.RunCommand(): g.proj -t proj4=+proj=longlat datum=nad83 datum_trans=-1
D1/4: grass.script.core.start_command(): g.proj -t proj4=+proj=longlat datum=nad83 datum_trans=-1
GUI D1/4: gcmd.RunCommand(): get return code 1 (0.070078 sec)
GUI D2/4: gcmd.RunCommand(): error D1/4: G_set_program_name(): g.proj
ERROR 1: PROJ: proj_create: Error -13: major axis or radius = 0 or not given

GRASS_INFO_ERROR(21047,1): Can't parse PROJ.4-style parameter string
GRASS_INFO_END(21047,1)

GUI D3/4: gcmd.RunCommand(): return stdout = None
GUI D1/4: gcmd.RunCommand(): g.proj -jf proj4=+proj=longlat  +lat_0=0.0 +lon_0=0.0 +ellps=grs80  +a=6378137.000  +rf=298.257222101 +no_defs location=location datum=nad83
D1/4: grass.script.core.start_command(): g.proj -jf proj4=+proj=longlat  +lat_0=0.0 +lon_0=0.0 +ellps=grs80  +a=6378137.000  +rf=298.257222101 +no_defs location=location datum=nad83
GUI D1/4: gcmd.RunCommand(): get return code 1 (0.036904 sec)
GUI D2/4: gcmd.RunCommand(): error D1/4: G_set_program_name(): g.proj
ERROR 1: PROJ: proj_create: Error -9: unknown elliptical parameter name

GRASS_INFO_ERROR(21048,1): Can't parse PROJ.4-style parameter string
GRASS_INFO_END(21048,1)

GUI D3/4: gcmd.RunCommand(): return stdout = None

[bal-agates]

More debugging info. Grass 7.8.3 new location sequence with advanced method.

1) Advanced methods - select coordinate system from parameter list 2) Projection code. Initial "Next" gray. Type "ll". Now "Next" enabled. 3) Datum with associated ellipsoid. 4) Datum code. Initial "Next" gray. Type "nad83". Now "Next" enabled. 5) Error Dialog displayed.

With debugging info turned on the first error I see is:

GUI D1/1: gcmd.RunCommand(): g.proj -t proj4=+proj=longlat datum=nad83 datum_trans=-1
GUI D1/1: gcmd.RunCommand(): get return code 1 (0.074414 sec)

Comparing this to Grass 7.6.0 the sequence 1-4 is the same but is followed by a dialog to specify the NAD83 region (1=whole nad83 region, 2=used in Florida,...)

In Grass 7.6.0 the debug output for the equivalent command is:

GUI D1/1: gcmd.RunCommand(): g.proj -t datum_trans=-1 proj4=+proj=longlat datum=nad83
GUI D1/1: gcmd.RunCommand(): get return code 0 (0.632159 sec)

The order of parameters is different but I do not think that should matter.

When I manually execute the above command in Grass 7.6.0 I see

GRASS 7.6.0 (co_dem):~ > g.proj -t datum_trans=-1 proj4=+proj=longlat datum=nad83
WARNING: Datum <unknown> not recognised by GRASS and no parameters found
Datum set to <nad83>
Ellipsoid set to <grs80>
---
1
Used in whole nad83 region
towgs84=0.000,0.000,0.000
Default 3-Parameter Transformation (May not be optimum for older datums; use this only if no more appropriate options are available.)
...

With Grass 7.8.3 started when I execute the above command I see

GRASS 7.8.3 (co_dem):~ > g.proj -t datum_trans=-1 proj4=+proj=longlat datum=nad83
ERROR 1: PROJ: proj_create: Error -13 (major axis or radius = 0 or not given)
ERROR: Can't parse PROJ.4-style parameter string
GRASS 7.8.3 (co_dem):~ >

It will take more digging to determine the root cause of failure. Need to determine where the region list was coming from. Possible earlier PROJ4 returned a region list that Grass knew how to handle and it no longer does? Or there was more functionality in g.proj?

I looked at the change suggested by petrasovaa for the string/bytes problem but for this debugging I did not apply any changes. I think that comes into effect after this first error.

[bal-agates]

I have a working theory as to the cause of this problem but do not have the build environment to prove my theory. I suspect either 1) something has changed in the GDAL library related to PROJ4 strings or 2) the PROJ4 string passed to GDAL has changed. I suspect the failure is at:

./general/g.proj/input.c:165 if (OSRImportFromProj4(hSRS, proj4string) != OGRERR_NONE)

I suspect proj4string = "+proj=longlat +no_defs". With the enclosed program I have verified this string causes the observed error using GDAL 3.1.0. The enclosed program also verifies that proj4string = "+proj=longlat +datum=NAD83 +no_defs" works. For this case in ./general/g.proj/main.c the datum should available in datum->answer. The function input_proj4 is declared

int input_proj4(char *proj4params)

and from "proj4params" one might guess the intention was to have multiple proj4 parameters (i.e. +proj and +datum) but I have found no mechanism for inproj4->answer to have "+datum=NAD83" concatenated. I suspect inproj4->answer is "+proj=longlat". I have found no recent changes in Grass that could account for a change in this area. The parser (./lib/gis/parser.c) is complex and I easily could have missed something. I think the calling chain is:

./general/g.proj/main.c:266 input_proj4(inproj4->answer);
./general/g.proj/input.c:165 if (OSRImportFromProj4(hSRS, proj4string) != OGRERR_NONE) 

I ran into problems (there were a myriad of library dependencies) trying to install an earlier version of GDAL so I could use my test program to see if anything changed in GDAL. If something has changed in GDAL a possible solution is to change input_proj4 passing 2 parameters, something like:

int input_proj4(char *proj, char *datum)

and then when building proj4string concatenate these values. I could provide proposed code changes if someone can verify my assumptions and rule out other causes.

#include <stdlib.h>
#include <stdio.h>
#include <ogr_api.h>
#include <ogr_srs_api.h>

int main(int argc, char** argv) 
{
    char* proj4string1 = "+proj=longlat +no_defs";
    char* proj4string2 = "+proj=longlat +datum=NAD83 +no_defs";
    char* proj4string3 = "+proj=longlat +datum_NAD83";
    char* proj4string4 = "+no_defs";
    char* proj4strings[] = { proj4string1, proj4string2, proj4string3, proj4string4, 0 };
    OGRSpatialReferenceH hSRS;
    int i;

        printf("%s was compiled against GDAL %s and is running against GDAL %s\n",
           argv[0], GDAL_RELEASE_NAME, GDALVersionInfo("RELEASE_NAME"));

    for (i = 0; proj4strings[i] != 0; i++) {
        printf("Testing string: '%s'\n", proj4strings[i]);
        hSRS = OSRNewSpatialReference(NULL);
        if (OSRImportFromProj4(hSRS, proj4strings[i]) != OGRERR_NONE) {
            printf("  ERROR: Can't parse PROJ.4-style parameter string\n");
        } else {
            printf("  SUCCESS\n");
        }
        OSRDestroySpatialReference(hSRS);
    }
}

/* Captured output macOS 10.13.6 on 2020/7/9
$ gdal_test1
gdal_test1 was compiled against GDAL 3.1.0 and is running against GDAL 3.1.0
Testing string: '+proj=longlat +no_defs'
ERROR 1: PROJ: proj_create: Error -13 (major axis or radius = 0 or not given)
  ERROR: Can't parse PROJ.4-style parameter string
Testing string: '+proj=longlat +datum=NAD83 +no_defs'
  SUCCESS
Testing string: '+proj=longlat +datum_NAD83'
  SUCCESS
Testing string: '+no_defs'
ERROR 1: PROJ: proj_create: unrecognized format / unknown name
  ERROR: Can't parse PROJ.4-style parameter string
*/

[metzm]

The reason is that the behaviour of PROJ has changed drastically with PROJ6+. The problem is that GRASS is not yet fully adapted to the new behaviour. Previously +proj=longlat +no_defs worked because the default datum was WGS84 and a proj string was supposed to define a crs. Now, there is no default datum, and a proj string can mean different things, usually a transformation pipeline and not a crs. Now, a datum transformation parameter no longer makes sense because this datum transformation parameter describes transformation from some other datum to WGS84 datum. This datum transformation parameter is only relevant for reprojection where now PROJ has a list of all sorts of transformations from one datum to another and selects an appropriate transformation automatically.

That means that the location wizard needs to be updated for PROJ6+ and g.proj needs to be updated to provide a valid proj4-style string, i.e. not +proj=longlat +no_defs but +proj=longlat datum=wg84 +no_defs +type=crs or another datum provided with the datum option for PROJ6+.

[mrael2]

I believe the error I'm seeing is related to what's already been reported, but here are additional details about how I duplicated the problem:

macOS 11.0.1 (Big Sur, Intel) GRASS GIS 7.8.4 (stable) as downloaded from http://grassmac.wikidot.com/downloads

1. Launch GRASS GIS 7.8.4 using its icon.
2. In the startup window, enter a valid database directory.
3. Click "New" to create new location.
4. Enter appropriate database, project, and location title.
5. Leave both "Set default region..." and "Create user mapset" boxes unchecked.
6. Click "Next."
7. Under "Advanced," click the button labeled "Select coordinate system parameters from a list."
8. Click "Next."
9. Select "utm - Universal Traverse Mercator (UTM)" from the list.
10. Click "Next."
11. Select the button labeled "Datum with associated ellipsoid."
12. Enter "13" for the projection zone.
13. For Southern Hemisphere, select "No".
14. Click "Next."
15. Select datum of "nad83" from the list.
16. Click "Next."
17. Error appears as follows:

ERROR 1: PROJ: proj_create: Error -9: unknown elliptical parameter name ERROR: Can't parse PROJ.4-style parameter string.

It sounds like changes are needed in both the GUI and g.proj to correct this issue. Are there any workarounds in the meantime?

[bal-agates]

The only short term workaround I know of is to use an EPSG code rather than using advanced settings. For your example

NAD83 UTM ZONE 13N = EPSG:26913

I used QGIS to find the EPSG code. I am not sure of the easiest way to find an EPSG code.

[cmbarton]

This is a real problem. It affects all current versions of GRASS AFAICT, from 7.8.x through 7.9

You cannot create a location from known projection information unless you parse it into a PROJ format or search for it online at EPSG.org. And neither of these is very satisfactory.

On top of this the search tool doesn't work. Put something into it that I know is in the list and nothing happens.

[cmbarton]

One possibility is that it is somehow looking for PROJ4 style info and we have migrated to PROJ5 and 6. But this is a shot in the dark.

[cmbarton]

A minor update. The search works in GRASS 7.9 for the EPSG code list but not for projection and datum. In GRASS 7.8, the search does not work in any place.

[nilason]

On top of this the search tool doesn't work. Put something into it that I know is in the list and nothing happens.

Just to illustrate. The search field is not functional in 7.8.5:

whereas it works nicely in 7.9:

[petrasovaa]

The search looks like a different issue, please create a new one to keep things organized.

[nilason]

The search looks like a different issue, please create a new one to keep things organized.

You're absolutely right. Done with #1319.

[metzm]

Regarding the original issue, the location wizard produces something like +proj=utm +zone=13 +ellps=grs80 +a=6378137.0 +rf=298.257222101 +no_defs The unknown elliptical parameter name is grs80, should be upper-case GRS80. The names recognized by PROJ are a mix of lower-case and upper-case, e.g. "GRS80" and "clrk80", and names provided to PROJ are case-sensitive. If using a proj definition, it is safer to provide only +a= +rf= without +ellps=. Needs to be fixed in the location wizard.

[cmbarton]

I've tested with utm and UTM, wgs84 and WGS84. And use the case for latlon as it is pulled from the PROJ-generated list (lower case 'll'). All generate the same error. So while case may be important, it doesn't seem to be what is causing this error.

Also, the names in the lists for projection and datum appear as they are generated by PROJ.

[petrasovaa]

Regarding the original issue, the location wizard produces something like +proj=utm +zone=13 +ellps=grs80 +a=6378137.0 +rf=298.257222101 +no_defs The unknown elliptical parameter name is grs80, should be upper-case GRS80. The names recognized by PROJ are a mix of lower-case and upper-case, e.g. "GRS80" and "clrk80", and names provided to PROJ are case-sensitive. If using a proj definition, it is safer to provide only +a= +rf= without +ellps=. Needs to be fixed in the location wizard.

Is this related to any of the new changes you did for new PROJ? Are any changes in location wizard needed (for master)? For example, there is still this transformation dialog, which seems to be obsolete. Perhaps you could summarize any necessary changes in a new issue. I would be happy to help with the GUI part, but I don't understand what needs to be done (for 78 or master).

[nilason]

@bal-agates

I have a working theory as to the cause of this problem but do not have the build environment to prove my theory.

I noticed you use MacPorts. If this is still the case you have everything you need for building from source (outside of MP).

I attached a file with configuration command I use: configure_grass_macports.txt that may need modification according to your setup.

There is a conflict between system and MP iconv, but if you're interested I'll help with that.

[nilason]

I suppose this is the issue @cmbarton references to in ML.

[cmbarton]

Yes. This is it.

[neteler]

Is this related to any of the new changes you did for new PROJ? Are any changes in location wizard needed (for master)? For example, there is still this transformation dialog, which seems to be obsolete. Perhaps you could summarize any necessary changes in a new issue. I would be happy to help with the GUI part, but I don't understand what needs to be done (for 78 or master).

I just tested with G79 (Fedora 33) and the same issue happens:

• Define custom CRS
• Projection code: ll
• [x] Datum with associated ellipsoid
• Datum code: nad83
• see error:

[cmbarton]

Would it be better to drop back to PROJ 5 for this upcoming release?

[neteler]

Would it be better to drop back to PROJ 5 for this upcoming release?

This is not possible - most distros have more recent versions.

I understood that "only" some uppercase vs lowercase is the issue? See https://github.com/OSGeo/grass/issues/758#issuecomment-776871542

[cmbarton]

Would it be better to drop back to PROJ 5 for this upcoming release?

This is not possible - most distros have more recent versions.

I understood that "only" some uppercase vs lowercase is the issue? See #758 (comment)

Too bad. M. Metz suggests other issues. See addtional comments by M. Metz

[nilason]

Projections is not my field of expertise, but seems @metzm is right with the issue of capital names (eg. wgs84 vs. WGS84), but they are defined in ellipse.table and datum.table (not in GUI). I'm still exploring this, but thought I'd mention this as a heads up or to pass the baton.

Eg. changing the files as follows, seem to solve the issue with datum wgs84 and nad83:

ellipse.table:

# from:
wgs84       "World Geodetic System 1984"    a=6378137.000   f=1/298.257223563
grs80       "Geodetic Reference System 1980" a=6378137.000  f=1/298.257222101

# to:
WGS84       "World Geodetic System 1984"    a=6378137.000   f=1/298.257223563
GRS80       "Geodetic Reference System 1980" a=6378137.000  f=1/298.257222101

datum.table:

# from:
wgs84  "WGS_1984"     wgs84     dx=0.0      dy=0.0      dz=0.0
nad83   "North_American_Datum_1983"       grs80     dx=0.0      dy=0.0  dz=0.0
# to:
wgs84  "WGS_1984"     WGS84     dx=0.0      dy=0.0      dz=0.0
nad83   "North_American_Datum_1983"       GRS80     dx=0.0      dy=0.0  dz=0.0

[neteler]

Isn't it that it is read twice, i.e. separately in libgis and location wizard? See these and following lines therein:

libgis: https://github.com/OSGeo/grass/blob/802b3c84beaad21ad15056257e50f19ff6cb951c/lib/gis/get_ellipse.c#L259

GUI location wizard: https://github.com/OSGeo/grass/blob/02bfae03d5db8c80a93a0e340ff3cdcd509a2785/gui/wxpython/location_wizard/wizard.py#L2599

[nilason]

Isn't it that it is read twice, i.e. separately in libgis and location wizard? See these and following lines therein:

I doubt it, they both read the same file and "wrong" spelling, but they are not connected directly as far as I can tell. Location wizard reads it, uses it in creating the proj4string which is passed to g.proj.

[nilason]

Should emphasise: the examples with eg. "projection code = ll" (lat/lon) and "datum code = nad83" works with above mentioned changes as well as similar combos using wsg84 as datum (which also gives the same error without changing the ellipse and datum files). Probably PROJ 5+ is less forgiving (case sensitive) than before.

Reference of ellipse codes: https://github.com/OSGeo/PROJ/blob/master/src/ellps.cpp

[neteler]

Do we still need our private copy of the ellipse file?

[nilason]

Do we still need our private copy of the ellipse file?

Very relevant question. And the same applies to some of the other files ending up in etc/proj like the datum.table and projections. I suspect there were valid reasons at some point in GRASS history to introduce them. Now with a very active PROJ project perhaps not anymore, or at least partially not.

[nilason]

This might however be fixed in GUI anyway, skipping the proj4string part +ellps=xxx as its values are inserted anyway eg. +a=yyyy +rf=xxxx as defined in ellipse.tables, which actually "does the job". Then there will be no need to change the files.

[nilason]

@cmbarton , @bal-agates , @mrael2 Could you please test if #1550 solves your issues. While I'm quite sure it solves certain combinations, I'm curious if there are cases which it doesn't.

[cmbarton]

@cmbarton , @bal-agates , @mrael2 Could you please test if #1550 solves your issues. While I'm quite sure it solves certain combinations, I'm curious if there are cases which it doesn't.

I will try this out tomorrow. Thanks much!

[cmbarton]

@cmbarton , @bal-agates , @mrael2 Could you please test if #1550 solves your issues. While I'm quite sure it solves certain combinations, I'm curious if there are cases which it doesn't.

I will try this out tomorrow. Thanks much!

I just build this for the Mac and it seems to work fine! I created a UTM/ETRS89 location for Spain and a LCC/NAD83 location for the US. Both were created without errors (though I can't test if they were created correctly from a projection perspective).

[nilason]

I just build this for the Mac and it seems to work fine! I created a UTM/ETRS89 location for Spain and a LCC/NAD83 location for the US. Both were created without errors (though I can't test if they were created correctly from a projection perspective).

What does g.proj -j give?

[cmbarton]

I just build this for the Mac and it seems to work fine! I created a UTM/ETRS89 location for Spain and a LCC/NAD83 location for the US. Both were created without errors (though I can't test if they were created correctly from a projection perspective).

What does g.proj -j give?

In one of the new locations or in any location?

[cmbarton]

This is for the Lambert CC NAD83 location. Looks good at a glance.

+proj=lcc +lat_0=23 +lon_0=-96 +lat_1=33 +lat_2=45 +x_0=0 +y_0=0 +no_defs +a=6378137 +rf=298.257222101 +towgs84=0.000,0.000,0.000 +type=crs +to_meter=1