Pseudo-geo camera extraction support

Former-commit-id: ffb842c85a

opendm/config.py

@@ -96,7 +96,7 @@ def config(argv=None):
     parser.add_argument('--end-with', '-e',
                         metavar='<string>',
                         action=StoreValue,
-                        default='odm_orthophoto',
+                        default='odm_report',
                         choices=processopts,
                         help=('Can be one of:' + ' | '.join(processopts)))
 

opendm/pseudogeo.py

@@ -7,7 +7,10 @@ from opendm import log
 def get_pseudogeo_utm():
     return '+proj=utm +zone=30 +ellps=WGS84 +datum=WGS84 +units=m +no_defs'
 
-def add_pseudo_georeferencing(geotiff, scale=1.0):
+def get_pseudogeo_scale():
+    return 0.1 # Arbitrarily chosen
+
+def add_pseudo_georeferencing(geotiff):
     if not io.file_exists(geotiff):
         log.ODM_WARNING("Cannot add pseudo georeferencing, %s does not exist" % geotiff)
         return
@@ -19,7 +22,7 @@ def add_pseudo_georeferencing(geotiff, scale=1.0):
         srs = osr.SpatialReference()
         srs.ImportFromProj4(get_pseudogeo_utm())
         dst_ds.SetProjection( srs.ExportToWkt() )
-        dst_ds.SetGeoTransform( [ 0.0, scale, 0.0, 0.0, 0.0, -scale ] )
+        dst_ds.SetGeoTransform( [ 0.0, get_pseudogeo_scale(), 0.0, 0.0, 0.0, -get_pseudogeo_scale() ] )
         dst_ds = None
 
     except Exception as e:

opendm/shots.py

@@ -1,8 +1,9 @@
 import os, json
 from opendm import log
-from opendm.pseudogeo import get_pseudogeo_utm
+from opendm.pseudogeo import get_pseudogeo_utm, get_pseudogeo_scale
 from opendm.location import transformer
 from pyproj import CRS
+import gdal
 import numpy as np
 import cv2
 
@@ -14,7 +15,7 @@ def get_origin(shot):
     """The origin of the pose in world coordinates."""
     return -get_rotation_matrix(np.array(shot['rotation'])).T.dot(np.array(shot['translation']))
 
-def get_geojson_shots_from_opensfm(reconstruction_file, geocoords_transformation_file=None, utm_srs=None):
+def get_geojson_shots_from_opensfm(reconstruction_file, geocoords_transformation_file=None, utm_srs=None, pseudo_geotiff=None):
     """
     Extract shots from OpenSfM's reconstruction.json
     """
@@ -22,10 +23,25 @@ def get_geojson_shots_from_opensfm(reconstruction_file, geocoords_transformation
     # Read transform (if available)
     if geocoords_transformation_file is not None and utm_srs is not None and os.path.exists(geocoords_transformation_file):
         geocoords = np.loadtxt(geocoords_transformation_file, usecols=range(4))
-    else:
+    elif pseudo_geotiff is not None and os.path.exists(pseudo_geotiff):
         # pseudogeo transform
         utm_srs = get_pseudogeo_utm()
-        geocoords = np.identity(4)
+
+        # the pseudo-georeferencing CRS UL corner is at 0,0
+        # but our shot coordinates aren't, so we need to offset them
+        raster = gdal.Open(pseudo_geotiff)
+        ulx, xres, _, uly, _, yres  = raster.GetGeoTransform()
+        lrx = ulx + (raster.RasterXSize * xres)
+        lry = uly + (raster.RasterYSize * yres)
+
+        geocoords = np.array([[1.0 / get_pseudogeo_scale() ** 2, 0, 0, ulx + lrx / 2.0],
+                              [0, 1.0 / get_pseudogeo_scale() ** 2, 0, uly + lry / 2.0],
+                              [0, 0, 1, 0],
+                              [0, 0, 0, 1]])
+        raster = None
+    else:
+        # Can't deal with this
+        return
 
     crstrans = transformer(CRS.from_proj4(utm_srs), CRS.from_epsg("4326"))
 
@@ -47,10 +63,10 @@ def get_geojson_shots_from_opensfm(reconstruction_file, geocoords_transformation
                     continue
                 
                 cam = cameras[cam]
-                R, T = geocoords[:3, :3], geocoords[:3, 3]
+                Rs, T = geocoords[:3, :3], geocoords[:3, 3]
                 origin = get_origin(shot)
 
-                utm_coords = np.dot(R, origin) + T
+                utm_coords = np.dot(Rs, origin) + T
                 trans_coords = crstrans.TransformPoint(utm_coords[0], utm_coords[1], utm_coords[2])
 
                 feats.append({

stages/odm_dem.py

@@ -130,8 +130,7 @@ class ODMDEMStage(types.ODM_Stage):
                                             overwrite=True)
 
                     if pseudo_georeference:
-                        # 0.1 is arbitrary
-                        pseudogeo.add_pseudo_georeferencing(dem_geotiff_path, 0.1)
+                        pseudogeo.add_pseudo_georeferencing(dem_geotiff_path)
                     
                     progress += 30
                     self.update_progress(progress)

stages/odm_orthophoto.py

@@ -160,8 +160,7 @@ class ODMOrthoPhotoStage(types.ODM_Stage):
                 geotiffcreated = True
             if not geotiffcreated:
                 if io.file_exists(tree.odm_orthophoto_render):
-                    # 0.1 is arbitrary
-                    pseudogeo.add_pseudo_georeferencing(tree.odm_orthophoto_render, 0.1)
+                    pseudogeo.add_pseudo_georeferencing(tree.odm_orthophoto_render)
                     log.ODM_INFO("Renaming %s --> %s" % (tree.odm_orthophoto_render, tree.odm_orthophoto_tif))
                     os.rename(tree.odm_orthophoto_render, tree.odm_orthophoto_tif)
                 else:

stages/odm_report.py

@@ -22,12 +22,15 @@ class ODMReport(types.ODM_Stage):
                 shots = get_geojson_shots_from_opensfm(tree.opensfm_reconstruction, tree.opensfm_transformation, reconstruction.get_proj_srs())    
             else:
                 # Psuedo geo
-                shots = get_geojson_shots_from_opensfm(tree.opensfm_reconstruction)
+                shots = get_geojson_shots_from_opensfm(tree.opensfm_reconstruction, pseudo_geotiff=tree.odm_orthophoto_tif)
 
-            with open(shots_geojson, "w") as fout:
-                fout.write(json.dumps(shots))
-                
-            logger.info("Wrote %s" % shots_geojson)
+            if shots:
+                with open(shots_geojson, "w") as fout:
+                    fout.write(json.dumps(shots))
+
+                log.ODM_INFO("Wrote %s" % shots_geojson)
+            else:
+                log.ODM_WARNING("Cannot extract shots")
         else:
             log.ODM_WARNING('Found a valid shots file in: %s' %
                             tree.shots_geojson)