diff --git a/Dockerfile b/Dockerfile index c3080707..ef8d6131 100644 --- a/Dockerfile +++ b/Dockerfile @@ -9,6 +9,7 @@ ENV PROJ_LIB=/usr/share/proj RUN mkdir /webodm WORKDIR /webodm +# Install Node.js RUN apt-get -qq update && apt-get -qq install -y --no-install-recommends wget RUN wget --no-check-certificate https://deb.nodesource.com/setup_12.x -O /tmp/node.sh && bash /tmp/node.sh RUN apt-get -qq update && apt-get -qq install -y nodejs @@ -19,7 +20,7 @@ RUN printf "Package: *\nPin: release a=testing\nPin-Priority: 750\n" > /etc/apt/ RUN printf "deb http://ftp.us.debian.org/debian/ stable main contrib non-free\ndeb-src http://ftp.us.debian.org/debian/ stable main contrib non-free" > /etc/apt/sources.list.d/stable.list RUN printf "deb http://ftp.us.debian.org/debian/ testing main contrib non-free\ndeb-src http://ftp.us.debian.org/debian/ testing main contrib non-free" > /etc/apt/sources.list.d/testing.list -# Install Python3, Node.js GDAL, nginx, letsencrypt, psql +# Install Python3, GDAL, nginx, letsencrypt, psql RUN apt-get -qq update && apt-get -qq install -t testing -y --no-install-recommends python3 python3-pip git g++ python3-dev libpq-dev binutils libproj-dev gdal-bin python3-gdal nginx certbot grass-core && apt-get -qq install -y --no-install-recommends gettext-base cron postgresql-client-9.6 RUN update-alternatives --install /usr/bin/python python /usr/bin/python2.7 1 && update-alternatives --install /usr/bin/python python /usr/bin/python3.8 2 diff --git a/app/cogeo.py b/app/cogeo.py index 28be25b8..19cacbb9 100644 --- a/app/cogeo.py +++ b/app/cogeo.py @@ -5,6 +5,7 @@ import shutil import rasterio import re import subprocess +from pipes import quote from rio_cogeo.cogeo import cog_validate, cog_translate from rio_tiler.utils import has_alpha_band from webodm import settings @@ -17,16 +18,23 @@ def valid_cogeo(src_path): :param src_path: path to GeoTIFF :return: true if the GeoTIFF is a cogeo, false otherwise """ - return cog_validate(src_path, strict=True) + try: + from app.vendor.validate_cloud_optimized_geotiff import validate + warnings, errors, details = validate(src_path, full_check=True) + return not errors and not warnings + except ModuleNotFoundError: + logger.warning("Using legacy cog_validate (osgeo.gdal package not found)") + # Legacy + return cog_validate(src_path, strict=True) -def assure_cogeo(src_path, use_legacy=False): +def assure_cogeo(src_path): """ Guarantee that the .tif passed as an argument is a Cloud Optimized GeoTIFF (cogeo) If the path is not a cogeo, it is destructively converted into a cogeo. If the file cannot be converted, the function does not change the file :param src_path: path to GeoTIFF (cogeo or not) - :param use_legacy: whether to force the use of legacy implementation. By default the best implementation is used. + :param force_use_legacy: whether to force the use of legacy implementation. By default the best implementation is used. :return: None """ @@ -57,11 +65,11 @@ def assure_cogeo(src_path, use_legacy=False): # This shouldn't happen use_legacy = True - if True or use_legacy: - logger.info("Using legacy implementation (GDAL >= 3.1 not found)") - make_cogeo_legacy(src_path) + if use_legacy: + logger.warning("Using legacy implementation (GDAL >= 3.1 not found)") + return make_cogeo_legacy(src_path) else: - make_cogeo_gdal(src_path) + return make_cogeo_gdal(src_path) def get_gdal_version(): # Bit of a hack without installing @@ -71,7 +79,7 @@ def get_gdal_version(): return None # Get version - version_output = subprocess.check_output("%s --version" % gdal_translate) + version_output = subprocess.check_output([gdal_translate, "--version"]).decode('utf-8') m = re.match(r"GDAL\s+([\d+])\.([\d+])\.([\d+]),\s+released", version_output) if not m: @@ -89,7 +97,33 @@ def make_cogeo_gdal(src_path): tmpfile = tempfile.mktemp('_cogeo.tif', dir=settings.MEDIA_TMP) swapfile = tempfile.mktemp('_cogeo_swap.tif', dir=settings.MEDIA_TMP) - # gdal_translate -of COG -co BLOCKSIZE=256 -co COMPRESS=deflate -co NUM_THREADS=4 -co BIGTIFF=IF_SAFER -co QUALITY=100 -co SPARSE_OK=ON --config GDAL_NUM_THREADS ALL_CPUS brighton.tif cog.tif + try: + subprocess.run(["gdal_translate", "-of", "COG", + "-co", "BLOCKSIZE=256", + "-co", "COMPRESS=deflate", + "-co", "NUM_THREADS=ALL_CPUS", + "-co", "BIGTIFF=IF_SAFER", + "--config", "GDAL_NUM_THREADS", "ALL_CPUS", + quote(src_path), quote(tmpfile)]) + except Exception as e: + logger.warning("Cannot create Cloud Optimized GeoTIFF: %s" % str(e)) + + if os.path.isfile(tmpfile): + shutil.move(src_path, swapfile) # Move to swap location + + try: + shutil.move(tmpfile, src_path) + except IOError as e: + logger.warning("Cannot move %s to %s: %s" % (tmpfile, src_path, str(e))) + shutil.move(swapfile, src_path) # Attempt to restore + raise e + + if os.path.isfile(swapfile): + os.remove(swapfile) + + return True + else: + return False def make_cogeo_legacy(src_path): """ diff --git a/app/vendor/__init__.py b/app/vendor/__init__.py new file mode 100644 index 00000000..e69de29b diff --git a/app/vendor/validate_cloud_optimized_geotiff.py b/app/vendor/validate_cloud_optimized_geotiff.py new file mode 100644 index 00000000..1b34f26e --- /dev/null +++ b/app/vendor/validate_cloud_optimized_geotiff.py @@ -0,0 +1,407 @@ +#!/usr/bin/env python3 +# -*- coding: utf-8 -*- +# ***************************************************************************** +# $Id$ +# +# Project: GDAL +# Purpose: Validate Cloud Optimized GeoTIFF file structure +# Author: Even Rouault, +# +# ***************************************************************************** +# Copyright (c) 2017, Even Rouault +# +# Permission is hereby granted, free of charge, to any person obtaining a +# copy of this software and associated documentation files (the "Software"), +# to deal in the Software without restriction, including without limitation +# the rights to use, copy, modify, merge, publish, distribute, sublicense, +# and/or sell copies of the Software, and to permit persons to whom the +# Software is furnished to do so, subject to the following conditions: +# +# The above copyright notice and this permission notice shall be included +# in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS +# OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, +# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL +# THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER +# LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +# FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER +# DEALINGS IN THE SOFTWARE. +# ***************************************************************************** + +import os.path +import struct +import sys +from osgeo import gdal + + +def Usage(): + print('Usage: validate_cloud_optimized_geotiff.py [-q] [--full-check=yes/no/auto] test.tif') + print('') + print('Options:') + print('-q: quiet mode') + print('--full-check=yes/no/auto: check tile/strip leader/trailer bytes. auto=yes for local files, and no for remote files') + return 1 + + +class ValidateCloudOptimizedGeoTIFFException(Exception): + pass + + +def full_check_band(f, band_name, band, errors, + block_order_row_major, + block_leader_size_as_uint4, + block_trailer_last_4_bytes_repeated, + mask_interleaved_with_imagery): + + block_size = band.GetBlockSize() + mask_band = None + if mask_interleaved_with_imagery: + mask_band = band.GetMaskBand() + mask_block_size = mask_band.GetBlockSize() + if block_size != mask_block_size: + errors += [band_name + ': mask block size is different from its imagery band'] + mask_band = None + + yblocks = (band.YSize + block_size[1] - 1) // block_size[1] + xblocks = (band.XSize + block_size[0] - 1) // block_size[0] + last_offset = 0 + for y in range(yblocks): + for x in range(xblocks): + + offset = band.GetMetadataItem('BLOCK_OFFSET_%d_%d' % (x, y), 'TIFF') + offset = int(offset) if offset is not None else 0 + bytecount = band.GetMetadataItem('BLOCK_SIZE_%d_%d' % (x, y), 'TIFF') + bytecount = int(bytecount) if bytecount is not None else 0 + + if offset > 0: + if block_order_row_major and offset < last_offset: + errors += [band_name + + ': offset of block (%d, %d) is smaller than previous block' % (x, y)] + + if block_leader_size_as_uint4: + gdal.VSIFSeekL(f, offset - 4, 0) + leader_size = struct.unpack('= 4: + gdal.VSIFSeekL(f, offset + bytecount - 4, 0) + last_bytes = gdal.VSIFReadL(8, 1, f) + if last_bytes[0:4] != last_bytes[4:8]: + errors += [band_name + + ': for block (%d, %d), trailer bytes are invalid' % (x, y)] + + if mask_band: + offset_mask = mask_band.GetMetadataItem('BLOCK_OFFSET_%d_%d' % (x, y), 'TIFF') + offset_mask = int(offset_mask) if offset_mask is not None else 0 + if offset > 0 and offset_mask > 0: + #bytecount_mask = int(mask_band.GetMetadataItem('BLOCK_SIZE_%d_%d' % (x,y), 'TIFF')) + expected_offset_mask = offset + bytecount + \ + (4 if block_leader_size_as_uint4 else 0) + \ + (4 if block_trailer_last_4_bytes_repeated else 0) + if offset_mask != expected_offset_mask: + errors += ['Mask of ' + band_name + ': for block (%d, %d), offset is %d, whereas %d was expected' % ( + x, y, offset_mask, expected_offset_mask)] + elif offset == 0 and offset_mask > 0: + if block_order_row_major and offset_mask < last_offset: + errors += ['Mask of ' + band_name + + ': offset of block (%d, %d) is smaller than previous block' % (x, y)] + + offset = offset_mask + + last_offset = offset + + +def validate(ds, check_tiled=True, full_check=False): + """Check if a file is a (Geo)TIFF with cloud optimized compatible structure. + + Args: + ds: GDAL Dataset for the file to inspect. + check_tiled: Set to False to ignore missing tiling. + full_check: Set to TRUe to check tile/strip leader/trailer bytes. Might be slow on remote files + + Returns: + A tuple, whose first element is an array of error messages + (empty if there is no error), and the second element, a dictionary + with the structure of the GeoTIFF file. + + Raises: + ValidateCloudOptimizedGeoTIFFException: Unable to open the file or the + file is not a Tiff. + """ + + if int(gdal.VersionInfo('VERSION_NUM')) < 2020000: + raise ValidateCloudOptimizedGeoTIFFException( + 'GDAL 2.2 or above required') + + unicode_type = type(''.encode('utf-8').decode('utf-8')) + if isinstance(ds, (str, unicode_type)): + gdal.PushErrorHandler() + ds = gdal.Open(ds) + gdal.PopErrorHandler() + if ds is None: + raise ValidateCloudOptimizedGeoTIFFException( + 'Invalid file : %s' % gdal.GetLastErrorMsg()) + if ds.GetDriver().ShortName != 'GTiff': + raise ValidateCloudOptimizedGeoTIFFException( + 'The file is not a GeoTIFF') + + details = {} + errors = [] + warnings = [] + filename = ds.GetDescription() + main_band = ds.GetRasterBand(1) + ovr_count = main_band.GetOverviewCount() + filelist = ds.GetFileList() + if filelist is not None and filename + '.ovr' in filelist: + errors += [ + 'Overviews found in external .ovr file. They should be internal'] + + if main_band.XSize > 512 or main_band.YSize > 512: + if check_tiled: + block_size = main_band.GetBlockSize() + if block_size[0] == main_band.XSize and block_size[0] > 1024: + errors += [ + 'The file is greater than 512xH or Wx512, but is not tiled'] + + if ovr_count == 0: + warnings += [ + 'The file is greater than 512xH or Wx512, it is recommended ' + 'to include internal overviews'] + + ifd_offset = int(main_band.GetMetadataItem('IFD_OFFSET', 'TIFF')) + ifd_offsets = [ifd_offset] + + block_order_row_major = False + block_leader_size_as_uint4 = False + block_trailer_last_4_bytes_repeated = False + mask_interleaved_with_imagery = False + + if ifd_offset not in (8, 16): + + # Check if there is GDAL hidden structural metadata + f = gdal.VSIFOpenL(filename, 'rb') + if not f: + raise ValidateCloudOptimizedGeoTIFFException("Cannot open file") + signature = struct.unpack('B' * 4, gdal.VSIFReadL(4, 1, f)) + bigtiff = signature in ((0x49, 0x49, 0x2B, 0x00), (0x4D, 0x4D, 0x00, 0x2B)) + if bigtiff: + expected_ifd_pos = 16 + else: + expected_ifd_pos = 8 + gdal.VSIFSeekL(f, expected_ifd_pos, 0) + pattern = "GDAL_STRUCTURAL_METADATA_SIZE=%06d bytes\n" % 0 + got = gdal.VSIFReadL(len(pattern), 1, f).decode('LATIN1') + if len(got) == len(pattern) and got.startswith('GDAL_STRUCTURAL_METADATA_SIZE='): + size = int(got[len('GDAL_STRUCTURAL_METADATA_SIZE='):][0:6]) + extra_md = gdal.VSIFReadL(size, 1, f).decode('LATIN1') + block_order_row_major = 'BLOCK_ORDER=ROW_MAJOR' in extra_md + block_leader_size_as_uint4 = 'BLOCK_LEADER=SIZE_AS_UINT4' in extra_md + block_trailer_last_4_bytes_repeated = 'BLOCK_TRAILER=LAST_4_BYTES_REPEATED' in extra_md + mask_interleaved_with_imagery = 'MASK_INTERLEAVED_WITH_IMAGERY=YES' in extra_md + if 'KNOWN_INCOMPATIBLE_EDITION=YES' in extra_md: + errors += ["KNOWN_INCOMPATIBLE_EDITION=YES is declared in the file"] + expected_ifd_pos += len(pattern) + size + expected_ifd_pos += expected_ifd_pos % 2 # IFD offset starts on a 2-byte boundary + gdal.VSIFCloseL(f) + + if expected_ifd_pos != ifd_offsets[0]: + errors += [ + 'The offset of the main IFD should be %d. It is %d instead' % (expected_ifd_pos, ifd_offsets[0])] + + details['ifd_offsets'] = {} + details['ifd_offsets']['main'] = ifd_offset + + for i in range(ovr_count): + # Check that overviews are by descending sizes + ovr_band = ds.GetRasterBand(1).GetOverview(i) + if i == 0: + if (ovr_band.XSize > main_band.XSize or + ovr_band.YSize > main_band.YSize): + errors += [ + 'First overview has larger dimension than main band'] + else: + prev_ovr_band = ds.GetRasterBand(1).GetOverview(i - 1) + if (ovr_band.XSize > prev_ovr_band.XSize or + ovr_band.YSize > prev_ovr_band.YSize): + errors += [ + 'Overview of index %d has larger dimension than ' + 'overview of index %d' % (i, i - 1)] + + if check_tiled: + block_size = ovr_band.GetBlockSize() + if block_size[0] == ovr_band.XSize and block_size[0] > 1024: + errors += [ + 'Overview of index %d is not tiled' % i] + + # Check that the IFD of descending overviews are sorted by increasing + # offsets + ifd_offset = int(ovr_band.GetMetadataItem('IFD_OFFSET', 'TIFF')) + ifd_offsets.append(ifd_offset) + details['ifd_offsets']['overview_%d' % i] = ifd_offset + if ifd_offsets[-1] < ifd_offsets[-2]: + if i == 0: + errors += [ + 'The offset of the IFD for overview of index %d is %d, ' + 'whereas it should be greater than the one of the main ' + 'image, which is at byte %d' % + (i, ifd_offsets[-1], ifd_offsets[-2])] + else: + errors += [ + 'The offset of the IFD for overview of index %d is %d, ' + 'whereas it should be greater than the one of index %d, ' + 'which is at byte %d' % + (i, ifd_offsets[-1], i - 1, ifd_offsets[-2])] + + # Check that the imagery starts by the smallest overview and ends with + # the main resolution dataset + + def get_block_offset(band): + blockxsize, blockysize = band.GetBlockSize() + for y in range(int((band.YSize + blockysize - 1) / blockysize)): + for x in range(int((band.XSize + blockxsize - 1) / blockxsize)): + block_offset = band.GetMetadataItem('BLOCK_OFFSET_%d_%d' % (x, y), 'TIFF') + if block_offset: + return int(block_offset) + return 0 + + block_offset = get_block_offset(main_band) + data_offsets = [block_offset] + details['data_offsets'] = {} + details['data_offsets']['main'] = block_offset + for i in range(ovr_count): + ovr_band = ds.GetRasterBand(1).GetOverview(i) + block_offset = get_block_offset(ovr_band) + data_offsets.append(block_offset) + details['data_offsets']['overview_%d' % i] = block_offset + + if data_offsets[-1] != 0 and data_offsets[-1] < ifd_offsets[-1]: + if ovr_count > 0: + errors += [ + 'The offset of the first block of the smallest overview ' + 'should be after its IFD'] + else: + errors += [ + 'The offset of the first block of the image should ' + 'be after its IFD'] + for i in range(len(data_offsets) - 2, 0, -1): + if data_offsets[i] != 0 and data_offsets[i] < data_offsets[i + 1]: + errors += [ + 'The offset of the first block of overview of index %d should ' + 'be after the one of the overview of index %d' % + (i - 1, i)] + if len(data_offsets) >= 2 and data_offsets[0] != 0 and data_offsets[0] < data_offsets[1]: + errors += [ + 'The offset of the first block of the main resolution image ' + 'should be after the one of the overview of index %d' % + (ovr_count - 1)] + + if full_check and (block_order_row_major or block_leader_size_as_uint4 or + block_trailer_last_4_bytes_repeated or + mask_interleaved_with_imagery): + f = gdal.VSIFOpenL(filename, 'rb') + if not f: + raise ValidateCloudOptimizedGeoTIFFException("Cannot open file") + + full_check_band(f, 'Main resolution image', main_band, errors, + block_order_row_major, + block_leader_size_as_uint4, + block_trailer_last_4_bytes_repeated, + mask_interleaved_with_imagery) + if main_band.GetMaskFlags() == gdal.GMF_PER_DATASET and \ + (filename + '.msk') not in ds.GetFileList(): + full_check_band(f, 'Mask band of main resolution image', + main_band.GetMaskBand(), errors, + block_order_row_major, + block_leader_size_as_uint4, + block_trailer_last_4_bytes_repeated, False) + for i in range(ovr_count): + ovr_band = ds.GetRasterBand(1).GetOverview(i) + full_check_band(f, 'Overview %d' % i, ovr_band, errors, + block_order_row_major, + block_leader_size_as_uint4, + block_trailer_last_4_bytes_repeated, + mask_interleaved_with_imagery) + if ovr_band.GetMaskFlags() == gdal.GMF_PER_DATASET and \ + (filename + '.msk') not in ds.GetFileList(): + full_check_band(f, 'Mask band of overview %d' % i, + ovr_band.GetMaskBand(), errors, + block_order_row_major, + block_leader_size_as_uint4, + block_trailer_last_4_bytes_repeated, False) + gdal.VSIFCloseL(f) + + return warnings, errors, details + + +def main(): + """Return 0 in case of success, 1 for failure.""" + + i = 1 + filename = None + quiet = False + full_check = None + while i < len(sys.argv): + if sys.argv[i] == '-q': + quiet = True + elif sys.argv[i] == '--full-check=yes': + full_check = True + elif sys.argv[i] == '--full-check=no': + full_check = False + elif sys.argv[i] == '--full-check=auto': + full_check = None + elif sys.argv[i][0] == '-': + return Usage() + elif filename is None: + filename = sys.argv[i] + else: + return Usage() + + i += 1 + + if filename is None: + return Usage() + + if full_check is None: + full_check = filename.startswith('/vsimem/') or os.path.exists(filename) + + try: + ret = 0 + warnings, errors, details = validate(filename, full_check=full_check) + if warnings: + if not quiet: + print('The following warnings were found:') + for warning in warnings: + print(' - ' + warning) + print('') + if errors: + if not quiet: + print('%s is NOT a valid cloud optimized GeoTIFF.' % filename) + print('The following errors were found:') + for error in errors: + print(' - ' + error) + print('') + ret = 1 + else: + if not quiet: + print('%s is a valid cloud optimized GeoTIFF' % filename) + + if not quiet and not warnings and not errors: + headers_size = min(details['data_offsets'][k] for k in details['data_offsets']) + if headers_size == 0: + headers_size = gdal.VSIStatL(filename).size + print('\nThe size of all IFD headers is %d bytes' % headers_size) + except ValidateCloudOptimizedGeoTIFFException as e: + if not quiet: + print('%s is NOT a valid cloud optimized GeoTIFF : %s' % + (filename, str(e))) + ret = 1 + + return ret + + +if __name__ == '__main__': + sys.exit(main())