Adapted elevation map plugin to new GRASS engine format

pull/829/head
Piero Toffanin 2020-03-11 22:31:37 -04:00
rodzic b542baca16
commit d6f0f47a0a
3 zmienionych plików z 97 dodań i 78 usunięć

Wyświetl plik

@ -37,21 +37,18 @@ class TaskElevationMapGenerate(TaskView):
noise_filter_size = float(request.data.get('noise_filter_size', 2))
current_dir = os.path.dirname(os.path.abspath(__file__))
context.add_param('dsm_file', dsm)
context.add_param('dsm', dsm)
context.add_param('interval', interval)
context.add_param('format', format)
context.add_param('noise_filter_size', noise_filter_size)
context.add_param('epsg', epsg)
context.add_param('python_script_path', os.path.join(current_dir, "elevationmap.py"))
context.add_param('python_path', plugin.get_python_packages_path())
if dtm != None:
context.add_param('dtm', '--dtm {}'.format(dtm))
else:
context.add_param('dtm', '')
context.add_param('dtm', dtm)
context.set_location(dsm)
celery_task_id = execute_grass_script.delay(os.path.join(current_dir, "calc_elevation_map.grass"), context.serialize()).task_id
celery_task_id = execute_grass_script.delay(os.path.join(current_dir, "elevationmap.py"), context.serialize()).task_id
return Response({'celery_task_id': celery_task_id}, status=status.HTTP_200_OK)
except GrassEngineException as e:

Wyświetl plik

@ -1,45 +0,0 @@
# dsm_file: GeoTIFF DSM containing the surface to calculate the levels of elevation
# interval: Interval that will split the different levels
# format: OGR output format
# noise_filter_size: Area in meters where we will clean up noise in the contours
# epsg: target EPSG code
# python_script_path: Path of the python script
# python_path: Path to python modules
# dtm: Optional text to include the GeoTIFF DTM
#
# ------
# output: If successful, prints the full path to the elevation map file. If an error.txt file is found, print its content. Otherwise, print 'error'"
ext=""
if [ "${format}" = "GeoJSON" ]; then
ext="json"
elif [ "${format}" = "GPKG" ]; then
ext="gpkg"
elif [ "${format}" = "DXF" ]; then
ext="dxf"
elif [ "${format}" = "ESRI Shapefile" ]; then
ext="shp"
fi
PYTHONPATH="${python_path}" "${python_script_path}" "${dsm_file}" ${interval} --epsg ${epsg} --noise_filter_size ${noise_filter_size} ${dtm} -o output.json
if [ $$ext != "json" ]; then
ogr2ogr -f "${format}" output.$$ext output.json > /dev/null
fi
if [ -e "output.$$ext" ]; then
# ESRI ShapeFile extra steps to compress into a zip archive
# we leverage Python's shutil in this case
if [ "${format}" = "ESRI Shapefile" ]; then
ext="zip"
mkdir elevationmap/
mv output* elevationmap/
echo "import shutil;shutil.make_archive('output', 'zip', 'elevationmap/')" | python
fi
echo "$$(pwd)/output.$$ext"
elif [ -e "error.txt" ]; then
cat "error.txt"
else
echo "error"
fi

Wyświetl plik

@ -1,23 +1,84 @@
#!/usr/bin/env python3
#%module
#% description: This script takes a GeoTIFF file, calculates its heighmap, and outputs it as a GeoJSON
#%end
#%option
#% key: dsm
#% type: string
#% required: yes
#% multiple: no
#% description: The path for the dsm file
#%end
#%option
#% key: intervals
#% type: double
#% required: yes
#% multiple: no
#% description: The intervals used to generate the diferent elevation levels
#%end
#%option
#% key: format
#% type: string
#% required: yes
#% multiple: no
#% description: OGR output format
#%end
#%option
#% key: dtm
#% type: string
#% required: no
#% multiple: no
#% description: The path for the dtm file
#%end
#%option
#% key: epsg
#% type: string
#% required: yes
#% multiple: no
#% description: The epsg code that will be used for output
#%end
#%option
#% key: noise_filter_size
#% type: double
#% required: yes
#% multiple: no
#% description: Area in meters where we will clean up noise in the contours
#%end
import cv2, math, argparse
import numpy as np
import rasterio as rio
from rasterio import warp, transform
from geojson import Feature, FeatureCollection, MultiPolygon, dumps
from geojson import Feature, FeatureCollection, MultiPolygon, dumps
import subprocess
import os
import glob
import shutil
import sys
import grass.script as grass
def main():
ext = ""
if opts['format'] == "GeoJSON":
ext = "json"
elif opts['format'] == "GPKG":
ext = "gpkg"
elif opts['format'] == "DXF":
ext = "dxf"
elif opts['format'] == "ESRI Shapefile":
ext = "shp"
def main(args):
# Open dsm
dsm = rio.open(args.dsm)
dsm = rio.open(opts['dsm'])
# Read the tiff as an numpy masked array
dsm_array = dsm.read(1, masked = True)
# Create a kernel based on the parameter 'noise_filter_size' and the tiff resolution
kernel = get_kernel(args.noise_filter_size, dsm)
kernel = get_kernel(float(opts['noise_filter_size']), dsm)
# Check if we want to use the dtm also
if args.dtm != None:
if opts['dtm'] != '':
# Open the dtm
dtm = rio.open(args.dtm)
dtm = rio.open(opts['dtm'])
# Assert that the dtm and dsm have the same bounds and resolution
assert_same_bounds_and_resolution(dsm, dtm)
# Calculate the different between the dsm and dtm
@ -26,7 +87,7 @@ def main(args):
array = dsm_array
# Calculate the ranges based on the parameter 'intervals' and the elevation array
ranges = calculate_ranges(args.intervals, array)
ranges = calculate_ranges(opts['intervals'], array)
features = []
@ -40,15 +101,32 @@ def main(args):
# Check if we found something
if len(contours) > 0:
# Transform contours from pixels to coordinates
mapped_contours = [map_pixels_to_coordinates(dsm, args.epsg, to_pixel_format(contour)) for contour in contours]
mapped_contours = [map_pixels_to_coordinates(dsm, opts['epsg'], to_pixel_format(contour)) for contour in contours]
# Build the MultiPolygon for based on the contours and their hierarchy
built_multi_polygon = LevelBuilder(bottom, top, mapped_contours, hierarchy[0]).build_multi_polygon()
features.append(built_multi_polygon)
# Write the GeoJSON to a file
dump = dumps(FeatureCollection(features))
with open(args.output, 'w+') as output:
output.write(dump)
with open("output.json", 'w+') as output:
output.write(dump)
if ext != "json":
subprocess.check_call(["ogr2ogr", "-f", opts['format'], "output.%s" % ext, "output.json"], stdout=subprocess.DEVNULL)
if os.path.isfile("output.%s" % ext):
if opts['format'] == "ESRI Shapefile":
ext="zip"
os.makedirs("contours")
contour_files = glob.glob("output.*")
for cf in contour_files:
shutil.move(cf, os.path.join("contours", os.path.basename(cf)))
shutil.make_archive('output', 'zip', 'contours/')
print(os.path.join(os.getcwd(), "output.%s" % ext))
else:
print("error")
def get_kernel(noise_filter_size, dsm):
"""Generate a kernel for noise filtering. Will return none if the noise_filter_size isn't positive"""
@ -151,18 +229,7 @@ class LevelBuilder:
multi_polygon = MultiPolygon(polygons)
return Feature(geometry = multi_polygon, properties = { 'bottom': int(self.bottom), 'top': int(self.top) })
if __name__ == '__main__':
parser = argparse.ArgumentParser(description = 'This script takes a GeoTIFF file, calculates its heighmap, and outputs it as a GeoJSON')
parser.add_argument('dsm', type = str, help = 'The path for the dsm file')
parser.add_argument('intervals', type = str, help = 'The intervals used to generate the diferent elevation levels')
parser.add_argument('-d', '--dtm', type = str, help = 'The path for the dtm file')
parser.add_argument('-e', '--epsg', type = int, help = 'The epsg code that will be used for output', default = 4326)
parser.add_argument('-k', '--noise_filter_size', type = float, help = 'Area in meters where we will clean up noise in the contours', default = 2)
parser.add_argument('-o', '--output', type = str, help = 'The path for the output file', default = "output.json")
args = parser.parse_args()
try:
main(args)
except Exception as e:
with open('error.txt', 'w+') as output:
output.write(str(e))
if __name__ == "__main__":
opts, _ = grass.parser()
sys.exit(main())