Effects of Precision in Viewshed Analysis

I have added high-precision calculation option to my viewshed analysis plugin for QGIS as promised. After spending the last couple of days sick, I am back here at last :)

Adding precision support was easier in Python than I've thought. The hard part was deciding on where to use the Decimal objects exactly. In the end, I have completely separated the precision calculation methods with unprecise ones. The last shared method is 'doVAnalysis'. The paths are then separated as in:


if (self.precision):
    p_visible = self.calculatePrecisionViewshedForPoint(targetX,targetY)
else:
    p_visible = self.calculateViewshedForPoint(targetX,targetY)


The rest of the method chain is like:

doVAnalysis -> calculateViewshedForPoint -> findLineCoeff -> isVisible -> readInterpolatedValue

doVAnalysis -> calculatePrecisionViewshedForPoint -> findPrecisionLineCoeff -> isPrecisionVisible -> readPrecisionInterpolatedValue

Here is an example showing the differences in the code:

Unprecise Arithmetic:

                xvw = 0
                yvw = float(self.dv)
                xpw = xp-xv
                ypw = float(self.readValue(xp,yp))
                c = self.findLineCoeff(xvw, yvw, xpw, ypw)
                aw = c[0]
                bw = c[1]


Precise Arithmetic:

                xvwD = D(0)
                yvwD = self.dvD
                xpwD = D(xp-xv)
                ypwD = D(repr(self.readValue(xp,yp)))
                cD = self.findPrecisionLineCoeff(xvwD, yvwD, xpwD, ypwD)
                awD = cD[0]
                bwD = cD[1]


The class 'D' is an alias for the Decimal class as shown below:


import decimal
D = decimal.Decimal # We use Python's function assignment


Upon reflection, I've realized the strict separation was not necessary and many operations on Decimals were not really any different than float operations. The differences come into play when converting a float value to Decimal. As you can see above, we have to use the 'repr()' output of the float for that.

So, what was the difference when using Decimal class? Well, the answer is in three parts: 1. nothing 2. something 3. too much

1. nothing:

In terms of results, precise arithmetics didn't produce a different viewshed for the limited test cases I've tried.

In the above map, the vector area for precision viewshed and floating-point viewshed for the same settings are overlapped. The red horizontal lines indicate floating-point viewshed area while the green vertical lines indicate precision viewshed area. As you can see they are perfectly identical. Since this is a ON or OFF type of binary result, it is to be expected. The errors in the floating-point calculations just don't accumulate enough to be reflected in the results. But there are real differences, like...

2. something:

When we look at the numerical results themselves, we see the differences. The heart of viewshed calculation is the comparison of a point's (interpolated) height with the height of the passing line-of-sight (L2) at that point. I have written the parameters of all such comparisons for a small scale map in a spreadsheet file and then calculated the differences between Ph (point height) and Lh (line-of-sight height). Here were the results:

Max diff. btw. precise/unprecise Ph:         0.0000487
Max diff. btw. precise/unprecise Lh:         0.0050000
Max diff. btw. precise/unprecise Ph-Lh:    0.0050367

It looks like whenever the height of Ph is close to the Lh within 0.005 meters the precision will start to affect the results. For viewshed analysis this is for many applications way below the sensitivity limit of the input data; so it is not essential. But for other types of operations, like the ones I've listed in the previous post, it may be significant.


3. too much:

The real difference in my test came in the amount of time each version has taken.

In a tiny raster sized 58x52, the times were 1.71 versus 11.42 seconds. This is almost a 7-fold increase..oops!  For a larger map with 332x287 pixels, the times became 106 versus 12395 seconds which is a 117-fold increase!

In the end, it seems like 'precision arithmetic' is not the thing for viewshed analysis and better left for financial applications and such. In an other topic, 106 seconds for a 332x287 map is still too long to be useful in the real world. So, I am planning to implement some faster (but less accurate) methods and make them selectable from the interface.

That's all folks; now back to my literature survey study...

summary: 'this was much ado about nothing...but was fun!'
mood: I feel this is juuuuust the beginning...