2010-04-20

Visibility Analysis: one way of doing it

It has been 18 days since I've started working on my Visibility Analysis Plugin for QGIS mini-project. Between my day job and many other activities I am more or less satisfied with my progress.


Currently VAP is able to calculate the viewshed for any point on the terrain with a relatively high precision. Of course, there are many improvements waiting to be made, but in this post I will just focus on the mathematical method I used to do that analysis.

An example analysis for visibility. Which areas how often are seen while kayaking through the river (or in other words, which areas see the most of the river)

I had some initial ideas about how to do VA using line formulas and some computer graphics techniques before starting. Some research in academic circles revealed that the same method was the standard approach used everywhere, to my surprise. This method, which I will call as "Discrete Interpolated Brute Force" or DIBF, also gives out the most accurate results by using the full information in the DEM RSG data.

I will try to describe the method very quickly for the interested:

1. Ultra Quick Introduction

As you can see on the left there are two types of elevation data:
  • RSGs: Regular Square Grids
  • TINs: Triangulated Irregular Networks
VAP can only calculate viewsheds for the RSG data (at least for now) by calculating the visibility of every point on terrain (Pi) from a known selected viewpoint (V)

2. XY Line and Z Line Equations

From basic euclidean line equation y=ax+b  and two points on that line, we can find the coefficients a and b.

There are two such lines we have to calculate.

1. L1: V-Pi line on the X-Y coordinate system
2. L2: V-Pi line on the L1-Z coordinate system
3. Checking for Visibility

Now we trace the L1 for every x starting from V point. We then interpolate the height values from the two squares that our y lies in between using a proximity weight.

For every intermediate I point, before going for the next, we check if the I point blocks the view. To do this, we just check if the height value of the I point is higher or lower than the z value of the L2 line at that x coordinate.

If none of the I points are blocking, then Pi is visible from V. Then we move on the next Pi point and repeat everything.

one more tip: if the incline of L1 leans to y axis as opposed to x axis, there will be some loss of resolution and accuracy in the analysis. To avoid this, just swap x and y values and we're done! :)  All the other code can remain the same except when reading height values from a point where original coordinates must be used obviously (well maybe not so obvious, as I've made that mistake and created a nasty bug hard to follow)

For swapping values, you can use Python's multi-variable assignment which evaluates right-hand side first, as in: x,y = y,x  . In fact, I'm in love with that feature :)

One disadvantage of DIBF algorithm is that it is as slow as a turtle. For every candidate point, we calculate two lines and operate on every single intermediate point as well as their immediate neighbors. I will return to that matter in a later post.


MadChuckle

summary: 'there are sometimes just one way to do it right.....and that one way is usually slow as hell !'
mood: serene

No comments:

Post a Comment