Merge DEM step finished

pull/979/head
Piero Toffanin 2019-04-29 14:01:55 -04:00
rodzic c7a0ca9dea
commit 73bf8ca2a6
4 zmienionych plików z 54 dodań i 27 usunięć

Wyświetl plik

@ -9,7 +9,7 @@ import os
from rasterio import logging
def euclidean_merge_dems(input_dems, output_dem):
def euclidean_merge_dems(input_dems, output_dem, creation_options={}):
"""
Based on https://github.com/mapbox/rio-merge-rgba
and ideas from Anna Petrasova
@ -22,13 +22,7 @@ def euclidean_merge_dems(input_dems, output_dem):
inputs = []
bounds=None
precision=7
creation_options={}
log_ = logging.getLogger()
debug_enabled = log_.isEnabledFor(logging.DEBUG)
if debug_enabled:
logging.disable(logging.DEBUG)
existing_dems = []
for dem in input_dems:
if not io.file_exists(dem):
@ -40,6 +34,12 @@ def euclidean_merge_dems(input_dems, output_dem):
log.ODM_WARNING("No input DEMs, skipping euclidean merge.")
return
# Suppress DEBUG logging for rasterio operations
log_ = logging.getLogger()
debug_enabled = log_.isEnabledFor(logging.DEBUG)
if debug_enabled:
logging.disable(logging.DEBUG)
with rasterio.open(existing_dems[0]) as first:
src_nodata = 0
nodatavals = first.get_nodatavals()
@ -110,10 +110,10 @@ def euclidean_merge_dems(input_dems, output_dem):
profile["transform"] = output_transform
profile["height"] = output_height
profile["width"] = output_width
profile["tiled"] = True
profile["blockxsize"] = 512
profile["blockysize"] = 512
profile["compress"] = 'lzw'
profile["tiled"] = creation_options.get('TILED', 'YES') == 'YES'
profile["blockxsize"] = creation_options.get('BLOCKXSIZE', 512)
profile["blockysize"] = creation_options.get('BLOCKYSIZE', 512)
profile["compress"] = creation_options.get('COMPRESS', 'LZW')
profile["nodata"] = src_nodata
# Creation opts
@ -133,14 +133,10 @@ def euclidean_merge_dems(input_dems, output_dem):
dst_count = first.count
dst_shape = (dst_count, dst_rows, dst_cols)
log.ODM_DEBUG("Temp shape: %r", dst_shape)
# dstarr = np.full(dst_shape, src_nodata, dtype=dtype)
dstarr = np.zeros(dst_shape, dtype=dtype)
distsum = np.zeros(dst_shape, dtype=dtype)
# Read up srcs until
# a. everything is data; i.e. no nodata
# b. no sources left
for src_d, src_e in sources:
# The full_cover behavior is problematic here as it includes
# extra pixels along the bottom right when the sources are
@ -187,5 +183,9 @@ def euclidean_merge_dems(input_dems, output_dem):
for _, euclidean_geotiff in inputs:
if io.file_exists(euclidean_geotiff):
os.remove(euclidean_geotiff)
# Restore logging
if debug_enabled:
logging.disable(logging.NOTSET)
return output_dem

Wyświetl plik

@ -0,0 +1,9 @@
def get_dem_vars(args):
return {
'TILED': 'YES',
'COMPRESS': 'LZW',
'BLOCKXSIZE': 512,
'BLOCKYSIZE': 512,
'NUM_THREADS': args.max_concurrency,
}

Wyświetl plik

@ -7,7 +7,7 @@ from opendm import system
from opendm import context
from opendm import types
from opendm import gsd
from opendm.dem import commands
from opendm.dem import commands, utils
from opendm.cropper import Cropper
class ODMDEMStage(types.ODM_Stage):
@ -80,13 +80,7 @@ class ODMDEMStage(types.ODM_Stage):
if args.crop > 0:
bounds_file_path = os.path.join(tree.odm_georeferencing, 'odm_georeferenced_model.bounds.gpkg')
if os.path.exists(bounds_file_path):
Cropper.crop(bounds_file_path, os.path.join(odm_dem_root, "{}.tif".format(product)), {
'TILED': 'YES',
'COMPRESS': 'LZW',
'BLOCKXSIZE': 512,
'BLOCKYSIZE': 512,
'NUM_THREADS': self.params.get('max_concurrency')
})
Cropper.crop(bounds_file_path, os.path.join(odm_dem_root, "{}.tif".format(product)), utils.get_dem_vars(args))
else:
log.ODM_WARNING('Found existing outputs in: %s' % odm_dem_root)
else:

Wyświetl plik

@ -6,7 +6,8 @@ from opendm import types
from opendm import io
from opendm import system
from opendm import orthophoto
from opendm.dem import pdal
from opendm.dem import pdal, utils
from opendm.dem.merge import euclidean_merge_dems
from opensfm.large import metadataset
from opendm.cropper import Cropper
from opendm.concurrency import get_max_memory, get_max_memory_mb
@ -237,10 +238,33 @@ class ODMMergeStage(types.ODM_Stage):
else:
log.ODM_WARNING("Found merged orthophoto in %s" % tree.odm_orthophoto_tif)
# Merge DEM
# Merge DEMs
# TODO: crop DEM if necessary
def merge_dems(dem_filename, human_name):
dem_file = os.path.join(tree.path("odm_dem", dem_filename))
if not io.file_exists(dem_file) or self.rerun():
all_dems = get_submodel_paths(tree.submodels_path, "odm_dem", dem_filename)
log.ODM_INFO("Merging %ss" % human_name)
# Merge
dem_vars = utils.get_dem_vars(args)
euclidean_merge_dems(all_dems, dem_file, dem_vars)
if io.file_exists(dem_file):
# Crop
if args.crop > 0:
Cropper.crop(merged_bounds_file, dem_file, dem_vars)
log.ODM_INFO("Created %s" % dem_file)
else:
log.ODM_WARNING("Cannot merge %s, %s was not created" % (human_name, dem_file))
else:
log.ODM_WARNING("Found merged %s in %s" % (human_name, dsm_file))
if args.dsm:
merge_dems("dsm.tif", "DSM")
if args.dtm:
merge_dems("dtm.tif", "DTM")
# Stop the pipeline short! We're done.
self.next_stage = None