[GRASS5] datum transforms...

On Tue, Jan 28, 2003 at 03:32:53PM -0500, Helena wrote:

I would vote for m.proj to be removed because it does not work with GRASS
data, but it may be useful to put somewhere information about cs2cs (which
I did not know about).

So, now we have 2 votes to remove it, no other votes.
Other opinions?

Maybe keep m.proj in the manual for some time with a note that it was
retired and people should use cs2cs (where do you get it?)

You get it here:

http://www.remotesensing.org/proj/
PROJ4 software

Paul Kelly added recently a hint to m.proj2 html page that it supports
datum transformation.
But m.proj doesn't.

The best place of course for this (and other essential non-GRASS tools)
would be in an on-line GRASS tutorial if it ever happens.

This list is already there:
http://grass.itc.it/
Menu: Related projects

[ok, the cs2cs hint I have added today]

Markus

But this is just my opinion (I have never used m.proj except for the
book)

Helena

Markus Neteler wrote:
>
> Hello Paul,
>
> (cc to 'grass5' as it is of general interest)
>
> On Mon, Jan 27, 2003 at 10:40:34AM +0000, Paul Kelly wrote:
> > Hello Markus
> > >From what I can see m.proj does the same thing as the proj or cs2cs
> > programs included with the PROJ.4 library. Does it pre-date its general
> > availability or is there something else it does that I have missed?
>
> The m.proj is there for a long time, maybe as a sort of convenient
> tool for 'proj'.
> The problem is that results calculated with m.proj/m.proj2 (this is the cmd
> version) are shifted by around 100m.
>
> > Do you have any examples of what people use it for?
>
> Yes: E.G. they have a paper map in a different projection (for
> example in Germany UTM and Gauss-Krueger are used in parallel).
> to check points etc some users use m.proj (maybe because they
> don't know that nowadays the nice 'cs2cs' exists).
>
> Then people may also use it to reproject point lists from GPS devices etc.
> before importing them with s.proj.
>
> > It may be simpler to re-write it using most of the code from cs2cs, but
> > that seems a bit pointless.
>
> Sounds right. But then we should think about removing m.proj[2]
> from CVS to avoid that some people continue to use it and receive
> shifted point data.
>
> > If m.proj doesn't actually operate on files in the GRASS database but uses
> > external text files, should it be there at all.....
>
> This might be discussed on 'grass5' if we want to support that
> in general in future. The problem is that many users take results for
> truth when working with GRASS tools. So I feel that it is dangerous
> when the results from m.proj[2] and cs2cs differ.
>
> An example: UTM Zone 32-> Gauss-Krueger (Transverse Mercator)
> # proj/nad/epsg
> # DHDN / Germany zone 3
> # datum from:
> # http://www.colorado.Edu/geography/gcraft/notes/datum/dlist.html
> cs2cs +init=epsg:32632 +to +init=epsg:31493 +towgs84=606,23,413
> 486000 5912000 0
> 3486066.49 5913928.42 6.189
>
> (this result is rather identical to other proprietary software, it differs
> by some tens of centimeters)
>
> m.proj delivers:
>
> Initializing INPUT PROJECTION:
> Please specify projection name
> >utm
> Initializing INPUT projection ELLIPSOID:
> >wgs84
> Enter INPUT Projection Zone [zone] : 32
> INPUT Projection: Would you like to use South Hemisphere (y/[n]) ? n
> Enter plural for INPUT units [meters]:
> Parameters used in INPUT projection:
> +proj=utm +a=6378137.0000000000 +es=0.0066943800 +zone=32 +unfact=1.0000000000
>
> Initializing OUTPUT PROJECTION:
> Please specify projection name
> >tmerc
> Initializing OUTPUT projection ELLIPSOID:
> Please specify ellipsoid name
> >bessel
> Enter OUTPUT Central Parallel [lat_0] (23N) : 0N
> Enter OUTPUT Central Meridian [lon_0] (96W) : 9E
> Enter OUTPUT Scale Factor at the Central Meridian [k_0] (1.000000) :
> Enter OUTPUT False Easting [x_0] (0.000000) : 3500000
> Enter OUTPUT False Northing [y_0] (0.000000) :
> Enter plural for OUTPUT units [meters]:
> Parameters used in OUTPUT projection:
> +proj=tmerc +a=6377397.1550000003 +es=0.0066743722
> +lat_0=0.0000000000 +lon_0=9.0000000000 +k=1.0000000000
> +x_0=3500000.0000000000 +y_0=0.0000000000 +unfact=1.0000000000
>
> Universe Transverse Mercator -> Transverse Mercator Conversion:
> X_in (meters) Y_in (meters) X_out (meters) Y_out (meters)
> ------------ ----------- ------------- -------------
> 486000.00 5912000.00 3485996.11 5913755.52
>
> -----------
> Differences cs2cs/m.proj:
> 3486066.49-3485996.11=70.38
> 5913928.42-5913755.52=172.9
>
> So the results of m.proj/m.proj2 are shifted as expected.
>
> As conclusion:
> Is it right that either m.proj/m.proj2 should ask for the
> datum (and pass to proj) or both commands be removed?
>
> Markus
>
> _______________________________________________
> grass5 mailing list
> grass5@grass.itc.it
> http://grass.itc.it/mailman/listinfo/grass5

_______________________________________________
grass5 mailing list
grass5@grass.itc.it
http://grass.itc.it/mailman/listinfo/grass5

On Thu, 13 Mar 2003, Markus Neteler wrote:

On Thu, Mar 13, 2003 at 03:31:43PM +0000, Paul Kelly wrote:
>
> On Thu, 13 Mar 2003, Markus Neteler wrote:
>
> > Suggestion: the [srv].proj could print the *used* proj info as well,
> > along with the transformed boundary coordinates (the latter is already
> > in r.proj). Then the user would know if datum/ellipsoid were picked up
> > as expected.
>
> I have done this already for debugging it, but it involved printing out
> the parameters from inside the library function pj_get_kv(). Printing to
> stderr from a library function probably isn't good as it mightn't suit all
> the programs that would be calling it. But it would be more awkward to put
> it into the user programs as well as extra work having to change all 3 of
> them. I'm not sure what would be the best way but will think about it when
> I have time.

I've worked out a better way of doing this using the new PROJ.4
pj_get_def() function, which expands the projection paramters into a
string as would be used for cs2cs input. So today I have implemented this
in the three programs [rsv].proj. But something unusual that this has
turned up is that v.proj appears to be using the same parameters for the
input and output projection, and a quick glance at the source code shows
that it sets both the output *and* input paramters after changing location
to the input location (r.proj and s.proj set the output parameters right
at the start before changing location).

I must have missed something as v.proj can't have been broken all this
time; does anybody use it and can confirm what I am missing? I have only
ever used r.proj and s.proj for real work before.

Paul

Dear developers,

based on help and suggestion from Paul Kelly I have
updated the PROJ_INFO file in the Spearfish sample data set.

The ellipsoid and datum was wrong (was WGS84, but it is Clarke66).

Paul suggested:

On Fri, May 09, 2003 at 05:01:35PM +0100, Paul Kelly wrote:

Hello again Markus
I think I have found a way of writing the PROJ_INFO so it is backward
compatible; it relies on the fact that if GRASS doesn't recognise an
option it will pass it on unaltered to PROJ so we just put the nadgrids
datum parameter in at the end:

name: UTM
proj: utm
ellps: clark66
a: 6378206.4000000004
es: 0.0067686580
f: 294.9786982000
zone: 13
nadgrids: conus

This won't be exactly the standard way of writing a PROJ_INFO file but it
will work and there is no problem with it so I think it should be fine to
use it to avoid having different spearfish versions like you say.

Paul

The updated package is here:

http://mpa.itc.it/markus/tmp/spearfish_grass50data.tar.gz

please
- try with any 5.0.x
- try on various platforms

You may simply run
g.region -l
which uses PROJ4 to report the current region in Lat/Long coordinates.

If reports are ok, I'll move the file to
http://grass.itc.it/data.html

Thanks,
(especially to Paul)

Markus

On Fri, 9 May 2003, Markus Neteler wrote:

The updated package is here:

http://mpa.itc.it/markus/tmp/spearfish_grass50data.tar.gz

please
- try with any 5.0.x
- try on various platforms

You may simply run
g.region -l
which uses PROJ4 to report the current region in Lat/Long coordinates.

If reports are ok, I'll move the file to
http://grass.itc.it/data.html

I think it should be safe to add
datum: nad27
to the parameters. This makes it clear what is going on and will be
ignored in GRASS versions before 5.0.0pre4 when there wasn't any datum
support at all. And people using 5.0.0pre4 to 5.0.2 will be able to take
advantage of the nad2nad datum conversions that were operational in those
versions.

So the new spearfish PROJ_INFO is:

name: UTM
datum: nad27
nadgrids: conus
proj: utm
ellps: clark66
a: 6378206.4000000004
es: 0.0067686580
f: 294.9786982000
zone: 13

I have also updated g.setproj so the PROJ_INFO files it writes are this
style, with either 'nadgrids: xxx' or 'towgs84: xxx' below the datum: line.
The old format 'datumparams: towgs84=xxx' is still accepted by pj_get_kv()
etc. and will be useful if additional datum shift parameters are added to
PROJ.4 in the future (which was my original intention) but I feel
compatibility of the PROJ_INFO files with older GRASS versions is more
important.

GRASS versions 5.0.0pre4 to 5.0.2 should now only have a problem if a datum
name that wasn't present in those versions is used. (It can't be fixed by
updating the src/libes/gis/datum.table file (i.e. $(GISBASE)/etc/datum.table)
as in 5.0.0pre4 to 5.0.2 the datums were hard coded into
src/libes/proj/pj_datums.c; this has now been changed). Versions older than
5.0.0pre4 totally ignored the datum so it will cause no problems there.

This is a probably good time to remind again that for 5.0.0pre4 to 5.0.2,
when 'datum: nad27' was present in the PROJ_INFO file a datum shift parameter
of 'nadgrids: conus' was automatically assumed and any other datum
parameters in the PROJ_INFO file were ignored. For the latest CVS version
the 'nadgrids: conus' datum parameter line must be added explicitly instead of
the dx: dy: dz 3-parameter transformation parameters. But this can be fixed by
running the latest g.setproj

Paul

On Mon, May 12, 2003 at 09:33:24AM +0100, Paul Kelly wrote:

On Fri, 9 May 2003, Markus Neteler wrote:

>
> The updated package is here:
>
> http://mpa.itc.it/markus/tmp/spearfish_grass50data.tar.gz

[...]

Thanks for your efforts, Paul. I have updated the PROJ_INFO file
in above package again according to your suggestion.

After receiving confirmation that above PROJ_INFO is valid also
for older GRASS versions than CVS HEAD, I'll update the
package with wrong PROJ_INFO on the GRASS web site.

Markus

On Mon, May 12, 2003 at 09:33:24AM +0100, Paul Kelly wrote:

On Fri, 9 May 2003, Markus Neteler wrote:
>
> The updated package is here:
>
> http://mpa.itc.it/markus/tmp/spearfish_grass50data.tar.gz
>
> please
> - try with any 5.0.x
> - try on various platforms
>
> You may simply run
> g.region -l
> which uses PROJ4 to report the current region in Lat/Long coordinates.
>
> If reports are ok, I'll move the file to
> http://grass.itc.it/data.html

Are there objections to update the spearfish_grass50data.tar.gz file?

If not, I'll replace the old version with wrong datum.

Markus

On Wed, 4 Jun 2003, Markus Neteler wrote:

Are there objections to update the spearfish_grass50data.tar.gz file?

If not, I'll replace the old version with wrong datum.

I say go ahead; I think it will be fine, and who knows: there may already
be people experimenting with importing the WGS84-referenced data
mentioned on the website and wondering why it does not align properly.

Paul