[GRASS5] [bug #3876] (grass) g.proj: proj4= fails, wkt= works partialy

GRASS GIS GRASS GIS developer
1

this bug's URL: http://intevation.de/rt/webrt?serial_num=3876
-------------------------------------------------------------------------

Subject: g.proj: proj4= fails, wkt= works partialy

Platform: GNU/Linux/i386
GRASS Version: CVS 2005-11-16

Here is my initial projection:

GRASS 6.1.cvs (huha2):~ > g.proj -p
-PROJ_INFO-------------------------------------------------
name : Latitude-Longitude
datum : S-42
towgs84 : 33.4,-146.6,-76.3,-0.359,-0.053,0.844,-0.84
proj : ll
ellps : krassovsky
-PROJ_UNITS------------------------------------------------
unit : degree
units : degrees
meters : 1.0

I wanted to change it with "g.proj wkt=wkt_file". It kindoff works, but:
1. the PERMANENT/WIND and PERMANENT/DEFAULT_WIND are left unchanged (still refer to latlong).
2. only the WIND in working mapset is changed properly
3. surprisingly, a correct DEFAULT_WIND appears in current working mapset - but it should go to PERMANENT.

Trying to change a projection to a following PROJ.4 description:

GRASS 6.1.cvs (huha2):~ > cat puwg1965_IV.prf
+proj=sterea +lat_0=50.625 +lon_0=21.08333333333333 +k=0.9998 +x_0=4637000 +y_0=5467000 +ellps=krass +towgs84=33.4297,-146.5746,-76.2865,-0.35867,-0.05283,0.84354,-0.84077 +units=m

It fails, no change is done:

GRASS 6.1.cvs (huha2):~ > g.proj -c proj4=puwg1965_IV.prf

WARNING! A projection file already exists for this location

This file contains all the parameters for the
location's projection: ll

    Overriding this information implies that the old projection parameters
    were incorrect. If you change the parameters, all existing data will be
    interpreted differently by the projection software.
    GRASS will not re-project your data automatically

Would you still like to overwrite the current projection information (y/n) [n] yProjection information updated!

GRASS 6.1.cvs (huha2):~ > g.proj -p
-PROJ_INFO-------------------------------------------------
name : Latitude-Longitude
datum : S-42
towgs84 : 33.4,-146.6,-76.3,-0.359,-0.053,0.844,-0.84
proj : ll
ellps : krassovsky
-PROJ_UNITS------------------------------------------------
unit : degree
units : degrees
meters : 1.0

Maciek

-------------------------------------------- Managed by Request Tracker

2

Request Tracker wrote:

I wanted to change it with "g.proj wkt=wkt_file". It kindoff works, but:
1. the PERMANENT/WIND and PERMANENT/DEFAULT_WIND are left unchanged (still refer to latlong).
2. only the WIND in working mapset is changed properly
3. surprisingly, a correct DEFAULT_WIND appears in current working mapset - but it should go to PERMANENT.

I take it you mean "g.proj -c wkt=wkt_file". You are right in that it won't work properly if you are not in the PERMANENT mapset. You need to have write access to the PERMANENT mapset anyway if you are going to change the projection files. I think changing the projection is a dramatic-enough thing to do to require the user to have selected the PERMANENT mapset before doing this, so I have changed g.proj to exit with an error if the user is not in PERMANENT.

This functionality (overwrite existing projection files) was not well tested; I have noticed also that the region was overwritten every time regardless of whether it had changed or not (if projection changed between UTM, XY, latlong or other, or if the zone changed {i.e. to or from UTM or changing UTM zone}). So if none of those change then the region file doesn't need to be changed and it won't change it. (Ideally of course the region file should not hold projection information but it is a legacy "feature" that it does and it is not really feasible to change it at this stage.)

If you have more multiple mapsets and the default region changes, you will need to run "g.region -d" in each mapset to update the region. This is another "feature" and also holds for g.setproj. I've changed g.proj so it will output a warning if you're likely to need to run g.region -d.

Trying to change a projection to a following PROJ.4 description:

GRASS 6.1.cvs (huha2):~ > cat puwg1965_IV.prf
+proj=sterea +lat_0=50.625 +lon_0=21.08333333333333 +k=0.9998 +x_0=4637000 +y_0=5467000 +ellps=krass +towgs84=33.4297,-146.5746,-76.2865,-0.35867,-0.05283,0.84354,-0.84077 +units=m

It fails, no change is done:

GRASS 6.1.cvs (huha2):~ > g.proj -c proj4=puwg1965_IV.prf

That should not work and should give an error (but would not have). I've fixed the return code checking so that it does. You need to run it like:
cat puwg1965_IV.prf | g.proj -c proj4=-
it is like that because PROJ.4-style descriptions are likely to be much shorter and more feasible to type on the command line than to read in from a file, and means you can do neat things with init codes and EPSG numbers easily. With the input wkt option that would of course work (it reads from a file there).
See the man page for more examples:
http://grass.itc.it/grass61/manuals/html61_user/g.proj.html

Paul

3

Paul,

Thanks for your fixes.

There are more possible bugs:

1.
GRASS 6.1.cvs (poloniny_ll):~ > g.proj location=poloniny_65
wkt=puwg1965_IV.prj

New location is not created, no error output, no information. What
happens actually (if anything)?

2.
GRASS 6.1.cvs (poloniny_ll):~ > g.proj --help

Description:
Prints and manipulates GRASS projection information files.

Usage:
g.proj [-pdjwefc] [georef=name] [wkt=name] [proj4=name]
   [location=name]
            ^
                  |
        That's propably why I assumed I
        should specify a file name here
        (part of my initial bug report).

Maciek

--------------------
W polskim Internecie s± setki milionów stron. My przekazujemy Tobie tylko najlepsze z nich!
http://katalog.epf.pl/

4

Another one:

3.
when the current location is xy, "g.proj -c wkt=" doesn't make any
change and no error is printed:

GRASS 6.1.cvs (poloniny_65):~ > g.proj -p
XY location (unprojected)
GRASS 6.1.cvs (poloniny_65):~ > g.proj -c wkt=puwg1965_IV.prj
GRASS 6.1.cvs (poloniny_65):~ > g.proj -p
XY location (unprojected)

Maciek

--------------------
W polskim Internecie s± setki milionów stron. My przekazujemy Tobie tylko najlepsze z nich!
http://katalog.epf.pl/

5

Maciek Sieczka wrote:

Another one:

3.
when the current location is xy, "g.proj -c wkt=" doesn't make any
change and no error is printed:

I can see how that would have happened before but I think the second change I made this afternoon should have fixed that. Have you tried with the latest version?

Paul

6

Maciek Sieczka wrote:

Paul,

Thanks for your fixes.

There are more possible bugs:

1.
GRASS 6.1.cvs (poloniny_ll):~ > g.proj location=poloniny_65
wkt=puwg1965_IV.prj

New location is not created, no error output, no information. What
happens actually (if anything)?

In theory it should verify the input file (puwg1965_IV.prj) is correct/acceptable and then not do anything because you haven't specified any flags. None of the flags are required, but I suppose it is required to have at least one or else nothing happens. I don't know if there's an easy way to program it in the parser though. If there isn't it could be done in the module with just a little bit more work.

2.
GRASS 6.1.cvs (poloniny_ll):~ > g.proj --help

Description:
Prints and manipulates GRASS projection information files.

Usage:
g.proj [-pdjwefc] [georef=name] [wkt=name] [proj4=name]
   [location=name]
            ^
                  |
        That's propably why I assumed I
        should specify a file name here
        (part of my initial bug report).

Hmm I didn't realise it was possible to change the text shown there. But it is (using the key_desc option) and I have.

Paul

7

> GRASS 6.1.cvs (poloniny_ll):~ > g.proj --help
>
> Description:
> Prints and manipulates GRASS projection information files.
>
> Usage:
> g.proj [-pdjwefc] [georef=name] [wkt=name] [proj4=name]
> [location=name]
> ^
> |
> That's propably why I assumed I
> should specify a file name here
> (part of my initial bug report).

Hmm I didn't realise it was possible to change the text shown there.
But it is (using the key_desc option) and I have.

changed "name" to "string" in lib/gis/parser.c. The auto-gen help pages
and TCL dialogs already use that.

Hamish