The aim is to generate total AGL maps for all conceivable viewpoints in a given area (i.e. DEM) - this is frequently in the order of 10s to 100s of thousands of points. To do this I have written a Python script but I am hitting a problem with the chunk of code that iteratively combines the individual AGL rasters into a single output.
The script takes as input a vector viewpoint layer and a DEM and it seems to work insofar as it loops as it should and generates a final output. However comparison with a manual 'testbed' dataset generated for a small sample of points reveals that the resultant values do not match. Trials with single points worked fine (i.e. for a single point the resultant AGL layer generated by the script was idennitcal to one manually generated) but in the case of two points the output seemed to be the sum of the 2nd AGL raster twice (rather than AGL1 + AGL2). I created a debug version of the code using elif for a small number of points and saving the temporary raster each time with a unique filename and it worked fine - the AGL layers generated at each iteration matched the manual versions as did the final AGL output. However, this is not really practical for 100K+ viewpoints.
If anyone can advise I would be deeply grateful - I am a very amateur python coder (as will become painfully clear & apologies for in advance for any violence done to programming language) but I'm giving it my best shot.
Attached is the python script as intended (hideAGLv8.py) and the debug version for 3 points (using elif) that does actually work for those 3 points insofar as it generates the expected output (hideAGLv4.py).
Cheers - Mark Gillings (UK Archaeologist)
