7

I want to write an algorithm that performs viewshed analysis in Python programming language using digital elevation model. I generated points from pixels from the raster data in QGIS. I calculated the latitude and longitude of these points. I determined one of the points as the observer point.

I exported the points as CSV format. I defined the CSV file using Pandas and GeoPandas libraries in Jupyter notebook. Then I created a buffer on the raster data. I got the remaining points in the buffer.

I'm having trouble doing the rest of the steps. How can I write an algorithm for viewshed or is there is an easier way to do viewshed analysis?

import gdal
import pandas as pd
import geopandas as gpd
import matplotlib.pyplot as plt
from rasterio.plot import show
from shapely.geometry import Point

ds = pd.read_csv("aster_latlon.csv")
geometry = [Point(xy) for xy in zip(ds.longitude, ds.latitude)]
crs = {'init': 'epsg:32620'}
gdf = gpd.GeoDataFrame(ds, geometry = geometry, crs = crs)

obs = pd.read_csv("observerpoint.csv")
geometry = [Point(xy) for xy in zip(obs.longitude, obs.latitude)]
crs = {'init': 'epsg:32620'}
obsGdf = gpd.GeoDataFrame(obs, geometry = geometry, crs = crs)

buffer = gdf.geometry[72477].buffer(0.08)
selected = gdf.intersects(buffer)
ds["Selected"] = selected 

x,y = buffer.exterior.xy 
l = x,y

lat = list(l[1])
long = list(l[0])

1 Answer 1

10

The easiest way to calculate a visibility raster in python is to call the ad hoc functions of a GIS from your code. This can be done in several ways (non-exhaustive list):

  1. Using the r.viewshed module from grass
  2. Using the visibility module of saga GIS
  3. Using the gdal_viewshed function from gdal.

Depending on the tool chosen, the method of calling will be slightly different. You should therefore look at the api documentation to know how to do this in detail. Here are some links that may help to do this:

An other solution is to write directly an algorithm of visibility. The easiest way to do this is to implement the algorithms that have been used in the tools included in the GIS. For example r.viewshed use an algorithm write by Haverkort, Toma and Zhuan, presented in the documentation and in this research paper :

Haverkort, H.,Toma L. and Zhuang Y. Computing Visibility on Terrains in External Memory. In the Proceedings of the 9th Workshop on Algorithm Engineering and Experiments / Workshop on Analytic Algorithms and Combinatoric. http://citeseerx.ist.psu.edu/viewdoc/download?doi=10.1.1.76.4282&rep=rep1&type=pdf.

gdal_viewshed use a Wang, Robinson and White algorithm, presented in this research paper :

Wang, Jianjun, Robinson, Gary J., and White, Kevin, Generating Viewsheds without Using Sightlines. Photogrammetric Engineering and Remote Sensing. https://www.asprs.org/wp-content/uploads/pers/2000journal/january/2000_jan_87-90.pdf

You can also use the Pixscape algorithm, prensented here:

Sahraoui, Y, Vuidel, G, Joly, D, Foltête, J-C. Integrated GIS software for computing landscape visibility metrics. Transactions in GIS. 2018; 22: 1309– 1323. https://doi.org/10.1111/tgis.12457

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.