cars.core.projection ==================== .. py:module:: cars.core.projection .. autoapi-nested-parse:: Projection module: contains some general purpose functions using polygons and data projections Functions --------- .. autoapisummary:: cars.core.projection.compute_dem_intersection_with_poly cars.core.projection.polygon_projection cars.core.projection.polygon_projection_crs cars.core.projection.geo_to_ecef cars.core.projection.ecef_to_enu cars.core.projection.geo_to_enu cars.core.projection.enu_to_aer cars.core.projection.geo_to_aer cars.core.projection.point_cloud_conversion cars.core.projection.point_cloud_conversion_crs cars.core.projection.get_xyz_np_array_from_dataset cars.core.projection.get_converted_xy_np_arrays_from_dataset cars.core.projection.point_cloud_conversion_dataset cars.core.projection.point_cloud_conversion_dataframe cars.core.projection.ground_polygon_from_envelopes cars.core.projection.ground_intersection_envelopes cars.core.projection.get_ground_direction cars.core.projection.get_ground_angles cars.core.projection.get_output_crs cars.core.projection.guess_vcrs_from_file_name Module Contents --------------- .. py:function:: compute_dem_intersection_with_poly(srtm_file: str, ref_poly: shapely.geometry.Polygon, ref_epsg: int) -> shapely.geometry.Polygon Compute the intersection polygon between the defined dem regions and the reference polygon in input :raise Exception: when the input dem doesn't intersect the reference polygon :param srtm_file: srtm file :param ref_poly: reference polygon :param ref_epsg: reference epsg code :return: The intersection polygon between the defined dem regions and the reference polygon in input .. py:function:: polygon_projection(poly: shapely.geometry.Polygon, epsg_in: int, epsg_out: int) -> shapely.geometry.Polygon Projects a polygon from an initial epsg code to another :param poly: poly to project :param epsg_in: initial epsg code :param epsg_out: final epsg code :return: The polygon in the final projection .. py:function:: polygon_projection_crs(poly: shapely.geometry.Polygon, crs_in: pyproj.CRS, crs_out: pyproj.CRS) -> shapely.geometry.Polygon Projects a polygon from an initial crs to another :param poly: poly to project :param crs_in: initial crs :param crs_out: final crs :return: The polygon in the final projection .. py:function:: geo_to_ecef(lat: numpy.ndarray, lon: numpy.ndarray, alt: numpy.ndarray) -> Tuple[numpy.ndarray, numpy.ndarray, numpy.ndarray] Point transformation from Geodetic of ellipsoid WGS-84) to ECEF ECEF: Earth-centered, Earth-fixed :param lat: input geodetic latitude (angle in degree) :param lon: input geodetic longitude (angle in degree) :param alt: input altitude above geodetic ellipsoid (meters) :return: ECEF (Earth centered, Earth fixed) x, y, z coordinates tuple (in meters) .. py:function:: ecef_to_enu(x_ecef: numpy.ndarray, y_ecef: numpy.ndarray, z_ecef: numpy.ndarray, lat0: numpy.ndarray, lon0: numpy.ndarray, alt0: numpy.ndarray) -> Tuple[numpy.ndarray, numpy.ndarray, numpy.ndarray] Coordinates conversion from ECEF Earth Centered to East North Up Coordinate from a reference point (lat0, lon0, alt0) See Wikipedia page for details: https://en.wikipedia.org/wiki/Geographic_coordinate_conversion :param x_ecef: target x ECEF coordinate (meters) :param y_ecef: target y ECEF coordinate (meters) :param z_ecef: target z ECEF coordinate (meters) :param lat0: Reference geodetic latitude :param lon0: Reference geodetic longitude :param alt0: Reference altitude above geodetic ellipsoid (meters) :return: ENU (xEast, yNorth zUp) target coordinates tuple (meters) .. py:function:: geo_to_enu(lat: numpy.ndarray, lon: numpy.ndarray, alt: numpy.ndarray, lat0: numpy.ndarray, lon0: numpy.ndarray, alt0: numpy.ndarray) -> Tuple[numpy.ndarray, numpy.ndarray, numpy.ndarray] Point transformation from WGS-84 Geodetic coordinates to to ENU. Use geo_to_ecef and ecef_to_enu functions. :param lat: input geodetic latitude (angle in degree) :param lon: input geodetic longitude (angle in degree) :param alt: input altitude above geodetic ellipsoid (meters) :param lat0: Reference geodetic latitude :param lon0: Reference geodetic longitude :param alt0: Reference altitude above geodetic ellipsoid (meters) :return: ENU (xEast, yNorth zUp) target coordinates tuple (meters) .. py:function:: enu_to_aer(x_east: numpy.ndarray, y_north: numpy.ndarray, z_up: numpy.ndarray) -> Tuple[numpy.ndarray, numpy.ndarray, numpy.ndarray] ENU coordinates to Azimuth, Elevation angle, Range from ENU origin Beware: Elevation angle is not the altitude. :param x_east: ENU East coordinate (meters) :param y_north: ENU North coordinate (meters) :param z_up: ENU Up coordinate (meters) :return: Azimuth, Elevation Angle, Slant Range (degrees, degrees, meters) .. py:function:: geo_to_aer(lat: numpy.ndarray, lon: numpy.ndarray, alt: numpy.ndarray, lat0: numpy.ndarray, lon0: numpy.ndarray, alt0: numpy.ndarray) -> Tuple[numpy.ndarray, numpy.ndarray, numpy.ndarray] Gives Azimuth, Elevation angle and Slant Range from a Reference to a Point with geodetic coordinates. :param lat: input geodetic latitude (angle in degree) :param lon: input geodetic longitude (angle in degree) :param alt: input altitude above geodetic ellipsoid (meters) :param lat0: Reference geodetic latitude :param lon0: Reference geodetic longitude :param alt0: Reference altitude above geodetic ellipsoid (meters) :return: Azimuth, Elevation Angle, Slant Range (degrees, degrees, meters) .. py:function:: point_cloud_conversion(cloud_in: numpy.ndarray, epsg_in: int, epsg_out: int) -> numpy.ndarray Convert a point cloud from a SRS to another one. :param cloud_in: cloud to project :param epsg_in: EPSG code of the input SRS :param epsg_out: EPSG code of the output SRS :return: Projected point cloud .. py:function:: point_cloud_conversion_crs(cloud_in: numpy.ndarray, crs_in: int, crs_out: int) -> numpy.ndarray Convert a point cloud from a SRS to another one. :param cloud_in: cloud to project :param crs_in: crs of the input SRS :param crs_out: crs of the output SRS :return: Projected point cloud .. py:function:: get_xyz_np_array_from_dataset(cloud_in: xarray.Dataset) -> Tuple[numpy.array, List[int]] Get a numpy array of size (nb_points, 3) with the columns being the x, y and z coordinates from a dataset as given in output of the triangulation. The original epipolar geometry shape is also given in output in order to reshape the output numpy array in its original geometry if necessary. :param cloud_in: input xarray dataset :return: a tuple composed of the xyz numàpy array and its original shape .. py:function:: get_converted_xy_np_arrays_from_dataset(cloud_in: xarray.Dataset, epsg_out: int) -> Tuple[numpy.array, numpy.array] Get the x and y coordinates as numpy array in the new referential indicated by epsg_out. TODO: add test :param cloud_in: input xarray dataset :param epsg_out: target epsg code :return: a tuple composed of the x and y numpy arrays .. py:function:: point_cloud_conversion_dataset(cloud: xarray.Dataset, epsg_out: int) Convert a point cloud as an xarray.Dataset to another epsg (inplace) TODO: add test :param cloud: cloud to project :param epsg_out: EPSG code of the output SRS .. py:function:: point_cloud_conversion_dataframe(cloud: pandas.DataFrame, epsg_in: int, epsg_out: int) Convert a point cloud as a panda.DataFrame to another epsg (inplace) :param cloud: cloud to project :param epsg_in: EPSG code of the input SRS :param epsg_out: EPSG code of the output SRS .. py:function:: ground_polygon_from_envelopes(poly_envelope1: shapely.geometry.Polygon, poly_envelope2: shapely.geometry.Polygon, epsg1: int, epsg2: int, tgt_epsg: int = 4326) -> Tuple[shapely.geometry.Polygon, Tuple[int, int, int, int]] compute the ground polygon of the intersection of two envelopes :raise: Exception when the envelopes don't intersect one to each other :param poly_envelope1: path to the first envelope :param poly_envelope2: path to the second envelope :param epsg1: EPSG code of poly_envelope1 :param epsg2: EPSG code of poly_envelope2 :param tgt_epsg: EPSG code of the new projection (default value is set to 4326) :return: a tuple with the shapely polygon of the intersection and the intersection's bounding box (described by a tuple (minx, miny, maxx, maxy)) .. py:function:: ground_intersection_envelopes(sensor1, sensor2, geomodel1, geomodel2, geometry_plugin: str, envelope_1_path: str, envelope_2_path: str, out_intersect_path: str, envelope_file_driver: str = 'GeoJSON', intersect_file_driver: str = 'GPKG') -> Tuple[shapely.geometry.Polygon, Tuple[int, int, int, int]] Compute ground intersection of two images with envelopes: 1/ Create envelopes for left, right images 2/ Read vectors and polygons with adequate EPSG codes. 3/ compute the ground polygon of the intersection of two envelopes 4/ Write the GPKG vector Returns a shapely polygon and intersection bounding box :raise: Exception when the envelopes don't intersect one to each other :param sensor1: path to left sensor image :param sensor2: path to right sensor image :param geomodel1: path and attributes for left geomodel :param geomodel2: path and attributes for right geomodel :param geometry_plugin: name of geometry plugin to use :param envelope_1_path: Path to the output shapefile left :param envelope_2_path: Path to the output shapefile right :param dem_dir: Directory containing DEM tiles :param default_alt: Default altitude above ellipsoid :param out_intersect_path: out vector file path to create :param envelope_file_driver: driver used to write envelope files :param intersect_file_driver: driver used to write the intersection file :return: a tuple with the shapely polygon of the intersection and the intersection's bounding box (described by a tuple (minx, miny, maxx, maxy)) .. py:function:: get_ground_direction(sensor, geomodel, geometry_plugin: str, x_coord: float = None, y_coord: float = None, z0_coord: float = None, z_coord: float = None) -> numpy.ndarray For a given image (x,y) point, compute the direction vector to ground The function uses the direct localization operation and makes a z variation to get a ground direction vector. By default, (x,y) is put at image center and z0, z at RPC geometric model limits. :param TODO :param geometry_plugin: name of geometry plugin to use :param x_coord: X Coordinate in input image sensor :param y_coord: Y Coordinate in input image sensor :param z0_coord: Z altitude reference coordinate :param z_coord: Z Altitude coordinate to take the image :param dem: path to the dem directory :param geoid: path to the geoid file :return: (lat0,lon0,alt0, lat,lon,alt) origin and end vector coordinates .. py:function:: get_ground_angles(sensor1, sensor2, geomodel1, geomodel2, geometry_plugin, x1_coord: float = None, y1_coord: float = None, z1_0_coord: float = None, z1_coord: float = None, x2_coord: float = None, y2_coord: float = None, z2_0_coord: float = None, z2_coord: float = None) -> Tuple[numpy.ndarray, numpy.ndarray, numpy.ndarray, numpy.ndarray, numpy.ndarray] For a given image (x,y) point, compute the Azimuth angle, Elevation angle (not the altitude !) and Range from Ground z0 perspective for both stereo image (img1: left and img2: right) Calculate also the convergence angle between the two satellites from ground. The function use get_ground_direction function to have coordinates of ground direction vector and compute angles and range. Ref: Jeong, Jaehoon. (2017). IMAGING GEOMETRY AND POSITIONING ACCURACY OF DUAL SATELLITE STEREO IMAGES: A REVIEW. ISPRS Annals of Photogrammetry, Remote Sensing and Spatial Information Sciences. IV-2/W4. 235-242. 10.5194/isprs-annals-IV-2-W4-235-2017. Perspectives: get bisector elevation (BIE), and asymmetry angle :param TODO :param geometry_plugin: name of geometry plugin to use :param x1_coord: X Coordinate in input left image1 sensor :param y1_coord: Y Coordinate in input left image1 sensor :param z1_0_coord: Left image1 Z altitude origin coordinate for ground direction vector :param z1_coord: Left image1 Z altitude end coordinate for ground direction vector :param x2_coord: X Coordinate in input right image2 sensor :param y2_coord: Y Coordinate in input right image2 sensor :param z2_0_coord: Right image2 Z altitude origin coordinate for ground direction vector :param z2_coord: Right image2 Z altitude end coordinate for ground direction vector :return: Left Azimuth, Left Elevation Angle, Right Azimuth, Right Elevation Angle, Convergence Angle .. py:function:: get_output_crs(epsg, out_conf) Détermine le CRS de sortie en fonction de la config. .. py:function:: guess_vcrs_from_file_name(filepath) Tries to detect the geoid's EPSG from the file name