cars.core.projection

Projection module: contains some general purpose functions using polygons and data projections

Module Contents

Functions

compute_dem_intersection_with_poly(...)

Compute the intersection polygon between the defined dem regions

polygon_projection(→ shapely.geometry.Polygon)

Projects a polygon from an initial epsg code to another

geo_to_ecef(→ Tuple[numpy.ndarray, numpy.ndarray, ...)

Point transformation from Geodetic of ellipsoid WGS-84) to ECEF

ecef_to_enu(→ Tuple[numpy.ndarray, numpy.ndarray, ...)

Coordinates conversion from ECEF Earth Centered to

geo_to_enu(→ Tuple[numpy.ndarray, numpy.ndarray, ...)

Point transformation from WGS-84 Geodetic coordinates to to ENU.

enu_to_aer(→ Tuple[numpy.ndarray, numpy.ndarray, ...)

ENU coordinates to Azimuth, Elevation angle, Range from ENU origin

geo_to_aer(→ Tuple[numpy.ndarray, numpy.ndarray, ...)

Gives Azimuth, Elevation angle and Slant Range

points_cloud_conversion(→ numpy.ndarray)

Convert a point cloud from a SRS to another one.

get_xyz_np_array_from_dataset(→ Tuple[numpy.array, ...)

Get a numpy array of size (nb_points, 3) with the columns

get_converted_xy_np_arrays_from_dataset(...)

Get the x and y coordinates as numpy array

points_cloud_conversion_dataset(cloud, epsg_out)

Convert a point cloud as an xarray.Dataset to another epsg (inplace)

points_cloud_conversion_dataframe(cloud, epsg_in, epsg_out)

Convert a point cloud as a panda.DataFrame to another epsg (inplace)

ground_polygon_from_envelopes(...)

compute the ground polygon of the intersection of two envelopes

ground_intersection_envelopes(...)

Compute ground intersection of two images with envelopes:

project_coordinates_on_line(→ Union[float, numpy.ndarray])

Project coordinates (x,y) on a line starting from origin with a

get_ground_direction(→ numpy.ndarray)

For a given image (x,y) point, compute the direction vector to ground

get_ground_angles(→ Tuple[numpy.ndarray, ...)

For a given image (x,y) point, compute the Azimuth angle,

cars.core.projection.compute_dem_intersection_with_poly(srtm_dir: 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

Raises

Exception – when the input dem doesn’t intersect the reference polygon

Parameters
  • srtm_dir – srtm directory

  • ref_poly – reference polygon

  • ref_epsg – reference epsg code

Returns

The intersection polygon between the defined dem regions and the reference polygon in input

cars.core.projection.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

Parameters
  • poly – poly to project

  • epsg_in – initial epsg code

  • epsg_out – final epsg code

Returns

The polygon in the final projection

cars.core.projection.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

Parameters
  • lat – input geodetic latitude (angle in degree)

  • lon – input geodetic longitude (angle in degree)

  • alt – input altitude above geodetic ellipsoid (meters)

Returns

ECEF (Earth centered, Earth fixed) x, y, z coordinates tuple (in meters)

cars.core.projection.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

Parameters
  • x_ecef – target x ECEF coordinate (meters)

  • y_ecef – target y ECEF coordinate (meters)

  • z_ecef – target z ECEF coordinate (meters)

  • lat0 – Reference geodetic latitude

  • lon0 – Reference geodetic longitude

  • alt0 – Reference altitude above geodetic ellipsoid (meters)

Returns

ENU (xEast, yNorth zUp) target coordinates tuple (meters)

cars.core.projection.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.

Parameters
  • lat – input geodetic latitude (angle in degree)

  • lon – input geodetic longitude (angle in degree)

  • alt – input altitude above geodetic ellipsoid (meters)

  • lat0 – Reference geodetic latitude

  • lon0 – Reference geodetic longitude

  • alt0 – Reference altitude above geodetic ellipsoid (meters)

Returns

ENU (xEast, yNorth zUp) target coordinates tuple (meters)

cars.core.projection.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.

Parameters
  • x_east – ENU East coordinate (meters)

  • y_north – ENU North coordinate (meters)

  • z_up – ENU Up coordinate (meters)

Returns

Azimuth, Elevation Angle, Slant Range (degrees, degrees, meters)

cars.core.projection.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.

Parameters
  • lat – input geodetic latitude (angle in degree)

  • lon – input geodetic longitude (angle in degree)

  • alt – input altitude above geodetic ellipsoid (meters)

  • lat0 – Reference geodetic latitude

  • lon0 – Reference geodetic longitude

  • alt0 – Reference altitude above geodetic ellipsoid (meters)

Returns

Azimuth, Elevation Angle, Slant Range (degrees, degrees, meters)

cars.core.projection.points_cloud_conversion(cloud_in: numpy.ndarray, epsg_in: int, epsg_out: int) numpy.ndarray

Convert a point cloud from a SRS to another one.

Parameters
  • cloud_in – cloud to project

  • epsg_in – EPSG code of the input SRS

  • epsg_out – EPSG code of the output SRS

Returns

Projected point cloud

cars.core.projection.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.

Parameters

cloud_in – input xarray dataset

Returns

a tuple composed of the xyz numàpy array and its original shape

cars.core.projection.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

Parameters
  • cloud_in – input xarray dataset

  • epsg_out – target epsg code

Returns

a tuple composed of the x and y numpy arrays

cars.core.projection.points_cloud_conversion_dataset(cloud: xarray.Dataset, epsg_out: int)

Convert a point cloud as an xarray.Dataset to another epsg (inplace) TODO: add test

Parameters
  • cloud – cloud to project

  • epsg_out – EPSG code of the output SRS

cars.core.projection.points_cloud_conversion_dataframe(cloud: pandas.DataFrame, epsg_in: int, epsg_out: int)

Convert a point cloud as a panda.DataFrame to another epsg (inplace)

Parameters
  • cloud – cloud to project

  • epsg_in – EPSG code of the input SRS

  • epsg_out – EPSG code of the output SRS

cars.core.projection.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 TODO: refacto with externals (OTB) and steps.

Raise

Exception when the envelopes don’t intersect one to each other

Parameters
  • poly_envelope1 – path to the first envelope

  • poly_envelope2 – path to the second envelope

  • epsg1 – EPSG code of poly_envelope1

  • epsg2 – EPSG code of poly_envelope2

  • tgt_epsg – EPSG code of the new projection (default value is set to 4326)

Returns

a tuple with the shapely polygon of the intersection and the intersection’s bounding box (described by a tuple (minx, miny, maxx, maxy))

cars.core.projection.ground_intersection_envelopes(sensor1, sensor2, geomodel1, geomodel2, geometry_plugin: str, shp1_path: str, shp2_path: str, out_intersect_path: str) 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

Parameters
  • sensor1 – path to left sensor image

  • sensor2 – path to right sensor image

  • geomodel1 – path and attributes for left geomodel

  • geomodel2 – path and attributes for right geomodel

  • geometry_plugin – name of geometry plugin to use

  • shp1_path – Path to the output shapefile left

  • shp2_path – Path to the output shapefile right

  • dem_dir – Directory containing DEM tiles

  • default_alt – Default altitude above ellipsoid

  • out_intersect_path – out vector file path to create

Returns

a tuple with the shapely polygon of the intersection and the intersection’s bounding box (described by a tuple (minx, miny, maxx, maxy))

cars.core.projection.project_coordinates_on_line(x_coord: Union[float, numpy.ndarray], y_coord: Union[float, numpy.ndarray], origin: Union[List[float], numpy.ndarray], vec: Union[List[float], numpy.ndarray]) Union[float, numpy.ndarray]

Project coordinates (x,y) on a line starting from origin with a direction vector vec, and return the euclidean distances between projected points and origin.

Parameters
  • x_coord – scalar or vector of coordinates x

  • y_coord – scalar or vector of coordinates x

  • origin – coordinates of origin point for line

  • vec – direction vector of line

Returns

vector of distances of projected points to origin

cars.core.projection.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

cars.core.projection.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