shareloc.geomodels.rpc ====================== .. py:module:: shareloc.geomodels.rpc .. autoapi-nested-parse:: 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. Classes ------- .. autoapisummary:: shareloc.geomodels.rpc.RPC Functions --------- .. autoapisummary:: shareloc.geomodels.rpc.polynomial_equation shareloc.geomodels.rpc.compute_rational_function_polynomial shareloc.geomodels.rpc.derivative_polynomial_latitude shareloc.geomodels.rpc.derivative_polynomial_longitude shareloc.geomodels.rpc.compute_loc_inverse_derivates_numba Module Contents --------------- .. py:class:: RPC(rpc_params) Bases: :py:obj:`shareloc.geomodels.geomodel_template.GeoModelTemplate` RPC class including direct and inverse localization instance methods .. py:attribute:: offset_alt :value: None .. py:attribute:: scale_alt :value: None .. py:attribute:: offset_col :value: None .. py:attribute:: scale_col :value: None .. py:attribute:: offset_row :value: None .. py:attribute:: scale_row :value: None .. py:attribute:: offset_x :value: None .. py:attribute:: scale_x :value: None .. py:attribute:: offset_y :value: None .. py:attribute:: scale_y :value: None .. py:attribute:: datum :value: None .. py:attribute:: type :value: 'RPC' .. py:attribute:: lim_extrapol :value: 1.0001 .. py:attribute:: monomes .. py:attribute:: monomes_deriv_1 .. py:attribute:: monomes_deriv_2 .. py:attribute:: inverse_coefficient :value: False .. py:attribute:: direct_coefficient :value: False .. py:attribute:: alt_minmax .. py:attribute:: col0 .. py:attribute:: colmax .. py:attribute:: row0 .. py:attribute:: rowmax .. py:method:: load(geomodel_path) :classmethod: 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) .. py:method:: direct_loc_h(row, col, alt, fill_nan=False, using_direct_coef=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 :param using_direct_coef: equals True if you want to use direct coefficients :type using_direct_coef: boolean :return: ground position (lon,lat,h) :rtype: numpy.ndarray 2D dimension with (N,3) shape, where N is number of input coordinates .. py:method:: direct_loc_grid_h(row0, col0, steprow, stepcol, nbrow, nbcol, alt) calculates a direct loc grid (lat, lon) from the direct RPCs at constant altitude TODO: not tested. :param row0: grid origin (row) :type row0: int :param col0: grid origin (col) :type col0: int :param steprow: grid step (row) :type steprow: int :param stepcol: grid step (col) :type stepcol: int :param nbrow: grid nb row :type nbrow: int :param nbcol: grid nb col :type nbcol: int :param alt: altitude of the grid :type alt: float :return: direct localization grid longitude and latitude :rtype: Tuple(numpy.array, numpy.array) .. py:method:: direct_loc_dtm(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 .. py:method:: inverse_loc(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) .. py:method:: filter_coordinates(first_coord, second_coord, fill_nan=False, direction='direct') Filter nan input values :param first_coord: first coordinate :type first_coord: 1D numpy.ndarray dtype=float64 :param second_coord: second coordinate :type second_coord: 1D numpy.ndarray dtype=float64 :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 :param direction: direct or inverse localisation :type direction: str in ('direct', 'inverse') :return: filtered coordinates :rtype: list of numpy.array (index of nan, first filtered, second filtered) .. py:method:: compute_loc_inverse_derivates(lon, lat, alt) Inverse loc partial derivatives analytical compute :param lon: longitude coordinate :param lat: latitude coordinate :param alt: altitude coordinate :return: partials derivatives of inverse localization :rtype: Tuple(dcol_dlon np.array, dcol_dlat np.array, drow_dlon np.array, drow_dlat np.array) .. py:method:: direct_loc_inverse_iterative(row, col, alt, nb_iter_max=10, fill_nan=False) Iterative direct localization using inverse RPC :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 :type alt: float :param nb_iter_max: max number of iteration :type alt: int :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: list of numpy.array .. py:method:: get_alt_min_max() returns altitudes min and max layers :return: alt_min,lat_max :rtype: list .. py:method:: get_dtm_alt_offset(corners, dtm) returns min/max altitude offset between dtm coordinates system and RPC one :param corners: corners of the DTM's footprint :type corners: np.ndarray (4x2) :param dtm: DTM to get alt offset from :type dtm: DTMIntersection or bindings_cpp.DTMIntersection :return: min/max altimetric difference between RPC's epsg minus dtm alti expressed in dtm epsg :rtype: list of float (1x2) .. py:method:: los_extrema(row, col, alt_min=None, alt_max=None, fill_nan=False, epsg=None) compute los extrema :param row: line sensor position :type row: float :param col: column sensor position :type col: float :param alt_min: los alt min :type alt_min: float :param alt_max: los alt max :type alt_max: float :param epsg: epsg code of the dtm :type epsg: int :return: los extrema :rtype: numpy.array (2x3) .. py:function:: polynomial_equation(xnorm, ynorm, znorm, coeff) Compute polynomial equation :param xnorm: Normalized longitude (for inverse) or column (for direct) position :type xnorm: float 64 :param ynorm: Normalized latitude (for inverse) or line (for direct) position :type ynorm: float 64 :param znorm: Normalized altitude position :type znorm: float 64 :param coeff: coefficients :type coeff: 1D np.array dtype np.float 64 :return: rational :rtype: float 64 .. py:function:: compute_rational_function_polynomial(lon_col_norm, lat_row_norm, alt_norm, num_col, den_col, num_lin, den_lin, scale_col, offset_col, scale_lin, offset_lin) Compute rational function polynomial using numba to reduce calculation time on multiple points. useful to compute direct and inverse localization using direct or inverse RPC. :param lon_col_norm: Normalized longitude (for inverse) or column (for direct) position :type lon_col_norm: 1D np.array dtype np.float 64 :param lat_row_norm: Normalized latitude (for inverse) or line (for direct) position :type lat_row_norm: 1D np.array dtype np.float 64 :param alt_norm: Normalized altitude position :type alt_norm: 1D np.array dtype np.float 64 :param num_col: Column numerator coefficients :type num_col: 1D np.array dtype np.float 64 :param den_col: Column denominator coefficients :type den_col: 1D np.array dtype np.float 64 :param num_lin: Line numerator coefficients :type num_lin: 1D np.array dtype np.float 64 :param den_lin: Line denominator coefficients :type den_lin: 1D np.array dtype np.float 64 :param scale_col: Column scale :type scale_col: float 64 :param offset_col: Column offset :type offset_col: float 64 :param scale_lin: Line scale :type scale_lin: float 64 :param offset_lin: Line offset :type offset_lin: float 64 :return: for inverse localization : sensor position (row, col). for direct localization : ground position (lon, lat) :rtype: Tuple(np.ndarray, np.ndarray) .. py:function:: derivative_polynomial_latitude(lon_norm, lat_norm, alt_norm, coeff) Compute latitude derivative polynomial equation :param lon_norm: Normalized longitude position :type lon_norm: float 64 :param lat_norm: Normalized latitude position :type lat_norm: float 64 :param alt_norm: Normalized altitude position :type alt_norm: float 64 :param coeff: coefficients :type coeff: 1D np.array dtype np.float 64 :return: rational derivative :rtype: float 64 .. py:function:: derivative_polynomial_longitude(lon_norm, lat_norm, alt_norm, coeff) Compute longitude derivative polynomial equation :param lon_norm: Normalized longitude position :type lon_norm: float 64 :param lat_norm: Normalized latitude position :type lat_norm: float 64 :param alt_norm: Normalized altitude position :type alt_norm: float 64 :param coeff: coefficients :type coeff: 1D np.array dtype np.float 64 :return: rational derivative :rtype: float 64 .. py:function:: compute_loc_inverse_derivates_numba(lon_norm, lat_norm, alt_norm, num_col, den_col, num_lin, den_lin, scale_col, scale_lon, scale_lin, scale_lat) Analytically compute the partials derivatives of inverse localization using numba to reduce calculation time on multiple points :param lon_norm: Normalized longitude position :type lon_norm: 1D np.array dtype np.float 64 :param lat_norm: Normalized latitude position :type lat_norm: 1D np.array dtype np.float 64 :param alt_norm: Normalized altitude position :type alt_norm: 1D np.array dtype np.float 64 :param num_col: Column numerator coefficients :type num_col: 1D np.array dtype np.float 64 :param den_col: Column denominator coefficients :type den_col: 1D np.array dtype np.float 64 :param num_lin: Line numerator coefficients :type num_lin: 1D np.array dtype np.float 64 :param den_lin: Line denominator coefficients :type den_lin: 1D np.array dtype np.float 64 :param scale_col: Column scale :type scale_col: float 64 :param scale_lon: Geodetic longitude scale :type scale_lon: float 64 :param scale_lin: Line scale :type scale_lin: float 64 :param scale_lat: Geodetic latitude scale :type scale_lat: float 64 :return: partials derivatives of inverse localization :rtype: Tuples(dcol_dlon np.array, dcol_dlat np.array, drow_dlon np.array, drow_dlat np.array)