[GRASS5] Re: [GRASSLIST:4019] d.zoom in lat/lon?

[Hamish]

Hamish:
> > I'm noticing when using d.zoom in a lat/lon location that the
> > resolution is not preserved, and I need to run 'g.region
> > res=$ORIG_RES -a' to fix it afterwards.
> > Is this a bug or a feature?
> 
Markus:
> To me it looks like a bug 
...
> I briefly looked into this a few weeks ago but didn't it.
> Maybe you are able to fix that?


Bug #1:

Right then. Increasing the significant digits on the rounding
calculation in set.c seems to have fixed it for resolutions up to 
0:00:00.0001 (3mm)

The "13" below is arbitrarily large, it was just an easy value to 
insert in vi. ;)


Index: src/display/d.zoom/cmd/set.c
===================================================================
RCS file:
/home/grass/grassrepository/grass/src/display/d.zoom/cmd/set.c,v
retrieving revision 1.3
diff -u -r1.3 set.c
--- src/display/d.zoom/cmd/set.c        25 Sep 2002 04:46:44 -0000     
1.3
+++ src/display/d.zoom/cmd/set.c        27 Jul 2004 07:59:45 -0000
@@ -50,8 +50,8 @@

     resetres = 1;
     while ( resetres ) {
-       nsr = round_to(window->ns_res, 3);
-       ewr = round_to(window->ew_res, 3);
+       nsr = round_to(window->ns_res, 13);
+       ewr = round_to(window->ew_res, 13);

        td = ceil (north / nsr);
        tnorth = td * nsr;



 *Note there are several more calls to round_to() using sd=3 in the 
code. It's all a bit too convoluded for me to follow so I've left 
them as-is for now.


The change makes the round_to() function in set.c do:

+while( rint(foo) < pow(10,12) )
-while( rint(foo) < pow(10,3) )

Will that break anything on wierd platforms/compilers?

Anyone have ideas why 3 was used aside from it being about right 
for res=meters?

But I'm thinking I'm missing the whole point of the "resetres" while-loop 
by making the rounding so exact? ie, why not just 'nsr=window->nsres' ...



Markus:
>(at least it's annoying because when zooming on continental scale you 
> often get an error that you left the max latitude etc).

Bug #1.5:

I think this happens regardless, e.g.

G:~> g.region n=90N s=90S e=180E w=180W res=0:05
G:~> d.redraw
G:~> d.zoom
[middle click for unzoom]
> 360 -> reset

north:90:02:50.28N south:90:02:50.28S east:179:59:19.32W west:179:59:19.32E
North, south, east limit of default region reached.
WARNING: G_set_window(): Illegal latitude for North
ERROR: region for current mapset is invalid
       line 3: <north: 90:02:50.28N>
       run "g.region"


i.e. 90N + res/2,  90S - res/2

The sig. digits patch attached makes the above error into:
north: 90:05N  south: 90:05S  east: 180E  west: 180W

this is because of the following code on lines 56-63:
	td = ceil (north / nsr);
	tnorth = td * nsr;
	... etc

when north is 90. (nsr is nsres=0:05 here)


cheap fix is on line 64 add something like
if (window->proj == PROJECTION_LL ) {
  if(tnorth>90)
    north=90;
  if(tsouth<-90)
    tsouth=-90;
  //e,w?
}

but that's not a proper fix, I think a lot of the set.c code assumes 
res>1 and might need to be rewritten.
??

So that's the where & how it's screwing up anyway ... hopefully that's
enough for it to be easily fixed by someone more familiar with the code.



Hamish