r.horizon generates incorrect results in single point mode

[rmyorston]

I wished to use r.horizon to calculate horizon values in single point mode. Initially I used a 50m DEM and the results looked reasonable. Later I obtained a 10m DEM for the region and performed the same calculation. The results were very different. When I used map mode with the 10m DEM the result for the chosen location and angle agreed with that of the 50m DEM.

The problem seems to be in test_low_res, which is intended to switch to a low resolution grid at a certain point in the calculation. I've built r.horizon without the call to test_low_res (just setting succes2 to 1) and the results with the 10m DEM agree pretty well with those for the 50m one.

In the code with the comment /*skip to the next lowres cell */ delta values for x and y are calculated. These depend on cosangle and sinangle. It appears that these aren't being set in single point mode, so both are zero. This results in mindel being 32000 and hence the point calculated is off the grid. In my case this stopped the horizon calculation after 2200m, far too early.

I think there are other problems with test_low_res, but I didn't investigate further.

[neteler]

See also (unrelated?) ticket #1575

[annakrat]

A test case for UTM32 location, on artificial surface:

g.region n=5118356 s=5088351 e=666798 w=636793 res=10
r.mapcalc --overwrite expression="test_surface = 100*sin( 0.2 *row()) + 100*sin( 0.2 *col())"
r.horizon elevation=test_surface direction=0   step=1 coordinates=656096,5102603  -d fi=horizon.csv maxd=800 --o

import numpy as np
import matplotlib.pyplot as plt

my_data = np.genfromtxt('horizon.csv', delimiter=',')
my_data = my_data[1:, :]
ax = plt.subplot(111, polar=True)
bars = ax.plot(my_data[:, 0] / 180 * np.pi, ( my_data[:, 1]) / 180 * np.pi)
plt.show()