From d4102e6d25250f26d7ae250d0685dce983573dd9 Mon Sep 17 00:00:00 2001 From: Emmanuel Dubois Date: Mon, 2 Oct 2023 12:05:32 +0200 Subject: [PATCH 01/12] feat: first geomodel structure to aggregate rpc and grid geomodel structure --- shareloc/geomodels/__init__.py | 8 +- shareloc/geomodels/geomodel.py | 112 ++++++++++++++++++++++++ shareloc/geomodels/geomodel_template.py | 58 ++++++++++++ shareloc/geomodels/grid.py | 8 +- shareloc/geomodels/rpc.py | 5 +- 5 files changed, 187 insertions(+), 4 deletions(-) create mode 100644 shareloc/geomodels/geomodel.py create mode 100644 shareloc/geomodels/geomodel_template.py diff --git a/shareloc/geomodels/__init__.py b/shareloc/geomodels/__init__.py index 8bea709..00772b2 100644 --- a/shareloc/geomodels/__init__.py +++ b/shareloc/geomodels/__init__.py @@ -1,7 +1,7 @@ #!/usr/bin/env python # coding: utf8 # -# Copyright (c) 2022 Centre National d'Etudes Spatiales (CNES). +# Copyright (c) 2023 Centre National d'Etudes Spatiales (CNES). # # This file is part of Shareloc # (see https://github.com/CNES/shareloc). @@ -20,4 +20,10 @@ # """ Shareloc geomodels module +Imports are used to simplify calls to module API Coregistration. """ +# Demcompare imports +from . import grid, rpc +from .geomodel import GeoModel + +__all__ = ["rpc", "grid", "GeoModel"] # To avoid flake8 F401 diff --git a/shareloc/geomodels/geomodel.py b/shareloc/geomodels/geomodel.py new file mode 100644 index 0000000..450167c --- /dev/null +++ b/shareloc/geomodels/geomodel.py @@ -0,0 +1,112 @@ +#!/usr/bin/env python +# coding: utf8 +# +# Copyright (c) 2023 Centre National d'Etudes Spatiales (CNES). +# +# This file is part of shareloc +# (see https://github.com/CNES/shareloc). +# +# Licensed under the Apache License, Version 2.0 (the "License"); +# you may not use this file except in compliance with the License. +# You may obtain a copy of the License at +# +# http://www.apache.org/licenses/LICENSE-2.0 +# +# Unless required by applicable law or agreed to in writing, software +# distributed under the License is distributed on an "AS IS" BASIS, +# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +# See the License for the specific language governing permissions and +# limitations under the License. +# +""" +This module contains the GeoModel class factory. +This main API GeoModel class generates an object linked +with geomodel configuration "geomodel_name" +(from registered GeoModelTemplate class) + +** Experimental and not used as API for now** +""" + +# Standard imports +import logging +from typing import Any, Dict + + +class GeoModel: + """ + GeoModel factory: + A class designed for registered all available geomodels + and instantiate them when needed. + """ + + # Dict (geomodel_name: str, class: object) containing registered geomodels + available_geomodels: Dict[str, Any] = {} + + def __new__(cls, geomodel_path: str, geomodel_type: str = "RPC"): + """ + Return a GeoModelTemplate child instance + associated with the "geomodel_name" given in the configuration + through create_geomodel local method for clarity. + + :param geomodel_path: Path of geomodel file + :type geomodel_path: string + :param geomodel_type: Type of the geomodel, default "rpc", used as key for specific geomodel subclass instance + :type geomodel_type: string + """ + return cls.create_geomodel(geomodel_path, geomodel_type) + + @classmethod + def create_geomodel(cls, geomodel_path: str, geomodel_type: str = "RPC"): + """ + Factory command to create the geomodel from geomodel_type + Return a GeoModelTemplate child instance + associated with the "geomodel_type" + + :param geomodel_path: Path of geomodel file + :type geomodel_path: string + :param geomodel_type: Type of the geomodel, default "rpc", used as key for specific geomodel subclass instance + :type geomodel_type: string + """ + + # If no type is given, the default is "rpc" + + # Create geomodel object with geomodel_path parameter from geomodel_type if exists + try: + geomodel_class = cls.available_geomodels[geomodel_type] + geomodel = geomodel_class(geomodel_path) + logging.debug("GeoModel type name: %s", geomodel_type) + except KeyError: + logging.error("Geomodel type named %s is not supported", geomodel_type) + raise + + return geomodel + + @classmethod + def print_avalaible_geomodels(cls): + """ + Print all registered applications + """ + for geomodel_type in cls.available_geomodels: + print(geomodel_type) + + @classmethod + def register(cls, geomodel_type: str): + """ + Allows to register the GeoModelTemplate subclass in the factory + with its geomodel geomodel_type through decorator + + :param geomodel_type: the subclass name to be registered + :type geomodel_type: string + """ + + def decorator(geomodel_subclass): + """ + Register the geomodel subclass in the available methods + + :param geomodel_subclass: the subclass to be registered + :type geomodel_subclass: object + """ + cls.available_geomodels[geomodel_type] = geomodel_subclass + return geomodel_subclass + + return decorator diff --git a/shareloc/geomodels/geomodel_template.py b/shareloc/geomodels/geomodel_template.py new file mode 100644 index 0000000..e1c6cdf --- /dev/null +++ b/shareloc/geomodels/geomodel_template.py @@ -0,0 +1,58 @@ +#!/usr/bin/env python +# coding: utf8 +# +# Copyright (c) 2023 Centre National d'Etudes Spatiales (CNES). +# +# This file is part of shareloc +# (see https://github.com/CNES/shareloc). +# +# Licensed under the Apache License, Version 2.0 (the "License"); +# you may not use this file except in compliance with the License. +# You may obtain a copy of the License at +# +# http://www.apache.org/licenses/LICENSE-2.0 +# +# Unless required by applicable law or agreed to in writing, software +# distributed under the License is distributed on an "AS IS" BASIS, +# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +# See the License for the specific language governing permissions and +# limitations under the License. +# +""" +This module contains the coregistration class template. +It contains the structure for all coregistration methods in subclasses and +generic coregistration code to avoid duplication. +""" + +# Standard imports +from abc import ABCMeta, abstractmethod +from typing import Dict + +# Third party imports + + +# pylint: disable=too-few-public-methods +class GeoModelTemplate(metaclass=ABCMeta): + """ + Class for general specification of a geometric model + declined in rpc.py and grid.py + """ + + @abstractmethod + def __init__(self, geomodel_file: str, geomodel_params: Dict = None): + """ + Return the geomodel object associated with the geomodel_type + given in the configuration + + :param geomodel_file: path of the geomodel file to instanciate. + :type geomodel_file: string + :param geomodel_params: Geomodels parameters as dict (unused, just to fit RPC) TODO: evolve with new API + :type geomodel_params: Dict + """ + # geomodel filename path + self.filename = geomodel_file + + # geomodel_params if exists (to evolve) + self.geomodels_params = geomodel_params + + # diff --git a/shareloc/geomodels/grid.py b/shareloc/geomodels/grid.py index 1bfa8e8..cb9fbcb 100755 --- a/shareloc/geomodels/grid.py +++ b/shareloc/geomodels/grid.py @@ -30,6 +30,8 @@ import numpy as np # Shareloc imports +from shareloc.geomodels.geomodel import GeoModel +from shareloc.geomodels.geomodel_template import GeoModelTemplate from shareloc.image import Image from shareloc.math_utils import interpol_bilin_grid, interpol_bilin_vectorized from shareloc.proj_utils import coordinates_conversion @@ -37,7 +39,8 @@ # gitlab issue #58 # pylint: disable=too-many-instance-attributes -class Grid: +@GeoModel.register("grid") +class Grid(GeoModelTemplate): """ multi H direct localization grid handling class. please refer to the main documentation grid format @@ -83,7 +86,8 @@ def __init__(self, grid_filename): :param grid_filename: grid filename (Geotiff) :type grid_filename: string """ - self.filename = grid_filename + # Instanciate GeoModelTemplate generic init with shared parameters + super().__init__(grid_filename) self.row0 = None self.col0 = None self.nbrow = None diff --git a/shareloc/geomodels/rpc.py b/shareloc/geomodels/rpc.py index 071a781..1b79abb 100755 --- a/shareloc/geomodels/rpc.py +++ b/shareloc/geomodels/rpc.py @@ -36,6 +36,8 @@ from numba import config, njit, prange # Shareloc imports +from shareloc.geomodels.geomodel import GeoModel +from shareloc.geomodels.geomodel_template import GeoModelTemplate from shareloc.proj_utils import coordinates_conversion # Set numba type of threading layer before parallel target compilation @@ -119,7 +121,8 @@ def identify_geotiff_rpc(image_filename): return None -class RPC: +@GeoModel.register("RPC") +class RPC(GeoModelTemplate): """ RPC class including direct and inverse localization instance methods """ From 1a4837cfe2661233b16e45c5ddb4e190ab2ec088 Mon Sep 17 00:00:00 2001 From: Emmanuel Dubois Date: Wed, 18 Oct 2023 17:18:27 +0200 Subject: [PATCH 02/12] feat: align geomodels Grid and RPC API, clean RPC.from_any, update doc, clean geomodels --- docs/source/getting_started.rst | 2 +- docs/source/user_manual_conventions.rst | 6 - docs/source/user_manual_geometric_models.rst | 2 +- shareloc/geofunctions/localization.py | 3 +- shareloc/geofunctions/rectification_grid.py | 2 +- shareloc/geomodels/geomodel.py | 6 +- shareloc/geomodels/geomodel_template.py | 26 +- shareloc/geomodels/grid.py | 24 +- shareloc/geomodels/rpc.py | 355 +++++++++---------- tests/geofunctions/test_localization.py | 28 +- tests/geofunctions/test_rectification.py | 72 +--- tests/geofunctions/test_triangulation.py | 12 +- tests/geomodels/test_los.py | 4 +- tests/geomodels/test_rpc.py | 38 +- 14 files changed, 271 insertions(+), 309 deletions(-) diff --git a/docs/source/getting_started.rst b/docs/source/getting_started.rst index a95ec74..57b5726 100644 --- a/docs/source/getting_started.rst +++ b/docs/source/getting_started.rst @@ -40,7 +40,7 @@ Quick Start >>> # Create RPC object from downloaded geometry file >>> rpc_geom_file = "left_image.geom" - >>> rpc = RPC.from_any(rpc_geom_file) + >>> rpc = RPC(rpc_geom_file) >>> # Create Localization object from created RPC >>> loc = Localization(rpc) diff --git a/docs/source/user_manual_conventions.rst b/docs/source/user_manual_conventions.rst index a91bffb..14c2f76 100644 --- a/docs/source/user_manual_conventions.rst +++ b/docs/source/user_manual_conventions.rst @@ -19,12 +19,6 @@ Geometric model convention for coordinates is, by default, [0.5,0.5] at the cent pixel convention -When dealing with shareloc RPC (``shareloc.geomodels.rpc.RPC``), the center at [0,0] can be changed by setting the ``topleftconvention`` option to ``False``. - -.. code-block:: bash - - @classmethod - def from_any(cls, primary_file, secondary_file=None, topleftconvention=True): Image convention ================ diff --git a/docs/source/user_manual_geometric_models.rst b/docs/source/user_manual_geometric_models.rst index 0caa518..fa70023 100644 --- a/docs/source/user_manual_geometric_models.rst +++ b/docs/source/user_manual_geometric_models.rst @@ -52,7 +52,7 @@ RPC class API Example $ python3 >>> from shareloc.geomodels.rpc import RPC >>> file_dimap = "RPC_PHR1B_P_201709281038393_SEN_PRG_FC_178609-001.XML") - >>> rpc_dimap = RPC.from_any(file_dimap) + >>> rpc_dimap = RPC(file_dimap) Direct location grids diff --git a/shareloc/geofunctions/localization.py b/shareloc/geofunctions/localization.py index 6739fb0..c0f58c8 100755 --- a/shareloc/geofunctions/localization.py +++ b/shareloc/geofunctions/localization.py @@ -53,7 +53,7 @@ def __init__(self, model, elevation=None, image=None, epsg=None): :param epsg: coordinate system of world points, if None model coordiante system will be used :type epsg: int """ - self.use_rpc = model.type == "rpc" + self.use_rpc = model.type == "RPC" self.model = model self.default_elevation = 0.0 self.dtm = None @@ -137,6 +137,7 @@ def inverse(self, lon, lat, h=None, using_geotransform=False): """ if not self.use_rpc and not hasattr(self.model, "pred_ofset_scale_lon"): + # for grids only self.model.estimate_inverse_loc_predictor() if h is None: h = self.default_elevation diff --git a/shareloc/geofunctions/rectification_grid.py b/shareloc/geofunctions/rectification_grid.py index 53b02a5..c9a8aff 100755 --- a/shareloc/geofunctions/rectification_grid.py +++ b/shareloc/geofunctions/rectification_grid.py @@ -40,7 +40,7 @@ def __init__(self, grid_filename): :param grid_filename: grid filename :type filename: string """ - self.filename = grid_filename + self.grid_filename = grid_filename dataset = rio.open(grid_filename) diff --git a/shareloc/geomodels/geomodel.py b/shareloc/geomodels/geomodel.py index 450167c..e6ee48b 100644 --- a/shareloc/geomodels/geomodel.py +++ b/shareloc/geomodels/geomodel.py @@ -25,6 +25,8 @@ (from registered GeoModelTemplate class) ** Experimental and not used as API for now** +** Will be used with an automatic switcher between Grid file format and +RPC file format obj = GeoModel("file_geom") only without geomodel_type """ # Standard imports @@ -48,6 +50,8 @@ def __new__(cls, geomodel_path: str, geomodel_type: str = "RPC"): associated with the "geomodel_name" given in the configuration through create_geomodel local method for clarity. + TODO: optional geomodel_type would profit to have an automatic switcher between geomodels (RPC, grids, ...) + :param geomodel_path: Path of geomodel file :type geomodel_path: string :param geomodel_type: Type of the geomodel, default "rpc", used as key for specific geomodel subclass instance @@ -68,7 +72,7 @@ def create_geomodel(cls, geomodel_path: str, geomodel_type: str = "RPC"): :type geomodel_type: string """ - # If no type is given, the default is "rpc" + # If no type is given, the default is "RPC" # Create geomodel object with geomodel_path parameter from geomodel_type if exists try: diff --git a/shareloc/geomodels/geomodel_template.py b/shareloc/geomodels/geomodel_template.py index e1c6cdf..dca847c 100644 --- a/shareloc/geomodels/geomodel_template.py +++ b/shareloc/geomodels/geomodel_template.py @@ -26,12 +26,18 @@ # Standard imports from abc import ABCMeta, abstractmethod -from typing import Dict -# Third party imports +# Global variable for optimization mode (functions in C) +# SHARELOC_OPTIM_GEOMODEL = False +# TODO: Override functions depending on optimization or not + +# if(SHARELOC_OPTIM_GEOMODEL == True): +# GeoModelTemplate.direct_loc_dtm = GeoModelTemplate.direct_loc_dtm_optim # pylint: disable=too-few-public-methods + + class GeoModelTemplate(metaclass=ABCMeta): """ Class for general specification of a geometric model @@ -39,20 +45,16 @@ class GeoModelTemplate(metaclass=ABCMeta): """ @abstractmethod - def __init__(self, geomodel_file: str, geomodel_params: Dict = None): + def __init__(self, geomodel_path: str): """ Return the geomodel object associated with the geomodel_type given in the configuration - :param geomodel_file: path of the geomodel file to instanciate. - :type geomodel_file: string - :param geomodel_params: Geomodels parameters as dict (unused, just to fit RPC) TODO: evolve with new API - :type geomodel_params: Dict + :param geomodel_path: path of the geomodel file to instanciate. + :type geomodel_path: string """ # geomodel filename path - self.filename = geomodel_file - - # geomodel_params if exists (to evolve) - self.geomodels_params = geomodel_params + self.geomodel_path: str = geomodel_path - # + # geomodel type. Set by the subclass + self.type: str diff --git a/shareloc/geomodels/grid.py b/shareloc/geomodels/grid.py index cb9fbcb..819df8a 100755 --- a/shareloc/geomodels/grid.py +++ b/shareloc/geomodels/grid.py @@ -19,9 +19,8 @@ # limitations under the License. # """ -localisation functions from multi h direct grids. +Localisation functions from multi h direct grids. """ -# pylint: disable=no-member # Standard imports import logging @@ -39,14 +38,15 @@ # gitlab issue #58 # pylint: disable=too-many-instance-attributes +# pylint: disable=no-member @GeoModel.register("grid") class Grid(GeoModelTemplate): """ multi H direct localization grid handling class. please refer to the main documentation grid format - :param filename: grid path - :type filename: str + Derives from GeoModelTemplate + :param row0: grid first pixel center along Y axis (row). :type row0: float :param col0: grid first pixel center along X axis (column). @@ -79,15 +79,18 @@ class Grid(GeoModelTemplate): :type type: str """ - def __init__(self, grid_filename): + def __init__(self, geomodel_path: str): """ Grid Constructor - :param grid_filename: grid filename (Geotiff) - :type grid_filename: string + :param geomodel_path: grid filename (Geotiff) + :type geomodel_path: string """ # Instanciate GeoModelTemplate generic init with shared parameters - super().__init__(grid_filename) + super().__init__(geomodel_path) + self.type = "multi H grid" + + # GeoModel Grid parameters definition (see documentation) self.row0 = None self.col0 = None self.nbrow = None @@ -102,8 +105,9 @@ def __init__(self, grid_filename): self.rowmax = None self.colmax = None self.epsg = 0 + + # Load function of grid parameters from grid file to grid object self.load() - self.type = "multi H grid" def load(self): """ @@ -116,7 +120,7 @@ def load(self): - lat_data : [alt,row,col] """ - grid_image = Image(self.filename, read_data=True) + grid_image = Image(self.geomodel_path, read_data=True) if grid_image.dataset.driver != "GTiff": raise TypeError( "Only Geotiff grids are accepted. Please refer to the documentation for grid supported format." diff --git a/shareloc/geomodels/rpc.py b/shareloc/geomodels/rpc.py index 1b79abb..e881e82 100755 --- a/shareloc/geomodels/rpc.py +++ b/shareloc/geomodels/rpc.py @@ -23,11 +23,11 @@ This module contains the RPC class corresponding to the RPC models. RPC models covered are : DIMAP V1, DIMAP V2, DIMAP V3, ossim (geom file), geotiff. """ -# pylint: disable=no-member # Standard imports import logging from os.path import basename +from typing import Dict from xml.dom import minidom # Third party imports @@ -121,6 +121,9 @@ def identify_geotiff_rpc(image_filename): return None +# pylint: disable=no-member + + @GeoModel.register("RPC") class RPC(GeoModelTemplate): """ @@ -129,13 +132,25 @@ class RPC(GeoModelTemplate): # gitlab issue #61 # pylint: disable=too-many-instance-attributes - def __init__(self, rpc_params): + def __init__(self, geomodel_path: str): + # Instanciate GeoModelTemplate generic init with shared parameters + super().__init__(geomodel_path) + self.type = "RPC" + + # initiate epsg and datum, can overriden by rpc_params self.epsg = None self.datum = None - for key, value in rpc_params.items(): + + # RPC parameters as Dict + self.rpc_params: Dict = {} + + # RPC parameters are load from geomodel_path to rpc params + self.load() + + # set class parameters from rpc_params (epsg and datum can be overriden) + for key, value in self.rpc_params.items(): setattr(self, key, value) - self.type = "rpc" if self.epsg is None: self.epsg = 4326 if self.datum is None: @@ -143,31 +158,32 @@ def __init__(self, rpc_params): self.lim_extrapol = 1.0001 + # Monomes seems not used in shareloc code: Clean ? # Each monome: c[0]*X**c[1]*Y**c[2]*Z**c[3] - monomes_order = [ - [1, 0, 0, 0], - [1, 1, 0, 0], - [1, 0, 1, 0], - [1, 0, 0, 1], - [1, 1, 1, 0], - [1, 1, 0, 1], - [1, 0, 1, 1], - [1, 2, 0, 0], - [1, 0, 2, 0], - [1, 0, 0, 2], - [1, 1, 1, 1], - [1, 3, 0, 0], - [1, 1, 2, 0], - [1, 1, 0, 2], - [1, 2, 1, 0], - [1, 0, 3, 0], - [1, 0, 1, 2], - [1, 2, 0, 1], - [1, 0, 2, 1], - [1, 0, 0, 3], - ] - - self.monomes = np.array(monomes_order) + self.monomes = np.array( + [ + [1, 0, 0, 0], + [1, 1, 0, 0], + [1, 0, 1, 0], + [1, 0, 0, 1], + [1, 1, 1, 0], + [1, 1, 0, 1], + [1, 0, 1, 1], + [1, 2, 0, 0], + [1, 0, 2, 0], + [1, 0, 0, 2], + [1, 1, 1, 1], + [1, 3, 0, 0], + [1, 1, 2, 0], + [1, 1, 0, 2], + [1, 2, 1, 0], + [1, 0, 3, 0], + [1, 0, 1, 2], + [1, 2, 0, 1], + [1, 0, 2, 1], + [1, 0, 0, 3], + ] + ) # monomial coefficients of 1st variable derivative self.monomes_deriv_1 = np.array( @@ -246,153 +262,164 @@ def __init__(self, rpc_params): self.row0 = self.offset_row - self.scale_row self.rowmax = self.offset_row + self.scale_row - @classmethod - def from_dimap(cls, dimap_filepath, topleftconvention=True): + def load(self): """ - Load from Dimap + Load from any RPC (auto identify driver) + from filename (dimap, ossim kwl, geotiff) - param dimap_filepath: dimap xml file - :type dimap_filepath: str - :param topleftconvention: [0,0] position - :type topleftconvention: boolean - If False : [0,0] is at the center of the Top Left pixel - If True : [0,0] is at the top left of the Top Left pixel (OSSIM) - """ - dimap_version = identify_dimap(dimap_filepath) - if identify_dimap(dimap_filepath) is not None: - if float(dimap_version) < 2.0: - return cls.from_dimap_v1(dimap_filepath, topleftconvention) - if float(dimap_version) >= 2.0: - return cls.read_dimap_coeff(dimap_filepath, topleftconvention) - else: - raise ValueError("can''t read dimap file") + TODO: topleftconvention always to True, set a standard and remove the option - return None + topleftconvention boolean: [0,0] position + If False : [0,0] is at the center of the Top Left pixel + If True : [0,0] is at the top left of the Top Left pixel (OSSIM) + """ + # Set topleftconvention (keeping historic option): to clean + topleftconvention = True - @classmethod - def read_dimap_coeff(cls, dimap_filepath, topleftconvention=True): + # If ends with XML --> DIMAP + if basename(self.geomodel_path.upper()).endswith("XML"): + dimap_version = identify_dimap(self.geomodel_path) + if dimap_version is not None: + if float(dimap_version) < 2.0: + self.load_dimap_v1(topleftconvention) + if float(dimap_version) >= 2.0: + self.load_dimap_coeff(topleftconvention) + else: + # If not DIMAP, identify ossim + ossim_model = identify_ossim_kwl(self.geomodel_path) + if ossim_model is not None: + self.load_ossim_kwl(topleftconvention) + else: + # Otherwise, check if RPC is in geotif + geotiff_rpc_dict = identify_geotiff_rpc(self.geomodel_path) + if geotiff_rpc_dict is not None: + self.load_geotiff(topleftconvention) + else: + # Raise error if no file recognized. + raise ValueError("can not read rpc file") + + def load_dimap_coeff(self, topleftconvention=True): """ Load from Dimap v2 and V3 - :param dimap_filepath: dimap xml file - :type dimap_filepath: str :param topleftconvention: [0,0] position :type topleftconvention: boolean If False : [0,0] is at the center of the Top Left pixel If True : [0,0] is at the top left of the Top Left pixel (OSSIM) """ - rpc_params = {} - if not basename(dimap_filepath).upper().endswith("XML"): + if not basename(self.geomodel_path).upper().endswith("XML"): raise ValueError("dimap must ends with .xml") - xmldoc = minidom.parse(dimap_filepath) + xmldoc = minidom.parse(self.geomodel_path) mtd = xmldoc.getElementsByTagName("Metadata_Identification") version = mtd[0].getElementsByTagName("METADATA_FORMAT")[0].attributes.items()[0][1] - rpc_params["driver_type"] = "dimap_v" + version + self.rpc_params["driver_type"] = "dimap_v" + version global_rfm = xmldoc.getElementsByTagName("Global_RFM")[0] normalisation_coeffs = global_rfm.getElementsByTagName("RFM_Validity")[0] - rpc_params["offset_row"] = float(normalisation_coeffs.getElementsByTagName("LINE_OFF")[0].firstChild.data) - rpc_params["offset_col"] = float(normalisation_coeffs.getElementsByTagName("SAMP_OFF")[0].firstChild.data) + self.rpc_params["offset_row"] = float(normalisation_coeffs.getElementsByTagName("LINE_OFF")[0].firstChild.data) + self.rpc_params["offset_col"] = float(normalisation_coeffs.getElementsByTagName("SAMP_OFF")[0].firstChild.data) if float(version) >= 3: direct_coeffs = global_rfm.getElementsByTagName("ImagetoGround_Values")[0] inverse_coeffs = global_rfm.getElementsByTagName("GroundtoImage_Values")[0] - rpc_params["num_x"] = [ + self.rpc_params["num_x"] = [ float(direct_coeffs.getElementsByTagName(f"LON_NUM_COEFF_{index}")[0].firstChild.data) for index in range(1, 21) ] - rpc_params["den_x"] = [ + self.rpc_params["den_x"] = [ float(direct_coeffs.getElementsByTagName(f"LON_DEN_COEFF_{index}")[0].firstChild.data) for index in range(1, 21) ] - rpc_params["num_y"] = [ + self.rpc_params["num_y"] = [ float(direct_coeffs.getElementsByTagName(f"LAT_NUM_COEFF_{index}")[0].firstChild.data) for index in range(1, 21) ] - rpc_params["den_y"] = [ + self.rpc_params["den_y"] = [ float(direct_coeffs.getElementsByTagName(f"LAT_DEN_COEFF_{index}")[0].firstChild.data) for index in range(1, 21) ] - rpc_params["offset_col"] -= 0.5 - rpc_params["offset_row"] -= 0.5 + self.rpc_params["offset_col"] -= 0.5 + self.rpc_params["offset_row"] -= 0.5 else: direct_coeffs = global_rfm.getElementsByTagName("Direct_Model")[0] inverse_coeffs = global_rfm.getElementsByTagName("Inverse_Model")[0] - rpc_params["num_x"] = [ + self.rpc_params["num_x"] = [ float(direct_coeffs.getElementsByTagName(f"SAMP_NUM_COEFF_{index}")[0].firstChild.data) for index in range(1, 21) ] - rpc_params["den_x"] = [ + self.rpc_params["den_x"] = [ float(direct_coeffs.getElementsByTagName(f"SAMP_DEN_COEFF_{index}")[0].firstChild.data) for index in range(1, 21) ] - rpc_params["num_y"] = [ + self.rpc_params["num_y"] = [ float(direct_coeffs.getElementsByTagName(f"LINE_NUM_COEFF_{index}")[0].firstChild.data) for index in range(1, 21) ] - rpc_params["den_y"] = [ + self.rpc_params["den_y"] = [ float(direct_coeffs.getElementsByTagName(f"LINE_DEN_COEFF_{index}")[0].firstChild.data) for index in range(1, 21) ] - rpc_params["offset_col"] -= 1.0 - rpc_params["offset_row"] -= 1.0 + self.rpc_params["offset_col"] -= 1.0 + self.rpc_params["offset_row"] -= 1.0 - rpc_params["num_col"] = [ + self.rpc_params["num_col"] = [ float(inverse_coeffs.getElementsByTagName(f"SAMP_NUM_COEFF_{index}")[0].firstChild.data) for index in range(1, 21) ] - rpc_params["den_col"] = [ + self.rpc_params["den_col"] = [ float(inverse_coeffs.getElementsByTagName(f"SAMP_DEN_COEFF_{index}")[0].firstChild.data) for index in range(1, 21) ] - rpc_params["num_row"] = [ + self.rpc_params["num_row"] = [ float(inverse_coeffs.getElementsByTagName(f"LINE_NUM_COEFF_{index}")[0].firstChild.data) for index in range(1, 21) ] - rpc_params["den_row"] = [ + self.rpc_params["den_row"] = [ float(inverse_coeffs.getElementsByTagName(f"LINE_DEN_COEFF_{index}")[0].firstChild.data) for index in range(1, 21) ] - rpc_params["scale_col"] = float(normalisation_coeffs.getElementsByTagName("SAMP_SCALE")[0].firstChild.data) - rpc_params["scale_row"] = float(normalisation_coeffs.getElementsByTagName("LINE_SCALE")[0].firstChild.data) - rpc_params["offset_alt"] = float(normalisation_coeffs.getElementsByTagName("HEIGHT_OFF")[0].firstChild.data) - rpc_params["scale_alt"] = float(normalisation_coeffs.getElementsByTagName("HEIGHT_SCALE")[0].firstChild.data) - rpc_params["offset_x"] = float(normalisation_coeffs.getElementsByTagName("LONG_OFF")[0].firstChild.data) - rpc_params["scale_x"] = float(normalisation_coeffs.getElementsByTagName("LONG_SCALE")[0].firstChild.data) - rpc_params["offset_y"] = float(normalisation_coeffs.getElementsByTagName("LAT_OFF")[0].firstChild.data) - rpc_params["scale_y"] = float(normalisation_coeffs.getElementsByTagName("LAT_SCALE")[0].firstChild.data) + self.rpc_params["scale_col"] = float(normalisation_coeffs.getElementsByTagName("SAMP_SCALE")[0].firstChild.data) + self.rpc_params["scale_row"] = float(normalisation_coeffs.getElementsByTagName("LINE_SCALE")[0].firstChild.data) + self.rpc_params["offset_alt"] = float( + normalisation_coeffs.getElementsByTagName("HEIGHT_OFF")[0].firstChild.data + ) + self.rpc_params["scale_alt"] = float( + normalisation_coeffs.getElementsByTagName("HEIGHT_SCALE")[0].firstChild.data + ) + self.rpc_params["offset_x"] = float(normalisation_coeffs.getElementsByTagName("LONG_OFF")[0].firstChild.data) + self.rpc_params["scale_x"] = float(normalisation_coeffs.getElementsByTagName("LONG_SCALE")[0].firstChild.data) + self.rpc_params["offset_y"] = float(normalisation_coeffs.getElementsByTagName("LAT_OFF")[0].firstChild.data) + self.rpc_params["scale_y"] = float(normalisation_coeffs.getElementsByTagName("LAT_SCALE")[0].firstChild.data) # If top left convention, 0.5 pixel shift added on col/row offsets if topleftconvention: - rpc_params["offset_col"] += 0.5 - rpc_params["offset_row"] += 0.5 - return cls(rpc_params) + self.rpc_params["offset_col"] += 0.5 + self.rpc_params["offset_row"] += 0.5 - @classmethod - def from_dimap_v1(cls, dimap_filepath, topleftconvention=True): + def load_dimap_v1(self, topleftconvention=True): """ Load from dimap v1 - :param dimap_filepath: dimap xml file - :type dimap_filepath: str + ** Deprecated, to clean ? ** + :param topleftconvention: [0,0] position :type topleftconvention: boolean If False : [0,0] is at the center of the Top Left pixel If True : [0,0] is at the top left of the Top Left pixel (OSSIM) """ - if not basename(dimap_filepath).upper().endswith("XML"): + if not basename(self.geomodel_path).upper().endswith("XML"): raise ValueError("dimap must ends with .xml") - xmldoc = minidom.parse(dimap_filepath) + xmldoc = minidom.parse(self.geomodel_path) mtd = xmldoc.getElementsByTagName("Metadata_Identification") version = mtd[0].getElementsByTagName("METADATA_PROFILE")[0].attributes.items()[0][1] - rpc_params = {"driver_type": "dimap_v" + version} + self.rpc_params = {"driver_type": "dimap_v" + version} global_rfm = xmldoc.getElementsByTagName("Global_RFM") rfm_validity = xmldoc.getElementsByTagName("RFM_Validity") @@ -412,50 +439,47 @@ def from_dimap_v1(cls, dimap_filepath, topleftconvention=True): scale_row = float(rfm_validity[0].getElementsByTagName("Row")[0].getElementsByTagName("A")[0].firstChild.data) offset_row = float(rfm_validity[0].getElementsByTagName("Row")[0].getElementsByTagName("B")[0].firstChild.data) - rpc_params["offset_col"] = offset_col - rpc_params["scale_col"] = scale_col - rpc_params["offset_row"] = offset_row - rpc_params["scale_row"] = scale_row - rpc_params["offset_alt"] = offset_alt - rpc_params["scale_alt"] = scale_alt - rpc_params["offset_x"] = offset_lon - rpc_params["scale_x"] = scale_lon - rpc_params["offset_y"] = offset_lat - rpc_params["scale_y"] = scale_lat - rpc_params["num_x"] = coeff_lon[0:20] - rpc_params["den_x"] = coeff_lon[20::] - rpc_params["num_y"] = coeff_lat[0:20] - rpc_params["den_y"] = coeff_lat[20::] - rpc_params["num_col"] = coeff_col[0:20] - rpc_params["den_col"] = coeff_col[20::] - rpc_params["num_row"] = coeff_lig[0:20] - rpc_params["den_row"] = coeff_lig[20::] - rpc_params["offset_col"] -= 1.0 - rpc_params["offset_row"] -= 1.0 + self.rpc_params["offset_col"] = offset_col + self.rpc_params["scale_col"] = scale_col + self.rpc_params["offset_row"] = offset_row + self.rpc_params["scale_row"] = scale_row + self.rpc_params["offset_alt"] = offset_alt + self.rpc_params["scale_alt"] = scale_alt + self.rpc_params["offset_x"] = offset_lon + self.rpc_params["scale_x"] = scale_lon + self.rpc_params["offset_y"] = offset_lat + self.rpc_params["scale_y"] = scale_lat + self.rpc_params["num_x"] = coeff_lon[0:20] + self.rpc_params["den_x"] = coeff_lon[20::] + self.rpc_params["num_y"] = coeff_lat[0:20] + self.rpc_params["den_y"] = coeff_lat[20::] + self.rpc_params["num_col"] = coeff_col[0:20] + self.rpc_params["den_col"] = coeff_col[20::] + self.rpc_params["num_row"] = coeff_lig[0:20] + self.rpc_params["den_row"] = coeff_lig[20::] + self.rpc_params["offset_col"] -= 1.0 + self.rpc_params["offset_row"] -= 1.0 # If top left convention, 0.5 pixel shift added on col/row offsets if topleftconvention: - rpc_params["offset_col"] += 0.5 - rpc_params["offset_row"] += 0.5 - return cls(rpc_params) + self.rpc_params["offset_col"] += 0.5 + self.rpc_params["offset_row"] += 0.5 - @classmethod - def from_geotiff(cls, image_filename, topleftconvention=True): + def load_geotiff(self, topleftconvention=True): """ Load from a geotiff image file - :param image_filename: image filename - :type image_filename: str + :param topleftconvention: [0,0] position :type topleftconvention: boolean If False : [0,0] is at the center of the Top Left pixel If True : [0,0] is at the top left of the Top Left pixel (OSSIM) """ - dataset = rio.open(image_filename) + dataset = rio.open(self.geomodel_path) rpc_dict = dataset.tags(ns="RPC") if not rpc_dict: - logging.error("%s does not contains RPCs ", image_filename) + logging.error("%s does not contains RPCs ", self.geomodel_path) raise ValueError - rpc_params = { + self.rpc_params = { "den_row": parse_coeff_line(rpc_dict["LINE_DEN_COEFF"]), "num_row": parse_coeff_line(rpc_dict["LINE_NUM_COEFF"]), "num_col": parse_coeff_line(rpc_dict["SAMP_NUM_COEFF"]), @@ -478,12 +502,10 @@ def from_geotiff(cls, image_filename, topleftconvention=True): # inverse coeff are not defined # If top left convention, 0.5 pixel shift added on col/row offsets if topleftconvention: - rpc_params["offset_col"] += 0.5 - rpc_params["offset_row"] += 0.5 - return cls(rpc_params) + self.rpc_params["offset_col"] += 0.5 + self.rpc_params["offset_row"] += 0.5 - @classmethod - def from_ossim_kwl(cls, ossim_kwl_filename, topleftconvention=True): + def load_ossim_kwl(self, topleftconvention=True): """ Load from a geom file @@ -492,10 +514,9 @@ def from_ossim_kwl(cls, ossim_kwl_filename, topleftconvention=True): If False : [0,0] is at the center of the Top Left pixel If True : [0,0] is at the top left of the Top Left pixel (OSSIM) """ - rpc_params = {} # OSSIM keyword list - rpc_params["driver_type"] = "ossim_kwl" - with open(ossim_kwl_filename, "r", encoding="utf-8") as ossim_file: + self.rpc_params["driver_type"] = "ossim_kwl" + with open(self.geomodel_path, "r", encoding="utf-8") as ossim_file: content = ossim_file.readlines() geom_dict = {} @@ -503,71 +524,43 @@ def from_ossim_kwl(cls, ossim_kwl_filename, topleftconvention=True): (key, val) = line.split(": ") geom_dict[key] = val.rstrip() - rpc_params["den_row"] = [np.nan] * 20 - rpc_params["num_row"] = [np.nan] * 20 - rpc_params["den_col"] = [np.nan] * 20 - rpc_params["num_col"] = [np.nan] * 20 + self.rpc_params["den_row"] = [np.nan] * 20 + self.rpc_params["num_row"] = [np.nan] * 20 + self.rpc_params["den_col"] = [np.nan] * 20 + self.rpc_params["num_col"] = [np.nan] * 20 for index in range(0, 20): axis = "line" num_den = "den" key = f"{axis}_{num_den}_coeff_{index:02d}" - rpc_params["den_row"][index] = float(geom_dict[key]) + self.rpc_params["den_row"][index] = float(geom_dict[key]) num_den = "num" key = f"{axis}_{num_den}_coeff_{index:02d}" - rpc_params["num_row"][index] = float(geom_dict[key]) + self.rpc_params["num_row"][index] = float(geom_dict[key]) axis = "samp" key = f"{axis}_{num_den}_coeff_{index:02d}" - rpc_params["num_col"][index] = float(geom_dict[key]) + self.rpc_params["num_col"][index] = float(geom_dict[key]) num_den = "den" key = f"{axis}_{num_den}_coeff_{index:02d}" - rpc_params["den_col"][index] = float(geom_dict[key]) - rpc_params["offset_col"] = float(geom_dict["samp_off"]) - rpc_params["scale_col"] = float(geom_dict["samp_scale"]) - rpc_params["offset_row"] = float(geom_dict["line_off"]) - rpc_params["scale_row"] = float(geom_dict["line_scale"]) - rpc_params["offset_alt"] = float(geom_dict["height_off"]) - rpc_params["scale_alt"] = float(geom_dict["height_scale"]) - rpc_params["offset_x"] = float(geom_dict["long_off"]) - rpc_params["scale_x"] = float(geom_dict["long_scale"]) - rpc_params["offset_y"] = float(geom_dict["lat_off"]) - rpc_params["scale_y"] = float(geom_dict["lat_scale"]) + self.rpc_params["den_col"][index] = float(geom_dict[key]) + self.rpc_params["offset_col"] = float(geom_dict["samp_off"]) + self.rpc_params["scale_col"] = float(geom_dict["samp_scale"]) + self.rpc_params["offset_row"] = float(geom_dict["line_off"]) + self.rpc_params["scale_row"] = float(geom_dict["line_scale"]) + self.rpc_params["offset_alt"] = float(geom_dict["height_off"]) + self.rpc_params["scale_alt"] = float(geom_dict["height_scale"]) + self.rpc_params["offset_x"] = float(geom_dict["long_off"]) + self.rpc_params["scale_x"] = float(geom_dict["long_scale"]) + self.rpc_params["offset_y"] = float(geom_dict["lat_off"]) + self.rpc_params["scale_y"] = float(geom_dict["lat_scale"]) # inverse coeff are not defined - rpc_params["num_x"] = None - rpc_params["den_x"] = None - rpc_params["num_y"] = None - rpc_params["den_y"] = None + self.rpc_params["num_x"] = None + self.rpc_params["den_x"] = None + self.rpc_params["num_y"] = None + self.rpc_params["den_y"] = None # If top left convention, 0.5 pixel shift added on col/row offsets if topleftconvention: - rpc_params["offset_col"] += 0.5 - rpc_params["offset_row"] += 0.5 - return cls(rpc_params) - - @classmethod - def from_any(cls, primary_file, topleftconvention=True): - """ - Load from any RPC (auto identify driver) - - :param primary_file: rpc filename (dimap, ossim kwl, geotiff) - :type primary_file: str - :param topleftconvention: [0,0] position - :type topleftconvention: boolean - If False : [0,0] is at the center of the Top Left pixel - If True : [0,0] is at the top left of the Top Left pixel (OSSIM) - """ - if basename(primary_file.upper()).endswith("XML"): - dimap_version = identify_dimap(primary_file) - if dimap_version is not None: - if float(dimap_version) < 2.0: - return cls.from_dimap_v1(primary_file, topleftconvention) - if float(dimap_version) >= 2.0: - return cls.read_dimap_coeff(primary_file, topleftconvention) - ossim_model = identify_ossim_kwl(primary_file) - if ossim_model is not None: - return cls.from_ossim_kwl(primary_file, topleftconvention) - geotiff_rpc_dict = identify_geotiff_rpc(primary_file) - if geotiff_rpc_dict is not None: - return cls.from_geotiff(primary_file, topleftconvention) - raise ValueError("can not read rpc file") + self.rpc_params["offset_col"] += 0.5 + self.rpc_params["offset_row"] += 0.5 def direct_loc_h(self, row, col, alt, fill_nan=False): """ diff --git a/tests/geofunctions/test_localization.py b/tests/geofunctions/test_localization.py index 3a45592..e6367fc 100644 --- a/tests/geofunctions/test_localization.py +++ b/tests/geofunctions/test_localization.py @@ -57,9 +57,9 @@ def test_localize_direct_rpc(): # first instanciate the RPC geometric model # data = os.path.join(data_path(), "rpc/phr_ventoux/", "left_image.geom") - # geom_model = RPC.from_any(data) + # geom_model = RPC(data) data = os.path.join(data_path(), "rpc/phr_ventoux/", "RPC_PHR1B_P_201308051042194_SEN_690908101-001.XML") - geom_model = RPC.from_any(data) + geom_model = RPC(data) # then read the Image to retrieve its geotransform image_filename = os.path.join(data_path(), "image/phr_ventoux/", "left_image.tif") image_left = Image(image_filename) @@ -94,7 +94,7 @@ def test_localize_direct_grid(): # first instanciate the Grid geometric model # data = os.path.join(data_path(), "rpc/phr_ventoux/", "left_image.geom") - # geom_model_1 = RPC.from_any(data) + # geom_model_1 = RPC(data) data = os.path.join(data_path(), "grid/phr_ventoux/", "GRID_PHR1B_P_201308051042194_SEN_690908101-001.tif") geom_model = Grid(data) # then read the Image to retrieve its geotransform @@ -230,7 +230,7 @@ def test_extent(): Test extent """ data_left = os.path.join(data_path(), "rectification", "left_image") - geom_model = RPC.from_any(data_left + ".geom", topleftconvention=True) + geom_model = RPC(data_left + ".geom") image_filename = os.path.join(data_path(), "image/phr_ventoux/", "left_image_pixsize_0_5.tif") image = Image(image_filename) loc_rpc_image = Localization(geom_model, elevation=None, image=image) @@ -246,7 +246,7 @@ def test_sensor_loc_dir_dtm_geoid(col, row, valid_coord): Test direct localization using image geotransform """ data = os.path.join(data_path(), "rectification", "left_image") - geom_model_left = RPC.from_any(data + ".geom", topleftconvention=True) + geom_model_left = RPC(data + ".geom") image_filename = os.path.join(data_path(), "image/phr_ventoux/", "left_image.tif") image_left = Image(image_filename) @@ -267,7 +267,7 @@ def test_sensor_loc_dir_dtm_geoid_utm(col, row, valid_coord): Test direct localization using image geotransform """ data = os.path.join(data_path(), "rectification", "left_image") - geom_model_left = RPC.from_any(data + ".geom", topleftconvention=True) + geom_model_left = RPC(data + ".geom") image_filename = os.path.join(data_path(), "image/phr_ventoux/", "left_image.tif") image_left = Image(image_filename) @@ -297,7 +297,7 @@ def test_sensor_loc_dir_vs_loc_rpc(row, col, h): data_folder = data_path() fichier_dimap = os.path.join(data_folder, f"rpc/PHRDIMAP_{id_scene}.XML") - fctrat = RPC.from_any(fichier_dimap) + fctrat = RPC(fichier_dimap) loc_rpc = Localization(fctrat) lonlatalt_rpc = loc_rpc.direct(row, col, h) @@ -392,7 +392,7 @@ def test_sensor_loc_inv_vs_loc_rpc(lon, lat, alt): data_folder = data_path() fichier_dimap = os.path.join(data_folder, f"rpc/PHRDIMAP_{id_scene}.XML") - fctrat = RPC.from_any(fichier_dimap, topleftconvention=True) + fctrat = RPC(fichier_dimap) loc_rpc = Localization(fctrat) [row_rpc, col_rpc, __] = loc_rpc.inverse(lon, lat, alt) @@ -508,7 +508,7 @@ def test_loc_dir_loc_inv_rpc(id_scene, rpc, row, col, h): data_folder = data_path() fichier_dimap = os.path.join(data_folder, "rpc", rpc) - fctrat = RPC.from_any(fichier_dimap, topleftconvention=True) + fctrat = RPC(fichier_dimap) (inv_row, inv_col, __) = fctrat.inverse_loc(lonlatalt[0][0], lonlatalt[0][1], lonlatalt[0][2]) assert row == pytest.approx(inv_row, abs=1e-2) assert col == pytest.approx(inv_col, abs=1e-2) @@ -560,7 +560,7 @@ def test_colocalization(col, row, alt): data_folder = data_path() id_scene = "P1BP--2018122638935449CP" file_dimap = os.path.join(data_folder, f"rpc/PHRDIMAP_{id_scene}.XML") - fctrat = RPC.from_dimap_v1(file_dimap) + fctrat = RPC(file_dimap) row_coloc, col_coloc, _ = coloc_rpc(fctrat, fctrat, row, col, alt) @@ -575,12 +575,12 @@ def test_sensor_coloc_using_geotransform(col, row, h): Test direct localization using image geotransform """ data_left = os.path.join(data_path(), "rectification", "left_image") - geom_model_left = RPC.from_any(data_left + ".geom", topleftconvention=True) + geom_model_left = RPC(data_left + ".geom") image_filename_left = os.path.join(data_path(), "image/phr_ventoux/", "left_image_pixsize_0_5.tif") image_left = Image(image_filename_left) data_right = os.path.join(data_path(), "rectification", "right_image") - geom_model_right = RPC.from_any(data_right + ".geom", topleftconvention=True) + geom_model_right = RPC(data_right + ".geom") image_filename_right = os.path.join(data_path(), "image/phr_ventoux/", "right_image_pixsize_0_5.tif") image_right = Image(image_filename_right) @@ -610,7 +610,7 @@ def test_sensor_loc_utm(col, row): Test direct localization using image geotransform """ data_left = os.path.join(data_path(), "rectification", "left_image") - geom_model = RPC.from_any(data_left + ".geom", topleftconvention=True) + geom_model = RPC(data_left + ".geom") epsg = 32631 loc_wgs = Localization(geom_model) loc_utm = Localization(geom_model, epsg=epsg) @@ -633,7 +633,7 @@ def test_sensor_loc_dir_dtm_multi_points(): left_im = Image(os.path.join(data_path(), "rectification", "left_image.tif")) - geom_model = RPC.from_any(os.path.join(data_path(), "rectification", "left_image.geom"), topleftconvention=True) + geom_model = RPC(os.path.join(data_path(), "rectification", "left_image.geom")) dtm_file = os.path.join(data_path(), "dtm", "srtm_ventoux", "srtm90_non_void_filled", "N44E005.hgt") geoid_file = os.path.join(data_path(), "dtm", "geoid", "egm96_15.gtx") diff --git a/tests/geofunctions/test_rectification.py b/tests/geofunctions/test_rectification.py index 26a9a28..a6e9556 100644 --- a/tests/geofunctions/test_rectification.py +++ b/tests/geofunctions/test_rectification.py @@ -62,12 +62,8 @@ def test_compute_stereorectification_epipolar_grids_geomodel_rpc(): left_im = Image(os.path.join(data_path(), "rectification", "left_image.tif")) right_im = Image(os.path.join(data_path(), "rectification", "right_image.tif")) - geom_model_left = RPC.from_any( - os.path.join(data_path(), "rectification", "left_image.geom"), topleftconvention=True - ) - geom_model_right = RPC.from_any( - os.path.join(data_path(), "rectification", "right_image.geom"), topleftconvention=True - ) + geom_model_left = RPC(os.path.join(data_path(), "rectification", "left_image.geom")) + geom_model_right = RPC(os.path.join(data_path(), "rectification", "right_image.geom")) epi_step = 30 elevation_offset = 50 @@ -113,12 +109,8 @@ def test_compute_stereorectification_epipolar_grids_geomodel_rpc_alti(): left_im = Image(os.path.join(data_path(), "rectification", "left_image.tif")) right_im = Image(os.path.join(data_path(), "rectification", "right_image.tif")) - geom_model_left = RPC.from_any( - os.path.join(data_path(), "rectification", "left_image.geom"), topleftconvention=True - ) - geom_model_right = RPC.from_any( - os.path.join(data_path(), "rectification", "right_image.geom"), topleftconvention=True - ) + geom_model_left = RPC(os.path.join(data_path(), "rectification", "left_image.geom")) + geom_model_right = RPC(os.path.join(data_path(), "rectification", "right_image.geom")) epi_step = 30 elevation_offset = 50 @@ -158,12 +150,8 @@ def test_compute_stereorectification_epipolar_grids_geomodel_rpc_dtm_geoid(): """ # first instantiate geometric models left and right (here RPC geometrics model) - geom_model_left = RPC.from_any( - os.path.join(data_path(), "rectification", "left_image.geom"), topleftconvention=True - ) - geom_model_right = RPC.from_any( - os.path.join(data_path(), "rectification", "right_image.geom"), topleftconvention=True - ) + geom_model_left = RPC(os.path.join(data_path(), "rectification", "left_image.geom")) + geom_model_right = RPC(os.path.join(data_path(), "rectification", "right_image.geom")) # read the images left_im = Image(os.path.join(data_path(), "rectification", "left_image.tif")) @@ -218,12 +206,8 @@ def test_compute_stereorectification_epipolar_grids_geomodel_rpc_dtm_geoid_roi() left_im = Image(os.path.join(data_path(), "rectification", "left_image.tif")) right_im = Image(os.path.join(data_path(), "rectification", "right_image.tif")) - geom_model_left = RPC.from_any( - os.path.join(data_path(), "rectification", "left_image.geom"), topleftconvention=True - ) - geom_model_right = RPC.from_any( - os.path.join(data_path(), "rectification", "right_image.geom"), topleftconvention=True - ) + geom_model_left = RPC(os.path.join(data_path(), "rectification", "left_image.geom")) + geom_model_right = RPC(os.path.join(data_path(), "rectification", "right_image.geom")) dtm_file = os.path.join(data_path(), "dtm", "srtm_ventoux", "srtm90_non_void_filled", "N44E005.hgt") geoid_file = os.path.join(data_path(), "dtm", "geoid", "egm96_15.gtx") @@ -438,12 +422,8 @@ def test_prepare_rectification(): """ left_im = Image(os.path.join(data_path(), "rectification", "left_image.tif")) - geom_model_left = RPC.from_any( - os.path.join(data_path(), "rectification", "left_image.geom"), topleftconvention=True - ) - geom_model_right = RPC.from_any( - os.path.join(data_path(), "rectification", "right_image.geom"), topleftconvention=True - ) + geom_model_left = RPC(os.path.join(data_path(), "rectification", "left_image.geom")) + geom_model_right = RPC(os.path.join(data_path(), "rectification", "right_image.geom")) epi_step = 30 elevation_offset = 50 @@ -478,12 +458,8 @@ def test_prepare_rectification_footprint(): """ left_im = Image(os.path.join(data_path(), "rectification", "left_image.tif")) - geom_model_left = RPC.from_any( - os.path.join(data_path(), "rectification", "left_image.geom"), topleftconvention=True - ) - geom_model_right = RPC.from_any( - os.path.join(data_path(), "rectification", "right_image.geom"), topleftconvention=True - ) + geom_model_left = RPC(os.path.join(data_path(), "rectification", "left_image.geom")) + geom_model_right = RPC(os.path.join(data_path(), "rectification", "right_image.geom")) epi_step = 30 elevation_offset = 50 @@ -509,12 +485,8 @@ def test_rectification_moving_along_lines(): """ Test moving along line in epipolar geometry """ - geom_model_left = RPC.from_any( - os.path.join(data_path(), "rectification", "left_image.geom"), topleftconvention=True - ) - geom_model_right = RPC.from_any( - os.path.join(data_path(), "rectification", "right_image.geom"), topleftconvention=True - ) + geom_model_left = RPC(os.path.join(data_path(), "rectification", "left_image.geom")) + geom_model_right = RPC(os.path.join(data_path(), "rectification", "right_image.geom")) current_coords = np.array([[5000.5, 5000.5, 0.0]], dtype=np.float64) mean_spacing = 1 @@ -545,12 +517,8 @@ def test_rectification_moving_to_next_line(): """ Test moving to next line in epipolar geometry """ - geom_model_left = RPC.from_any( - os.path.join(data_path(), "rectification", "left_image.geom"), topleftconvention=True - ) - geom_model_right = RPC.from_any( - os.path.join(data_path(), "rectification", "right_image.geom"), topleftconvention=True - ) + geom_model_left = RPC(os.path.join(data_path(), "rectification", "left_image.geom")) + geom_model_right = RPC(os.path.join(data_path(), "rectification", "right_image.geom")) current_coords = np.array([[5000.5, 5000.5, 0.0]], dtype=np.float64) mean_spacing = 1 @@ -697,12 +665,8 @@ def test_rectification_grid_pos_inside_prepare_footprint_bounding_box(): # Compute shareloc epipolar footprint left_im = Image(os.path.join(data_path(), "rectification", "left_image.tif")) - geom_model_left = RPC.from_any( - os.path.join(data_path(), "rectification", "left_image.geom"), topleftconvention=True - ) - geom_model_right = RPC.from_any( - os.path.join(data_path(), "rectification", "right_image.geom"), topleftconvention=True - ) + geom_model_left = RPC(os.path.join(data_path(), "rectification", "left_image.geom")) + geom_model_right = RPC(os.path.join(data_path(), "rectification", "right_image.geom")) epi_step = 30 elevation_offset = 50 default_elev = 0.0 diff --git a/tests/geofunctions/test_triangulation.py b/tests/geofunctions/test_triangulation.py index 2b6cc2b..87a929d 100644 --- a/tests/geofunctions/test_triangulation.py +++ b/tests/geofunctions/test_triangulation.py @@ -156,10 +156,10 @@ def test_epi_triangulation_sift_rpc(): data_folder = data_path() id_scene = "PHR1B_P_201709281038045_SEN_PRG_FC_178608-001" file_geom = os.path.join(data_folder, f"rpc/{id_scene}.geom") - geom_model_left = RPC.from_any(file_geom, topleftconvention=True) + geom_model_left = RPC(file_geom) id_scene = "PHR1B_P_201709281038393_SEN_PRG_FC_178609-001" file_geom = os.path.join(data_folder, f"rpc/{id_scene}.geom") - geom_model_right = RPC.from_any(file_geom, topleftconvention=True) + geom_model_right = RPC(file_geom) grid_left_filename = os.path.join(data_path(), "rectification_grids", "left_epipolar_grid.tif") grid_right_filename = os.path.join(data_path(), "rectification_grids", "right_epipolar_grid.tif") @@ -217,10 +217,10 @@ def test_epi_triangulation_disp_rpc(): data_folder = data_path() id_scene = "PHR1B_P_201709281038045_SEN_PRG_FC_178608-001" file_geom = os.path.join(data_folder, f"rpc/{id_scene}.geom") - geom_model_left = RPC.from_any(file_geom, topleftconvention=True) + geom_model_left = RPC(file_geom) id_scene = "PHR1B_P_201709281038393_SEN_PRG_FC_178609-001" file_geom = os.path.join(data_folder, f"rpc/{id_scene}.geom") - geom_model_right = RPC.from_any(file_geom, topleftconvention=True) + geom_model_right = RPC(file_geom) # grid_left_filename = os.path.join(data_path(), "rectification_grids", # "grid_{}.tif".format(id_scene_left)) @@ -258,9 +258,9 @@ def test_epi_triangulation_disp_rpc_roi(): """ data_folder = data_path() file_geom = os.path.join(data_folder, "rpc/phr_ventoux/left_image.geom") - geom_model_left = RPC.from_any(file_geom, topleftconvention=True) + geom_model_left = RPC(file_geom) file_geom = os.path.join(data_folder, "rpc/phr_ventoux/right_image.geom") - geom_model_right = RPC.from_any(file_geom, topleftconvention=True) + geom_model_right = RPC(file_geom) grid_left_filename = os.path.join(data_path(), "rectification_grids", "left_epipolar_grid_ventoux.tif") grid_right_filename = os.path.join(data_path(), "rectification_grids", "right_epipolar_grid_ventoux.tif") diff --git a/tests/geomodels/test_los.py b/tests/geomodels/test_los.py index 733b1f5..71fc6c6 100644 --- a/tests/geomodels/test_los.py +++ b/tests/geomodels/test_los.py @@ -51,7 +51,7 @@ def test_los_creation_with_different_alt(alt): # Load geometrical model id_scene = "P1BP--2017092838284574CP" file_dimap = os.path.join(data_folder, f"rpc/RPC_{id_scene}.XML") - geometrical_model = RPC.from_any(file_dimap) + geometrical_model = RPC(file_dimap) # Get alt max if alt is None: alt_min, alt_max = geometrical_model.get_alt_min_max() @@ -103,7 +103,7 @@ def test_compare_two_altitude(): # Load geometrical model id_scene = "P1BP--2017092838284574CP" file_dimap = os.path.join(data_folder, f"rpc/RPC_{id_scene}.XML") - geometrical_model = RPC.from_any(file_dimap) + geometrical_model = RPC(file_dimap) # Create los model_los = LOS(matches_left, geometrical_model, [310, 850]) diff --git a/tests/geomodels/test_rpc.py b/tests/geomodels/test_rpc.py index 6b128db..cdb707a 100755 --- a/tests/geomodels/test_rpc.py +++ b/tests/geomodels/test_rpc.py @@ -50,25 +50,25 @@ def test_rpc_drivers(): # Test DIMAP RPC id_scene = "P1BP--2018122638935449CP" file_dimap = os.path.join(data_folder, f"rpc/PHRDIMAP_{id_scene}.XML") - fctrat_dimap = RPC.from_any(file_dimap) + fctrat_dimap = RPC(file_dimap) assert fctrat_dimap.driver_type == "dimap_v1.4" # Test OSSIM KWL GEOM RPC id_scene = "PHR1B_P_201709281038393_SEN_PRG_FC_178609-001" file_geom = os.path.join(data_folder, f"rpc/{id_scene}.geom") - fctrat_geom = RPC.from_any(file_geom) + fctrat_geom = RPC(file_geom) assert fctrat_geom.driver_type == "ossim_kwl" # Test DIMAPv2 RPC id_scene = "PHR1B_P_201709281038393_SEN_PRG_FC_178609-001" file_dimap_v2 = os.path.join(data_folder, f"rpc/RPC_{id_scene}.XML") - fctrat_dimap_v2 = RPC.from_any(file_dimap_v2, topleftconvention=True) + fctrat_dimap_v2 = RPC(file_dimap_v2) assert fctrat_dimap_v2.driver_type == "dimap_v2.15" # Test fake RPC fake_rpc = os.path.join(data_folder, "rpc/fake_rpc.txt") try: - RPC.from_any(fake_rpc, topleftconvention=True) # Raise ValueError->True + RPC(fake_rpc) # Raise ValueError->True raise AssertionError() # Assert false if no exception raised except ValueError: assert True @@ -123,7 +123,7 @@ def test_rpc_ossim_kwl(id_scene, lon, lat, alt, row_vt, col_vt): """ data_folder = data_path() file_geom = os.path.join(data_folder, f"rpc/{id_scene}.geom") - fctrat_geom = RPC.from_any(file_geom, topleftconvention=True) + fctrat_geom = RPC(file_geom) (row, col, __) = fctrat_geom.inverse_loc(lon, lat, alt) assert col == pytest.approx(col_vt, abs=1e-2) @@ -165,7 +165,7 @@ def test_rpc_from_any(id_scene, lon, lat, alt, row_vt, col_vt): """ data_folder = data_path() rpc_file = os.path.join(data_folder, "rpc", id_scene) - fctrat = RPC.from_any(rpc_file, topleftconvention=True) + fctrat = RPC(rpc_file) (row, col, __) = fctrat.inverse_loc(lon, lat, alt) assert fctrat.epsg == 4326 assert fctrat.datum == "ellipsoid" @@ -183,7 +183,7 @@ def test_rpc_direct_iterative(id_scene, lon, lat, alt): """ data_folder = data_path() rpc_file = os.path.join(data_folder, "rpc", id_scene) - fctrat = RPC.from_any(rpc_file, topleftconvention=True) + fctrat = RPC(rpc_file) (row, col, __) = fctrat.inverse_loc(lon, lat, alt) (lon2, lat2, __) = fctrat.direct_loc_inverse_iterative(row, col, alt) assert lon == pytest.approx(lon2, abs=1e-2) @@ -205,7 +205,7 @@ def test_rpc_from_geotiff_without_rpc(prod, can_read): data_folder = data_path() rpc_file = os.path.join(data_folder, "rpc", prod) try: - RPC.from_geotiff(rpc_file, topleftconvention=True) + RPC(rpc_file) assert can_read except ValueError: assert not can_read @@ -221,14 +221,14 @@ def test_rpc_dimap_v2(id_scene, lon, lat, alt, row_vt, col_vt): """ data_folder = data_path() file_dimap = os.path.join(data_folder, f"rpc/RPC_{id_scene}.XML") - fctrat_dimap = RPC.from_dimap(file_dimap, topleftconvention=True) + fctrat_dimap = RPC(file_dimap) (row, col, __) = fctrat_dimap.inverse_loc(lon, lat, alt) assert col == pytest.approx(col_vt, abs=1e-2) assert row == pytest.approx(row_vt, abs=1e-2) file_geom = os.path.join(data_folder, f"rpc/{id_scene}.geom") - fctrat_geom = RPC.from_any(file_geom, topleftconvention=True) + fctrat_geom = RPC(file_geom) (row, col, __) = fctrat_geom.inverse_loc(lon, lat, alt) assert col == pytest.approx(col_vt, abs=1e-2) assert row == pytest.approx(row_vt, abs=1e-2) @@ -243,7 +243,7 @@ def test_rpc_phrdimap(col, row, alt): id_scene = "P1BP--2018122638935449CP" file_dimap = os.path.join(data_folder, f"rpc/PHRDIMAP_{id_scene}.XML") - fctrat = RPC.from_dimap(file_dimap) + fctrat = RPC(file_dimap) (lonlatalt) = fctrat.direct_loc_h(row, col, alt) @@ -261,7 +261,7 @@ def test_rpc_direct_inverse_iterative_vs_direct(col, row, alt): id_scene = "P1BP--2018122638935449CP" file_dimap = os.path.join(data_folder, f"rpc/PHRDIMAP_{id_scene}.XML") - fctrat = RPC.from_dimap_v1(file_dimap) + fctrat = RPC(file_dimap) # (col,lig,alt)=(100,1000,400) lonlatalt = fctrat.direct_loc_h(row, col, alt) @@ -281,7 +281,7 @@ def test_rpc_direct_inverse_iterative_vs_direct_multiple_points(): id_scene = "P1BP--2018122638935449CP" file_dimap = os.path.join(data_folder, f"rpc/PHRDIMAP_{id_scene}.XML") - fctrat = RPC.from_dimap_v1(file_dimap) + fctrat = RPC(file_dimap) (col, row, alt) = (np.array([600, 610]), np.array([200, 210]), np.array([125])) p_direct = fctrat.direct_loc_h(row, col, alt) @@ -308,7 +308,7 @@ def test_rpc_direct_iterative_nan(): id_scene = "P1BP--2018122638935449CP" file_dimap = os.path.join(data_folder, f"rpc/PHRDIMAP_{id_scene}.XML") - fctrat = RPC.from_dimap_v1(file_dimap) + fctrat = RPC(file_dimap) (col, row, alt) = (np.array([600, np.nan]), np.array([200, 210]), np.array([125])) p_direct0 = fctrat.direct_loc_inverse_iterative(row, col, alt, fill_nan=True) @@ -327,7 +327,7 @@ def test_rpc_direct_iterative_all_nan(): id_scene = "P1BP--2018122638935449CP" file_dimap = os.path.join(data_folder, f"rpc/PHRDIMAP_{id_scene}.XML") - fctrat = RPC.from_dimap_v1(file_dimap) + fctrat = RPC(file_dimap) (col, row, alt) = (np.array([600, np.nan]), np.array([200, 210]), np.array([125])) direct_loc_tab = fctrat.direct_loc_inverse_iterative(row, col, alt, fill_nan=True) @@ -353,7 +353,7 @@ def test_rpc_direct_inverse_iterative(col, row, alt): id_scene = "P1BP--2018122638935449CP" file_dimap = os.path.join(data_folder, f"rpc/PHRDIMAP_{id_scene}.XML") - fctrat = RPC.from_dimap_v1(file_dimap) + fctrat = RPC(file_dimap) lonlatalt = fctrat.direct_loc_h(row, col, alt) (row_inv, col_inv, __) = fctrat.inverse_loc(lonlatalt[0][0], lonlatalt[0][1], lonlatalt[0][2]) @@ -374,7 +374,7 @@ def test_rpc_direct_dtm(id_scene, index_x, index_y): data_folder = data_path() rpc_file = os.path.join(data_folder, "rpc", id_scene) - fctrat = RPC.from_any(rpc_file, topleftconvention=True) + fctrat = RPC(rpc_file) id_scene = "P1BP--2017092838284574CP" data_folder_mnt = data_path("ellipsoide", id_scene) @@ -400,7 +400,7 @@ def test_rpc_los_extrapolation(id_scene, row, col): """ data_folder = data_path() file_dimap = os.path.join(data_folder, f"rpc/RPC_{id_scene}.XML") - fctrat = RPC.from_dimap(file_dimap, topleftconvention=True) + fctrat = RPC(file_dimap) los_edges = fctrat.los_extrema(row, col) altmin = -10 altmax = 2000 @@ -417,7 +417,7 @@ def test_rpc_minmax(): data_folder = data_path() id_scene = "P1BP--2018122638935449CP" fichier_dimap = os.path.join(data_folder, f"rpc/PHRDIMAP_{id_scene}.XML") - fctrat = RPC.from_any(fichier_dimap) + fctrat = RPC(fichier_dimap) (h_min, h_max) = fctrat.get_alt_min_max() assert h_min == 532.5 assert h_max == 617.5 From 5e27a329a697ed407375f374494992f6af42720b Mon Sep 17 00:00:00 2001 From: Emmanuel Dubois Date: Wed, 18 Oct 2023 17:57:41 +0200 Subject: [PATCH 03/12] feat: add geomodel template interface loc functions --- shareloc/geomodels/geomodel_template.py | 51 ++++++++++++++++++++++++- 1 file changed, 49 insertions(+), 2 deletions(-) diff --git a/shareloc/geomodels/geomodel_template.py b/shareloc/geomodels/geomodel_template.py index dca847c..5ccf790 100644 --- a/shareloc/geomodels/geomodel_template.py +++ b/shareloc/geomodels/geomodel_template.py @@ -35,8 +35,6 @@ # if(SHARELOC_OPTIM_GEOMODEL == True): # GeoModelTemplate.direct_loc_dtm = GeoModelTemplate.direct_loc_dtm_optim -# pylint: disable=too-few-public-methods - class GeoModelTemplate(metaclass=ABCMeta): """ @@ -58,3 +56,52 @@ def __init__(self, geomodel_path: str): # geomodel type. Set by the subclass self.type: str + + # Define GeoModelTemplate functions interface + + @abstractmethod + def direct_loc_h(self, row, col, alt, fill_nan=False): + """ + direct localization at constant altitude + + :param row: line sensor position + :type row: float or 1D numpy.ndarray dtype=float64 + :param col: column sensor position + :type col: float or 1D numpy.ndarray dtype=float64 + :param alt: altitude + :param fill_nan: fill numpy.nan values with lon and lat offset if true (same as OTB/OSSIM), nan is returned + otherwise + :type fill_nan: boolean + :return: ground position (lon,lat,h) + :rtype: numpy.ndarray 2D dimension with (N,3) shape, where N is number of input coordinates + """ + + @abstractmethod + def direct_loc_dtm(self, row, col, dtm): + """ + direct localization on dtm + + :param row: line sensor position + :type row: float + :param col: column sensor position + :type col: float + :param dtm: dtm intersection model + :type dtm: shareloc.geofunctions.dtm_intersection + :return: ground position (lon,lat,h) in dtm coordinates system + :rtype: numpy.ndarray 2D dimension with (N,3) shape, where N is number of input coordinates + """ + + @abstractmethod + def inverse_loc(self, lon, lat, alt): + """ + Inverse localization + + :param lon: longitude position + :type lon: float or 1D numpy.ndarray dtype=float64 + :param lat: latitude position + :type lat: float or 1D numpy.ndarray dtype=float64 + :param alt: altitude + :type alt: float + :return: sensor position (row, col, alt) + :rtype: tuple(1D np.array row position, 1D np.array col position, 1D np.array alt) + """ From 7cbd963082f5876b6668aa48b632558a923a1db9 Mon Sep 17 00:00:00 2001 From: Emmanuel Dubois Date: Thu, 16 Nov 2023 13:45:16 +0100 Subject: [PATCH 04/12] feat: replace RPC and Grid by GeoModel --- .pylintrc | 2 +- shareloc/geomodels/geomodel.py | 2 +- tests/geofunctions/test_localization.py | 36 +++++++++--------- tests/geofunctions/test_rectification.py | 47 ++++++++++++------------ tests/geofunctions/test_triangulation.py | 17 ++++----- tests/geomodels/test_grid.py | 4 +- tests/geomodels/test_los.py | 6 +-- tests/geomodels/test_rpc.py | 41 +++++++++++---------- 8 files changed, 77 insertions(+), 78 deletions(-) diff --git a/.pylintrc b/.pylintrc index 16d5e65..73c1b42 100644 --- a/.pylintrc +++ b/.pylintrc @@ -250,7 +250,7 @@ ignore-on-opaque-inference=yes # List of class names for which member attributes should not be checked (useful # for classes with dynamically set attributes). This supports the use of # qualified names. -ignored-classes=optparse.Values,thread._local,_thread._local +ignored-classes=optparse.Values,thread._local,_thread._local, GeoModel # List of module names for which member attributes should not be checked # (useful for modules/projects where namespaces are manipulated during runtime diff --git a/shareloc/geomodels/geomodel.py b/shareloc/geomodels/geomodel.py index e6ee48b..36dcb81 100644 --- a/shareloc/geomodels/geomodel.py +++ b/shareloc/geomodels/geomodel.py @@ -77,7 +77,7 @@ def create_geomodel(cls, geomodel_path: str, geomodel_type: str = "RPC"): # Create geomodel object with geomodel_path parameter from geomodel_type if exists try: geomodel_class = cls.available_geomodels[geomodel_type] - geomodel = geomodel_class(geomodel_path) + geomodel = geomodel_class(geomodel_path).load() logging.debug("GeoModel type name: %s", geomodel_type) except KeyError: logging.error("Geomodel type named %s is not supported", geomodel_type) diff --git a/tests/geofunctions/test_localization.py b/tests/geofunctions/test_localization.py index e6367fc..8fe2882 100644 --- a/tests/geofunctions/test_localization.py +++ b/tests/geofunctions/test_localization.py @@ -34,10 +34,10 @@ from shareloc.geofunctions.dtm_intersection import DTMIntersection from shareloc.geofunctions.localization import Localization from shareloc.geofunctions.localization import coloc as coloc_rpc +from shareloc.geomodels import GeoModel # Shareloc imports -from shareloc.geomodels.grid import Grid, coloc -from shareloc.geomodels.rpc import RPC +from shareloc.geomodels.grid import coloc from shareloc.image import Image from shareloc.proj_utils import coordinates_conversion @@ -57,9 +57,9 @@ def test_localize_direct_rpc(): # first instanciate the RPC geometric model # data = os.path.join(data_path(), "rpc/phr_ventoux/", "left_image.geom") - # geom_model = RPC(data) + # geom_model = GeoModel(data) data = os.path.join(data_path(), "rpc/phr_ventoux/", "RPC_PHR1B_P_201308051042194_SEN_690908101-001.XML") - geom_model = RPC(data) + geom_model = GeoModel(data) # then read the Image to retrieve its geotransform image_filename = os.path.join(data_path(), "image/phr_ventoux/", "left_image.tif") image_left = Image(image_filename) @@ -94,9 +94,9 @@ def test_localize_direct_grid(): # first instanciate the Grid geometric model # data = os.path.join(data_path(), "rpc/phr_ventoux/", "left_image.geom") - # geom_model_1 = RPC(data) + # geom_model_1 = GeoModel(data) data = os.path.join(data_path(), "grid/phr_ventoux/", "GRID_PHR1B_P_201308051042194_SEN_690908101-001.tif") - geom_model = Grid(data) + geom_model = GeoModel(data) # then read the Image to retrieve its geotransform image_filename = os.path.join(data_path(), "image/phr_ventoux/", "left_image.tif") image_left = Image(image_filename) @@ -136,7 +136,7 @@ def prepare_loc(alti="geoide", id_scene="P1BP--2017030824934340CP"): dtm = DTMIntersection(fic) # load grid model gld = os.path.join(data_folder, grid_name) - gri = Grid(gld) + gri = GeoModel(gld) return dtm, gri @@ -230,7 +230,7 @@ def test_extent(): Test extent """ data_left = os.path.join(data_path(), "rectification", "left_image") - geom_model = RPC(data_left + ".geom") + geom_model = GeoModel(data_left + ".geom") image_filename = os.path.join(data_path(), "image/phr_ventoux/", "left_image_pixsize_0_5.tif") image = Image(image_filename) loc_rpc_image = Localization(geom_model, elevation=None, image=image) @@ -246,7 +246,7 @@ def test_sensor_loc_dir_dtm_geoid(col, row, valid_coord): Test direct localization using image geotransform """ data = os.path.join(data_path(), "rectification", "left_image") - geom_model_left = RPC(data + ".geom") + geom_model_left = GeoModel(data + ".geom") image_filename = os.path.join(data_path(), "image/phr_ventoux/", "left_image.tif") image_left = Image(image_filename) @@ -267,7 +267,7 @@ def test_sensor_loc_dir_dtm_geoid_utm(col, row, valid_coord): Test direct localization using image geotransform """ data = os.path.join(data_path(), "rectification", "left_image") - geom_model_left = RPC(data + ".geom") + geom_model_left = GeoModel(data + ".geom") image_filename = os.path.join(data_path(), "image/phr_ventoux/", "left_image.tif") image_left = Image(image_filename) @@ -297,7 +297,7 @@ def test_sensor_loc_dir_vs_loc_rpc(row, col, h): data_folder = data_path() fichier_dimap = os.path.join(data_folder, f"rpc/PHRDIMAP_{id_scene}.XML") - fctrat = RPC(fichier_dimap) + fctrat = GeoModel(fichier_dimap) loc_rpc = Localization(fctrat) lonlatalt_rpc = loc_rpc.direct(row, col, h) @@ -392,7 +392,7 @@ def test_sensor_loc_inv_vs_loc_rpc(lon, lat, alt): data_folder = data_path() fichier_dimap = os.path.join(data_folder, f"rpc/PHRDIMAP_{id_scene}.XML") - fctrat = RPC(fichier_dimap) + fctrat = GeoModel(fichier_dimap) loc_rpc = Localization(fctrat) [row_rpc, col_rpc, __] = loc_rpc.inverse(lon, lat, alt) @@ -508,7 +508,7 @@ def test_loc_dir_loc_inv_rpc(id_scene, rpc, row, col, h): data_folder = data_path() fichier_dimap = os.path.join(data_folder, "rpc", rpc) - fctrat = RPC(fichier_dimap) + fctrat = GeoModel(fichier_dimap) (inv_row, inv_col, __) = fctrat.inverse_loc(lonlatalt[0][0], lonlatalt[0][1], lonlatalt[0][2]) assert row == pytest.approx(inv_row, abs=1e-2) assert col == pytest.approx(inv_col, abs=1e-2) @@ -560,7 +560,7 @@ def test_colocalization(col, row, alt): data_folder = data_path() id_scene = "P1BP--2018122638935449CP" file_dimap = os.path.join(data_folder, f"rpc/PHRDIMAP_{id_scene}.XML") - fctrat = RPC(file_dimap) + fctrat = GeoModel(file_dimap) row_coloc, col_coloc, _ = coloc_rpc(fctrat, fctrat, row, col, alt) @@ -575,12 +575,12 @@ def test_sensor_coloc_using_geotransform(col, row, h): Test direct localization using image geotransform """ data_left = os.path.join(data_path(), "rectification", "left_image") - geom_model_left = RPC(data_left + ".geom") + geom_model_left = GeoModel(data_left + ".geom") image_filename_left = os.path.join(data_path(), "image/phr_ventoux/", "left_image_pixsize_0_5.tif") image_left = Image(image_filename_left) data_right = os.path.join(data_path(), "rectification", "right_image") - geom_model_right = RPC(data_right + ".geom") + geom_model_right = GeoModel(data_right + ".geom") image_filename_right = os.path.join(data_path(), "image/phr_ventoux/", "right_image_pixsize_0_5.tif") image_right = Image(image_filename_right) @@ -610,7 +610,7 @@ def test_sensor_loc_utm(col, row): Test direct localization using image geotransform """ data_left = os.path.join(data_path(), "rectification", "left_image") - geom_model = RPC(data_left + ".geom") + geom_model = GeoModel(data_left + ".geom") epsg = 32631 loc_wgs = Localization(geom_model) loc_utm = Localization(geom_model, epsg=epsg) @@ -633,7 +633,7 @@ def test_sensor_loc_dir_dtm_multi_points(): left_im = Image(os.path.join(data_path(), "rectification", "left_image.tif")) - geom_model = RPC(os.path.join(data_path(), "rectification", "left_image.geom")) + geom_model = GeoModel(os.path.join(data_path(), "rectification", "left_image.geom")) dtm_file = os.path.join(data_path(), "dtm", "srtm_ventoux", "srtm90_non_void_filled", "N44E005.hgt") geoid_file = os.path.join(data_path(), "dtm", "geoid", "egm96_15.gtx") diff --git a/tests/geofunctions/test_rectification.py b/tests/geofunctions/test_rectification.py index a6e9556..31e7887 100644 --- a/tests/geofunctions/test_rectification.py +++ b/tests/geofunctions/test_rectification.py @@ -43,8 +43,7 @@ prepare_rectification, ) from shareloc.geofunctions.rectification_grid import RectificationGrid -from shareloc.geomodels.grid import Grid -from shareloc.geomodels.rpc import RPC +from shareloc.geomodels import GeoModel from shareloc.image import Image # Shareloc test imports @@ -62,8 +61,8 @@ def test_compute_stereorectification_epipolar_grids_geomodel_rpc(): left_im = Image(os.path.join(data_path(), "rectification", "left_image.tif")) right_im = Image(os.path.join(data_path(), "rectification", "right_image.tif")) - geom_model_left = RPC(os.path.join(data_path(), "rectification", "left_image.geom")) - geom_model_right = RPC(os.path.join(data_path(), "rectification", "right_image.geom")) + geom_model_left = GeoModel(os.path.join(data_path(), "rectification", "left_image.geom")) + geom_model_right = GeoModel(os.path.join(data_path(), "rectification", "right_image.geom")) epi_step = 30 elevation_offset = 50 @@ -109,8 +108,8 @@ def test_compute_stereorectification_epipolar_grids_geomodel_rpc_alti(): left_im = Image(os.path.join(data_path(), "rectification", "left_image.tif")) right_im = Image(os.path.join(data_path(), "rectification", "right_image.tif")) - geom_model_left = RPC(os.path.join(data_path(), "rectification", "left_image.geom")) - geom_model_right = RPC(os.path.join(data_path(), "rectification", "right_image.geom")) + geom_model_left = GeoModel(os.path.join(data_path(), "rectification", "left_image.geom")) + geom_model_right = GeoModel(os.path.join(data_path(), "rectification", "right_image.geom")) epi_step = 30 elevation_offset = 50 @@ -150,8 +149,8 @@ def test_compute_stereorectification_epipolar_grids_geomodel_rpc_dtm_geoid(): """ # first instantiate geometric models left and right (here RPC geometrics model) - geom_model_left = RPC(os.path.join(data_path(), "rectification", "left_image.geom")) - geom_model_right = RPC(os.path.join(data_path(), "rectification", "right_image.geom")) + geom_model_left = GeoModel(os.path.join(data_path(), "rectification", "left_image.geom")) + geom_model_right = GeoModel(os.path.join(data_path(), "rectification", "right_image.geom")) # read the images left_im = Image(os.path.join(data_path(), "rectification", "left_image.tif")) @@ -206,8 +205,8 @@ def test_compute_stereorectification_epipolar_grids_geomodel_rpc_dtm_geoid_roi() left_im = Image(os.path.join(data_path(), "rectification", "left_image.tif")) right_im = Image(os.path.join(data_path(), "rectification", "right_image.tif")) - geom_model_left = RPC(os.path.join(data_path(), "rectification", "left_image.geom")) - geom_model_right = RPC(os.path.join(data_path(), "rectification", "right_image.geom")) + geom_model_left = GeoModel(os.path.join(data_path(), "rectification", "left_image.geom")) + geom_model_right = GeoModel(os.path.join(data_path(), "rectification", "right_image.geom")) dtm_file = os.path.join(data_path(), "dtm", "srtm_ventoux", "srtm90_non_void_filled", "N44E005.hgt") geoid_file = os.path.join(data_path(), "dtm", "geoid", "egm96_15.gtx") @@ -256,10 +255,10 @@ def test_compute_stereorectification_epipolar_grids_geomodel_grid(): """ # first instantiate geometric models left and right (here Grid geometric model) - geom_model_left = Grid( + geom_model_left = GeoModel( os.path.join(data_path(), "grid/phr_ventoux/GRID_PHR1B_P_201308051042194_SEN_690908101-001.tif") ) - geom_model_right = Grid( + geom_model_right = GeoModel( os.path.join(data_path(), "grid/phr_ventoux/GRID_PHR1B_P_201308051042523_SEN_690908101-002.tif") ) @@ -312,10 +311,10 @@ def test_compute_stereorectification_epipolar_grids_geomodel_grid_dtm_geoid(): """ # first instantiate geometric models left and right (here Grid geometrics model) - geom_model_left = Grid( + geom_model_left = GeoModel( os.path.join(data_path(), "grid/phr_ventoux/GRID_PHR1B_P_201308051042194_SEN_690908101-001.tif") ) - geom_model_right = Grid( + geom_model_right = GeoModel( os.path.join(data_path(), "grid/phr_ventoux/GRID_PHR1B_P_201308051042523_SEN_690908101-002.tif") ) @@ -422,8 +421,8 @@ def test_prepare_rectification(): """ left_im = Image(os.path.join(data_path(), "rectification", "left_image.tif")) - geom_model_left = RPC(os.path.join(data_path(), "rectification", "left_image.geom")) - geom_model_right = RPC(os.path.join(data_path(), "rectification", "right_image.geom")) + geom_model_left = GeoModel(os.path.join(data_path(), "rectification", "left_image.geom")) + geom_model_right = GeoModel(os.path.join(data_path(), "rectification", "right_image.geom")) epi_step = 30 elevation_offset = 50 @@ -458,8 +457,8 @@ def test_prepare_rectification_footprint(): """ left_im = Image(os.path.join(data_path(), "rectification", "left_image.tif")) - geom_model_left = RPC(os.path.join(data_path(), "rectification", "left_image.geom")) - geom_model_right = RPC(os.path.join(data_path(), "rectification", "right_image.geom")) + geom_model_left = GeoModel(os.path.join(data_path(), "rectification", "left_image.geom")) + geom_model_right = GeoModel(os.path.join(data_path(), "rectification", "right_image.geom")) epi_step = 30 elevation_offset = 50 @@ -485,8 +484,8 @@ def test_rectification_moving_along_lines(): """ Test moving along line in epipolar geometry """ - geom_model_left = RPC(os.path.join(data_path(), "rectification", "left_image.geom")) - geom_model_right = RPC(os.path.join(data_path(), "rectification", "right_image.geom")) + geom_model_left = GeoModel(os.path.join(data_path(), "rectification", "left_image.geom")) + geom_model_right = GeoModel(os.path.join(data_path(), "rectification", "right_image.geom")) current_coords = np.array([[5000.5, 5000.5, 0.0]], dtype=np.float64) mean_spacing = 1 @@ -517,8 +516,8 @@ def test_rectification_moving_to_next_line(): """ Test moving to next line in epipolar geometry """ - geom_model_left = RPC(os.path.join(data_path(), "rectification", "left_image.geom")) - geom_model_right = RPC(os.path.join(data_path(), "rectification", "right_image.geom")) + geom_model_left = GeoModel(os.path.join(data_path(), "rectification", "left_image.geom")) + geom_model_right = GeoModel(os.path.join(data_path(), "rectification", "right_image.geom")) current_coords = np.array([[5000.5, 5000.5, 0.0]], dtype=np.float64) mean_spacing = 1 @@ -665,8 +664,8 @@ def test_rectification_grid_pos_inside_prepare_footprint_bounding_box(): # Compute shareloc epipolar footprint left_im = Image(os.path.join(data_path(), "rectification", "left_image.tif")) - geom_model_left = RPC(os.path.join(data_path(), "rectification", "left_image.geom")) - geom_model_right = RPC(os.path.join(data_path(), "rectification", "right_image.geom")) + geom_model_left = GeoModel(os.path.join(data_path(), "rectification", "left_image.geom")) + geom_model_right = GeoModel(os.path.join(data_path(), "rectification", "right_image.geom")) epi_step = 30 elevation_offset = 50 default_elev = 0.0 diff --git a/tests/geofunctions/test_triangulation.py b/tests/geofunctions/test_triangulation.py index 87a929d..07b9800 100644 --- a/tests/geofunctions/test_triangulation.py +++ b/tests/geofunctions/test_triangulation.py @@ -36,8 +36,7 @@ # Shareloc imports from shareloc.geofunctions.triangulation import distance_point_los, epipolar_triangulation, sensor_triangulation -from shareloc.geomodels.grid import Grid -from shareloc.geomodels.rpc import RPC +from shareloc.geomodels import GeoModel # Shareloc test imports from ..helpers import data_path @@ -86,7 +85,7 @@ def prepare_loc(alti="geoide", id_scene="P1BP--2017030824934340CP"): data_folder = data_path(alti, id_scene) # load grid gld = os.path.join(data_folder, f"GRID_{id_scene}.tif") - gri = Grid(gld) + gri = GeoModel(gld) return gri @@ -156,10 +155,10 @@ def test_epi_triangulation_sift_rpc(): data_folder = data_path() id_scene = "PHR1B_P_201709281038045_SEN_PRG_FC_178608-001" file_geom = os.path.join(data_folder, f"rpc/{id_scene}.geom") - geom_model_left = RPC(file_geom) + geom_model_left = GeoModel(file_geom) id_scene = "PHR1B_P_201709281038393_SEN_PRG_FC_178609-001" file_geom = os.path.join(data_folder, f"rpc/{id_scene}.geom") - geom_model_right = RPC(file_geom) + geom_model_right = GeoModel(file_geom) grid_left_filename = os.path.join(data_path(), "rectification_grids", "left_epipolar_grid.tif") grid_right_filename = os.path.join(data_path(), "rectification_grids", "right_epipolar_grid.tif") @@ -217,10 +216,10 @@ def test_epi_triangulation_disp_rpc(): data_folder = data_path() id_scene = "PHR1B_P_201709281038045_SEN_PRG_FC_178608-001" file_geom = os.path.join(data_folder, f"rpc/{id_scene}.geom") - geom_model_left = RPC(file_geom) + geom_model_left = GeoModel(file_geom) id_scene = "PHR1B_P_201709281038393_SEN_PRG_FC_178609-001" file_geom = os.path.join(data_folder, f"rpc/{id_scene}.geom") - geom_model_right = RPC(file_geom) + geom_model_right = GeoModel(file_geom) # grid_left_filename = os.path.join(data_path(), "rectification_grids", # "grid_{}.tif".format(id_scene_left)) @@ -258,9 +257,9 @@ def test_epi_triangulation_disp_rpc_roi(): """ data_folder = data_path() file_geom = os.path.join(data_folder, "rpc/phr_ventoux/left_image.geom") - geom_model_left = RPC(file_geom) + geom_model_left = GeoModel(file_geom) file_geom = os.path.join(data_folder, "rpc/phr_ventoux/right_image.geom") - geom_model_right = RPC(file_geom) + geom_model_right = GeoModel(file_geom) grid_left_filename = os.path.join(data_path(), "rectification_grids", "left_epipolar_grid_ventoux.tif") grid_right_filename = os.path.join(data_path(), "rectification_grids", "right_epipolar_grid_ventoux.tif") diff --git a/tests/geomodels/test_grid.py b/tests/geomodels/test_grid.py index 90ecac1..4ab9cb6 100755 --- a/tests/geomodels/test_grid.py +++ b/tests/geomodels/test_grid.py @@ -28,7 +28,7 @@ import scipy # Shareloc imports -from shareloc.geomodels.grid import Grid +from shareloc.geomodels import GeoModel from shareloc.image import Image # Shareloc test imports @@ -41,7 +41,7 @@ def fixture_get_geotiff_grid(): get grid and associated path """ geotiff_grid_path = data_path("ellipsoide", "loc_direct_grid_PHR_2013072139303958CP.tif") - gri_geotiff = Grid(geotiff_grid_path) + gri_geotiff = GeoModel(geotiff_grid_path) grid_image = Image(geotiff_grid_path, read_data=True) return gri_geotiff, grid_image diff --git a/tests/geomodels/test_los.py b/tests/geomodels/test_los.py index 71fc6c6..f105025 100644 --- a/tests/geomodels/test_los.py +++ b/tests/geomodels/test_los.py @@ -28,8 +28,8 @@ import pytest # Shareloc imports +from shareloc.geomodels import GeoModel from shareloc.geomodels.los import LOS -from shareloc.geomodels.rpc import RPC from shareloc.proj_utils import coordinates_conversion # Shareloc test imports @@ -51,7 +51,7 @@ def test_los_creation_with_different_alt(alt): # Load geometrical model id_scene = "P1BP--2017092838284574CP" file_dimap = os.path.join(data_folder, f"rpc/RPC_{id_scene}.XML") - geometrical_model = RPC(file_dimap) + geometrical_model = GeoModel(file_dimap) # Get alt max if alt is None: alt_min, alt_max = geometrical_model.get_alt_min_max() @@ -103,7 +103,7 @@ def test_compare_two_altitude(): # Load geometrical model id_scene = "P1BP--2017092838284574CP" file_dimap = os.path.join(data_folder, f"rpc/RPC_{id_scene}.XML") - geometrical_model = RPC(file_dimap) + geometrical_model = GeoModel(file_dimap) # Create los model_los = LOS(matches_left, geometrical_model, [310, 850]) diff --git a/tests/geomodels/test_rpc.py b/tests/geomodels/test_rpc.py index cdb707a..42ae96c 100755 --- a/tests/geomodels/test_rpc.py +++ b/tests/geomodels/test_rpc.py @@ -35,7 +35,8 @@ # Shareloc imports from shareloc.geofunctions.dtm_intersection import DTMIntersection -from shareloc.geomodels.rpc import RPC, identify_dimap, identify_ossim_kwl +from shareloc.geomodels import GeoModel +from shareloc.geomodels.rpc import identify_dimap, identify_ossim_kwl # Shareloc test imports from ..helpers import data_path @@ -50,25 +51,25 @@ def test_rpc_drivers(): # Test DIMAP RPC id_scene = "P1BP--2018122638935449CP" file_dimap = os.path.join(data_folder, f"rpc/PHRDIMAP_{id_scene}.XML") - fctrat_dimap = RPC(file_dimap) + fctrat_dimap = GeoModel(file_dimap) assert fctrat_dimap.driver_type == "dimap_v1.4" # Test OSSIM KWL GEOM RPC id_scene = "PHR1B_P_201709281038393_SEN_PRG_FC_178609-001" file_geom = os.path.join(data_folder, f"rpc/{id_scene}.geom") - fctrat_geom = RPC(file_geom) + fctrat_geom = GeoModel(file_geom) assert fctrat_geom.driver_type == "ossim_kwl" # Test DIMAPv2 RPC id_scene = "PHR1B_P_201709281038393_SEN_PRG_FC_178609-001" file_dimap_v2 = os.path.join(data_folder, f"rpc/RPC_{id_scene}.XML") - fctrat_dimap_v2 = RPC(file_dimap_v2) + fctrat_dimap_v2 = GeoModel(file_dimap_v2) assert fctrat_dimap_v2.driver_type == "dimap_v2.15" # Test fake RPC fake_rpc = os.path.join(data_folder, "rpc/fake_rpc.txt") try: - RPC(fake_rpc) # Raise ValueError->True + GeoModel(fake_rpc) # Raise ValueError->True raise AssertionError() # Assert false if no exception raised except ValueError: assert True @@ -123,7 +124,7 @@ def test_rpc_ossim_kwl(id_scene, lon, lat, alt, row_vt, col_vt): """ data_folder = data_path() file_geom = os.path.join(data_folder, f"rpc/{id_scene}.geom") - fctrat_geom = RPC(file_geom) + fctrat_geom = GeoModel(file_geom) (row, col, __) = fctrat_geom.inverse_loc(lon, lat, alt) assert col == pytest.approx(col_vt, abs=1e-2) @@ -165,7 +166,7 @@ def test_rpc_from_any(id_scene, lon, lat, alt, row_vt, col_vt): """ data_folder = data_path() rpc_file = os.path.join(data_folder, "rpc", id_scene) - fctrat = RPC(rpc_file) + fctrat = GeoModel(rpc_file) (row, col, __) = fctrat.inverse_loc(lon, lat, alt) assert fctrat.epsg == 4326 assert fctrat.datum == "ellipsoid" @@ -183,7 +184,7 @@ def test_rpc_direct_iterative(id_scene, lon, lat, alt): """ data_folder = data_path() rpc_file = os.path.join(data_folder, "rpc", id_scene) - fctrat = RPC(rpc_file) + fctrat = GeoModel(rpc_file) (row, col, __) = fctrat.inverse_loc(lon, lat, alt) (lon2, lat2, __) = fctrat.direct_loc_inverse_iterative(row, col, alt) assert lon == pytest.approx(lon2, abs=1e-2) @@ -205,7 +206,7 @@ def test_rpc_from_geotiff_without_rpc(prod, can_read): data_folder = data_path() rpc_file = os.path.join(data_folder, "rpc", prod) try: - RPC(rpc_file) + GeoModel(rpc_file) assert can_read except ValueError: assert not can_read @@ -221,14 +222,14 @@ def test_rpc_dimap_v2(id_scene, lon, lat, alt, row_vt, col_vt): """ data_folder = data_path() file_dimap = os.path.join(data_folder, f"rpc/RPC_{id_scene}.XML") - fctrat_dimap = RPC(file_dimap) + fctrat_dimap = GeoModel(file_dimap) (row, col, __) = fctrat_dimap.inverse_loc(lon, lat, alt) assert col == pytest.approx(col_vt, abs=1e-2) assert row == pytest.approx(row_vt, abs=1e-2) file_geom = os.path.join(data_folder, f"rpc/{id_scene}.geom") - fctrat_geom = RPC(file_geom) + fctrat_geom = GeoModel(file_geom) (row, col, __) = fctrat_geom.inverse_loc(lon, lat, alt) assert col == pytest.approx(col_vt, abs=1e-2) assert row == pytest.approx(row_vt, abs=1e-2) @@ -243,7 +244,7 @@ def test_rpc_phrdimap(col, row, alt): id_scene = "P1BP--2018122638935449CP" file_dimap = os.path.join(data_folder, f"rpc/PHRDIMAP_{id_scene}.XML") - fctrat = RPC(file_dimap) + fctrat = GeoModel(file_dimap) (lonlatalt) = fctrat.direct_loc_h(row, col, alt) @@ -261,7 +262,7 @@ def test_rpc_direct_inverse_iterative_vs_direct(col, row, alt): id_scene = "P1BP--2018122638935449CP" file_dimap = os.path.join(data_folder, f"rpc/PHRDIMAP_{id_scene}.XML") - fctrat = RPC(file_dimap) + fctrat = GeoModel(file_dimap) # (col,lig,alt)=(100,1000,400) lonlatalt = fctrat.direct_loc_h(row, col, alt) @@ -281,7 +282,7 @@ def test_rpc_direct_inverse_iterative_vs_direct_multiple_points(): id_scene = "P1BP--2018122638935449CP" file_dimap = os.path.join(data_folder, f"rpc/PHRDIMAP_{id_scene}.XML") - fctrat = RPC(file_dimap) + fctrat = GeoModel(file_dimap) (col, row, alt) = (np.array([600, 610]), np.array([200, 210]), np.array([125])) p_direct = fctrat.direct_loc_h(row, col, alt) @@ -308,7 +309,7 @@ def test_rpc_direct_iterative_nan(): id_scene = "P1BP--2018122638935449CP" file_dimap = os.path.join(data_folder, f"rpc/PHRDIMAP_{id_scene}.XML") - fctrat = RPC(file_dimap) + fctrat = GeoModel(file_dimap) (col, row, alt) = (np.array([600, np.nan]), np.array([200, 210]), np.array([125])) p_direct0 = fctrat.direct_loc_inverse_iterative(row, col, alt, fill_nan=True) @@ -327,7 +328,7 @@ def test_rpc_direct_iterative_all_nan(): id_scene = "P1BP--2018122638935449CP" file_dimap = os.path.join(data_folder, f"rpc/PHRDIMAP_{id_scene}.XML") - fctrat = RPC(file_dimap) + fctrat = GeoModel(file_dimap) (col, row, alt) = (np.array([600, np.nan]), np.array([200, 210]), np.array([125])) direct_loc_tab = fctrat.direct_loc_inverse_iterative(row, col, alt, fill_nan=True) @@ -353,7 +354,7 @@ def test_rpc_direct_inverse_iterative(col, row, alt): id_scene = "P1BP--2018122638935449CP" file_dimap = os.path.join(data_folder, f"rpc/PHRDIMAP_{id_scene}.XML") - fctrat = RPC(file_dimap) + fctrat = GeoModel(file_dimap) lonlatalt = fctrat.direct_loc_h(row, col, alt) (row_inv, col_inv, __) = fctrat.inverse_loc(lonlatalt[0][0], lonlatalt[0][1], lonlatalt[0][2]) @@ -374,7 +375,7 @@ def test_rpc_direct_dtm(id_scene, index_x, index_y): data_folder = data_path() rpc_file = os.path.join(data_folder, "rpc", id_scene) - fctrat = RPC(rpc_file) + fctrat = GeoModel(rpc_file) id_scene = "P1BP--2017092838284574CP" data_folder_mnt = data_path("ellipsoide", id_scene) @@ -400,7 +401,7 @@ def test_rpc_los_extrapolation(id_scene, row, col): """ data_folder = data_path() file_dimap = os.path.join(data_folder, f"rpc/RPC_{id_scene}.XML") - fctrat = RPC(file_dimap) + fctrat = GeoModel(file_dimap) los_edges = fctrat.los_extrema(row, col) altmin = -10 altmax = 2000 @@ -417,7 +418,7 @@ def test_rpc_minmax(): data_folder = data_path() id_scene = "P1BP--2018122638935449CP" fichier_dimap = os.path.join(data_folder, f"rpc/PHRDIMAP_{id_scene}.XML") - fctrat = RPC(fichier_dimap) + fctrat = GeoModel(fichier_dimap) (h_min, h_max) = fctrat.get_alt_min_max() assert h_min == 532.5 assert h_max == 617.5 From a7f08b70b2fba09bdb1c870f0017b7745b82285c Mon Sep 17 00:00:00 2001 From: GUINET Jonathan Date: Tue, 21 Nov 2023 09:42:27 +0000 Subject: [PATCH 05/12] update rpc and grid readers --- shareloc/geomodels/geomodel.py | 2 +- shareloc/geomodels/geomodel_template.py | 15 +- shareloc/geomodels/grid.py | 22 +- shareloc/geomodels/rpc.py | 445 ++--------------------- shareloc/geomodels/rpc_readers.py | 416 +++++++++++++++++++++ tests/geofunctions/test_localization.py | 4 +- tests/geofunctions/test_rectification.py | 8 +- tests/geofunctions/test_triangulation.py | 2 +- tests/geomodels/test_grid.py | 2 +- tests/geomodels/test_rpc.py | 2 +- 10 files changed, 485 insertions(+), 433 deletions(-) create mode 100755 shareloc/geomodels/rpc_readers.py diff --git a/shareloc/geomodels/geomodel.py b/shareloc/geomodels/geomodel.py index 36dcb81..9bb9f8a 100644 --- a/shareloc/geomodels/geomodel.py +++ b/shareloc/geomodels/geomodel.py @@ -77,7 +77,7 @@ def create_geomodel(cls, geomodel_path: str, geomodel_type: str = "RPC"): # Create geomodel object with geomodel_path parameter from geomodel_type if exists try: geomodel_class = cls.available_geomodels[geomodel_type] - geomodel = geomodel_class(geomodel_path).load() + geomodel = geomodel_class.load(geomodel_path) logging.debug("GeoModel type name: %s", geomodel_type) except KeyError: logging.error("Geomodel type named %s is not supported", geomodel_type) diff --git a/shareloc/geomodels/geomodel_template.py b/shareloc/geomodels/geomodel_template.py index 5ccf790..388c242 100644 --- a/shareloc/geomodels/geomodel_template.py +++ b/shareloc/geomodels/geomodel_template.py @@ -43,17 +43,11 @@ class GeoModelTemplate(metaclass=ABCMeta): """ @abstractmethod - def __init__(self, geomodel_path: str): + def __init__(self): """ Return the geomodel object associated with the geomodel_type given in the configuration - - :param geomodel_path: path of the geomodel file to instanciate. - :type geomodel_path: string """ - # geomodel filename path - self.geomodel_path: str = geomodel_path - # geomodel type. Set by the subclass self.type: str @@ -105,3 +99,10 @@ def inverse_loc(self, lon, lat, alt): :return: sensor position (row, col, alt) :rtype: tuple(1D np.array row position, 1D np.array col position, 1D np.array alt) """ + + @classmethod + @abstractmethod + def load(cls, geomodel_path, **kwargs): + """ + load function with class specific args + """ diff --git a/shareloc/geomodels/grid.py b/shareloc/geomodels/grid.py index 819df8a..ae1b53c 100755 --- a/shareloc/geomodels/grid.py +++ b/shareloc/geomodels/grid.py @@ -87,9 +87,9 @@ def __init__(self, geomodel_path: str): :type geomodel_path: string """ # Instanciate GeoModelTemplate generic init with shared parameters - super().__init__(geomodel_path) + super().__init__() self.type = "multi H grid" - + self.geomodel_path = geomodel_path # GeoModel Grid parameters definition (see documentation) self.row0 = None self.col0 = None @@ -105,11 +105,10 @@ def __init__(self, geomodel_path: str): self.rowmax = None self.colmax = None self.epsg = 0 + self.read() - # Load function of grid parameters from grid file to grid object - self.load() - - def load(self): + @classmethod + def load(cls, geomodel_path): """ Load grid and fill Class attributes. @@ -119,7 +118,18 @@ def load(self): - lon_data : [alt,row,col] - lat_data : [alt,row,col] """ + return cls(geomodel_path) + + def read(self): + """ + Load grid and fill Class attributes. + The grid is read as an shareloc.image. Image and class attributes are filled. + Shareloc geotiff grids are stored by increasing altitude H0 ... Hx + 2 data cubes are defined: + - lon_data : [alt,row,col] + - lat_data : [alt,row,col] + """ grid_image = Image(self.geomodel_path, read_data=True) if grid_image.dataset.driver != "GTiff": raise TypeError( diff --git a/shareloc/geomodels/rpc.py b/shareloc/geomodels/rpc.py index e881e82..7ed9190 100755 --- a/shareloc/geomodels/rpc.py +++ b/shareloc/geomodels/rpc.py @@ -26,101 +26,20 @@ # Standard imports import logging -from os.path import basename -from typing import Dict -from xml.dom import minidom # Third party imports import numpy as np -import rasterio as rio from numba import config, njit, prange # Shareloc imports from shareloc.geomodels.geomodel import GeoModel from shareloc.geomodels.geomodel_template import GeoModelTemplate +from shareloc.geomodels.rpc_readers import rpc_reader from shareloc.proj_utils import coordinates_conversion # Set numba type of threading layer before parallel target compilation config.THREADING_LAYER = "omp" - -def parse_coeff_line(coeff_str): - """ - split str coef to float list - - :param coeff_str: line coef - :type coeff_str: str - :return: coeff list - :rtype: list() - """ - return [float(el) for el in coeff_str.split()] - - -def identify_dimap(xml_file): - """ - parse xml file to identify dimap and its version - - :param xml_file: dimap rpc file - :type xml_file: str - :return: dimap info : dimap_version and None if not an dimap file - :rtype: str - """ - try: - xmldoc = minidom.parse(xml_file) - mtd = xmldoc.getElementsByTagName("Metadata_Identification") - mtd_format = mtd[0].getElementsByTagName("METADATA_FORMAT")[0].firstChild.data - if mtd_format == "DIMAP_PHR": - version_tag = "METADATA_PROFILE" - else: - version_tag = "METADATA_FORMAT" - version = mtd[0].getElementsByTagName(version_tag)[0].attributes.items()[0][1] - except Exception: # pylint: disable=broad-except - return None - - return version - - -def identify_ossim_kwl(ossim_kwl_file): - """ - parse geom file to identify if it is an ossim model - - :param ossim_kwl_fil : ossim keyword list file - :type ossim_kwl_file: str - :return: ossimmodel or None if not an ossim kwl file - :rtype: str - """ - try: - with open(ossim_kwl_file, encoding="utf-8") as ossim_file: - content = ossim_file.readlines() - geom_dict = {} - for line in content: - (key, val) = line.split(": ") - geom_dict[key] = val.rstrip() - if "type" in geom_dict: - if geom_dict["type"].strip().startswith("ossim"): - return geom_dict["type"].strip() - return None - except Exception: # pylint: disable=broad-except - return None - - -def identify_geotiff_rpc(image_filename): - """ - read image file to identify if it is a geotiff which contains RPCs - - :param image_filename: image_filename - :type image_filename: str - :return: rpc info, rpc dict or None if not a geotiff with rpc - :rtype: str - """ - try: - dataset = rio.open(image_filename) - rpc_dict = dataset.tags(ns="RPC") - return rpc_dict - except Exception: # pylint: disable=broad-except - return None - - # pylint: disable=no-member @@ -132,25 +51,16 @@ class RPC(GeoModelTemplate): # gitlab issue #61 # pylint: disable=too-many-instance-attributes - def __init__(self, geomodel_path: str): - # Instanciate GeoModelTemplate generic init with shared parameters - super().__init__(geomodel_path) - self.type = "RPC" - - # initiate epsg and datum, can overriden by rpc_params + # gitlab issue #61 + # pylint: disable=too-many-instance-attributes + def __init__(self, rpc_params): + super().__init__() self.epsg = None self.datum = None - - # RPC parameters as Dict - self.rpc_params: Dict = {} - - # RPC parameters are load from geomodel_path to rpc params - self.load() - - # set class parameters from rpc_params (epsg and datum can be overriden) - for key, value in self.rpc_params.items(): + for key, value in rpc_params.items(): setattr(self, key, value) + self.type = "RPC" if self.epsg is None: self.epsg = 4326 if self.datum is None: @@ -158,32 +68,31 @@ def __init__(self, geomodel_path: str): self.lim_extrapol = 1.0001 - # Monomes seems not used in shareloc code: Clean ? # Each monome: c[0]*X**c[1]*Y**c[2]*Z**c[3] - self.monomes = np.array( - [ - [1, 0, 0, 0], - [1, 1, 0, 0], - [1, 0, 1, 0], - [1, 0, 0, 1], - [1, 1, 1, 0], - [1, 1, 0, 1], - [1, 0, 1, 1], - [1, 2, 0, 0], - [1, 0, 2, 0], - [1, 0, 0, 2], - [1, 1, 1, 1], - [1, 3, 0, 0], - [1, 1, 2, 0], - [1, 1, 0, 2], - [1, 2, 1, 0], - [1, 0, 3, 0], - [1, 0, 1, 2], - [1, 2, 0, 1], - [1, 0, 2, 1], - [1, 0, 0, 3], - ] - ) + monomes_order = [ + [1, 0, 0, 0], + [1, 1, 0, 0], + [1, 0, 1, 0], + [1, 0, 0, 1], + [1, 1, 1, 0], + [1, 1, 0, 1], + [1, 0, 1, 1], + [1, 2, 0, 0], + [1, 0, 2, 0], + [1, 0, 0, 2], + [1, 1, 1, 1], + [1, 3, 0, 0], + [1, 1, 2, 0], + [1, 1, 0, 2], + [1, 2, 1, 0], + [1, 0, 3, 0], + [1, 0, 1, 2], + [1, 2, 0, 1], + [1, 0, 2, 1], + [1, 0, 0, 3], + ] + + self.monomes = np.array(monomes_order) # monomial coefficients of 1st variable derivative self.monomes_deriv_1 = np.array( @@ -262,7 +171,8 @@ def __init__(self, geomodel_path: str): self.row0 = self.offset_row - self.scale_row self.rowmax = self.offset_row + self.scale_row - def load(self): + @classmethod + def load(cls, geomodel_path, topleftconvention=True): """ Load from any RPC (auto identify driver) from filename (dimap, ossim kwl, geotiff) @@ -274,293 +184,8 @@ def load(self): If True : [0,0] is at the top left of the Top Left pixel (OSSIM) """ # Set topleftconvention (keeping historic option): to clean - topleftconvention = True - - # If ends with XML --> DIMAP - if basename(self.geomodel_path.upper()).endswith("XML"): - dimap_version = identify_dimap(self.geomodel_path) - if dimap_version is not None: - if float(dimap_version) < 2.0: - self.load_dimap_v1(topleftconvention) - if float(dimap_version) >= 2.0: - self.load_dimap_coeff(topleftconvention) - else: - # If not DIMAP, identify ossim - ossim_model = identify_ossim_kwl(self.geomodel_path) - if ossim_model is not None: - self.load_ossim_kwl(topleftconvention) - else: - # Otherwise, check if RPC is in geotif - geotiff_rpc_dict = identify_geotiff_rpc(self.geomodel_path) - if geotiff_rpc_dict is not None: - self.load_geotiff(topleftconvention) - else: - # Raise error if no file recognized. - raise ValueError("can not read rpc file") - - def load_dimap_coeff(self, topleftconvention=True): - """ - Load from Dimap v2 and V3 - - :param topleftconvention: [0,0] position - :type topleftconvention: boolean - If False : [0,0] is at the center of the Top Left pixel - If True : [0,0] is at the top left of the Top Left pixel (OSSIM) - """ - - if not basename(self.geomodel_path).upper().endswith("XML"): - raise ValueError("dimap must ends with .xml") - - xmldoc = minidom.parse(self.geomodel_path) - - mtd = xmldoc.getElementsByTagName("Metadata_Identification") - version = mtd[0].getElementsByTagName("METADATA_FORMAT")[0].attributes.items()[0][1] - self.rpc_params["driver_type"] = "dimap_v" + version - global_rfm = xmldoc.getElementsByTagName("Global_RFM")[0] - normalisation_coeffs = global_rfm.getElementsByTagName("RFM_Validity")[0] - self.rpc_params["offset_row"] = float(normalisation_coeffs.getElementsByTagName("LINE_OFF")[0].firstChild.data) - self.rpc_params["offset_col"] = float(normalisation_coeffs.getElementsByTagName("SAMP_OFF")[0].firstChild.data) - if float(version) >= 3: - direct_coeffs = global_rfm.getElementsByTagName("ImagetoGround_Values")[0] - inverse_coeffs = global_rfm.getElementsByTagName("GroundtoImage_Values")[0] - - self.rpc_params["num_x"] = [ - float(direct_coeffs.getElementsByTagName(f"LON_NUM_COEFF_{index}")[0].firstChild.data) - for index in range(1, 21) - ] - self.rpc_params["den_x"] = [ - float(direct_coeffs.getElementsByTagName(f"LON_DEN_COEFF_{index}")[0].firstChild.data) - for index in range(1, 21) - ] - self.rpc_params["num_y"] = [ - float(direct_coeffs.getElementsByTagName(f"LAT_NUM_COEFF_{index}")[0].firstChild.data) - for index in range(1, 21) - ] - self.rpc_params["den_y"] = [ - float(direct_coeffs.getElementsByTagName(f"LAT_DEN_COEFF_{index}")[0].firstChild.data) - for index in range(1, 21) - ] - self.rpc_params["offset_col"] -= 0.5 - self.rpc_params["offset_row"] -= 0.5 - - else: - direct_coeffs = global_rfm.getElementsByTagName("Direct_Model")[0] - inverse_coeffs = global_rfm.getElementsByTagName("Inverse_Model")[0] - - self.rpc_params["num_x"] = [ - float(direct_coeffs.getElementsByTagName(f"SAMP_NUM_COEFF_{index}")[0].firstChild.data) - for index in range(1, 21) - ] - self.rpc_params["den_x"] = [ - float(direct_coeffs.getElementsByTagName(f"SAMP_DEN_COEFF_{index}")[0].firstChild.data) - for index in range(1, 21) - ] - self.rpc_params["num_y"] = [ - float(direct_coeffs.getElementsByTagName(f"LINE_NUM_COEFF_{index}")[0].firstChild.data) - for index in range(1, 21) - ] - self.rpc_params["den_y"] = [ - float(direct_coeffs.getElementsByTagName(f"LINE_DEN_COEFF_{index}")[0].firstChild.data) - for index in range(1, 21) - ] - self.rpc_params["offset_col"] -= 1.0 - self.rpc_params["offset_row"] -= 1.0 - - self.rpc_params["num_col"] = [ - float(inverse_coeffs.getElementsByTagName(f"SAMP_NUM_COEFF_{index}")[0].firstChild.data) - for index in range(1, 21) - ] - self.rpc_params["den_col"] = [ - float(inverse_coeffs.getElementsByTagName(f"SAMP_DEN_COEFF_{index}")[0].firstChild.data) - for index in range(1, 21) - ] - self.rpc_params["num_row"] = [ - float(inverse_coeffs.getElementsByTagName(f"LINE_NUM_COEFF_{index}")[0].firstChild.data) - for index in range(1, 21) - ] - self.rpc_params["den_row"] = [ - float(inverse_coeffs.getElementsByTagName(f"LINE_DEN_COEFF_{index}")[0].firstChild.data) - for index in range(1, 21) - ] - - self.rpc_params["scale_col"] = float(normalisation_coeffs.getElementsByTagName("SAMP_SCALE")[0].firstChild.data) - self.rpc_params["scale_row"] = float(normalisation_coeffs.getElementsByTagName("LINE_SCALE")[0].firstChild.data) - self.rpc_params["offset_alt"] = float( - normalisation_coeffs.getElementsByTagName("HEIGHT_OFF")[0].firstChild.data - ) - self.rpc_params["scale_alt"] = float( - normalisation_coeffs.getElementsByTagName("HEIGHT_SCALE")[0].firstChild.data - ) - self.rpc_params["offset_x"] = float(normalisation_coeffs.getElementsByTagName("LONG_OFF")[0].firstChild.data) - self.rpc_params["scale_x"] = float(normalisation_coeffs.getElementsByTagName("LONG_SCALE")[0].firstChild.data) - self.rpc_params["offset_y"] = float(normalisation_coeffs.getElementsByTagName("LAT_OFF")[0].firstChild.data) - self.rpc_params["scale_y"] = float(normalisation_coeffs.getElementsByTagName("LAT_SCALE")[0].firstChild.data) - # If top left convention, 0.5 pixel shift added on col/row offsets - if topleftconvention: - self.rpc_params["offset_col"] += 0.5 - self.rpc_params["offset_row"] += 0.5 - - def load_dimap_v1(self, topleftconvention=True): - """ - Load from dimap v1 - - ** Deprecated, to clean ? ** - - :param topleftconvention: [0,0] position - :type topleftconvention: boolean - If False : [0,0] is at the center of the Top Left pixel - If True : [0,0] is at the top left of the Top Left pixel (OSSIM) - """ - - if not basename(self.geomodel_path).upper().endswith("XML"): - raise ValueError("dimap must ends with .xml") - - xmldoc = minidom.parse(self.geomodel_path) - - mtd = xmldoc.getElementsByTagName("Metadata_Identification") - version = mtd[0].getElementsByTagName("METADATA_PROFILE")[0].attributes.items()[0][1] - self.rpc_params = {"driver_type": "dimap_v" + version} - - global_rfm = xmldoc.getElementsByTagName("Global_RFM") - rfm_validity = xmldoc.getElementsByTagName("RFM_Validity") - coeff_lon = [float(el) for el in global_rfm[0].getElementsByTagName("F_LON")[0].firstChild.data.split()] - coeff_lat = [float(el) for el in global_rfm[0].getElementsByTagName("F_LAT")[0].firstChild.data.split()] - coeff_col = [float(el) for el in global_rfm[0].getElementsByTagName("F_COL")[0].firstChild.data.split()] - coeff_lig = [float(el) for el in global_rfm[0].getElementsByTagName("F_ROW")[0].firstChild.data.split()] - - scale_lon = float(rfm_validity[0].getElementsByTagName("Lon")[0].getElementsByTagName("A")[0].firstChild.data) - offset_lon = float(rfm_validity[0].getElementsByTagName("Lon")[0].getElementsByTagName("B")[0].firstChild.data) - scale_lat = float(rfm_validity[0].getElementsByTagName("Lat")[0].getElementsByTagName("A")[0].firstChild.data) - offset_lat = float(rfm_validity[0].getElementsByTagName("Lat")[0].getElementsByTagName("B")[0].firstChild.data) - scale_alt = float(rfm_validity[0].getElementsByTagName("Alt")[0].getElementsByTagName("A")[0].firstChild.data) - offset_alt = float(rfm_validity[0].getElementsByTagName("Alt")[0].getElementsByTagName("B")[0].firstChild.data) - scale_col = float(rfm_validity[0].getElementsByTagName("Col")[0].getElementsByTagName("A")[0].firstChild.data) - offset_col = float(rfm_validity[0].getElementsByTagName("Col")[0].getElementsByTagName("B")[0].firstChild.data) - scale_row = float(rfm_validity[0].getElementsByTagName("Row")[0].getElementsByTagName("A")[0].firstChild.data) - offset_row = float(rfm_validity[0].getElementsByTagName("Row")[0].getElementsByTagName("B")[0].firstChild.data) - - self.rpc_params["offset_col"] = offset_col - self.rpc_params["scale_col"] = scale_col - self.rpc_params["offset_row"] = offset_row - self.rpc_params["scale_row"] = scale_row - self.rpc_params["offset_alt"] = offset_alt - self.rpc_params["scale_alt"] = scale_alt - self.rpc_params["offset_x"] = offset_lon - self.rpc_params["scale_x"] = scale_lon - self.rpc_params["offset_y"] = offset_lat - self.rpc_params["scale_y"] = scale_lat - self.rpc_params["num_x"] = coeff_lon[0:20] - self.rpc_params["den_x"] = coeff_lon[20::] - self.rpc_params["num_y"] = coeff_lat[0:20] - self.rpc_params["den_y"] = coeff_lat[20::] - self.rpc_params["num_col"] = coeff_col[0:20] - self.rpc_params["den_col"] = coeff_col[20::] - self.rpc_params["num_row"] = coeff_lig[0:20] - self.rpc_params["den_row"] = coeff_lig[20::] - self.rpc_params["offset_col"] -= 1.0 - self.rpc_params["offset_row"] -= 1.0 - # If top left convention, 0.5 pixel shift added on col/row offsets - if topleftconvention: - self.rpc_params["offset_col"] += 0.5 - self.rpc_params["offset_row"] += 0.5 - - def load_geotiff(self, topleftconvention=True): - """ - Load from a geotiff image file - - - :param topleftconvention: [0,0] position - :type topleftconvention: boolean - If False : [0,0] is at the center of the Top Left pixel - If True : [0,0] is at the top left of the Top Left pixel (OSSIM) - """ - dataset = rio.open(self.geomodel_path) - rpc_dict = dataset.tags(ns="RPC") - if not rpc_dict: - logging.error("%s does not contains RPCs ", self.geomodel_path) - raise ValueError - self.rpc_params = { - "den_row": parse_coeff_line(rpc_dict["LINE_DEN_COEFF"]), - "num_row": parse_coeff_line(rpc_dict["LINE_NUM_COEFF"]), - "num_col": parse_coeff_line(rpc_dict["SAMP_NUM_COEFF"]), - "den_col": parse_coeff_line(rpc_dict["SAMP_DEN_COEFF"]), - "offset_col": float(rpc_dict["SAMP_OFF"]), - "scale_col": float(rpc_dict["SAMP_SCALE"]), - "offset_row": float(rpc_dict["LINE_OFF"]), - "scale_row": float(rpc_dict["LINE_SCALE"]), - "offset_alt": float(rpc_dict["HEIGHT_OFF"]), - "scale_alt": float(rpc_dict["HEIGHT_SCALE"]), - "offset_x": float(rpc_dict["LONG_OFF"]), - "scale_x": float(rpc_dict["LONG_SCALE"]), - "offset_y": float(rpc_dict["LAT_OFF"]), - "scale_y": float(rpc_dict["LAT_SCALE"]), - "num_x": None, - "den_x": None, - "num_y": None, - "den_y": None, - } - # inverse coeff are not defined - # If top left convention, 0.5 pixel shift added on col/row offsets - if topleftconvention: - self.rpc_params["offset_col"] += 0.5 - self.rpc_params["offset_row"] += 0.5 - - def load_ossim_kwl(self, topleftconvention=True): - """ - Load from a geom file - - :param topleftconvention: [0,0] position - :type topleftconvention: boolean - If False : [0,0] is at the center of the Top Left pixel - If True : [0,0] is at the top left of the Top Left pixel (OSSIM) - """ - # OSSIM keyword list - self.rpc_params["driver_type"] = "ossim_kwl" - with open(self.geomodel_path, "r", encoding="utf-8") as ossim_file: - content = ossim_file.readlines() - - geom_dict = {} - for line in content: - (key, val) = line.split(": ") - geom_dict[key] = val.rstrip() - - self.rpc_params["den_row"] = [np.nan] * 20 - self.rpc_params["num_row"] = [np.nan] * 20 - self.rpc_params["den_col"] = [np.nan] * 20 - self.rpc_params["num_col"] = [np.nan] * 20 - for index in range(0, 20): - axis = "line" - num_den = "den" - key = f"{axis}_{num_den}_coeff_{index:02d}" - self.rpc_params["den_row"][index] = float(geom_dict[key]) - num_den = "num" - key = f"{axis}_{num_den}_coeff_{index:02d}" - self.rpc_params["num_row"][index] = float(geom_dict[key]) - axis = "samp" - key = f"{axis}_{num_den}_coeff_{index:02d}" - self.rpc_params["num_col"][index] = float(geom_dict[key]) - num_den = "den" - key = f"{axis}_{num_den}_coeff_{index:02d}" - self.rpc_params["den_col"][index] = float(geom_dict[key]) - self.rpc_params["offset_col"] = float(geom_dict["samp_off"]) - self.rpc_params["scale_col"] = float(geom_dict["samp_scale"]) - self.rpc_params["offset_row"] = float(geom_dict["line_off"]) - self.rpc_params["scale_row"] = float(geom_dict["line_scale"]) - self.rpc_params["offset_alt"] = float(geom_dict["height_off"]) - self.rpc_params["scale_alt"] = float(geom_dict["height_scale"]) - self.rpc_params["offset_x"] = float(geom_dict["long_off"]) - self.rpc_params["scale_x"] = float(geom_dict["long_scale"]) - self.rpc_params["offset_y"] = float(geom_dict["lat_off"]) - self.rpc_params["scale_y"] = float(geom_dict["lat_scale"]) - # inverse coeff are not defined - self.rpc_params["num_x"] = None - self.rpc_params["den_x"] = None - self.rpc_params["num_y"] = None - self.rpc_params["den_y"] = None - # If top left convention, 0.5 pixel shift added on col/row offsets - if topleftconvention: - self.rpc_params["offset_col"] += 0.5 - self.rpc_params["offset_row"] += 0.5 + cls.geomodel_path = geomodel_path + return cls(rpc_reader(geomodel_path, topleftconvention)) def direct_loc_h(self, row, col, alt, fill_nan=False): """ diff --git a/shareloc/geomodels/rpc_readers.py b/shareloc/geomodels/rpc_readers.py new file mode 100755 index 0000000..d7bac09 --- /dev/null +++ b/shareloc/geomodels/rpc_readers.py @@ -0,0 +1,416 @@ +#!/usr/bin/env python +# coding: utf8 +# +# Copyright (c) 2023 Centre National d'Etudes Spatiales (CNES). +# +# This file is part of Shareloc +# (see https://github.com/CNES/shareloc). +# +# Licensed under the Apache License, Version 2.0 (the "License"); +# you may not use this file except in compliance with the License. +# You may obtain a copy of the License at +# +# http://www.apache.org/licenses/LICENSE-2.0 +# +# Unless required by applicable law or agreed to in writing, software +# distributed under the License is distributed on an "AS IS" BASIS, +# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +# See the License for the specific language governing permissions and +# limitations under the License. +# +""" +This module contains the RPC formats handling to instantiate RPC models. +RPC models covered are : DIMAP V1, DIMAP V2, DIMAP V3, ossim (geom file), geotiff. +""" + +# Standard imports +import logging +from os.path import basename +from typing import Dict +from xml.dom import minidom + +# Third party imports +import numpy as np +import rasterio as rio + + +def rpc_reader(geomodel_path: str, topleftconvention: bool = True) -> Dict: + """ + Load from any RPC (auto identify driver) + from filename (dimap, ossim kwl, geotiff) + + TODO: topleftconvention always to True, set a standard and remove the option + + topleftconvention boolean: [0,0] position + If False : [0,0] is at the center of the Top Left pixel + If True : [0,0] is at the top left of the Top Left pixel (OSSIM) + + :param geomodel_path geomodel filename + :return rpc dict filled with parameters + """ + + # If ends with XML --> DIMAP + if basename(geomodel_path.upper()).endswith("XML"): + dimap_version = identify_dimap(geomodel_path) + if dimap_version is not None: + if float(dimap_version) < 2.0: + return rpc_reader_dimap_v1(geomodel_path, topleftconvention) + if float(dimap_version) >= 2.0: + return rpc_reader_dimap_v23(geomodel_path, topleftconvention) + ossim_model = identify_ossim_kwl(geomodel_path) + if ossim_model is not None: + return rpc_reader_ossim_kwl(geomodel_path, topleftconvention) + geotiff_rpc_dict = identify_geotiff_rpc(geomodel_path) + if geotiff_rpc_dict is not None: + return rpc_reader_geotiff(geomodel_path, topleftconvention) + raise ValueError("can not read rpc file") + + +def rpc_reader_dimap_v23(geomodel_path: str, topleftconvention: bool = True) -> Dict: + """ + Load from Dimap v2 and V3 + + :param geomodel_path: path to geomodel + :param topleftconvention: [0,0] position + :type topleftconvention: boolean + If False : [0,0] is at the center of the Top Left pixel + If True : [0,0] is at the top left of the Top Left pixel (OSSIM) + + :returns rpc_dict + """ + + if not basename(geomodel_path).upper().endswith("XML"): + raise ValueError("dimap must ends with .xml") + + xmldoc = minidom.parse(geomodel_path) + + mtd = xmldoc.getElementsByTagName("Metadata_Identification") + version = mtd[0].getElementsByTagName("METADATA_FORMAT")[0].attributes.items()[0][1] + rpc_params = {"driver_type": "dimap_v" + version} + global_rfm = xmldoc.getElementsByTagName("Global_RFM")[0] + normalisation_coeffs = global_rfm.getElementsByTagName("RFM_Validity")[0] + rpc_params["offset_row"] = float(normalisation_coeffs.getElementsByTagName("LINE_OFF")[0].firstChild.data) + rpc_params["offset_col"] = float(normalisation_coeffs.getElementsByTagName("SAMP_OFF")[0].firstChild.data) + if float(version) >= 3: + direct_coeffs = global_rfm.getElementsByTagName("ImagetoGround_Values")[0] + inverse_coeffs = global_rfm.getElementsByTagName("GroundtoImage_Values")[0] + + rpc_params["num_x"] = [ + float(direct_coeffs.getElementsByTagName(f"LON_NUM_COEFF_{index}")[0].firstChild.data) + for index in range(1, 21) + ] + rpc_params["den_x"] = [ + float(direct_coeffs.getElementsByTagName(f"LON_DEN_COEFF_{index}")[0].firstChild.data) + for index in range(1, 21) + ] + rpc_params["num_y"] = [ + float(direct_coeffs.getElementsByTagName(f"LAT_NUM_COEFF_{index}")[0].firstChild.data) + for index in range(1, 21) + ] + rpc_params["den_y"] = [ + float(direct_coeffs.getElementsByTagName(f"LAT_DEN_COEFF_{index}")[0].firstChild.data) + for index in range(1, 21) + ] + rpc_params["offset_col"] -= 0.5 + rpc_params["offset_row"] -= 0.5 + + else: + direct_coeffs = global_rfm.getElementsByTagName("Direct_Model")[0] + inverse_coeffs = global_rfm.getElementsByTagName("Inverse_Model")[0] + + rpc_params["num_x"] = [ + float(direct_coeffs.getElementsByTagName(f"SAMP_NUM_COEFF_{index}")[0].firstChild.data) + for index in range(1, 21) + ] + rpc_params["den_x"] = [ + float(direct_coeffs.getElementsByTagName(f"SAMP_DEN_COEFF_{index}")[0].firstChild.data) + for index in range(1, 21) + ] + rpc_params["num_y"] = [ + float(direct_coeffs.getElementsByTagName(f"LINE_NUM_COEFF_{index}")[0].firstChild.data) + for index in range(1, 21) + ] + rpc_params["den_y"] = [ + float(direct_coeffs.getElementsByTagName(f"LINE_DEN_COEFF_{index}")[0].firstChild.data) + for index in range(1, 21) + ] + rpc_params["offset_col"] -= 1.0 + rpc_params["offset_row"] -= 1.0 + + rpc_params["num_col"] = [ + float(inverse_coeffs.getElementsByTagName(f"SAMP_NUM_COEFF_{index}")[0].firstChild.data) + for index in range(1, 21) + ] + rpc_params["den_col"] = [ + float(inverse_coeffs.getElementsByTagName(f"SAMP_DEN_COEFF_{index}")[0].firstChild.data) + for index in range(1, 21) + ] + rpc_params["num_row"] = [ + float(inverse_coeffs.getElementsByTagName(f"LINE_NUM_COEFF_{index}")[0].firstChild.data) + for index in range(1, 21) + ] + rpc_params["den_row"] = [ + float(inverse_coeffs.getElementsByTagName(f"LINE_DEN_COEFF_{index}")[0].firstChild.data) + for index in range(1, 21) + ] + + rpc_params["scale_col"] = float(normalisation_coeffs.getElementsByTagName("SAMP_SCALE")[0].firstChild.data) + rpc_params["scale_row"] = float(normalisation_coeffs.getElementsByTagName("LINE_SCALE")[0].firstChild.data) + rpc_params["offset_alt"] = float(normalisation_coeffs.getElementsByTagName("HEIGHT_OFF")[0].firstChild.data) + rpc_params["scale_alt"] = float(normalisation_coeffs.getElementsByTagName("HEIGHT_SCALE")[0].firstChild.data) + rpc_params["offset_x"] = float(normalisation_coeffs.getElementsByTagName("LONG_OFF")[0].firstChild.data) + rpc_params["scale_x"] = float(normalisation_coeffs.getElementsByTagName("LONG_SCALE")[0].firstChild.data) + rpc_params["offset_y"] = float(normalisation_coeffs.getElementsByTagName("LAT_OFF")[0].firstChild.data) + rpc_params["scale_y"] = float(normalisation_coeffs.getElementsByTagName("LAT_SCALE")[0].firstChild.data) + # If top left convention, 0.5 pixel shift added on col/row offsets + if topleftconvention: + rpc_params["offset_col"] += 0.5 + rpc_params["offset_row"] += 0.5 + return rpc_params + + +def rpc_reader_dimap_v1(geomodel_path: str, topleftconvention: bool = True) -> Dict: + """ + Load from dimap v1 + + ** Deprecated, to clean ? ** + + :param geomodel_path: path to geomodel + :param topleftconvention: [0,0] position + :type topleftconvention: boolean + If False : [0,0] is at the center of the Top Left pixel + If True : [0,0] is at the top left of the Top Left pixel (OSSIM) + """ + + if not basename(geomodel_path).upper().endswith("XML"): + raise ValueError("dimap must ends with .xml") + + xmldoc = minidom.parse(geomodel_path) + + mtd = xmldoc.getElementsByTagName("Metadata_Identification") + version = mtd[0].getElementsByTagName("METADATA_PROFILE")[0].attributes.items()[0][1] + rpc_params = {"driver_type": "dimap_v" + version} + + global_rfm = xmldoc.getElementsByTagName("Global_RFM") + rfm_validity = xmldoc.getElementsByTagName("RFM_Validity") + coeff_lon = [float(el) for el in global_rfm[0].getElementsByTagName("F_LON")[0].firstChild.data.split()] + coeff_lat = [float(el) for el in global_rfm[0].getElementsByTagName("F_LAT")[0].firstChild.data.split()] + coeff_col = [float(el) for el in global_rfm[0].getElementsByTagName("F_COL")[0].firstChild.data.split()] + coeff_lig = [float(el) for el in global_rfm[0].getElementsByTagName("F_ROW")[0].firstChild.data.split()] + + scale_lon = float(rfm_validity[0].getElementsByTagName("Lon")[0].getElementsByTagName("A")[0].firstChild.data) + offset_lon = float(rfm_validity[0].getElementsByTagName("Lon")[0].getElementsByTagName("B")[0].firstChild.data) + scale_lat = float(rfm_validity[0].getElementsByTagName("Lat")[0].getElementsByTagName("A")[0].firstChild.data) + offset_lat = float(rfm_validity[0].getElementsByTagName("Lat")[0].getElementsByTagName("B")[0].firstChild.data) + scale_alt = float(rfm_validity[0].getElementsByTagName("Alt")[0].getElementsByTagName("A")[0].firstChild.data) + offset_alt = float(rfm_validity[0].getElementsByTagName("Alt")[0].getElementsByTagName("B")[0].firstChild.data) + scale_col = float(rfm_validity[0].getElementsByTagName("Col")[0].getElementsByTagName("A")[0].firstChild.data) + offset_col = float(rfm_validity[0].getElementsByTagName("Col")[0].getElementsByTagName("B")[0].firstChild.data) + scale_row = float(rfm_validity[0].getElementsByTagName("Row")[0].getElementsByTagName("A")[0].firstChild.data) + offset_row = float(rfm_validity[0].getElementsByTagName("Row")[0].getElementsByTagName("B")[0].firstChild.data) + + rpc_params["offset_col"] = offset_col + rpc_params["scale_col"] = scale_col + rpc_params["offset_row"] = offset_row + rpc_params["scale_row"] = scale_row + rpc_params["offset_alt"] = offset_alt + rpc_params["scale_alt"] = scale_alt + rpc_params["offset_x"] = offset_lon + rpc_params["scale_x"] = scale_lon + rpc_params["offset_y"] = offset_lat + rpc_params["scale_y"] = scale_lat + rpc_params["num_x"] = coeff_lon[0:20] + rpc_params["den_x"] = coeff_lon[20::] + rpc_params["num_y"] = coeff_lat[0:20] + rpc_params["den_y"] = coeff_lat[20::] + rpc_params["num_col"] = coeff_col[0:20] + rpc_params["den_col"] = coeff_col[20::] + rpc_params["num_row"] = coeff_lig[0:20] + rpc_params["den_row"] = coeff_lig[20::] + rpc_params["offset_col"] -= 1.0 + rpc_params["offset_row"] -= 1.0 + # If top left convention, 0.5 pixel shift added on col/row offsets + if topleftconvention: + rpc_params["offset_col"] += 0.5 + rpc_params["offset_row"] += 0.5 + return rpc_params + + +def rpc_reader_geotiff(geomodel_path, topleftconvention=True) -> Dict: + """ + Load from a geotiff image file + + :param geomodel_path: path to geomodel + :param topleftconvention: [0,0] position + :type topleftconvention: boolean + If False : [0,0] is at the center of the Top Left pixel + If True : [0,0] is at the top left of the Top Left pixel (OSSIM) + """ + dataset = rio.open(geomodel_path) + rpc_dict = dataset.tags(ns="RPC") + if not rpc_dict: + logging.error("%s does not contains RPCs ", geomodel_path) + raise ValueError + rpc_params = { + "den_row": parse_coeff_line(rpc_dict["LINE_DEN_COEFF"]), + "num_row": parse_coeff_line(rpc_dict["LINE_NUM_COEFF"]), + "num_col": parse_coeff_line(rpc_dict["SAMP_NUM_COEFF"]), + "den_col": parse_coeff_line(rpc_dict["SAMP_DEN_COEFF"]), + "offset_col": float(rpc_dict["SAMP_OFF"]), + "scale_col": float(rpc_dict["SAMP_SCALE"]), + "offset_row": float(rpc_dict["LINE_OFF"]), + "scale_row": float(rpc_dict["LINE_SCALE"]), + "offset_alt": float(rpc_dict["HEIGHT_OFF"]), + "scale_alt": float(rpc_dict["HEIGHT_SCALE"]), + "offset_x": float(rpc_dict["LONG_OFF"]), + "scale_x": float(rpc_dict["LONG_SCALE"]), + "offset_y": float(rpc_dict["LAT_OFF"]), + "scale_y": float(rpc_dict["LAT_SCALE"]), + "num_x": None, + "den_x": None, + "num_y": None, + "den_y": None, + } + # inverse coeff are not defined + # If top left convention, 0.5 pixel shift added on col/row offsets + if topleftconvention: + rpc_params["offset_col"] += 0.5 + rpc_params["offset_row"] += 0.5 + return rpc_params + + +def rpc_reader_ossim_kwl(geomodel_path: str, topleftconvention: bool = True) -> Dict: + """ + Load from a geom file + + :param geomodel_path: path to geomodel + :param topleftconvention: [0,0] position + :type topleftconvention: boolean + If False : [0,0] is at the center of the Top Left pixel + If True : [0,0] is at the top left of the Top Left pixel (OSSIM) + """ + # OSSIM keyword list + rpc_params = {"driver_type": "ossim_kwl"} + with open(geomodel_path, "r", encoding="utf-8") as ossim_file: + content = ossim_file.readlines() + + geom_dict = {} + for line in content: + (key, val) = line.split(": ") + geom_dict[key] = val.rstrip() + + rpc_params["den_row"] = [np.nan] * 20 + rpc_params["num_row"] = [np.nan] * 20 + rpc_params["den_col"] = [np.nan] * 20 + rpc_params["num_col"] = [np.nan] * 20 + for index in range(0, 20): + axis = "line" + num_den = "den" + key = f"{axis}_{num_den}_coeff_{index:02d}" + rpc_params["den_row"][index] = float(geom_dict[key]) + num_den = "num" + key = f"{axis}_{num_den}_coeff_{index:02d}" + rpc_params["num_row"][index] = float(geom_dict[key]) + axis = "samp" + key = f"{axis}_{num_den}_coeff_{index:02d}" + rpc_params["num_col"][index] = float(geom_dict[key]) + num_den = "den" + key = f"{axis}_{num_den}_coeff_{index:02d}" + rpc_params["den_col"][index] = float(geom_dict[key]) + rpc_params["offset_col"] = float(geom_dict["samp_off"]) + rpc_params["scale_col"] = float(geom_dict["samp_scale"]) + rpc_params["offset_row"] = float(geom_dict["line_off"]) + rpc_params["scale_row"] = float(geom_dict["line_scale"]) + rpc_params["offset_alt"] = float(geom_dict["height_off"]) + rpc_params["scale_alt"] = float(geom_dict["height_scale"]) + rpc_params["offset_x"] = float(geom_dict["long_off"]) + rpc_params["scale_x"] = float(geom_dict["long_scale"]) + rpc_params["offset_y"] = float(geom_dict["lat_off"]) + rpc_params["scale_y"] = float(geom_dict["lat_scale"]) + # inverse coeff are not defined + rpc_params["num_x"] = None + rpc_params["den_x"] = None + rpc_params["num_y"] = None + rpc_params["den_y"] = None + # If top left convention, 0.5 pixel shift added on col/row offsets + if topleftconvention: + rpc_params["offset_col"] += 0.5 + rpc_params["offset_row"] += 0.5 + return rpc_params + + +def identify_dimap(xml_file): + """ + parse xml file to identify dimap and its version + + :param xml_file: dimap rpc file + :type xml_file: str + :return: dimap info : dimap_version and None if not an dimap file + :rtype: str + """ + try: + xmldoc = minidom.parse(xml_file) + mtd = xmldoc.getElementsByTagName("Metadata_Identification") + mtd_format = mtd[0].getElementsByTagName("METADATA_FORMAT")[0].firstChild.data + if mtd_format == "DIMAP_PHR": + version_tag = "METADATA_PROFILE" + else: + version_tag = "METADATA_FORMAT" + version = mtd[0].getElementsByTagName(version_tag)[0].attributes.items()[0][1] + except Exception: # pylint: disable=broad-except + return None + + return version + + +def parse_coeff_line(coeff_str): + """ + split str coef to float list + + :param coeff_str: line coef + :type coeff_str: str + :return: coeff list + :rtype: list() + """ + return [float(el) for el in coeff_str.split()] + + +def identify_geotiff_rpc(image_filename): + """ + read image file to identify if it is a geotiff which contains RPCs + + :param image_filename: image_filename + :type image_filename: str + :return: rpc info, rpc dict or None if not a geotiff with rpc + :rtype: str + """ + try: + dataset = rio.open(image_filename) + rpc_dict = dataset.tags(ns="RPC") + return rpc_dict + except Exception: # pylint: disable=broad-except + return None + + +def identify_ossim_kwl(ossim_kwl_file): + """ + parse geom file to identify if it is an ossim model + + :param ossim_kwl_fil : ossim keyword list file + :type ossim_kwl_file: str + :return: ossimmodel or None if not an ossim kwl file + :rtype: str + """ + try: + with open(ossim_kwl_file, encoding="utf-8") as ossim_file: + content = ossim_file.readlines() + geom_dict = {} + for line in content: + (key, val) = line.split(": ") + geom_dict[key] = val.rstrip() + if "type" in geom_dict: + if geom_dict["type"].strip().startswith("ossim"): + return geom_dict["type"].strip() + return None + except Exception: # pylint: disable=broad-except + return None diff --git a/tests/geofunctions/test_localization.py b/tests/geofunctions/test_localization.py index 8fe2882..93866b6 100644 --- a/tests/geofunctions/test_localization.py +++ b/tests/geofunctions/test_localization.py @@ -96,7 +96,7 @@ def test_localize_direct_grid(): # data = os.path.join(data_path(), "rpc/phr_ventoux/", "left_image.geom") # geom_model_1 = GeoModel(data) data = os.path.join(data_path(), "grid/phr_ventoux/", "GRID_PHR1B_P_201308051042194_SEN_690908101-001.tif") - geom_model = GeoModel(data) + geom_model = GeoModel(data, "grid") # then read the Image to retrieve its geotransform image_filename = os.path.join(data_path(), "image/phr_ventoux/", "left_image.tif") image_left = Image(image_filename) @@ -136,7 +136,7 @@ def prepare_loc(alti="geoide", id_scene="P1BP--2017030824934340CP"): dtm = DTMIntersection(fic) # load grid model gld = os.path.join(data_folder, grid_name) - gri = GeoModel(gld) + gri = GeoModel(gld, "grid") return dtm, gri diff --git a/tests/geofunctions/test_rectification.py b/tests/geofunctions/test_rectification.py index 31e7887..1a82b09 100644 --- a/tests/geofunctions/test_rectification.py +++ b/tests/geofunctions/test_rectification.py @@ -256,10 +256,10 @@ def test_compute_stereorectification_epipolar_grids_geomodel_grid(): # first instantiate geometric models left and right (here Grid geometric model) geom_model_left = GeoModel( - os.path.join(data_path(), "grid/phr_ventoux/GRID_PHR1B_P_201308051042194_SEN_690908101-001.tif") + os.path.join(data_path(), "grid/phr_ventoux/GRID_PHR1B_P_201308051042194_SEN_690908101-001.tif"), "grid" ) geom_model_right = GeoModel( - os.path.join(data_path(), "grid/phr_ventoux/GRID_PHR1B_P_201308051042523_SEN_690908101-002.tif") + os.path.join(data_path(), "grid/phr_ventoux/GRID_PHR1B_P_201308051042523_SEN_690908101-002.tif"), "grid" ) # read the images @@ -312,10 +312,10 @@ def test_compute_stereorectification_epipolar_grids_geomodel_grid_dtm_geoid(): # first instantiate geometric models left and right (here Grid geometrics model) geom_model_left = GeoModel( - os.path.join(data_path(), "grid/phr_ventoux/GRID_PHR1B_P_201308051042194_SEN_690908101-001.tif") + os.path.join(data_path(), "grid/phr_ventoux/GRID_PHR1B_P_201308051042194_SEN_690908101-001.tif"), "grid" ) geom_model_right = GeoModel( - os.path.join(data_path(), "grid/phr_ventoux/GRID_PHR1B_P_201308051042523_SEN_690908101-002.tif") + os.path.join(data_path(), "grid/phr_ventoux/GRID_PHR1B_P_201308051042523_SEN_690908101-002.tif"), "grid" ) # read the images diff --git a/tests/geofunctions/test_triangulation.py b/tests/geofunctions/test_triangulation.py index 07b9800..d6e048a 100644 --- a/tests/geofunctions/test_triangulation.py +++ b/tests/geofunctions/test_triangulation.py @@ -85,7 +85,7 @@ def prepare_loc(alti="geoide", id_scene="P1BP--2017030824934340CP"): data_folder = data_path(alti, id_scene) # load grid gld = os.path.join(data_folder, f"GRID_{id_scene}.tif") - gri = GeoModel(gld) + gri = GeoModel(gld, "grid") return gri diff --git a/tests/geomodels/test_grid.py b/tests/geomodels/test_grid.py index 4ab9cb6..f40e58d 100755 --- a/tests/geomodels/test_grid.py +++ b/tests/geomodels/test_grid.py @@ -41,7 +41,7 @@ def fixture_get_geotiff_grid(): get grid and associated path """ geotiff_grid_path = data_path("ellipsoide", "loc_direct_grid_PHR_2013072139303958CP.tif") - gri_geotiff = GeoModel(geotiff_grid_path) + gri_geotiff = GeoModel(geotiff_grid_path, "grid") grid_image = Image(geotiff_grid_path, read_data=True) return gri_geotiff, grid_image diff --git a/tests/geomodels/test_rpc.py b/tests/geomodels/test_rpc.py index 42ae96c..1829e4b 100755 --- a/tests/geomodels/test_rpc.py +++ b/tests/geomodels/test_rpc.py @@ -36,7 +36,7 @@ # Shareloc imports from shareloc.geofunctions.dtm_intersection import DTMIntersection from shareloc.geomodels import GeoModel -from shareloc.geomodels.rpc import identify_dimap, identify_ossim_kwl +from shareloc.geomodels.rpc_readers import identify_dimap, identify_ossim_kwl # Shareloc test imports from ..helpers import data_path From 3cb1e3d0643a68f2ec0f275ed6161f215c864b0b Mon Sep 17 00:00:00 2001 From: GUINET Jonathan Date: Tue, 21 Nov 2023 12:29:08 +0000 Subject: [PATCH 06/12] update load() interface --- shareloc/geomodels/geomodel_template.py | 2 +- shareloc/geomodels/rpc.py | 4 ++-- 2 files changed, 3 insertions(+), 3 deletions(-) diff --git a/shareloc/geomodels/geomodel_template.py b/shareloc/geomodels/geomodel_template.py index 388c242..9804ae7 100644 --- a/shareloc/geomodels/geomodel_template.py +++ b/shareloc/geomodels/geomodel_template.py @@ -102,7 +102,7 @@ def inverse_loc(self, lon, lat, alt): @classmethod @abstractmethod - def load(cls, geomodel_path, **kwargs): + def load(cls, geomodel_path): """ load function with class specific args """ diff --git a/shareloc/geomodels/rpc.py b/shareloc/geomodels/rpc.py index 7ed9190..2f37850 100755 --- a/shareloc/geomodels/rpc.py +++ b/shareloc/geomodels/rpc.py @@ -172,7 +172,7 @@ def __init__(self, rpc_params): self.rowmax = self.offset_row + self.scale_row @classmethod - def load(cls, geomodel_path, topleftconvention=True): + def load(cls, geomodel_path): """ Load from any RPC (auto identify driver) from filename (dimap, ossim kwl, geotiff) @@ -185,7 +185,7 @@ def load(cls, geomodel_path, topleftconvention=True): """ # Set topleftconvention (keeping historic option): to clean cls.geomodel_path = geomodel_path - return cls(rpc_reader(geomodel_path, topleftconvention)) + return cls(rpc_reader(geomodel_path, topleftconvention=True)) def direct_loc_h(self, row, col, alt, fill_nan=False): """ From 885bfac2e6be0df706f7d52ccff6ae20516a0139 Mon Sep 17 00:00:00 2001 From: Jonathan Guinet Date: Wed, 29 Nov 2023 17:00:43 +0100 Subject: [PATCH 07/12] doc: update doc in functions arguments --- docs/source/getting_started.rst | 4 ++-- docs/source/user_manual_functions.rst | 14 +++++++------- docs/source/user_manual_geometric_models.rst | 8 ++++---- shareloc/geofunctions/localization.py | 6 +++--- shareloc/geofunctions/rectification.py | 20 ++++++++++---------- shareloc/geofunctions/triangulation.py | 8 ++++---- 6 files changed, 30 insertions(+), 30 deletions(-) diff --git a/docs/source/getting_started.rst b/docs/source/getting_started.rst index 57b5726..71966b4 100644 --- a/docs/source/getting_started.rst +++ b/docs/source/getting_started.rst @@ -35,12 +35,12 @@ Quick Start $ python3 >>> # Import shareloc modules rpc and localization - >>> from shareloc.geomodels.rpc import RPC + >>> from shareloc.geomodels import GeoModel >>> from shareloc.geofunctions.localization import Localization >>> # Create RPC object from downloaded geometry file >>> rpc_geom_file = "left_image.geom" - >>> rpc = RPC(rpc_geom_file) + >>> rpc = GeoModel(rpc_geom_file) >>> # Create Localization object from created RPC >>> loc = Localization(rpc) diff --git a/docs/source/user_manual_functions.rst b/docs/source/user_manual_functions.rst index b9c8fb1..9bec2af 100644 --- a/docs/source/user_manual_functions.rst +++ b/docs/source/user_manual_functions.rst @@ -40,7 +40,7 @@ It is possible to use geometric model directly using ``shareloc.grid.Grid`` and """ constructor :param model : geometric model - :type model : shareloc.grid or shareloc.rpc + :type model : GeoModelTemplate :param elevation : dtm or default elevation over ellipsoid if None elevation is set to 0 :type elevation : shareloc.dtm or float or np.ndarray :param image : image class to handle geotransform @@ -108,9 +108,9 @@ colocalization returns image positions (row2,col2) in image 2 from (row1,col1) p Colocalization : direct localization with model1, then inverse localization with model2 :param model1: geometric model 1 - :type model1: shareloc.grid or shareloc.rpc + :type model1: GeomodelTemplate :param model2: geometric model 2 - :type model2: shareloc.grid or shareloc.rpc + :type model2: GeomodelTemplate :param row: sensor row :type row: int or 1D numpy array :param col: sensor col @@ -163,9 +163,9 @@ where :math:`v_i` is the orientation of the :term:`LOS` i and :math:`s_i` the ha :param matches : matches in sensor coordinates Nx[row (left), col (left), row (right), col (right)] :type matches : np.array :param geometrical_model_left : left image geometrical model - :type geometrical_model_left : shareloc.grid or shareloc.rpc + :type geometrical_model_left : GeomodelTemplate :param geometrical_model_right : right image geometrical model - :type geometrical_model_right : shareloc.grid or shareloc.rpc + :type geometrical_model_right : GeomodelTemplate :param left_min_max : left min/max for los creation, if None model min/max will be used :type left_min_max : list :param right_min_max : right min/max for los creation, if None model min/max will be used @@ -207,11 +207,11 @@ Algorithm details can be found in reference below. :param left_im: left image :type left_im: shareloc.image object :param geom_model_left: geometric model of the left image - :type geom_model_left: shareloc.grid or shareloc.rpc + :type geom_model_left: GeomodelTemplate :param right_im: right image :type right_im: shareloc.image object :param geom_model_right: geometric model of the right image - :type geom_model_right: shareloc.grid or shareloc.rpc + :type geom_model_right: GeomodelTemplate :param elevation: elevation :type elevation: shareloc.dtm or float :param epi_step: epipolar step diff --git a/docs/source/user_manual_geometric_models.rst b/docs/source/user_manual_geometric_models.rst index fa70023..921742e 100644 --- a/docs/source/user_manual_geometric_models.rst +++ b/docs/source/user_manual_geometric_models.rst @@ -50,9 +50,9 @@ RPC class API Example $ wget https://raw.githubusercontent.com/CNES/shareloc/tests/data/rpc/RPC_PHR1B_P_201709281038393_SEN_PRG_FC_178609-001.XML $ python3 - >>> from shareloc.geomodels.rpc import RPC + >>> from shareloc.geomodels import GeoModel >>> file_dimap = "RPC_PHR1B_P_201709281038393_SEN_PRG_FC_178609-001.XML") - >>> rpc_dimap = RPC(file_dimap) + >>> rpc_dimap = GeoModel(file_dimap) Direct location grids @@ -140,9 +140,9 @@ Grid API Example $ wget https://raw.githubusercontent.com/CNES/shareloc/tests/data/ellipsoide/loc_direct_grid_PHR_2013072139303958CP.tif $ python3 - >>> from shareloc.geomodels.grid import Grid + >>> from shareloc.geomodels import GeoModel >>> geotiff_grid_path = "loc_direct_grid_PHR_2013072139303958CP.tif" - >>> geotiff_grid = Grid(geotiff_grid_path) + >>> geotiff_grid = GeoModel(geotiff_grid_path, "grid") References __________ diff --git a/shareloc/geofunctions/localization.py b/shareloc/geofunctions/localization.py index c0f58c8..7a714e1 100755 --- a/shareloc/geofunctions/localization.py +++ b/shareloc/geofunctions/localization.py @@ -45,7 +45,7 @@ def __init__(self, model, elevation=None, image=None, epsg=None): Localization constructor :param model: geometric model - :type model: shareloc.grid or shareloc.rpc + :type model: GeomodelTemplate :param elevation: dtm or default elevation over ellipsoid if None elevation is set to 0 :type elevation: shareloc.dtm or float or np.ndarray :param image: image class to handle geotransform @@ -167,9 +167,9 @@ def coloc(model1, model2, row, col, elevation=None, image1=None, image2=None, us Colocalization : direct localization with model1, then inverse localization with model2 :param model1: geometric model 1 - :type model1: shareloc.grid or shareloc.rpc + :type model1: GeomodelTemplate :param model2: geometric model 2 - :type model2: shareloc.grid or shareloc.rpc + :type model2: GeomodelTemplate :param row: sensor row :type row: int or 1D numpy array :param col: sensor col diff --git a/shareloc/geofunctions/rectification.py b/shareloc/geofunctions/rectification.py index c4677fd..e2c4e59 100644 --- a/shareloc/geofunctions/rectification.py +++ b/shareloc/geofunctions/rectification.py @@ -108,9 +108,9 @@ def compute_local_epipolar_line(geom_model_left, geom_model_right, left_point, e Estimate the beginning and the ending of local epipolar line in left image :param geom_model_left: geometric model of the left image - :type geom_model_left: shareloc.grid or shareloc.rpc + :type geom_model_left: GeomodelTemplate :param geom_model_right: geometric model of the right image - :type geom_model_right: shareloc.grid or shareloc.rpc + :type geom_model_right: GeomodelTemplate :param left_point: georeferenced coordinates in the left image :type left_point: 1D numpy array : [row coord, col coord, altitude] or 2D numpy array : (number of points, [row coord, col coord, altitude]) @@ -159,9 +159,9 @@ def prepare_rectification(left_im, geom_model_left, geom_model_right, elevation, :param left_im: left image :type left_im: shareloc.image object :param geom_model_left: geometric model of the left image - :type geom_model_left: shareloc.grid or shareloc.rpc + :type geom_model_left: GeomodelTemplate :param geom_model_right: geometric model of the right image - :type geom_model_right: shareloc.grid or shareloc.rpc + :type geom_model_right: GeomodelTemplate :param elevation: elevation :type elevation: shareloc.dtm or float :param epi_step: epipolar step @@ -269,9 +269,9 @@ def get_epipolar_extent( :param left_im: left image :type left_im: shareloc.image object :param geom_model_left: geometric model of the left image - :type geom_model_left: shareloc.grid or shareloc.rpc + :type geom_model_left: GeomodelTemplate :param geom_model_right: geometric model of the right image - :type geom_model_right: shareloc.grid or shareloc.rpc + :type geom_model_right: GeomodelTemplate :param elevation: elevation :type elevation: shareloc.dtm or float :param epi_step: epipolar step @@ -333,9 +333,9 @@ def moving_along_axis( Moving to the next line in epipolar geometry :param geom_model_left: geometric model of the left image - :type geom_model_left: shareloc.grid or shareloc.rpc + :type geom_model_left: GeomodelTemplate :param geom_model_right: geometric model of the right image - :type geom_model_right: shareloc.grid or shareloc.rpc + :type geom_model_right: GeoModelTemplate :param current_coords: current line in the left epipolar geometry or current georeferenced coordinates in left epipolar line :type current_coords: 1D np.array [row, col, altitude] @@ -389,11 +389,11 @@ def compute_stereorectification_epipolar_grids( :param left_im: left image :type left_im: shareloc.image object :param geom_model_left: geometric model of the left image - :type geom_model_left: shareloc.grid or shareloc.rpc + :type geom_model_left: GeomodelTemplate :param right_im: right image :type right_im: shareloc.image object :param geom_model_right: geometric model of the right image - :type geom_model_right: shareloc.grid or shareloc.rpc + :type geom_model_right: GeomodelTemplate :param elevation: elevation :type elevation: shareloc.dtm or float :param epi_step: epipolar step diff --git a/shareloc/geofunctions/triangulation.py b/shareloc/geofunctions/triangulation.py index 075d765..2e4c5d5 100644 --- a/shareloc/geofunctions/triangulation.py +++ b/shareloc/geofunctions/triangulation.py @@ -54,9 +54,9 @@ def sensor_triangulation( :param matches: matches in sensor coordinates Nx[row (left), col (left), row (right), col (right)] :type matches: np.array :param geometrical_model_left: left image geometrical model - :type geometrical_model_left: shareloc.grid or shareloc.rpc + :type geometrical_model_left: GeomodelTemplate :param geometrical_model_right: right image geometrical model - :type geometrical_model_right: shareloc.grid or shareloc.rpc + :type geometrical_model_right: GeomodelTemplate :param left_min_max: left min/max for los creation, if None model min/max will be used :type left_min_max: list :param right_min_max: right min/max for los creation, if None model min/max will be used @@ -190,9 +190,9 @@ def epipolar_triangulation( :param matches_type: 'disp' or 'sift' :type matches_type: str :param geometrical_model_left: left image geometrical model - :type geometrical_model_left: shareloc.grid or shareloc.rpc + :type geometrical_model_left: GeomodelTemplate :param geometrical_model_right: right image geometrical model - :type geometrical_model_right: shareloc.grid or shareloc.rpc + :type geometrical_model_right: GeomodelTemplate :param grid_left: left rectification grid filename :type grid_left: str :param grid_right: right rectification grid filename From aa5708ea0a86c388ca6878a646dc975ada6495b6 Mon Sep 17 00:00:00 2001 From: Jonathan Guinet Date: Fri, 1 Dec 2023 08:26:04 +0100 Subject: [PATCH 08/12] doc: doc update explicit model type --- docs/source/getting_started.rst | 2 +- docs/source/user_manual_geometric_models.rst | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/docs/source/getting_started.rst b/docs/source/getting_started.rst index 71966b4..db7fe58 100644 --- a/docs/source/getting_started.rst +++ b/docs/source/getting_started.rst @@ -40,7 +40,7 @@ Quick Start >>> # Create RPC object from downloaded geometry file >>> rpc_geom_file = "left_image.geom" - >>> rpc = GeoModel(rpc_geom_file) + >>> rpc = GeoModel(rpc_geom_file, "rpc") >>> # Create Localization object from created RPC >>> loc = Localization(rpc) diff --git a/docs/source/user_manual_geometric_models.rst b/docs/source/user_manual_geometric_models.rst index 921742e..6331ad5 100644 --- a/docs/source/user_manual_geometric_models.rst +++ b/docs/source/user_manual_geometric_models.rst @@ -52,7 +52,7 @@ RPC class API Example $ python3 >>> from shareloc.geomodels import GeoModel >>> file_dimap = "RPC_PHR1B_P_201709281038393_SEN_PRG_FC_178609-001.XML") - >>> rpc_dimap = GeoModel(file_dimap) + >>> rpc_dimap = GeoModel(file_dimap, "rpc") Direct location grids From 3f83471839616065d078dd07634dcf2f2e24a378 Mon Sep 17 00:00:00 2001 From: Jonathan Guinet Date: Fri, 1 Dec 2023 08:57:29 +0100 Subject: [PATCH 09/12] refac: moved class member init --- shareloc/geomodels/geomodel_template.py | 12 +++------ shareloc/geomodels/grid.py | 31 +++++++++++++++--------- shareloc/geomodels/rpc.py | 2 -- tests/geofunctions/test_rectification.py | 8 ++---- tests/geomodels/test_los.py | 2 +- 5 files changed, 26 insertions(+), 29 deletions(-) diff --git a/shareloc/geomodels/geomodel_template.py b/shareloc/geomodels/geomodel_template.py index 9804ae7..599b406 100644 --- a/shareloc/geomodels/geomodel_template.py +++ b/shareloc/geomodels/geomodel_template.py @@ -27,14 +27,6 @@ # Standard imports from abc import ABCMeta, abstractmethod -# Global variable for optimization mode (functions in C) -# SHARELOC_OPTIM_GEOMODEL = False - -# TODO: Override functions depending on optimization or not - -# if(SHARELOC_OPTIM_GEOMODEL == True): -# GeoModelTemplate.direct_loc_dtm = GeoModelTemplate.direct_loc_dtm_optim - class GeoModelTemplate(metaclass=ABCMeta): """ @@ -102,7 +94,9 @@ def inverse_loc(self, lon, lat, alt): @classmethod @abstractmethod - def load(cls, geomodel_path): + def load(cls, geomodel_path: str): """ load function with class specific args + + :param geomodel_path: filename of geomodel """ diff --git a/shareloc/geomodels/grid.py b/shareloc/geomodels/grid.py index ae1b53c..2c513c8 100755 --- a/shareloc/geomodels/grid.py +++ b/shareloc/geomodels/grid.py @@ -38,7 +38,6 @@ # gitlab issue #58 # pylint: disable=too-many-instance-attributes -# pylint: disable=no-member @GeoModel.register("grid") class Grid(GeoModelTemplate): """ @@ -105,6 +104,17 @@ def __init__(self, geomodel_path: str): self.rowmax = None self.colmax = None self.epsg = 0 + + # inverse loc predictor attributes + self.pred_col_min = None + self.pred_row_min = None + self.pred_col_max = None + self.pred_row_max = None + self.pred_ofset_scale_lon = None + self.pred_ofset_scale_lat = None + self.pred_ofset_scale_row = None + self.pred_ofset_scale_col = None + self.read() @classmethod @@ -545,7 +555,7 @@ def estimate_inverse_loc_predictor(self, nbrow_pred=3, nbcol_pred=3): a_max[imes, 5] = glonlat[0, irow, icol] imes += 1 - # Calcul des coeffcients + # Compute coefficients mat_a_min = np.array(a_min) mat_a_max = np.array(a_max) @@ -559,16 +569,15 @@ def estimate_inverse_loc_predictor(self, nbrow_pred=3, nbcol_pred=3): coef_col_max = t_aa_max_inv @ mat_a_max.T @ b_col coef_row_max = t_aa_max_inv @ mat_a_max.T @ b_row - # TODO: refactor to have only class parameters in __init__ # seems not clear to understand with inverse_loc_predictor function ... - setattr(self, "pred_col_min", coef_col_min.flatten()) # noqa: 502 - setattr(self, "pred_row_min", coef_row_min.flatten()) # noqa: 502 - setattr(self, "pred_col_max", coef_col_max.flatten()) # noqa: 502 - setattr(self, "pred_row_max", coef_row_max.flatten()) # noqa: 502 - setattr(self, "pred_ofset_scale_lon", [lon_ofset, lon_scale]) # noqa: 502 - setattr(self, "pred_ofset_scale_lat", [lat_ofset, lat_scale]) # noqa: 502 - setattr(self, "pred_ofset_scale_row", [row_ofset, row_scale]) # noqa: 502 - setattr(self, "pred_ofset_scale_col", [col_ofset, col_scale]) # noqa: 502 + self.pred_col_min = coef_col_min.flatten() + self.pred_row_min = coef_row_min.flatten() + self.pred_col_max = coef_col_max.flatten() + self.pred_row_max = coef_row_max.flatten() + self.pred_ofset_scale_lon = [lon_ofset, lon_scale] + self.pred_ofset_scale_lat = [lat_ofset, lat_scale] + self.pred_ofset_scale_row = [row_ofset, row_scale] + self.pred_ofset_scale_col = [col_ofset, col_scale] def inverse_loc_predictor(self, lon, lat, alt=0.0): """ diff --git a/shareloc/geomodels/rpc.py b/shareloc/geomodels/rpc.py index 2f37850..d9c8373 100755 --- a/shareloc/geomodels/rpc.py +++ b/shareloc/geomodels/rpc.py @@ -49,8 +49,6 @@ class RPC(GeoModelTemplate): RPC class including direct and inverse localization instance methods """ - # gitlab issue #61 - # pylint: disable=too-many-instance-attributes # gitlab issue #61 # pylint: disable=too-many-instance-attributes def __init__(self, rpc_params): diff --git a/tests/geofunctions/test_rectification.py b/tests/geofunctions/test_rectification.py index 1a82b09..e249cb6 100644 --- a/tests/geofunctions/test_rectification.py +++ b/tests/geofunctions/test_rectification.py @@ -548,12 +548,8 @@ def test_rectification_moving_to_axis_error(): """ Test moving along axis with wrong axis """ - geom_model_left = RPC.from_any( - os.path.join(data_path(), "rectification", "left_image.geom"), topleftconvention=True - ) - geom_model_right = RPC.from_any( - os.path.join(data_path(), "rectification", "right_image.geom"), topleftconvention=True - ) + geom_model_left = GeoModel(os.path.join(data_path(), "rectification", "left_image.geom")) + geom_model_right = GeoModel(os.path.join(data_path(), "rectification", "right_image.geom")) current_coords = np.array([5000.5, 5000.5, 0.0], dtype=np.float64) mean_spacing = 1 diff --git a/tests/geomodels/test_los.py b/tests/geomodels/test_los.py index f105025..3d516dc 100644 --- a/tests/geomodels/test_los.py +++ b/tests/geomodels/test_los.py @@ -127,7 +127,7 @@ def test_los_parameters(): # Load geometrical model id_scene = "P1BP--2017092838284574CP" file_dimap = os.path.join(data_folder, f"rpc/RPC_{id_scene}.XML") - geometrical_model = RPC.from_any(file_dimap) + geometrical_model = GeoModel(file_dimap) # Create los model_los = LOS(matches_left, geometrical_model, [310, 850]) From 26bbdba6133b31ab36f91185c242a7080139b95f Mon Sep 17 00:00:00 2001 From: Jonathan Guinet Date: Fri, 1 Dec 2023 09:03:22 +0100 Subject: [PATCH 10/12] doc: update shpinx --- docs/source/getting_started.rst | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/docs/source/getting_started.rst b/docs/source/getting_started.rst index db7fe58..be9b464 100644 --- a/docs/source/getting_started.rst +++ b/docs/source/getting_started.rst @@ -40,7 +40,7 @@ Quick Start >>> # Create RPC object from downloaded geometry file >>> rpc_geom_file = "left_image.geom" - >>> rpc = GeoModel(rpc_geom_file, "rpc") + >>> rpc = GeoModel(rpc_geom_file, "rpc") # "rpc" is the geomodel type in ("rpc", "grid") with default value "rpc" >>> # Create Localization object from created RPC >>> loc = Localization(rpc) From 213f20b07a53d371cd5c3f04393ec829d49a7b88 Mon Sep 17 00:00:00 2001 From: Jonathan Guinet Date: Fri, 1 Dec 2023 09:03:49 +0100 Subject: [PATCH 11/12] refac: change member of grid --- shareloc/geofunctions/localization.py | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/shareloc/geofunctions/localization.py b/shareloc/geofunctions/localization.py index 7a714e1..1b5037e 100755 --- a/shareloc/geofunctions/localization.py +++ b/shareloc/geofunctions/localization.py @@ -136,9 +136,10 @@ def inverse(self, lon, lat, h=None, using_geotransform=False): :rtype: Tuple(1D np.ndarray row position, 1D np.ndarray col position, 1D np.ndarray alt) """ - if not self.use_rpc and not hasattr(self.model, "pred_ofset_scale_lon"): + if not self.use_rpc: # for grids only - self.model.estimate_inverse_loc_predictor() + if self.model.pred_ofset_scale_lon is None: + self.model.estimate_inverse_loc_predictor() if h is None: h = self.default_elevation From bf9df23f312877750a345a9197fe7ca41dd88a58 Mon Sep 17 00:00:00 2001 From: Jonathan Guinet Date: Fri, 1 Dec 2023 11:31:16 +0100 Subject: [PATCH 12/12] refac: moved espg member to geomodeltemplate --- shareloc/geomodels/geomodel_template.py | 2 ++ shareloc/geomodels/grid.py | 1 - shareloc/geomodels/rpc.py | 15 ++++++++++++--- 3 files changed, 14 insertions(+), 4 deletions(-) diff --git a/shareloc/geomodels/geomodel_template.py b/shareloc/geomodels/geomodel_template.py index 599b406..93e4e7c 100644 --- a/shareloc/geomodels/geomodel_template.py +++ b/shareloc/geomodels/geomodel_template.py @@ -42,6 +42,8 @@ def __init__(self): """ # geomodel type. Set by the subclass self.type: str + # geomodel epsg projection code + self.epsg: int = None # Define GeoModelTemplate functions interface diff --git a/shareloc/geomodels/grid.py b/shareloc/geomodels/grid.py index 2c513c8..3aa704f 100755 --- a/shareloc/geomodels/grid.py +++ b/shareloc/geomodels/grid.py @@ -103,7 +103,6 @@ def __init__(self, geomodel_path: str): self.alts_down = None self.rowmax = None self.colmax = None - self.epsg = 0 # inverse loc predictor attributes self.pred_col_min = None diff --git a/shareloc/geomodels/rpc.py b/shareloc/geomodels/rpc.py index d9c8373..b5085c1 100755 --- a/shareloc/geomodels/rpc.py +++ b/shareloc/geomodels/rpc.py @@ -40,8 +40,6 @@ # Set numba type of threading layer before parallel target compilation config.THREADING_LAYER = "omp" -# pylint: disable=no-member - @GeoModel.register("RPC") class RPC(GeoModelTemplate): @@ -53,7 +51,18 @@ class RPC(GeoModelTemplate): # pylint: disable=too-many-instance-attributes def __init__(self, rpc_params): super().__init__() - self.epsg = None + + self.offset_alt = None + self.scale_alt = None + self.offset_col = None + self.scale_col = None + self.offset_row = None + self.scale_row = None + self.offset_x = None + self.scale_x = None + self.offset_y = None + self.scale_y = None + self.datum = None for key, value in rpc_params.items(): setattr(self, key, value)