[GRASSLIST:4910] Problems with importing polar data

Anton_B_Ivanov · November 8, 2002, 10:57pm

Hello GRASS gurus!

I am new to the GRASS program and it seems to be VERY powerful. However,
I have been a little frustrated trying to import my polar data into it.

Here's a problem : I have data (images) projected into a polar
stereographic projection (using some camera speicific software). I import
them using r.in.bin and set up proper lat/lon bounds for these images. The
projection setup is as ups (see PROJ_INFO file below).

However, when I'm trying display data, put a grid on it or make a
postscript map, the output projection look like a simple cylindrical or
latititude/longitude.

How do I tell GRASS to display my data in polar stereographic projection
and overlay grid on it?

Another problem is to project data from simple cylindrical projection
into polar stereographic (using r.proj) . I have no luck there either.

Any help will be highly appereciated!!!!

Thanks a lot in advance
Anton Ivanov

PS. here's PROJ_INFO and PROJ_UNITS pair for my location. I work with
Martian data, so geoid is that of Mars and no datum has been specified.

PROJ_INFO
----------------------------------
name: Universal Polar
Stereographic
proj: ups
ellps: mars_iau2000
a: 33961900.0000000000
es:0.0117373700
f: 169.8944472236
lat_0: -55.0000000000
lon_0: 0.0000000000
south: defined
-------------------------------
PROJ_UNITS
-------------------------------
unit: degree
units: degrees
meters: 59272.9493611000

Glynn_Clements · November 9, 2002, 2:53am

Anton B. Ivanov wrote:

  I am new to the GRASS program and it seems to be VERY powerful. However,
I have been a little frustrated trying to import my polar data into it.

  Here's a problem : I have data (images) projected into a polar
stereographic projection (using some camera speicific software). I import
them using r.in.bin and set up proper lat/lon bounds for these images. The
projection setup is as ups (see PROJ_INFO file below).

  However, when I'm trying display data, put a grid on it or make a
postscript map, the output projection look like a simple cylindrical or
latititude/longitude.

Is the data which you imported really in UPS, and not Lat/Lon?

How do I tell GRASS to display my data in polar stereographic projection
and overlay grid on it?

GRASS displays data "as is"; none of the display commands will attempt
to project it. The *only* commands which project data are r.proj,
v.proj and s.proj. Most of the other commands simply work with
arbitrary X/Y coordinates, without any consideration as to what those
coordinates actually represent (apart from a few commands which take
account of the units, e.g. feet/metres).

Another problem is to project data from simple cylindrical projection
into polar stereographic (using r.proj) . I have no luck there either.

If you have Lat/Lon data, import it into a Lat/Lon location. Then
create a UPS location, and run r.proj from that location.

  Any help will be highly appereciated!!!!

  Thanks a lot in advance
  Anton Ivanov

PS. here's PROJ_INFO and PROJ_UNITS pair for my location. I work with
Martian data, so geoid is that of Mars and no datum has been specified.

PROJ_INFO
----------------------------------
name: Universal Polar Stereographic
proj: ups

Again, is the data which you are importing really in UPS? The GRASS
import utilities don't perform any projection of the data; you have to
import it into a location whose projection matches the data.

--
Glynn Clements <glynn.clements@virgin.net>

Eric_G_Miller · November 9, 2002, 3:12am

On Fri, Nov 08, 2002 at 02:57:02PM -0800, Anton B. Ivanov wrote:

Hello GRASS gurus!

PS. here's PROJ_INFO and PROJ_UNITS pair for my location. I work with
Martian data, so geoid is that of Mars and no datum has been specified.

Mars?

PROJ_INFO
----------------------------------
name: Universal Polar
Stereographic
proj: ups
ellps: mars_iau2000
a: 33961900.0000000000
es:0.0117373700
f: 169.8944472236
lat_0: -55.0000000000
lon_0: 0.0000000000
south: defined
-------------------------------
PROJ_UNITS
-------------------------------
unit: degree
units: degrees
meters: 59272.9493611000

The above does not look correct, but I don't have experience with UPS.
Anyway, the proj doc OF-90-284.ps says UPS has the following default
parameters:

k: 0.9994 (scale factor)
lat_0: 0
x_0: 2000000 (meters)
y_0: 2000000 (meters)
lat_0: 90N (north) or 90S (south)
south: ?

Units should be meters AFAIK. So,

unit: meter
units: meters
meters: 1.0

g.setproj appears busted for this projection (can someone confirm?).

--
static const char signature =
"Copyright (c) 2002 Eric G. Miller <egm2@jps.net>";

Morten_Hulden1 · November 9, 2002, 7:44pm

On Fri, 8 Nov 2002, Eric G. Miller wrote:

On Fri, Nov 08, 2002 at 02:57:02PM -0800, Anton B. Ivanov wrote:
> Hello GRASS gurus!

> PS. here's PROJ_INFO and PROJ_UNITS pair for my location. I work with
> Martian data, so geoid is that of Mars and no datum has been specified.

Mars?

> PROJ_INFO
> ----------------------------------
> name: Universal Polar
> Stereographic
> proj: ups
> ellps: mars_iau2000
> a: 33961900.0000000000
> es:0.0117373700
> f: 169.8944472236
> lat_0: -55.0000000000
> lon_0: 0.0000000000
> south: defined
> -------------------------------
> PROJ_UNITS
> -------------------------------
> unit: degree
> units: degrees
> meters: 59272.9493611000

The above does not look correct, but I don't have experience with UPS.
Anyway, the proj doc OF-90-284.ps says UPS has the following default
parameters:

k: 0.9994 (scale factor)
lat_0: 0
x_0: 2000000 (meters)
y_0: 2000000 (meters)
lat_0: 90N (north) or 90S (south)
south: ?

Units should be meters AFAIK. So,

unit: meter
units: meters
meters: 1.0

g.setproj appears busted for this projection (can someone confirm?).

Yes.

If proj does not have built-in default values for k0, x0 and y0 for this
projection then src/libes/geo_init.c needs patching in order to add
k-factor, false easting and northing. Anyway, default values for lon_0 and
lat_0 should be fixed (though user can change the values when running
g.setproj). After patching g.setproj and m.proj should be recompiled.

Not sure if lat_0 really is necessary for this projection; the 'south'
flag really should tell proj what it needs to know.

--- geo_init.c 2002-06-16 17:29:25.000000000 +0200
+++ geo_init.c_new 2002-11-09 20:40:03.000000000 +0100
@@ -1150,11 +1150,23 @@

TABLE[UPS][LON0].ask = 1;
TABLE[UPS][LON0].def_exists = 1;
- TABLE[UPS][LON0].deflt = 20.0;
+ TABLE[UPS][LON0].deflt = 0.0;

TABLE[UPS][LAT0].ask = 1;
TABLE[UPS][LAT0].def_exists = 1;
- TABLE[UPS][LAT0].deflt = 55.0;
+ TABLE[UPS][LAT0].deflt = 90.0;
+
+ TABLE[UPS][KFACT].ask = 1;
+ TABLE[UPS][KFACT].def_exists = 1;
+ TABLE[UPS][KFACT].deflt = 0.9994;
+
+ TABLE[UPS][X0].ask = 1;
+ TABLE[UPS][X0].def_exists = 1;
+ TABLE[UPS][X0].deflt = 2000000.0;
+
+ TABLE[UPS][Y0].ask = 1;
+ TABLE[UPS][Y0].def_exists = 1;
+ TABLE[UPS][Y0].deflt = 2000000.0;

TABLE[UPS][SOUTH].ask = 1;

Morten_Hulden1 · November 9, 2002, 8:56pm

On Sat, 9 Nov 2002, Morten Hulden wrote:

If proj does not have built-in default values for k0, x0 and y0 for this
projection then src/libes/geo_init.c needs patching in order to add
k-factor, false easting and northing.

I checked by running external proj with -V: Looks like proj really uses
k_0=0.9994, x_0=2000000 and y_0=2000000 as default values for UPS
projection.

So patching geo_init.c is strictly not necessary and would only result in
a PROJ_INFO that reflects the values that proj would use anyway.

But when running g.setproj it's important to give latitude=90 or -90. 55N
for this projection does not make sense.

Morten Hulden

Eric_G_Miller · November 10, 2002, 4:03am

On Sat, Nov 09, 2002 at 09:56:25PM +0100, Morten Hulden wrote:

On Sat, 9 Nov 2002, Morten Hulden wrote:

> If proj does not have built-in default values for k0, x0 and y0 for this
> projection then src/libes/geo_init.c needs patching in order to add
> k-factor, false easting and northing.

I checked by running external proj with -V: Looks like proj really uses
k_0=0.9994, x_0=2000000 and y_0=2000000 as default values for UPS
projection.

So patching geo_init.c is strictly not necessary and would only result in
a PROJ_INFO that reflects the values that proj would use anyway.

But when running g.setproj it's important to give latitude=90 or -90. 55N
for this projection does not make sense.

I'd think patching geo_init.c to set all the defaults and not ask except for
lat_0 and south. I would just ask for south and then set lat_0
appropriately, but g.setproj is too scary for words... Still requires
user to know to use 90/-90 and n/y for south?, respectively...

Index: geo_init.c

RCS file: /grassrepository/grass/src/libes/gis/geo_init.c,v
retrieving revision 1.6
diff -r1.6 geo_init.c
1284c1284,1292
< TABLE[UPS][LON0].ask = 1;
---

        TABLE[UPS][KFACT].def_exists = 1;
        TABLE[UPS][KFACT].deflt = 0.9994;

        TABLE[UPS][X0].def_exists = 1;
        TABLE[UPS][X0].deflt = 2000000;

        TABLE[UPS][Y0].def_exists = 1;
        TABLE[UPS][Y0].deflt = 2000000;

1286c1294
< TABLE[UPS][LON0].deflt = 20.0;
---

TABLE[UPS][LON0].deflt = 0.0;

1290c1298
< TABLE[UPS][LAT0].deflt = 55.0;
---

TABLE[UPS][LAT0].deflt = 90.0;

--
static const char signature =
"Copyright (c) 2002 Eric G. Miller <egm2@jps.net>";

Morten_Hulden1 · November 10, 2002, 6:02am

On Sat, 9 Nov 2002, Eric G. Miller wrote:

On Sat, Nov 09, 2002 at 09:56:25PM +0100, Morten Hulden wrote:
>
> On Sat, 9 Nov 2002, Morten Hulden wrote:
>
> > If proj does not have built-in default values for k0, x0 and y0 for this
> > projection then src/libes/geo_init.c needs patching in order to add
> > k-factor, false easting and northing.
>
> I checked by running external proj with -V: Looks like proj really uses
> k_0=0.9994, x_0=2000000 and y_0=2000000 as default values for UPS
> projection.

...

I'd think patching geo_init.c to set all the defaults and not ask except for
lat_0 and south. I would just ask for south and then set lat_0
appropriately, but g.setproj is too scary for words... Still requires
user to know to use 90/-90 and n/y for south?, respectively...

also found out proj ignores lat_0 value in ups projection. so simplified
patch below. no need to ask for latitude at all and no need to touch the
too-ugly-for-words g.setproj.

the patch is purely cosmetic to make PROJ_INFO display real values used by
proj - ups projection should work correctly as is, even if values
PROJ_INFO are misleading.

better take this discussion over to developers' list.

--- geo_init.c 2002-06-16 17:29:25.000000000 +0200
+++ geo_init.c_new 2002-11-10 06:57:09.000000000 +0100
@@ -1150,11 +1150,16 @@

TABLE[UPS][LON0].ask = 1;
TABLE[UPS][LON0].def_exists = 1;
- TABLE[UPS][LON0].deflt = 20.0;
+ TABLE[UPS][LON0].deflt = 0.0;

- TABLE[UPS][LAT0].ask = 1;
- TABLE[UPS][LAT0].def_exists = 1;
- TABLE[UPS][LAT0].deflt = 55.0;
+ TABLE[UPS][KFACT].def_exists = 1;
+ TABLE[UPS][KFACT].deflt = 0.9994;
+
+ TABLE[UPS][X0].def_exists = 1;
+ TABLE[UPS][X0].deflt = 2000000.0;
+
+ TABLE[UPS][Y0].def_exists = 1;
+ TABLE[UPS][Y0].deflt = 2000000.0;

TABLE[UPS][SOUTH].ask = 1;