[GRASS-user] GRASS & earth curvature correction (viewsheds, LOS)

Dear list --

I just recently had the opportunity to revisit
the issue of earth curvature correction in r.los
(and, by extension, r.cva). The curvature of the
Earth's surface means that visibility will be
over-estimated if the line-of-sight algorithm
does not correct for its effect. I initially
implemented a tentative correction in r.cva and
it was with some anxiety that I saw that (untested)
method spread into the source code of r.los.

However, after some testing and comparison with
the output of ArcGIS' Viewshed module, I can now
confirm that the correction method works just fine.

So that's the good news. The remaining itch is that
atmospheric refraction is not corrected for in any
GRASS line-of-sight based module. That effect is
smaller, but it works against the effect of earth
curvature. This means that a module which corrects
for curvature, but not refraction, will tend to
over-compensate (by about 6/7).

I have not tested r.viewshed's curvature correction.
But if it uses the same method as r.los, then the
same observations apply.

As an alternative to all this mess, I have attached
a little script which allows the user to pre-process
the digital elevation model, in effect producing
a "fake" planar model that "fools" LOS algorithms into
operating in a plane instead of on a curved surface.
This little trick works as long as the region is
not more than a few hundred kilometers across. Kudos
to Bill Huber of Quantitative Decisions, who made the
details of the computations available via the user
forums of a well-known GIS producer. More details
can be found in the HTML man page included in the
attached ZIP file.

Pre-processing the elevation data and then turning
OFF any built-in correction methods of any LOS
module used will make the computations transparent
and more comparable.

Btw.: When comparing the output of ArcGIS' Viewshed
with that of r.los (or any other LOS modules I have
tried, including some from SEXTANTE), you will notice
that the ArcGIS viewsheds are always more contiguous
and less "noisy". That software seems to use some
internal tolerance or smoothing to achieve this effect.
I am not sure it's a desirable thing to emulate, though.
Better smooth the result manually and in a controlled
manner using e.g. a median filter, if desired.

------
Files attached to this email may be in ISO 26300 format (OASIS Open Document Format). If you have difficulty opening them, please visit http://iso26300.info for more information.

(attachments)

r.ecurv.zip (4.1 KB)

curvature. This means that a module which corrects
for curvature, but not refraction, will tend to
over-compensate (by about 6/7).

That's 1/7, of course :wink:

Ben

------
Files attached to this email may be in ISO 26300 format (OASIS Open Document Format). If you have difficulty opening them, please visit http://iso26300.info for more information.

Benjamin wrote:

I just recently had the opportunity to revisit
the issue of earth curvature correction in r.los
(and, by extension, r.cva). The curvature of the
Earth's surface means that visibility will be
over-estimated if the line-of-sight algorithm
does not correct for its effect. I initially
implemented a tentative correction in r.cva and
it was with some anxiety that I saw that (untested)
method spread into the source code of r.los.

However, after some testing and comparison with
the output of ArcGIS' Viewshed module, I can now
confirm that the correction method works just fine.

It gives confidence sure, but there's no way of knowing if their
black box is actually correct or not, maybe we make the same
mistakes? :slight_smile: All the better to follow/cite published articles.

So that's the good news. The remaining itch is that
atmospheric refraction is not corrected for in any
GRASS line-of-sight based module. That effect is
smaller, but it works against the effect of earth
curvature. This means that a module which corrects
for curvature, but not refraction, will tend to
over-compensate (by about 6/7).

I have not tested r.viewshed's curvature correction.
But if it uses the same method as r.los, then the
same observations apply.

I'm not sure, but I suspect it might. MarkusM?

As an alternative to all this mess, I have attached
a little script which allows the user to pre-process
the digital elevation model, in effect producing
a "fake" planar model that "fools" LOS algorithms into
operating in a plane instead of on a curved surface.
This little trick works as long as the region is
not more than a few hundred kilometers across. Kudos
to Bill Huber of Quantitative Decisions, who made the
details of the computations available via the user
forums of a well-known GIS producer. More details
can be found in the HTML man page included in the
attached ZIP file.

Pre-processing the elevation data and then turning
OFF any built-in correction methods of any LOS
module used will make the computations transparent
and more comparable.

certainly the ideal solution is to have both curvature and
refection corrections implemented as flags in the new & improved
r.viewshed. (say, is that now ready for moving into grass7? *)

Does refraction only work in the vertical or would you be able
to slightly see around some distant volcano horizontally?
..perhaps the phenomenon is related to the atm pressure
gradient, in which case purely a vertical-only effect..?

[*] @MarkusM: re. changing out-of-view to NULL in r.viewshed, I
now notice in the help page that NULL was previously reserved for
NULLs in the original DEM. But in that case I suspect that
everything in the shadow of the DEM's NULLs should also be NULL,
as the visibility will be unknowable. I'm mildly in favour of
keeping the code as-is, and changing the help page to match.

Btw.: When comparing the output of ArcGIS' Viewshed
with that of r.los (or any other LOS modules I have
tried, including some from SEXTANTE), you will notice
that the ArcGIS viewsheds are always more contiguous
and less "noisy". That software seems to use some
internal tolerance or smoothing to achieve this effect.

Do you think that noise happens because of a translation of the
coarse horizontal grid cell size when it is translated into the
vertical plane?

I am not sure it's a desirable thing to emulate, though.
Better smooth the result manually and in a controlled
manner using e.g. a median filter, if desired.

absolutely. If it's a common task maybe offer a flag to run
that automatically, but certainly don't mess with the results
by default.

thanks,
Hamish

On Wed, May 25, 2011 at 12:34 PM, Hamish <hamish_b@yahoo.com> wrote:

Benjamin wrote:

[...]

I have not tested r.viewshed's curvature correction.
But if it uses the same method as r.los, then the
same observations apply.

I'm not sure, but I suspect it might. MarkusM?

r.viewshed's correction for earth curvature is copied as is from (recent) r.los.

[...]

[*] @MarkusM: re. changing out-of-view to NULL in r.viewshed, I
now notice in the help page that NULL was previously reserved for
NULLs in the original DEM. But in that case I suspect that
everything in the shadow of the DEM's NULLs should also be NULL,
as the visibility will be unknowable. I'm mildly in favour of
keeping the code as-is, and changing the help page to match.

The original behaviour of r.viewshed was to distinguish between
visible and invisible cells and cells for which visibility could not
be determined because input elevation is NULL. r.los on the contrary
sets both invisible cells and cells for which visibility could not be
determined to NULL. I synchronized r.viewshed to r.los because r.los
is in the main branches whereas r.viewshed is in addons.
Distinguishing between invisible cells and cells for which visibility
could not be determined might provide useful information, though.

Markus M

Markus Metz wrote:

The original behaviour of r.viewshed was to distinguish
between visible and invisible cells and cells for which
visibility could not be determined because input elevation is
NULL. r.los on the contrary sets both invisible cells and cells
for which visibility could not be determined to NULL. I
synchronized r.viewshed to r.los because r.los is in the main
branches whereas r.viewshed is in addons.
Distinguishing between invisible cells and cells for which
visibility could not be determined might provide useful
information, though.

well, AFAIR r.los uses NULL for that simply because that's what I
decided to do when I added NULL support to the module some years
ago.
So to distinguish between uncalculated (because beyond the
requested search radius?), and in-search-region but non-visible?

Well, it's not to me, but I guess potentially that's useful to
someone somewhere. I hate to needlessly lose information, but
in general I think NULL for non-visible is what you'd want the
output to show most of the time. So again I'm mildly in favour
of keeping out-of-sight = NULL, but not of too strong an opinion
about it.

Hamish