[GRASS5] LatLong E/W problem in southern hemisphere

GRASS GIS GRASS GIS developer
1

Hi,

we have a map in LatLong, located somewhere in Argentina:

g.region rast=noaa -p
projection: 3 (Latitude-Longitude)
zone: 0
datum: wgs84
ellipsoid: a=6378137 es=0.00669438
north: 30:24S
south: 35:42:00.0036S
west: 58:17:59.9964W
east: 65:11:59.9892W <- note this correct value
nsres: 0:03:47.1429
ewres: 4:42:28.800096
rows: 84
cols: 75
                                                                                                          Mapset <PERMANENT> in Location <santafell> GRASS 5.7.cvs > g.region -g
n=-30.4
s=-35.700001
w=-58.299999
e=294.800003 <- note the wrong value
nsres=0.06309525
ewres=4.70800003

The bug seems to be somewhere in libgis, because:

r.in.bin ~/pub/avhrrpf.ndvi.1ntfsa.840101.11705583 output=noaa bytes=1 north=-30.400000 south=-35.700001 east=-65.199997 west=-58.299999 rows=84 cols=75
D0/0: E: 294.800003 N:-58.299999
D0/0: E: 294.800003 N:-58.299999
D0/0: E: 294.800003 N:-58.299999
Percent Complete: 100%
CREATING SUPPORT FILES FOR noaa

I have added the debug output to
lib/gis/adj_cellhd.c
(line 103)
   G_debug(0,"E: %f N:%F",cellhd->east,cellhd->west);

Anyone working in the southern hemisphere and observing similar
problems?

Markus

2

On Fri, 15 Oct 2004, Markus Neteler wrote:

Hi,

we have a map in LatLong, located somewhere in Argentina:

g.region rast=noaa -p
projection: 3 (Latitude-Longitude)
zone: 0
datum: wgs84
ellipsoid: a=6378137 es=0.00669438
north: 30:24S
south: 35:42:00.0036S
west: 58:17:59.9964W
east: 65:11:59.9892W <- note this correct value

Just asking: does Southern hemisphere mean we look at things upside down?
Shouldn't east be (yes) east of west, this east is almost 7 degrees west
of west, and I think the algorithm is (logically) assuming that you want
what is asked for, 353 degrees? Or am I missing something (I wouldn't be
surprised!)? Maybe a warning could be useful when the nsres/ewres ratio is
so extreme?

Roger

nsres: 0:03:47.1429
ewres: 4:42:28.800096
rows: 84
cols: 75
                                                                                                          Mapset <PERMANENT> in Location <santafell> GRASS 5.7.cvs > g.region -g
n=-30.4
s=-35.700001
w=-58.299999
e=294.800003 <- note the wrong value
nsres=0.06309525
ewres=4.70800003

The bug seems to be somewhere in libgis, because:

r.in.bin ~/pub/avhrrpf.ndvi.1ntfsa.840101.11705583 output=noaa bytes=1 north=-30.400000 south=-35.700001 east=-65.199997 west=-58.299999 rows=84 cols=75
D0/0: E: 294.800003 N:-58.299999
D0/0: E: 294.800003 N:-58.299999
D0/0: E: 294.800003 N:-58.299999
Percent Complete: 100%
CREATING SUPPORT FILES FOR noaa

I have added the debug output to
lib/gis/adj_cellhd.c
(line 103)
   G_debug(0,"E: %f N:%F",cellhd->east,cellhd->west);

Anyone working in the southern hemisphere and observing similar
problems?

Markus

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

--
Roger Bivand
Economic Geography Section, Department of Economics, Norwegian School of
Economics and Business Administration, Breiviksveien 40, N-5045 Bergen,
Norway. voice: +47 55 95 93 55; fax +47 55 95 93 93
e-mail: Roger.Bivand@nhh.no

3

On Fri, Oct 15, 2004 at 06:30:40PM +0200, Roger Bivand wrote:

On Fri, 15 Oct 2004, Markus Neteler wrote:

> Hi,
>
> we have a map in LatLong, located somewhere in Argentina:
>
> g.region rast=noaa -p
> projection: 3 (Latitude-Longitude)
> zone: 0
> datum: wgs84
> ellipsoid: a=6378137 es=0.00669438
> north: 30:24S
> south: 35:42:00.0036S
> west: 58:17:59.9964W
> east: 65:11:59.9892W <- note this correct value

Just asking: does Southern hemisphere mean we look at things upside down?
Shouldn't east be (yes) east of west, this east is almost 7 degrees west
of west, and I think the algorithm is (logically) assuming that you want
what is asked for, 353 degrees? Or am I missing something (I wouldn't be
surprised!)? Maybe a warning could be useful when the nsres/ewres ratio is
so extreme?

Roger, the problem (to me) is that above we have the expected value
65:11:59.9892W, which becomes later 294.800003
(360 - 294.800003 = 65.2 !). So the wrap around doesn't work (or whatever).

GRASS 5.7.cvs > g.region -g

> n=-30.4
> s=-35.700001
> w=-58.299999
> e=294.800003 <- note the wrong value
> nsres=0.06309525
> ewres=4.70800003
>
> The bug seems to be somewhere in libgis, because:
>
> r.in.bin ~/pub/avhrrpf.ndvi.1ntfsa.840101.11705583 output=noaa bytes=1 north=-30.400000 south=-35.700001 east=-65.199997 west=-58.299999 rows=84 cols=75
> D0/0: E: 294.800003 N:-58.299999
> D0/0: E: 294.800003 N:-58.299999
> D0/0: E: 294.800003 N:-58.299999
> Percent Complete: 100%
> CREATING SUPPORT FILES FOR noaa
>
> I have added the debug output to
> lib/gis/adj_cellhd.c
> (line 103)
> G_debug(0,"E: %f N:%F",cellhd->east,cellhd->west);
>
> Anyone working in the southern hemisphere and observing similar
> problems?

If east remained -65.199997 instead of becoming 294.800003, everything would
be ok and the test for East > West would work as well.

Markus

4

On Fri, 15 Oct 2004, Markus Neteler wrote:

On Fri, Oct 15, 2004 at 06:30:40PM +0200, Roger Bivand wrote:
> On Fri, 15 Oct 2004, Markus Neteler wrote:
>
> > Hi,
> >
> > we have a map in LatLong, located somewhere in Argentina:
> >
> > g.region rast=noaa -p
> > projection: 3 (Latitude-Longitude)
> > zone: 0
> > datum: wgs84
> > ellipsoid: a=6378137 es=0.00669438
> > north: 30:24S
> > south: 35:42:00.0036S
> > west: 58:17:59.9964W
> > east: 65:11:59.9892W <- note this correct value
>
> Just asking: does Southern hemisphere mean we look at things upside down?
> Shouldn't east be (yes) east of west, this east is almost 7 degrees west
> of west, and I think the algorithm is (logically) assuming that you want
> what is asked for, 353 degrees? Or am I missing something (I wouldn't be
> surprised!)? Maybe a warning could be useful when the nsres/ewres ratio is
> so extreme?

Roger, the problem (to me) is that above we have the expected value
65:11:59.9892W, which becomes later 294.800003
(360 - 294.800003 = 65.2 !). So the wrap around doesn't work (or whatever).

Why is the declared east west of west? Isn't the program just doing what
it was asked to do? It seems to say that, well, if the stated east is west
of west *in lat-long* this can be resolved by going round the other way,
which it then does. Isn't this an input error?

east=-65.199997 west=-58.299999

west= east=

seems more likely? I don't think it's the code, I think it's doing what it
was told?

Roger

GRASS 5.7.cvs > g.region -g
> > n=-30.4
> > s=-35.700001
> > w=-58.299999
> > e=294.800003 <- note the wrong value
> > nsres=0.06309525
> > ewres=4.70800003
> >
> > The bug seems to be somewhere in libgis, because:
> >
> > r.in.bin ~/pub/avhrrpf.ndvi.1ntfsa.840101.11705583 output=noaa bytes=1 north=-30.400000 south=-35.700001 east=-65.199997 west=-58.299999 rows=84 cols=75
> > D0/0: E: 294.800003 N:-58.299999
> > D0/0: E: 294.800003 N:-58.299999
> > D0/0: E: 294.800003 N:-58.299999
> > Percent Complete: 100%
> > CREATING SUPPORT FILES FOR noaa
> >
> > I have added the debug output to
> > lib/gis/adj_cellhd.c
> > (line 103)
> > G_debug(0,"E: %f N:%F",cellhd->east,cellhd->west);
> >
> > Anyone working in the southern hemisphere and observing similar
> > problems?

If east remained -65.199997 instead of becoming 294.800003, everything would
be ok and the test for East > West would work as well.

Markus

--
Roger Bivand
Economic Geography Section, Department of Economics, Norwegian School of
Economics and Business Administration, Breiviksveien 40, N-5045 Bergen,
Norway. voice: +47 55 95 93 55; fax +47 55 95 93 93
e-mail: Roger.Bivand@nhh.no

5

On Fri, Oct 15, 2004 at 06:53:39PM +0200, Roger Bivand wrote:

On Fri, 15 Oct 2004, Markus Neteler wrote:

> On Fri, Oct 15, 2004 at 06:30:40PM +0200, Roger Bivand wrote:
> > On Fri, 15 Oct 2004, Markus Neteler wrote:
> >
> > > Hi,
> > >
> > > we have a map in LatLong, located somewhere in Argentina:
> > >
> > > g.region rast=noaa -p
> > > projection: 3 (Latitude-Longitude)
> > > zone: 0
> > > datum: wgs84
> > > ellipsoid: a=6378137 es=0.00669438
> > > north: 30:24S
> > > south: 35:42:00.0036S
> > > west: 58:17:59.9964W
> > > east: 65:11:59.9892W <- note this correct value
> >
> > Just asking: does Southern hemisphere mean we look at things upside down?
> > Shouldn't east be (yes) east of west, this east is almost 7 degrees west
> > of west, and I think the algorithm is (logically) assuming that you want
> > what is asked for, 353 degrees? Or am I missing something (I wouldn't be
> > surprised!)? Maybe a warning could be useful when the nsres/ewres ratio is
> > so extreme?
>
> Roger, the problem (to me) is that above we have the expected value
> 65:11:59.9892W, which becomes later 294.800003
> (360 - 294.800003 = 65.2 !). So the wrap around doesn't work (or whatever).
>

Why is the declared east west of west? Isn't the program just doing what
it was asked to do? It seems to say that, well, if the stated east is west
of west *in lat-long* this can be resolved by going round the other way,
which it then does. Isn't this an input error?

east=-65.199997 west=-58.299999

west= east=

seems more likely? I don't think it's the code, I think it's doing what it
was told?

Roger,

let's see it from a user's perspective (my colleague):
He fetched the NOAA image from the web, file in plain binary
format. The metadata told him how to select the parameters
for r.in.bin.
First he mixed east and west and it gave *no* error (while it
does in case of north/south confusion). This delivered a messy
map. Eventually I found out that he mixed east/west, changing
that the map was imported correctly. But of course I was
wondering why r.in.bin didn't complain.

You may be right that the behaviour is correct. But in fact
not helpful for this particular case. I darkly remember that
this was discussed once here, but don't recall details.
It was probably the same trap of 'how to go around the globe...'.
But shouldn't this then raise alarm due to a wrong number
of columns? Ah, in fact the resolution was wrong then. So
it's the way how r.in.bin determines rows/cols/resolution.

Now it seems clear to me: no bug, but unfortunate case.

Markus

6

On Fri, 15 Oct 2004, Markus Neteler wrote:

On Fri, Oct 15, 2004 at 06:53:39PM +0200, Roger Bivand wrote:
> On Fri, 15 Oct 2004, Markus Neteler wrote:
>
> > On Fri, Oct 15, 2004 at 06:30:40PM +0200, Roger Bivand wrote:
> > > On Fri, 15 Oct 2004, Markus Neteler wrote:
> > >
> > > > Hi,
> > > >
> > > > we have a map in LatLong, located somewhere in Argentina:
> > > >
> > > > g.region rast=noaa -p
> > > > projection: 3 (Latitude-Longitude)
> > > > zone: 0
> > > > datum: wgs84
> > > > ellipsoid: a=6378137 es=0.00669438
> > > > north: 30:24S
> > > > south: 35:42:00.0036S
> > > > west: 58:17:59.9964W
> > > > east: 65:11:59.9892W <- note this correct value
> > >
> > > Just asking: does Southern hemisphere mean we look at things upside down?
> > > Shouldn't east be (yes) east of west, this east is almost 7 degrees west
> > > of west, and I think the algorithm is (logically) assuming that you want
> > > what is asked for, 353 degrees? Or am I missing something (I wouldn't be
> > > surprised!)? Maybe a warning could be useful when the nsres/ewres ratio is
> > > so extreme?
> >
> > Roger, the problem (to me) is that above we have the expected value
> > 65:11:59.9892W, which becomes later 294.800003
> > (360 - 294.800003 = 65.2 !). So the wrap around doesn't work (or whatever).
> >
>
> Why is the declared east west of west? Isn't the program just doing what
> it was asked to do? It seems to say that, well, if the stated east is west
> of west *in lat-long* this can be resolved by going round the other way,
> which it then does. Isn't this an input error?
>
> east=-65.199997 west=-58.299999
>
> west= east=
>
> seems more likely? I don't think it's the code, I think it's doing what it
> was told?
>
Roger,

let's see it from a user's perspective (my colleague):
He fetched the NOAA image from the web, file in plain binary
format. The metadata told him how to select the parameters
for r.in.bin.
First he mixed east and west and it gave *no* error (while it
does in case of north/south confusion). This delivered a messy
map. Eventually I found out that he mixed east/west, changing
that the map was imported correctly. But of course I was
wondering why r.in.bin didn't complain.

You may be right that the behaviour is correct. But in fact
not helpful for this particular case. I darkly remember that
this was discussed once here, but don't recall details.
It was probably the same trap of 'how to go around the globe...'.
But shouldn't this then raise alarm due to a wrong number
of columns? Ah, in fact the resolution was wrong then. So
it's the way how r.in.bin determines rows/cols/resolution.

Now it seems clear to me: no bug, but unfortunate case.

Right - I think a warning when the ratio of resolutions seems extreme
might be justified, or a note on the manpage to check the ewres for
sanity. With a projected data set, the E > W test would have worked.

Roger

Markus

--
Roger Bivand
Economic Geography Section, Department of Economics, Norwegian School of
Economics and Business Administration, Breiviksveien 40, N-5045 Bergen,
Norway. voice: +47 55 95 93 55; fax +47 55 95 93 93
e-mail: Roger.Bivand@nhh.no