diff --git a/opendm/orthophoto.py b/opendm/orthophoto.py index d5b14a37..c76ca9b9 100644 --- a/opendm/orthophoto.py +++ b/opendm/orthophoto.py @@ -19,6 +19,38 @@ from osgeo import ogr +def clamp_negative_pixels(raster_path): + """ + Checks a raster (georeferenced or not, any number of bands) for negative + pixel values. If any are found, they are zeroed in place on disk and a + warning is logged. + + Uses rasterio (rather than osgeo.gdal's Band.ReadAsArray/WriteArray) since + the latter lazily imports osgeo.gdal_array, a separate compiled extension + that can be ABI-incompatible with the installed numpy in some environments. + :return the total number of negative pixels found (and zeroed) across all bands + """ + if not io.file_exists(raster_path): + return 0 + + total_negative = 0 + try: + with rasterio.open(raster_path, 'r+') as ds: + arr = ds.read() + negative = arr < 0 + total_negative = int(np.count_nonzero(negative)) + if total_negative > 0: + arr[negative] = 0 + ds.write(arr) + except Exception as e: + log.WARNING("Cannot open %s to check for negative pixel values: %s" % (raster_path, str(e))) + return 0 + + if total_negative > 0: + log.WARNING("%s: found and zeroed %s negative pixel values" % (raster_path, total_negative)) + + return total_negative + def get_orthophoto_vars(args): return { 'TILED': 'NO' if args.orthophoto_no_tiled else 'YES', diff --git a/stages/mvstex.py b/stages/mvstex.py index de2045e3..0763bcf3 100644 --- a/stages/mvstex.py +++ b/stages/mvstex.py @@ -1,10 +1,11 @@ -import os, shutil +import os, shutil, glob from opendm import log from opendm import io from opendm import system from opendm import context from opendm import types +from opendm import orthophoto from opendm.multispectral import get_primary_band_name from opendm.photo import find_largest_photo_dim from opendm.objpacker import obj_pack @@ -80,11 +81,18 @@ def add_run(nvm_file, primary=True, band=None): # Format arguments to fit Mvs-Texturing app skipGlobalSeamLeveling = "" + skipLocalSeamLeveling = "" keepUnseenFaces = "" nadir = "" - if args.texturing_skip_global_seam_leveling: + # for multispectral data, force disabling of global seam leveling to avoid negative texture and reflectance values, + # or if user requested disabling for RGB/thermal data + if reconstruction.multi_camera or args.texturing_skip_global_seam_leveling: skipGlobalSeamLeveling = "--skip_global_seam_leveling" + # for multispectral data, force disabling of local seam leveling to avoid negative texture and reflectance values + if reconstruction.multi_camera: + skipLocalSeamLeveling = "--skip_local_seam_leveling" + # TODO: allow user-requested --skip_local_seam_leveling for RGB/thermal data? if args.texturing_keep_unseen_faces: keepUnseenFaces = "--keep_unseen_faces" if (r['nadir']): @@ -98,6 +106,7 @@ def add_run(nvm_file, primary=True, band=None): 'dataTerm': 'gmi', 'outlierRemovalType': 'gauss_clamping', 'skipGlobalSeamLeveling': skipGlobalSeamLeveling, + 'skipLocalSeamLeveling': skipLocalSeamLeveling, 'keepUnseenFaces': keepUnseenFaces, 'toneMapping': 'none', 'nadirMode': nadir, @@ -121,12 +130,20 @@ def add_run(nvm_file, primary=True, band=None): '-t {toneMapping} ' '{intermediate} ' '{skipGlobalSeamLeveling} ' + '{skipLocalSeamLeveling} ' '{keepUnseenFaces} ' '{nadirMode} ' '{labelingFile} ' '{numThreads} ' '{maxTextureSize} '.format(**kwargs)) + # For multispectral data, check the resulting texture atlases for negative pixel + # values (reflectances should be >= 0). Warn and zero if any found. Typically, disabling global and local + # seam leveling avoids negative pixel values, but force an explicit clamp in case something slips by + if reconstruction.multi_camera: + for texfile in glob.glob(kwargs['out_dir'] + "*.tif"): + orthophoto.clamp_negative_pixels(texfile) + if r['primary'] and (not r['nadir'] or args.skip_3dmodel): # Single material? if args.texturing_single_material: diff --git a/stages/odm_orthophoto.py b/stages/odm_orthophoto.py index b44c85a1..de835dec 100644 --- a/stages/odm_orthophoto.py +++ b/stages/odm_orthophoto.py @@ -109,6 +109,13 @@ def process(self, args, outputs): '-outputCornerFile "{corners}" {bands} {depth_idx} {inpaint} ' '{utm_offsets} {a_srs} {vars} {gdal_configs} '.format(**kwargs), env_vars={'OMP_NUM_THREADS': args.max_concurrency}) + # Multispectral reflectance values should never be negative, so catch, warn and zero + # any that slipped through before any further cutline/feathering/tiling steps + # operate on this raster. Typically, these will have been avoided in previous stages. However, if they + # show up in this stage, then we should figure out why and fix the problem (not just the symptom) + if reconstruction.multi_camera: + orthophoto.clamp_negative_pixels(kwargs['ortho']) + # Create georeferenced GeoTiff if reconstruction.is_georeferenced(): bounds_file_path = os.path.join(tree.odm_georeferencing, 'odm_georeferenced_model.bounds.gpkg')