diff --git a/SuperBuild/cmake/External-OpenSfM.cmake b/SuperBuild/cmake/External-OpenSfM.cmake index 5b798f8d..875bb343 100644 --- a/SuperBuild/cmake/External-OpenSfM.cmake +++ b/SuperBuild/cmake/External-OpenSfM.cmake @@ -19,7 +19,7 @@ ExternalProject_Add(${_proj_name} #--Download step-------------- DOWNLOAD_DIR ${SB_DOWNLOAD_DIR} GIT_REPOSITORY https://github.com/OpenDroneMap/OpenSfM/ - GIT_TAG 263 + GIT_TAG 264 #--Update/Patch step---------- UPDATE_COMMAND git submodule update --init --recursive #--Configure step------------- diff --git a/VERSION b/VERSION index ec1cf33c..2714f531 100644 --- a/VERSION +++ b/VERSION @@ -1 +1 @@ -2.6.3 +2.6.4 diff --git a/opendm/osfm.py b/opendm/osfm.py index 3fdec0ab..fc4752f7 100644 --- a/opendm/osfm.py +++ b/opendm/osfm.py @@ -23,7 +23,7 @@ from opensfm import report from opendm.multispectral import get_photos_by_band from opendm.gpu import has_gpus from opensfm import multiview -from opensfm.actions.export_geocoords import _get_transformation +from opensfm.actions.export_geocoords import _transform class OSFMContext: def __init__(self, opensfm_project_path): @@ -458,47 +458,33 @@ class OSFMContext: """ Load ground control point information. """ - ds = DataSet(self.opensfm_project_path) - gcps = ds.load_ground_control_points() - if not gcps: + gcp_stats_file = self.path("stats", "ground_control_points.json") + + if not io.file_exists(gcp_stats_file): return [] - reconstructions = ds.load_reconstruction() - reference = ds.load_reference() + gcps_stats = {} + try: + with open(gcp_stats_file) as f: + gcps_stats = json.loads(f.read()) + except: + log.ODM_INFO("Cannot parse %s" % gcp_stats_file) + if not gcps_stats: + return [] + + ds = DataSet(self.opensfm_project_path) + reference = ds.load_reference() projection = pyproj.Proj(proj4) - t = _get_transformation(reference, projection, (0, 0)) - A, b = t[:3, :3], t[:3, 3] result = [] - - for gcp in gcps: - if not gcp.coordinates.has_value: - continue - triangulated = None - - for rec in reconstructions: - triangulated = multiview.triangulate_gcp(gcp, rec.shots, 1.0, 0.1) - if triangulated is None: - continue - else: - break - - if triangulated is None: - continue - - triangulated_topocentric = np.dot(A.T, triangulated) - - coordinates_topocentric = np.array(gcp.coordinates.value) - coordinates = np.dot(A, coordinates_topocentric) + b - triangulated = triangulated + b - + for gcp in gcps_stats: + geocoords = _transform(gcp['coordinates'], reference, projection) result.append({ - 'id': gcp.id, - 'observations': [obs.shot_id for obs in gcp.observations], - 'triangulated': triangulated, - 'coordinates': coordinates, - 'error': np.abs(triangulated_topocentric - coordinates_topocentric) + 'id': gcp['id'], + 'observations': gcp['observations'], + 'coordinates': geocoords, + 'error': gcp['error'] }) return result diff --git a/opendm/types.py b/opendm/types.py index 568d101d..38fccb8d 100644 --- a/opendm/types.py +++ b/opendm/types.py @@ -235,6 +235,7 @@ class ODM_Tree(object): self.opensfm_reconstruction = os.path.join(self.opensfm, 'reconstruction.json') self.opensfm_reconstruction_nvm = os.path.join(self.opensfm, 'undistorted/reconstruction.nvm') self.opensfm_geocoords_reconstruction = os.path.join(self.opensfm, 'reconstruction.geocoords.json') + self.opensfm_topocentric_reconstruction = os.path.join(self.opensfm, 'reconstruction.topocentric.json') # OpenMVS self.openmvs_model = os.path.join(self.openmvs, 'scene_dense_dense_filtered.ply') diff --git a/stages/odm_georeferencing.py b/stages/odm_georeferencing.py index 6e908466..c9cbaff5 100644 --- a/stages/odm_georeferencing.py +++ b/stages/odm_georeferencing.py @@ -3,18 +3,22 @@ import struct import pipes import fiona import fiona.crs +import json from collections import OrderedDict +from pyproj import CRS from opendm import io from opendm import log from opendm import types from opendm import system from opendm import context +from opendm import location from opendm.cropper import Cropper from opendm import point_cloud from opendm.multispectral import get_primary_band_name from opendm.osfm import OSFMContext + class ODMGeoreferencingStage(types.ODM_Stage): def process(self, args, outputs): tree = outputs['tree'] @@ -24,6 +28,7 @@ class ODMGeoreferencingStage(types.ODM_Stage): gcp_export_file = tree.path("odm_georeferencing", "ground_control_points.gpkg") gcp_gml_export_file = tree.path("odm_georeferencing", "ground_control_points.gml") + gcp_geojson_export_file = tree.path("odm_georeferencing", "ground_control_points.geojson") if reconstruction.has_gcp() and (not io.file_exists(gcp_export_file) or self.rerun()): octx = OSFMContext(tree.opensfm) @@ -36,9 +41,6 @@ class ODMGeoreferencingStage(types.ODM_Stage): ('id', 'str'), ('observations_count', 'int'), ('observations_list', 'str'), - ('triangulated_x', 'float'), - ('triangulated_y', 'float'), - ('triangulated_z', 'float'), ('error_x', 'float'), ('error_y', 'float'), ('error_z', 'float'), @@ -58,10 +60,7 @@ class ODMGeoreferencingStage(types.ODM_Stage): 'properties': OrderedDict([ ('id', gcp['id']), ('observations_count', len(gcp['observations'])), - ('observations_list', ",".join(gcp['observations'])), - ('triangulated_x', gcp['triangulated'][0]), - ('triangulated_y', gcp['triangulated'][1]), - ('triangulated_z', gcp['triangulated'][2]), + ('observations_list', ",".join([obs['shot_id'] for obs in gcp['observations']])), ('error_x', gcp['error'][0]), ('error_y', gcp['error'][1]), ('error_z', gcp['error'][2]), @@ -73,10 +72,35 @@ class ODMGeoreferencingStage(types.ODM_Stage): system.run('ogr2ogr -of GML "{}" "{}"'.format(gcp_gml_export_file, gcp_export_file)) except Exception as e: log.ODM_WARNING("Cannot generate ground control points GML file: %s" % str(e)) + + # Write GeoJSON + geojson = { + 'type': 'FeatureCollection', + 'features': [] + } + + from_srs = CRS.from_proj4(reconstruction.georef.proj4()) + to_srs = CRS.from_epsg(4326) + transformer = location.transformer(from_srs, to_srs) + + for gcp in gcps: + properties = gcp.copy() + del properties['coordinates'] + + geojson['features'].append({ + 'geometry': { + 'type': 'Point', + 'coordinates': transformer.TransformPoint(*gcp['coordinates']), + }, + 'properties': properties + }) + + with open(gcp_geojson_export_file, 'w') as f: + f.write(json.dumps(geojson, indent=4)) + else: log.ODM_WARNING("GCPs could not be loaded for writing to %s" % gcp_export_file) - if not io.file_exists(tree.odm_georeferencing_model_laz) or self.rerun(): cmd = ('pdal translate -i "%s" -o \"%s\"' % (tree.filtered_point_cloud, tree.odm_georeferencing_model_laz)) stages = ["ferry"] diff --git a/stages/run_opensfm.py b/stages/run_opensfm.py index 7aa8411f..415778a4 100644 --- a/stages/run_opensfm.py +++ b/stages/run_opensfm.py @@ -68,14 +68,11 @@ class ODMOpenSfMStage(types.ODM_Stage): self.update_progress(75) # We now switch to a geographic CRS - geocoords_flag_file = octx.path("exported_geocoords.txt") - - if reconstruction.is_georeferenced() and (not io.file_exists(geocoords_flag_file) or self.rerun()): + if reconstruction.is_georeferenced() and (not io.file_exists(tree.opensfm_topocentric_reconstruction) or self.rerun()): octx.run('export_geocoords --reconstruction --proj "%s" --offset-x %s --offset-y %s' % (reconstruction.georef.proj4(), reconstruction.georef.utm_east_offset, reconstruction.georef.utm_north_offset)) - # Destructive + shutil.move(tree.opensfm_reconstruction, tree.opensfm_topocentric_reconstruction) shutil.move(tree.opensfm_geocoords_reconstruction, tree.opensfm_reconstruction) - octx.touch(geocoords_flag_file) else: log.ODM_WARNING("Will skip exporting %s" % tree.opensfm_geocoords_reconstruction)