[GRASS-dev] d.vect render= and corrupted vector display

Hi

I have tried different render methods with a vector box which fits the
region perfectly. For each render method the vector box is displayed
wrong on the X mointor. Details:

in spearfish:

$ g.region rast=elevation.10m -a
$ v.in.region output=frame type=line
$ d.mon x0

#1
$ d.vect frame rander=g
(right egde of the vector box invisible)

#2
$ d.erase; d.vect frame render=r
(right egde of the vector box invisible)

#3
$ d.erase; d.vect frame render=d
(right egde of the vector box invisible)

#4
$ d.erase; d.vect frame render=c
(nothing is visible)

#5
$ d.erase; d.vect frame render=c
(nothing is visible)

Then I repeated tests with width=2 and greater.

For cases 1-3 the right edge became visible, but the right and left
edges were thinner than top and bottom.

For case 4, no matter what width nothing is visible.

For case 5, width > 1 makes the right edge visible, but the other 3 are
never displayed.

Maciek

Maciej Sieczka wrote:

#5
$ d.erase; d.vect frame render=c

Typo. Should read: render=l.

(nothing is visible)

Maciek

Maciej Sieczka wrote:

I have tried different render methods with a vector box which fits the
region perfectly. For each render method the vector box is displayed
wrong on the X mointor. Details:

in spearfish:

$ g.region rast=elevation.10m -a
$ v.in.region output=frame type=line
$ d.mon x0

#1
$ d.vect frame rander=g
(right egde of the vector box invisible)

#2
$ d.erase; d.vect frame render=r
(right egde of the vector box invisible)

#3
$ d.erase; d.vect frame render=d
(right egde of the vector box invisible)

#4
$ d.erase; d.vect frame render=c
(nothing is visible)

#5
$ d.erase; d.vect frame render=c
(nothing is visible)

Then I repeated tests with width=2 and greater.

For cases 1-3 the right edge became visible, but the right and left
edges were thinner than top and bottom.

For case 4, no matter what width nothing is visible.

For case 5, width > 1 makes the right edge visible, but the other 3 are
never displayed.

I've changed the clipping/culling code to handle this particular case
(i.e. use <= 0 rather than < 0), so #4 and #5 now match the other
three.

Beyond that, the main issue is that the code which sets up the
coordinate mapping (D_setup/D_do_conversions) maps the current region
to the current frame *exactly*, without any margin. This is arguably
the right thing, even if the net result isn't always ideal (and the
result of v.in.region is the worst possible case).

This doesn't allow for thick lines, where the path itself (i.e. the
centre-line) might lie just outside the region but the "stroke" can
still be visible. It also doesn't allow for the fact that the tiniest
rounding error might push the path outside.

Note that the clipping/culling code doesn't make any allowance for the
line width. It can't; it doesn't know the line width (R_line_width()
does directly to the driver; the display library doesn't get to see
it).

I should probably add an option to set a "margin", i.e. the difference
between the clip window used for the actual rendering and the clip
window used for geometric clipping and culling.

There's also the issue that the actual line may vary by +/- half a
pixel from the "correct" position. That's an inevitable fact of
rasterisation; the driver can't fill half of a pixel. If the path lies
exactly along an edge of the "screen" (window/image), the line will
end up either just on-screen or just off-screen.

FWIW, if you try this with the PS driver and width=1 (which translates
to one point, not one device pixel), and render the PostScript at a
high resolution, you'll note that it really does draw half of each of
the two vertical lines (adding "[10] 0 setdash 1 setlinecap" makes it
more clear).

--
Glynn Clements <glynn@gclements.plus.com>

Glynn

Thank you for looking into this. Please read below.

Glynn Clements wrote:

I've changed the clipping/culling code to handle this particular case
(i.e. use <= 0 rather than < 0), so #4 and #5 now match the other
three.

Beyond that, the main issue is that the code which sets up the
coordinate mapping (D_setup/D_do_conversions) maps the current region
to the current frame *exactly*, without any margin. This is arguably
the right thing, even if the net result isn't always ideal (and the
result of v.in.region is the worst possible case).

Don't you think it is confussing for the user? If I do v.in.region, I'm
right to expect that the whole vector will be displayed when my working
region matches the v.in.region's output extent. And if this doesn't
happen (as it is currently), I start wondering whether there could be a
bug in v.in.region, or g.region, or d.vect an so on. Waste of time,
less trust in GRASS, bogus bug reports etc.

If it is possible to fix that without too much effort, I would be for it.

This doesn't allow for thick lines, where the path itself (i.e. the
centre-line) might lie just outside the region but the "stroke" can
still be visible. It also doesn't allow for the fact that the tiniest
rounding error might push the path outside.

Note that the clipping/culling code doesn't make any allowance for the
line width. It can't; it doesn't know the line width (R_line_width()
does directly to the driver; the display library doesn't get to see
it).

Propably I did not understand something, but as I can see render=c/l
work well together with width= in d.vect. Did you mean something else?

Maciek

P.S.

If this is usefull, Hamish made some guess about a possible nature of
the problem about a year ago; please read on
http://intevation.de/rt/webrt?serial_num=3613, message dated: Thu, Sep
8 2005 02:48:20.

Maciej Sieczka wrote:

> I've changed the clipping/culling code to handle this particular case
> (i.e. use <= 0 rather than < 0), so #4 and #5 now match the other
> three.

> Beyond that, the main issue is that the code which sets up the
> coordinate mapping (D_setup/D_do_conversions) maps the current region
> to the current frame *exactly*, without any margin. This is arguably
> the right thing, even if the net result isn't always ideal (and the
> result of v.in.region is the worst possible case).

Don't you think it is confussing for the user? If I do v.in.region, I'm
right to expect that the whole vector will be displayed when my working
region matches the v.in.region's output extent.

Hmm; not really. v.in.region generates a vector map corresponding to
the region's boundary. d.vect displays the portion of a vector map
which lies inside the region.

Is a region's boundary *inside* the region?

If you treat the boundary as an area, then the area which d.vect will
fill is exactly the region, i.e. inside the boundary == inside the
region. But is the boundary itself inside the region?

I'm aware that this may be counter-intuitive, but that isn't the same
thing as being wrong. The problem is that I don't see how to change
this without creating more significant problems.

And if this doesn't
happen (as it is currently), I start wondering whether there could be a
bug in v.in.region, or g.region, or d.vect an so on. Waste of time,
less trust in GRASS, bogus bug reports etc.

If it is possible to fix that without too much effort, I would be for it.

How do you propose "fixing" this case without breaking more sensible
cases. Bear in mind that this is (quite literally) a "boundary case".

> This doesn't allow for thick lines, where the path itself (i.e. the
> centre-line) might lie just outside the region but the "stroke" can
> still be visible. It also doesn't allow for the fact that the tiniest
> rounding error might push the path outside.
>
> Note that the clipping/culling code doesn't make any allowance for the
> line width. It can't; it doesn't know the line width (R_line_width()
> does directly to the driver; the display library doesn't get to see
> it).

Propably I did not understand something, but as I can see render=c/l
work well together with width= in d.vect. Did you mean something else?

They only work because:

a) I tweaked the clip/cull tests slightly (changed <0 to <=0), and

b) the coordinates of the region's boundary are integers, so
converting from cartographic coordinates to display coordinates
doesn't suffer from rounding error.

Had the right-hand edge been given a display X coordinate of
640.000001, render=c/l would have discarded the edge, as it lies
entirely outside the region.

The fact that part of the "stroke" (the filled area centred upon the
actual line) would still be visible isn't taken into account when
clipping/culling, as the display library doesn't know about the line
width.

That's why I say that there needs to be a D_line_width() or D_margin()
or similar, so that the clipping/culling code won't discard edges
which lie only just outside the region.

--
Glynn Clements <glynn@gclements.plus.com>

Glynn Clements wrote:

Maciej Sieczka wrote:

I've changed the clipping/culling code to handle this particular case
(i.e. use <= 0 rather than < 0), so #4 and #5 now match the other
three.
Beyond that, the main issue is that the code which sets up the
coordinate mapping (D_setup/D_do_conversions) maps the current region
to the current frame *exactly*, without any margin. This is arguably
the right thing, even if the net result isn't always ideal (and the
result of v.in.region is the worst possible case).

Don't you think it is confussing for the user? If I do v.in.region, I'm
right to expect that the whole vector will be displayed when my working
region matches the v.in.region's output extent.

Hmm; not really. v.in.region generates a vector map corresponding to
the region's boundary. d.vect displays the portion of a vector map
which lies inside the region.

Is a region's boundary *inside* the region?

1. This is not a boundary, but a line (v.in.region output=frame type=line).

2. Even if it was a boundary - either the *whole* boundary should be
displayed, or not. Currently, 75% of the boundary is displayed, while
the other 25% is not. This looks strange to the user.

Too bad if cannot be "fixed", though I'm 100% OK with whatever you
decide about it and thank you for taking your time to check it out and
explain it.

If you treat the boundary as an area, then the area which d.vect will
fill is exactly the region, i.e. inside the boundary == inside the
region. But is the boundary itself inside the region?

I'm aware that this may be counter-intuitive, but that isn't the same
thing as being wrong. The problem is that I don't see how to change
this without creating more significant problems.

And if this doesn't
happen (as it is currently), I start wondering whether there could be a
bug in v.in.region, or g.region, or d.vect an so on. Waste of time,
less trust in GRASS, bogus bug reports etc.

If it is possible to fix that without too much effort, I would be for it.

How do you propose "fixing" this case

I can't know that. I thought you might, so I asked.

without breaking more sensible cases.

What could be broken?

Bear in mind that this is
(quite literally) a "boundary case".

I don't agree. This issue might bias user's understanding of GRASS'
region, which is not trivial:

g.region vect=frame; d.vect frame

In a result, 25% of lines is missing on display, while region and
display were ordered to encompass the whole vector "frame". As the 3/4
of lines is displayed while 1/4 is not, it just looks "not good" to the
user, thus might incline him to think there is a defect somewhere
underneath.

Please note that if you decide this cannot be fixed/too intrusive to do
so/to much work compared to profit, please be sure I'm perfectly OK
with your decision. Yet, the issue would remain.

Maciek

Maciej Sieczka wrote:

>>> I've changed the clipping/culling code to handle this particular case
>>> (i.e. use <= 0 rather than < 0), so #4 and #5 now match the other
>>> three.
>>> Beyond that, the main issue is that the code which sets up the
>>> coordinate mapping (D_setup/D_do_conversions) maps the current region
>>> to the current frame *exactly*, without any margin. This is arguably
>>> the right thing, even if the net result isn't always ideal (and the
>>> result of v.in.region is the worst possible case).
>> Don't you think it is confussing for the user? If I do v.in.region, I'm
>> right to expect that the whole vector will be displayed when my working
>> region matches the v.in.region's output extent.

> Hmm; not really. v.in.region generates a vector map corresponding to
> the region's boundary. d.vect displays the portion of a vector map
> which lies inside the region.
>
> Is a region's boundary *inside* the region?

1. This is not a boundary, but a line (v.in.region output=frame type=line).

Note that I mean "boundary" in the general sense, not specifically
GV_BOUNDARY.

2. Even if it was a boundary - either the *whole* boundary should be
displayed, or not. Currently, 75% of the boundary is displayed, while
the other 25% is not. This looks strange to the user.

Typically, two of the region edges will be coincident with the
corresponding edges of the display frame, while the other two won't
(unless the region's aspect ratio exactly matches that of the display
frame).

For a vertical or horizontal line segment which is concident with an
edge of the display frame, it's essentially a 50/50 chance as to
whether that edge will lie inside or outside the frame; even the
slightest rounding error can decide it. Note that it is not guaranteed
that (x/y)*y == x when using floating point.

If you go out of your way to create lines which lie exactly on a
boundary, the results are bound to be unpredictable. Unfortunately,
that's exactly what v.in.region does.

Ultimately, the only solution is to stick to even line widths (i.e.
width=2, width=4 etc), or use the PS driver. Then you're guaranteed to
get exactly half of a line on each side (if you don't, that *is* a
bug).

> How do you propose "fixing" this case

I can't know that. I thought you might, so I asked.

> without breaking more sensible cases.

What could be broken?

Well, that depends upon what you change. AFAICT, changing anything
which is likely to affect the result of v.in.region type=line + d.vect
amounts to introducing an error.

One solution is to expand the region used by D_do_conversions()
(either pass in an expanded region or have D_do_conversions() expand
it) so that data which lies either within the current region or
exactly on its bounds always lies comfortably within the expanded
region.

But expanded by how much? A metre, a pixel, 1%, ...? It would have to
be consistent if you want to overlay maps.

Also, that would cause problems for any external code which attempts
to reproduce the conversion; ISTR that gis.m does this to convert
mouse coordinates to geographic coordinates.

> Bear in mind that this is
> (quite literally) a "boundary case".

I don't agree. This issue might bias user's understanding of GRASS'
region, which is not trivial:

That doesn't have any effect upon whether this is a boundary case.

g.region vect=frame; d.vect frame

In a result, 25% of lines is missing on display, while region and
display were ordered to encompass the whole vector "frame". As the 3/4
of lines is displayed while 1/4 is not, it just looks "not good" to the
user, thus might incline him to think there is a defect somewhere
underneath.

Please note that if you decide this cannot be fixed/too intrusive to do
so/to much work compared to profit, please be sure I'm perfectly OK
with your decision. Yet, the issue would remain.

If it was critical, it could be trivially fixed by prohibiting odd
line widths in the X and PNG drivers. Unfortunately, a line width of 2
may be much slower than a line width of 1 for the X driver, as it uses
a completely different (and more complex) algorithm.

But this is just a single (albeit rather striking) example of an
inevitable problem with raster graphics: aliasing. Even if you fix
this specific case by means of some fudge, there are an unlimited
number of similar cases.

--
Glynn Clements <glynn@gclements.plus.com>

Glynn Clements wrote:

Maciej Sieczka wrote:

I've changed the clipping/culling code to handle this particular case
(i.e. use <= 0 rather than < 0), so #4 and #5 now match the other
three.
Beyond that, the main issue is that the code which sets up the
coordinate mapping (D_setup/D_do_conversions) maps the current region
to the current frame *exactly*, without any margin. This is arguably
the right thing, even if the net result isn't always ideal (and the
result of v.in.region is the worst possible case).

Don't you think it is confussing for the user? If I do v.in.region, I'm
right to expect that the whole vector will be displayed when my working
region matches the v.in.region's output extent.

Hmm; not really. v.in.region generates a vector map corresponding to
the region's boundary. d.vect displays the portion of a vector map
which lies inside the region.

Is a region's boundary *inside* the region?

1. This is not a boundary, but a line (v.in.region output=frame type=line).

Note that I mean "boundary" in the general sense, not specifically
GV_BOUNDARY.

2. Even if it was a boundary - either the *whole* boundary should be
displayed, or not. Currently, 75% of the boundary is displayed, while
the other 25% is not. This looks strange to the user.

Typically, two of the region edges will be coincident with the
corresponding edges of the display frame, while the other two won't

Do you mean that either both horizontal or both vertical lines should
be not visible on display? I don't confirm. Alway a single line is
missing - either the bottom, or the right side edge, depending on the
with X monitor windows's aspect ratio. (Talking of d.vet width=1 here.)

To reproduce:

1. v.in.region output=frame type=line
2. d.vect frame
3. play with X monitor windows's aspect ratio

(unless the region's aspect ratio exactly matches that of the display
frame).

If you mean that in such case all edges should be displayed, I don't
confirm. To reproduce:

1. g.region n=4926750 s=4914430 w=594880 e=607200 res=10 -ag
(results is a square region)

2. d.mon x0; d.resize width=320 height=320
(square X display)

3. v.in.region square type=line; d.vect square
(only left and top edge are displayed)

For a vertical or horizontal line segment which is concident with an
edge of the display frame, it's essentially a 50/50 chance as to
whether that edge will lie inside or outside the frame; even the
slightest rounding error can decide it. Note that it is not guaranteed
that (x/y)*y == x when using floating point.

If you go out of your way to create lines which lie exactly on a
boundary, the results are bound to be unpredictable. Unfortunately,
that's exactly what v.in.region does.

Ultimately, the only solution is to stick to even line widths (i.e.
width=2, width=4 etc), or use the PS driver. Then you're guaranteed to
get exactly half of a line on each side (if you don't, that *is* a
bug).

Hmm. With line width=1 I get the effect as described above. With width
2 or 4 I get the effect that either both horizontal or both vertical
edges (depending on the with X monitor windows's aspect ratio) are
twice that thick as the other pair of lines. If it is not clear what I
mean, please look at the attached screenshot of d.vect frame width=2.
Do you mean this is a bug?

Bear in mind that this is
(quite literally) a "boundary case".

I don't agree. This issue might bias user's understanding of GRASS'
region, which is not trivial:

That doesn't have any effect upon whether this is a boundary case.

By "boundary case", do you mean to assume that v.in.region+d.vect is
unlikely to be used? Just curios.

Maciek

(attachments)

horizontal_edges_thinner.png
vertical_edges_thinner.png

Maciej Sieczka wrote:

> Typically, two of the region edges will be coincident with the
> corresponding edges of the display frame, while the other two won't

Do you mean that either both horizontal or both vertical lines should
be not visible on display?

No. I mean that one pair (both horizontal or both vertical) will
normally be visible, while the other two should each have a 50/50
chance of being visible.

In the latter case, various factors may combine to affect the actual
ratio. E.g. the fact that the transformation is rooted at the
north-west corner will typically mean that the top/left edges always
result in an exact zero, producing the same result every time.

As with anything related to floating-point, the architecture and any
compilation options may affect the outcome.

> For a vertical or horizontal line segment which is concident with an
> edge of the display frame, it's essentially a 50/50 chance as to
> whether that edge will lie inside or outside the frame; even the
> slightest rounding error can decide it. Note that it is not guaranteed
> that (x/y)*y == x when using floating point.
>
> If you go out of your way to create lines which lie exactly on a
> boundary, the results are bound to be unpredictable. Unfortunately,
> that's exactly what v.in.region does.
>
> Ultimately, the only solution is to stick to even line widths (i.e.
> width=2, width=4 etc), or use the PS driver. Then you're guaranteed to
> get exactly half of a line on each side (if you don't, that *is* a
> bug).

Hmm. With line width=1 I get the effect as described above. With width
2 or 4 I get the effect that either both horizontal or both vertical
edges (depending on the with X monitor windows's aspect ratio) are
twice that thick as the other pair of lines.

This is correct behaviour. One pair of region edges will always
exactly touch the corresponding edges of the frame, while the other
pair will usually have a margin.

For an edge which exactly touches the frame, any line which follows
that edge will get clipped down the middle.

If the aspect ratio of the region exactly matches that of the frame,
all four edges will coincide, and all four lines should be half-drawn.

If it is not clear what I
mean, please look at the attached screenshot of d.vect frame width=2.
Do you mean this is a bug?

No, I mean that anything else would be a bug.

>>> Bear in mind that this is
>>> (quite literally) a "boundary case".

>> I don't agree. This issue might bias user's understanding of GRASS'
>> region, which is not trivial:

> That doesn't have any effect upon whether this is a boundary case.

By "boundary case", do you mean to assume that v.in.region+d.vect is
unlikely to be used? Just curios.

No, I mean that it's a case where a "point" is neither inside nor
outside but on the boundary.

The term doesn't refer solely to geometry. In any context where some
continuous (non-discrete) space is partitioned into distinct
categories, you get values which fall on the boundary between
categories. E.g. if you partition the real numbers into positive and
negative, which is zero?

--
Glynn Clements <glynn@gclements.plus.com>

Glynn Clements wrote:

Maciej Sieczka wrote:

By "boundary case", do you mean to assume that v.in.region+d.vect is
unlikely to be used? Just curios.

No, I mean that it's a case where a "point" is neither inside nor
outside but on the boundary.

The term doesn't refer solely to geometry. In any context where some
continuous (non-discrete) space is partitioned into distinct
categories, you get values which fall on the boundary between
categories. E.g. if you partition the real numbers into positive and
negative, which is zero?

As to display, I guess I get it know. However, I'm not sure what to
think about this issue in regard to computational GRASS region.

Let's take a following eaxample; in spearfish:

$ g.region rast=landuse -a
$ v.in.region lu_border type=line
$ v.category lu_border out=lu_border_cat option=add
$ g.region vect=lu_border_cat
$ v.to.rast input=lu_border_cat output=lu_border_cat use=cat

The outcome is a raster which includes only the top and left lines,
while region was ordered to cover the whole "lu_border_cat" vector.

From user's point of view, half of the input vector was ignored by

v.to.rast, in error. Is this not fixable too?

Maciek

Maciej Sieczka wrote:

>> By "boundary case", do you mean to assume that v.in.region+d.vect is
>> unlikely to be used? Just curios.

> No, I mean that it's a case where a "point" is neither inside nor
> outside but on the boundary.
>
> The term doesn't refer solely to geometry. In any context where some
> continuous (non-discrete) space is partitioned into distinct
> categories, you get values which fall on the boundary between
> categories. E.g. if you partition the real numbers into positive and
> negative, which is zero?

As to display, I guess I get it know. However, I'm not sure what to
think about this issue in regard to computational GRASS region.

Let's take a following eaxample; in spearfish:

$ g.region rast=landuse -a
$ v.in.region lu_border type=line
$ v.category lu_border out=lu_border_cat option=add
$ g.region vect=lu_border_cat
$ v.to.rast input=lu_border_cat output=lu_border_cat use=cat

The outcome is a raster which includes only the top and left lines,
while region was ordered to cover the whole "lu_border_cat" vector.
From user's point of view, half of the input vector was ignored by
v.to.rast, in error. Is this not fixable too?

No. Exactly the same issues apply here.

--
Glynn Clements <glynn@gclements.plus.com>