OSGeo
Viewshed calculations with curvature coefficient are incorrect when SRS uses different units for geo transform and reference ellipsiod
Viewshed calculations do not correctly handle curvature coefficient when SRS defines reference ellipsoid in different units from the geo transform.
Expected behavior and actual behavior.
Would expect this to raise a warning / error that the results may not be valid. At the very least, there should be a warning in the viewshed function documentation.
Current documentation includes the warning "The algorithm as implemented currently will only output meaningful results if the georeferencing is in a projected coordinate reference system." However, this does not clearly address the issue described here. A warning "The curvature coefficient will only provide meaningful results if the SRS defines the reference ellipsoid in the same units used by the geo transform" could be added to the section on the curvature coefficient.
Steps to reproduce the problem.
Using a file in epsg:4326 calculate the viewshed at any point. Change the curvature coefficient from the default value of 1.0 to 0 (-cc 0) and re-run gdal_viewshed. Comparing the two viewsheds will show no difference between the files, where there should be a substantial difference between a cc of 0 (flat) and 1 (normal earth curvature).
This arises since the SRS defines the geo transform in degrees and the reference ellipsoid in meters. This causes the calculations to treat the earth's diameter as many times larger than it is, and the curvature coefficient to be ignored for all practical purposes.
Operating system
RHEL 9.3
GDAL version and provenance
3.8.0, built from source
@sa2329 Did you find a work around for this? A type to convert to before processing?
@temas Currently I'm working around this by projecting into a different CRS. Using UTM so units are meters, which matches ellipsoid units and seems to give correct results.
If I have time later, may experiment with leaving the data in epsg:4326 and changing the geotransform on the data to be in meters rather than degrees. Will add a comment here with results if I do.