Satellite–Derived Topography and Morphometry for VHR Coastal Habitat Mapping: The Pleiades–1 Tri–Stereo Enhancement

: The evolution of the coastal fringe is closely linked to the impact of climate change, speciﬁcally increases in sea level and storm intensity. The anthropic pressure that is inﬂicted on these fragile environments strengthens the risk. Therefore, numerous research projects look into the possibility of monitoring and understanding the coastal environment in order to better identify its dynamics and adaptation to the major changes that are currently taking place in the landscape. This new study aims to improve the habitat mapping/classiﬁcation at Very High Resolution (VHR) using Pleiades–1–derived topography, its morphometric by–products, and Pleiades–1–derived imageries. A tri–stereo dataset was acquired and processed by image pairing to obtain nine digital surface models (DSM) that were 0.50 m pixel size using the free software RSP (RPC Stereo Processor) and that were calibrated and validated with the 2018–LiDAR dataset that was available for the study area: the Emerald Coast in Brittany (France). Four morphometric predictors that were derived from the best of the nine generated DSMs were calculated via a freely available software (SAGA GIS): slope, aspect, topographic position index (TPI), and TPI–based landform classiﬁcation (TPILC). A maximum likelihood classiﬁcation of the area was calculated using nine classes: the salt marsh, dune, rock, urban, ﬁeld, forest, beach, road, and seawater classes. With an RMSE of 4 m, the DSM#2–3_1 (from images #2 and #3 with one ground control point) outperformed the other DSMs. The classiﬁcation results that were computed from the DSM#2–3_1 demonstrate the importance of the contribution of the morphometric predictors that were added to the reference Red–Green–Blue (RGB, 76.37% in overall accuracy, OA). The best combination of TPILC that was added to the RGB + DSM provided a gain of 13% in the OA, reaching 89.37%. These ﬁndings will help scientists and managers who are tasked with coastal risks at VHR.


Global Change
Coastal landscapes have faced significant changes over billions of years, and their evolution is concomitant with major climatic upheavals. In its sixth and most recent report, the Intergovernmental Panel on Climate Change (IPCC) indicates that climate change is occurring more rapidly than originally predicted, with unprecedented increases in sea levels, heat waves, and the faster melting of polar ice caps [1]. Currently, mankind is trying

Landuse/Landcover Observation Techniques
Observation techniques for tracking landscape changes are many and varied. At local-scale and very-high (VH) spatial resolution, unmanned aerial vehicles (UAV) are useful for the VH temporal resolution monitoring of coastal socio-ecosystems [3]. UAVs are cost-efficient and easily deployable for shoreline detection [4] and for the identification of seasonal variations in saltmarsh meadows [5]. They are, however, not well suited for monitoring areas at the landscape scale (several km 2 ) due to not only the restrictions that are imposed by legislation but also the technical limitations that are enforced by the number of flight times that are permitted by the battery capacity. In addition, in coastal areas, the meteorological and marine conditions require a maximum time of presence on the site due to the tides (±one hour after low water slack).
Manned aerial vehicles (MAV) serve as a robust alternative that leverages passive sensors with a basic red-green-blue (RGB) spectrum and sometimes the infrared spectrum [6] or an active light detection and ranging (LiDAR) sensor [7]. The complete solution makes it difficult to plan missions, and sensors such as LiDAR are rather expensive.

Spaceborne Acquisition and Stereoscopy
Yet, the analysis of the environment at the landscape scale is made possible by satellites which have a VH spatial, a multispectral, and even a hyperspectral resolution for the bestequipped satellites [8].
A spaceborne solution exists to obtain multispectral VH resolution images that are 0.50 m and 0.30 m pixels in size, which are provided by Pleiades-1 or WorldView-3 and 4, respectively. Satellite-based multispectral VH resolution mapping of the coastal fringe has been successfully performed in studies focusing on tropical [9] or temperate environments [10].
Since 2000, some remote sensing satellites that are specialized in stereo acquisition have been launched into orbit around the Earth, such as the Worldview-1, -2, -3, -4 constellations, GeoEye-1 and -2, the Pleiades-1 and, -1B constellations, and the 2021launched Pleiades Neo, whose images are not yet available for research at the time of this submission [11]. The operating principle of stereoscopy is to photograph an object or a landscape from two different angles in the same way as human vision is able to, with a specific overlap for determining the 3D information of the obtained images. Sometimes a third angle (nadir) of view can be available as a redundant observation to increase the accuracy when producing a digital elevation or surface model.
Satellite-based stereo topography has the capability of improving coastal mapping by improving the spectral discrimination of eco-geo-morphological objects [12]. However, when a tri-stereo product is used, do they augment/boost this coastal mapping? Do the morphometric parameters that are derived from the topography contribute to a better classification of coastal ecosystems than basic spectral information? This paper along with the experimental results that are presented in it will seek to answer these two questions.

The Study Site
The entire study site (76 km 2 terrestrial part) is located on the Emerald Coast in Brittany (France) along the Channel Sea (48.60 • N, 2.00 • W; Figure 1). It is characterized by a diversity of ecosystems that are shaped by the proximity of a megatidal sea and is one of the six areas with the highest tidal ranges in the world (up to 14 m) [13]. The Rance, a coastal river, ends its course in the bay of Saint Malo, dividing the area in two sub-sites. In terms of land cover, this study area is composed of temperate zone coastal vegetation, salt marshes, rocks, dunes, and fine sand beaches. The coastline is strongly indented, leaving multiple sandy beaches surrounded by rocky points and islets a little further offshore. Bays are also common on this coastal fringe and are featured by the presence of salt marsh meadows. In terms of urbanism, this Brittany coastal fringe is subject to strong anthropic pressures. The small fishing villages of the past have evolved into resort urban areas that now attract tourists in search of iodized air and marine landscapes. salt marshes, rocks, dunes, and fine sand beaches. The coastline is strongly indented, leaving multiple sandy beaches surrounded by rocky points and islets a little further offshore. Bays are also common on this coastal fringe and are featured by the presence of salt marsh meadows. In terms of urbanism, this Brittany coastal fringe is subject to strong anthropic pressures. The small fishing villages of the past have evolved into resort urban areas that now attract tourists in search of iodized air and marine landscapes. As with all coastal areas in the world, the Emerald Coast is not exempted from coastal risks. Its highly populated coastline increases the vulnerability of lowland populations.
Data collection is based on tri-stereo images from the Pleiades-1 satellite sensor (Table 1; Figure 2a). The satellite orbited over the study area on 28 November 2020 to collect three images at 11 h 26 min 14 s (UTC), then at 11 h 26 min 24 s (UTC), and finally at 11 h 26 min 32 s (UTC; Table 1). Each dataset of images contains panchromatic and multispectral images (R, G, B, NIR) that are 0.5 m and 2 m pixels in size, respectively (Table 1, Figure 2b,c). The images were delivered without initial geometric processing (primary level) and without radiometric processing.  As with all coastal areas in the world, the Emerald Coast is not exempted from coastal risks. Its highly populated coastline increases the vulnerability of lowland populations.
Data collection is based on tri-stereo images from the Pleiades-1 satellite sensor (Table 1; Figure 2a). The satellite orbited over the study area on 28 November 2020 to collect three images at 11 h 26 min 14 s (UTC), then at 11 h 26 min 24 s (UTC), and finally at 11 h 26 min 32 s (UTC; Table 1). Each dataset of images contains panchromatic and multi-spectral images (R, G, B, NIR) that are 0.5 m and 2 m pixels in size, respectively (Table 1, Figure 2b,c). The images were delivered without initial geometric processing (primary level) and without radiometric processing.

LiDAR Airborne Dataset
Thanks to multiple LiDAR data acquisition campaigns on the French territory by the French Navy's Hydrographic and Oceanographic Department (SHOM), a LiDAR point cloud was available for the study site. The land/sea continuum is guaranteed by the precision of the topo-bathymetric dataset (horizontal and vertical topographic accuracy of 0.20 m), which was acquired with a Leica HawkEye-3 sensor (Chiroptera + Deep channel).
The 2018 LiDAR point cloud was used to calibrate and validate the digital elevation model (DEM) using 36 validation points that were evenly distributed over the study area. The coordinates in XY (WGS84 UTM 30N) and Z (ellipsoidal height) were extracted for each point (Figure 3).

LiDAR Airborne Dataset
Thanks to multiple LiDAR data acquisition campaigns on the French territory by the French Navy's Hydrographic and Oceanographic Department (SHOM), a LiDAR point cloud was available for the study site. The land/sea continuum is guaranteed by the precision of the topo-bathymetric dataset (horizontal and vertical topographic accuracy of 0.20 m), which was acquired with a Leica HawkEye-3 sensor (Chiroptera + Deep channel).
The 2018 LiDAR point cloud was used to calibrate and validate the digital elevation model (DEM) using 36 validation points that were evenly distributed over the study area. The coordinates in XY (WGS84 UTM 30N) and Z (ellipsoidal height) were extracted for each point (Figure 3).

Coastal Landscape Classes
Nine classes that are representative of the coastal landscape were identified (Table 2, Figure 4): dune (white dune vegetation Ammophila arenaria), salt marsh (salt marsh vegetation composed mainly of Spartina, Salicornia, Suaeda, and Halimione portulacoides), rock, urban (building roof), forest (mix of deciduous and coniferous), field (cultivated and uncultivated), beach (wet and dry sand of grain of 0.06 to 2 mm), road (mainly asphalt), and seawater (shallow to deep salt water). A sub-study site was extracted for the classification tests (red rectangle in Figure 3). Remote Sens. 2022, 14, x FOR PEER REVIEW 5 of 19 Figure 3. 2018 LiDAR of Emerald Coast with evenly distributed calibration points (black circle) and validation points (black triangle); the red rectangle corresponds to the sub-study site.

Coastal Landscape Classes
Nine classes that are representative of the coastal landscape were identified (Table 2, Figure 4): dune (white dune vegetation Ammophila arenaria), salt marsh (salt marsh vegetation composed mainly of Spartina, Salicornia, Suaeda, and Halimione portulacoides), rock, urban (building roof), forest (mix of deciduous and coniferous), field (cultivated and uncultivated), beach (wet and dry sand of grain of 0.06 to 2 mm), road (mainly asphalt), and seawater (shallow to deep salt water). A sub-study site was extracted for the classification tests (red rectangle in Figure 3).

Coastal Landscape Classes
Nine classes that are representative of the coastal landscape were identified (Table 2, Figure 4): dune (white dune vegetation Ammophila arenaria), salt marsh (salt marsh vegetation composed mainly of Spartina, Salicornia, Suaeda, and Halimione portulacoides), rock, urban (building roof), forest (mix of deciduous and coniferous), field (cultivated and uncultivated), beach (wet and dry sand of grain of 0.06 to 2 mm), road (mainly asphalt), and seawater (shallow to deep salt water). A sub-study site was extracted for the classification tests (red rectangle in Figure 3).

Satellite-Derived Topography: Photogrammetry Reconstruction
Three couples of panchromatic images, which were characterized by three intersection angles, were processed using the opensource software RSP (RPC stereo processor, [15], Figure 5) from the tri-stereo images from Pleiades-1. RSP was used to process stereo satellite images so as to create point dense clouds and then to create a Digital Surface Model (DSM) using the photogrammetry technique. For each pair of images, RSP used the Rational Polynomial Coefficients (RPC) files, which were delivered with the Pleiades-1 images, which contain the geometric parameters for the same images to build the projection relationship between the 3D and 2D space. The additional Ground Control Points (GCP) increased the accuracy of the reconstruction horizontally and vertically.

Satellite-Derived Topography: Photogrammetry Reconstruction
Three couples of panchromatic images, which were characterized by three intersection angles, were processed using the opensource software RSP (RPC stereo processor, [15], Figure 5) from the tri-stereo images from Pleiades-1. RSP was used to process stereo satellite images so as to create point dense clouds and then to create a Digital Surface Model (DSM) using the photogrammetry technique. For each pair of images, RSP used the Rational Polynomial Coefficients (RPC) files, which were delivered with the Pleiades-1 images, which contain the geometric parameters for the same images to build the projection relationship between the 3D and 2D space. The additional Ground Control Points (GCP) increased the accuracy of the reconstruction horizontally and vertically.

Satellite-Derived Topography: Photogrammetry Reconstruction
Three couples of panchromatic images, which were characterized by three intersection angles, were processed using the opensource software RSP (RPC stereo processor, [15], Figure 5) from the tri-stereo images from Pleiades-1. RSP was used to process stereo satellite images so as to create point dense clouds and then to create a Digital Surface Model (DSM) using the photogrammetry technique. For each pair of images, RSP used the Rational Polynomial Coefficients (RPC) files, which were delivered with the Pleiades-1 images, which contain the geometric parameters for the same images to build the projection relationship between the 3D and 2D space. The additional Ground Control Points (GCP) increased the accuracy of the reconstruction horizontally and vertically.

Satellite-Derived Topography: Photogrammetry Reconstruction
Three couples of panchromatic images, which were characterized by three intersection angles, were processed using the opensource software RSP (RPC stereo processor, [15], Figure 5) from the tri-stereo images from Pleiades-1. RSP was used to process stereo satellite images so as to create point dense clouds and then to create a Digital Surface Model (DSM) using the photogrammetry technique. For each pair of images, RSP used the Rational Polynomial Coefficients (RPC) files, which were delivered with the Pleiades-1 images, which contain the geometric parameters for the same images to build the projection relationship between the 3D and 2D space. The additional Ground Control Points (GCP) increased the accuracy of the reconstruction horizontally and vertically.

Satellite-Derived Topography: Photogrammetry Reconstruction
Three couples of panchromatic images, which were characterized by three intersection angles, were processed using the opensource software RSP (RPC stereo processor, [15], Figure 5) from the tri-stereo images from Pleiades-1. RSP was used to process stereo satellite images so as to create point dense clouds and then to create a Digital Surface Model (DSM) using the photogrammetry technique. For each pair of images, RSP used the Rational Polynomial Coefficients (RPC) files, which were delivered with the Pleiades-1 images, which contain the geometric parameters for the same images to build the projection relationship between the 3D and 2D space. The additional Ground Control Points (GCP) increased the accuracy of the reconstruction horizontally and vertically.

Satellite-Derived Topography: Photogrammetry Reconstruction
Three couples of panchromatic images, which were characterized by three intersection angles, were processed using the opensource software RSP (RPC stereo processor, [15], Figure 5) from the tri-stereo images from Pleiades-1. RSP was used to process stereo satellite images so as to create point dense clouds and then to create a Digital Surface Model (DSM) using the photogrammetry technique. For each pair of images, RSP used the Rational Polynomial Coefficients (RPC) files, which were delivered with the Pleiades-1 images, which contain the geometric parameters for the same images to build the projection relationship between the 3D and 2D space. The additional Ground Control Points (GCP) increased the accuracy of the reconstruction horizontally and vertically.

Satellite-Derived Topography: Photogrammetry Reconstruction
Three couples of panchromatic images, which were characterized by three intersection angles, were processed using the opensource software RSP (RPC stereo processor, [15], Figure 5) from the tri-stereo images from Pleiades-1. RSP was used to process stereo satellite images so as to create point dense clouds and then to create a Digital Surface Model (DSM) using the photogrammetry technique. For each pair of images, RSP used the Rational Polynomial Coefficients (RPC) files, which were delivered with the Pleiades-1 images, which contain the geometric parameters for the same images to build the projection relationship between the 3D and 2D space. The additional Ground Control Points (GCP) increased the accuracy of the reconstruction horizontally and vertically.

Satellite-Derived Topography: Photogrammetry Reconstruction
Three couples of panchromatic images, which were characterized by three intersection angles, were processed using the opensource software RSP (RPC stereo processor, [15], Figure 5) from the tri-stereo images from Pleiades-1. RSP was used to process stereo satellite images so as to create point dense clouds and then to create a Digital Surface Model (DSM) using the photogrammetry technique. For each pair of images, RSP used the Rational Polynomial Coefficients (RPC) files, which were delivered with the Pleiades-1 images, which contain the geometric parameters for the same images to build the projection relationship between the 3D and 2D space. The additional Ground Control Points (GCP) increased the accuracy of the reconstruction horizontally and vertically.

Satellite-Derived Topography: Photogrammetry Reconstruction
Three couples of panchromatic images, which were characterized by three intersection angles, were processed using the opensource software RSP (RPC stereo processor, [15], Figure 5) from the tri-stereo images from Pleiades-1. RSP was used to process stereo satellite images so as to create point dense clouds and then to create a Digital Surface Model (DSM) using the photogrammetry technique. For each pair of images, RSP used the Rational Polynomial Coefficients (RPC) files, which were delivered with the Pleiades-1 images, which contain the geometric parameters for the same images to build the projection relationship between the 3D and 2D space. The additional Ground Control Points (GCP) increased the accuracy of the reconstruction horizontally and vertically.

Satellite-Derived Topography: Photogrammetry Reconstruction
Three couples of panchromatic images, which were characterized by three intersection angles, were processed using the opensource software RSP (RPC stereo processor, [15], Figure 5) from the tri-stereo images from Pleiades-1. RSP was used to process stereo satellite images so as to create point dense clouds and then to create a Digital Surface Model (DSM) using the photogrammetry technique. For each pair of images, RSP used the Rational Polynomial Coefficients (RPC) files, which were delivered with the Pleiades-1 images, which contain the geometric parameters for the same images to build the projection relationship between the 3D and 2D space. The additional Ground Control Points (GCP) increased the accuracy of the reconstruction horizontally and vertically.
where P is the value of the calculated DSM ellipsoidal height; O is the value of the LiDAR ellipsoidal height; and n is the number of observations. In addition, the DSMs were also evaluated at the ecosystem scale by referring to the nine identified classes. A mean of the ellipsoidal heights per class was calculated, and a dispersion parameter such as standard deviation was determined to show the homogeneity or heterogeneity of the statistical series.
The MS images underwent three pre-processing stages before they were integrated into the analyses: In addition, the DSMs were also evaluated at the ecosystem scale by referring to the nine identified classes. A mean of the ellipsoidal heights per class was calculated, and a dispersion parameter such as standard deviation was determined to show the homogeneity or heterogeneity of the statistical series.
The MS images underwent three pre-processing stages before they were integrated into the analyses: ENVI's FLAASH tool was used to correct the influence of the atmosphere (top of atmosphere, TOA) for each pixel in the images. The radiometric correction of the MS images consisted of converting the numerical radiance values into TOA reflectance values. The reflectance was calculated according to the radiance/irradiance (solar) ratio.
The MS images were also geometrically corrected using both the RPC files and the DSMs, which had been created previously from the panchromatic stereo images. This preprocessing stage corrected the distortions in the images that were related to the positioning of the satellite or the structure of the landform.
The MS images were pan-sharpened using the Gram-Schmidt algorithm. It resampled the 2 m pixel size of the MS image to be of the 0.5 m pixel size of the panchromatic images (see James et al., 2021).

Satellite-Derived Morphometry
The morphometric features such as slope, aspect, topographic position index (TPI), and TPI-based landform classification (TPILC) were calculated from the best DSMs using the SAGA GIS opensource software ( Figure 5) [16]. These indicators that were related to the topography of the study site were coupled with the basic RGB spectral information to reveal the contributions of each predictor. The near-infrared band of the image was also compared to the morphometric predictors.
The slope highlights the inclination of a pixel, and this aspect defines the orientation of the slope from a compass direction. The slope and aspect predictors generate raster images that are computed from the DSM. The percentage of slope is the ratio between the difference in altitude and the horizontal distance. A 3 × 3 pixel moving window compares the values of a pixel around its neighbors to define the slope percentage.
TPI computes the elevation or altitude of each pixel and subtracts it from the mean elevation or altitude of a neighborhood of that pixel of a grid raster [17]. Values that are lower than 0 correspond to valleys. Values higher than 0 are ridges, and those that are around 0 are flat areas. TPI-based landform classification was founded on the same principle as the basic TPI. Two different scales were combined to allow for the better identification of the topographic differences [18].

Classification Algorithm
A supervised machine learning classifier algorithm was tested: maximum likelihood (ML) with ENVI ® software ( Figure 5) [19]. ML is a probabilistic method that calculates the variance and covariance of each class by assuming that the statistics of each class in each band are normally distributed. A pixel is then assigned to the class with the most likely probability of membership.
An array of 500 calibration pixels and 500 validation pixels per class were extracted from the satellite-derived products ( Table 2). Overall accuracy (OA) is determined as the sum of the correctly classified pixels divided by the number of pixels. The producer accuracy (PA) corresponds to the accuracy of the map from the producer's point of view, i.e., from the algorithm. The result expressed in % indicates the fraction of correctly classified pixels of those that are known to belong to the class [20]. The OA of each by-product combination was evaluated through the calibration/validation pixels that were calculated with the confusion matrix.

Results
After the DSMs were derived and evaluated, the best one was further investigated to build topographic by-products. A ML algorithm was applied to this DSM using nine representative landscape classes from the study site. The contribution of each derived topographic band was evaluated at the landscape (OA) and class (PA) level.

Global Evaluation
The results of the overall DSM evaluation showed a slight increase in the overall accuracy of the DSMs with the addition of GCPs ( Figure 6). Thus, without GCP, with 1 GCP, and with 3 GCPs, the results were 4.03 m, 4 m, and 4.01 m for DSM#2-3 (Figures 6 and 7a). DSM#1-2 and DSM#1-3 increased their accuracy by 0.04 m and 0.02 m, respectively, as soon as a GCP was added (Figures 6 and 7b,c). However, image pairs #1-#2 and #1-#3 gave unsatisfactory results compared to image pair #2-#3.

Global Evaluation
The results of the overall DSM evaluation showed a slight increase in the overall accuracy of the DSMs with the addition of GCPs ( Figure 6). Thus, without GCP, with 1 GCP, and with 3 GCPs, the results were 4.03 m, 4 m, and 4.01 m for DSM#2-3 (Figures 6 and  7a). DSM#1-2 and DSM#1-3 increased their accuracy by 0.04 m and 0.02 m, respectively, as soon as a GCP was added (Figures 6 and 7b,c). However, image pairs #1-#2 and #1-#3 gave unsatisfactory results compared to image pair #2-#3. At the global scale, the analysis of the three best DSMs (DSM#2-3_1, DSM#1-2_3, and DSM#1-3_3) can be compared to the validation points that were extracted from the Li-DAR-2018 (Figure 7a-c). The distribution of the DSM#2-3_1 points shows a low dispersion and therefore a positive correlation between the LiDAR-2018 dataset and DSM#2-3_1 (Figure 7a). The DSM#1-3_3 appeared to be sparsely correlated to the LiDAR-2018 dataset.
(a) At the global scale, the analysis of the three best DSMs (DSM#2-3_1, DSM#1-2_3, and DSM#1-3_3) can be compared to the validation points that were extracted from the LiDAR-2018 (Figure 7a-c). The distribution of the DSM#2-3_1 points shows a low dispersion and therefore a positive correlation between the LiDAR-2018 dataset and DSM#2-3_1 (Figure 7a). The DSM#1-3_3 appeared to be sparsely correlated to the LiDAR-2018 dataset.

Class Level DSM Evaluation
The analysis of the DSM results from the Pleiades-1 images can also be examined at the class level (Figure 8) on the sub-study site. Based on the validation polygons from the classes, each DSM was evaluated.

Global Evaluation
The results of the overall DSM evaluation showed a slight increase in the overall accuracy of the DSMs with the addition of GCPs ( Figure 6). Thus, without GCP, with 1 GCP, and with 3 GCPs, the results were 4.03 m, 4 m, and 4.01 m for DSM#2-3 (Figures 6 and  7a). DSM#1-2 and DSM#1-3 increased their accuracy by 0.04 m and 0.02 m, respectively, as soon as a GCP was added (Figures 6 and 7b,c). However, image pairs #1-#2 and #1-#3 gave unsatisfactory results compared to image pair #2-#3.

Class Level DSM Evaluation
The analysis of the DSM results from the Pleiades-1 images can also be examined at the class level ( Figure 8) on the sub-study site. Based on the validation polygons from the classes, each DSM was evaluated.
Depending on the ecosystem, the mean heights between the DSMs within the same class can vary. This is the case for almost all of the classes if the three DSMs are compared to each other.
However, when DSM#2-3_1 and DSM#1-3_3 are compared to each other (due to their performance at the landscape scale), the mean differences in heights were: +10.

Morphometric Derivatives
Four main morphometric by-products were calculated: • The slope values ranged from 0 to 89 degrees, with 0° corresponding to a flat surface such as the seawater or flatland (in green) and 89° corresponding to a steep cliff (in Depending on the ecosystem, the mean heights between the DSMs within the same class can vary. This is the case for almost all of the classes if the three DSMs are compared to each other. However, when DSM#2-3_1 and DSM#1-3_3 are compared to each other (due to their performance at the landscape scale), the mean differences in heights were: +10.35 m, +10.22 m, +13.16 m, +15.24 m, +15.11 m, +13.95 m, +10.83 m, +7.71, and +6.15 m for the salt marsh, dune, rock, urban, forest, beach, road, and seawater classes, respectively.

Morphometric Derivatives
Four main morphometric by-products were calculated:

•
The slope values ranged from 0 to 89 • , with 0 • corresponding to a flat surface such as the seawater or flatland (in green) and 89 • corresponding to a steep cliff (in red in Figure 9a).

•
The aspect is categorized in 10 classes from 0 to 360 • , according to the main cardinal points (north, south, east, west; Figure 9b). A value of −1 corresponds to flat areas such as those for seawater. • TPI is the third morphometric contributor (Figure 9c).

Morphometric Derivatives
Four main morphometric by-products were calculated: • The slope values ranged from 0 to 89 degrees, with 0° corresponding to a flat surface such as the seawater or flatland (in green) and 89° corresponding to a steep cliff (in red in Figure 9a).

•
The aspect is categorized in 10 classes from 0 to 360 degrees, according to the main cardinal points (north, south, east, west; Figure 9b). A value of −1 corresponds to flat areas such as those for seawater. • TPI is the third morphometric contributor (Figure 9c).

Overall Accuracy at the Landscape Scale
At the landscape scale, each band combination was assessed. To test individual contribution, the RGB spectral composite was used as a basis, achieving an OA of 76.37%

Overall Accuracy at the Landscape Scale
At the landscape scale, each band combination was assessed. To test individual contribution, the RGB spectral composite was used as a basis, achieving an OA of 76.37% (Figure 10).
At the landscape scale, each band combination was assessed. To test individual contribution, the RGB spectral composite was used as a basis, achieving an OA of 76.37% ( Figure 10).
The NIR band was also tested and obtained a score of 84.34% combined with the RGB.
Thus, in addition to the RGB, the DSM significantly enhanced the OA by +12.51%, but RGB + DSM + slope decreased the score (−11.23%) compared to the RGB + DSM combination ( Figure 10).
The morphometric predictor combinations provided a strong contribution to the classification, with +9.05%, +11.45%, and +12.52% for the RGB + DSM + TPI and RGB + DSM + morphometric predictors and for RGB + DSM + aspect, respectively.
Finally, the best combination was the TPILC predictor combined with RGB + DSM, with a classification score of 89.37%, namely an augmentation of +13%.  The NIR band was also tested and obtained a score of 84.34% combined with the RGB. Thus, in addition to the RGB, the DSM significantly enhanced the OA by +12.51%, but RGB + DSM + slope decreased the score (−11.23%) compared to the RGB + DSM combination ( Figure 10).
The morphometric predictor combinations provided a strong contribution to the classification, with +9.05%, +11.45%, and +12.52% for the RGB + DSM + TPI and RGB + DSM + morphometric predictors and for RGB + DSM + aspect, respectively.
Finally, the best combination was the TPILC predictor combined with RGB + DSM, with a classification score of 89.37%, namely an augmentation of +13%.

Evaluation at the Class Level
The analysis of the confusion matrix of each morphometric predictor added to the basic RGB highlights a heterogeneity between the classes (Table 3, Figures 11 and 12). Table 3. Producer accuracy (in %) from confusion matrix using the maximum likelihood classifier for the salt marsh, dune, rock, urban, field, forest, beach, road, and seawater classes.

The Intersection Angle as a Key Determinant
At the global scale, nine DSM were computed via RSP from three Pleiades-1 satellite images captured at three different angles of incidence: 16.41° for image #1, 15.35° for image #2, and 16.05° for image #3 (Figure 2a and Table 1). DSM#2-3_1, which was derived from the stereo reconstruction of images #2 and #3 outperformed, the other DSMs. According to the results of the satellite photogrammetry, reconstructions with the closest intersection angles (5.13°) produced better point-measurement accuracy compared to the the 2018-LiDAR altimetric reference. Another answer can also be provided by focusing on the solar angle. [21]. Moreover, close or near-similar intersection angles increase the risk of "hidden sides" because of their proximity [22,23]. Thus, the urban class obtained the worst results with the RGB and increased in terms of classification performance when another predictor was added. The slope predictor increased by 2.46. The addition of a morphometric predictor allowed the 70.13% of classification performance to increase until the threshold of 84.13% as reached with the combination RGB + DSM + morphometric predictor. The NIR predictor holds up reasonably, achieving a score of with 71% (Table 3, Figures 11 and 12).
The trend for the forest class seems to be the same as the trend observed for the urban class. The slope provides a modest contribution of 2.47%, followed by the NIR predictor with a contribution of 16.67%, and then by the DSM contribution with a contribution of 32.27%. When added to the RGB + DSM, morphometric predictors aspect, TPILC, TPI, and the entire combination obtain values of +30.93%, +30.97%, +37.6%, and +40.2%, respectively.
The salt marsh, beach, and dune classes increased in classification accuracy when a morphometric variable was added to the reference RGB: +10.73%, +17.4%, and +36%, respectively with the RGB + DSM combination until the addition of TPI variable at 98.6% and 99.07% for the salt marsh and dune classes. The beach class obtained the best result with RGB + NIR, achieving results of 99.4%. In contrast to the salt marsh and beach classes, the rock class performed worst, with 13.27% and 26.33% with the RGB + TPI combinations and the complete combination of morphometric variables added to the RGB, respectively.
The road class obtained linear classification scores between 91.47% for RGB + DSM + TPI and 94.2% for the RGB + DSM + morphometric predictors.
The field class follows a slightly different pattern than the other classes since the worst result, which was still very acceptable, is achieved by the combination RGB + DSM + slope, with a score of 86.47%. The RGB combination increased the classification performance by only +0.2% and +0.93% for the RGB + DSM + aspect and +1.13% for the RGB + DSM combination. The RGB + DSM + TPILC, RGB + NIR, and RGB + DSM + morphometric predictors increased the classification performance by +2.2%, 2.8%, and +5.46%. The best combination for the field class was obtained by RGB + DSM + TPI, with a score of 93.6%.
The seawater class achieved high scores of almost 100% PA, regardless of the predictor ( Table 3, Figures 11 and 12).  (Figure 2a and Table 1). DSM#2-3_1, which was derived from the stereo reconstruction of images #2 and #3 outperformed, the other DSMs. According to the results of the satellite photogrammetry, reconstructions with the closest intersection angles (5.13 • ) produced better point-measurement accuracy compared to the the 2018-LiDAR altimetric reference. Another answer can also be provided by focusing on the solar angle [21]. Moreover, close or near-similar intersection angles increase the risk of "hidden sides" because of their proximity [22,23].

Discussion
However, many studies have shown interest in using tri-stereo satellite images to benefit, when possible, from a nadir view [24]. The benefit of such a tri-stereo enables a reduction in the shadows that are created by trees or buildings. This approach is highly valued by urban planners, as it limits the risk of shadows and hidden areas [25].
As for photogrammetric reconstructions from UAV images, a deficient RMSE could be explained by reconstruction artifacts related to the algorithm or to the images themselves [26].

Ground Control Point Effect
At a more local scale, three scenarios can explain the results that were obtained. The first is the number of GCP that were extracted from the 2018-LiDAR as calibration points for our topographic models. The tests that were performed without GCP, with one, or with three GCP, showed that when we added GCP points, the model is more accurate [27]. However, The RPC files that were delivered with the Pleiades-1 images could be more powerful. The RPC files integrate all of the parameters that are related to the photographs in order to correct for the satellite images. This approach could be investigated in a future study.

Information Reflected by the LiDAR Wavelengths
Two laser sensors were used for the whole LiDAR dataset: a mixed topo-bathymetric laser in the NIR and green spectra, and a stronger bathymetric laser in the green spectrum, which was specifically used for deeper areas. Depending on the nature (albedo) of the coastal habitat, the wavelength that is used by the LiDAR does not return the same information, being more or less reflected. Thus, the reflectance in the NIR is strong for eco-geo-systems with high chlorophyll dominance, and conversely, it is fully absorbed by habitats with a high water content, such as seawater, for example [28].
The density of the LiDAR points may also play a crucial role. The intertidal and coastal land area is denser in a number of points per square meter than marine area (deep water area) is. The detection of features such as trees is easy when the density of measured points is high [29].

Topographic Contribution to Habitat Classification
At the landscape scale, the contribution of the morphometric predictors surpassed the classification performance. Although the basic RGB performs well, as soon as morphometric predictors are added, the OA increases (from +1.28% to +13%). The positive effect of the slope raster and its by-products on the landscape classification were evidenced in several studies based on satellite sensors [30] as well as in UAV studies [31]. The classification test with the addition of the NIR band showed the relevance of a classification with morphometric variables. Indeed, when these variables were added to the RGB, the classification performance on the coastal fringe increased considerably contrary to the multispectral RGB + NIR combination without any other predictors.
At the ecosystem class level, the morphometric predictors produced different results depending on the coastal habitat type. For example, the topographic variables provided successful results for salt marshes. This particular ecosystem is located in a geographical area that is sheltered from strong prevailing marine currents. Coastal erosion has minimal or no impact on the meadows, allowing the different plant species to grow.

Conclusions
This research study on satellite photogrammetry with Pleiades-1 tri-stereo images is meaningful for the classification of coastal landscapes at VHR, especially in the context of climate change and increasing anthropic pressure on the coastal fringe.
Three pairs of Pleiades-1 panchromatic images at the 0.5 m pixel size were tested for DSM generation, and nine DSM were evaluated from 36 2018 LiDAR validation points.
The best DSM was derived from images #2 and #3 (DSM#2-3_1), which featured, respectively, with incidence angles of 15.35 • and 16.05 • and an intersection angle of 5.13 • . From this new DSM, four morphometric by-products were calculated: slope, aspect, topographic position index (TPI), and TPI-based landform classification (TPILC).
A pixel-based classifier, the probabilistic maximum likelihood, was applied to the 0.5 m pansharpened RGB images, which initially had a pixel side of 2 m. Nine classes (dune, salt marsh, rock, urban, field, forest, beach, road, and seawater) were examined to map the study site ( Figure 13). The best combination of morphometric predictors provided a gain of 13% in the OA, reaching 89.37%, when added to the RGB + DSM. These findings will help scientists and managers who tasked with the coastal risks at VHR. The best DSM was derived from images #2 and #3 (DSM#2-3_1), which featured, respectively, with incidence angles of 15.35° and 16.05° and an intersection angle of 5.13°. From this new DSM, four morphometric by-products were calculated: slope, aspect, topographic position index (TPI), and TPI-based landform classification (TPILC).
A pixel-based classifier, the probabilistic maximum likelihood, was applied to the 0.5 m pansharpened RGB images, which initially had a pixel side of 2 m. Nine classes (dune, salt marsh, rock, urban, field, forest, beach, road, and seawater) were examined to map the study site ( Figure 13). The best combination of morphometric predictors provided a gain of 13% in the OA, reaching 89.37%, when added to the RGB + DSM. These findings will help scientists and managers who tasked with the coastal risks at VHR.
In a study that will be published in the near future, photogrammetric reconstruction from a Pleiades-Neo panchromatic triplet at the 0.30 m spatial resolution will be evaluated, given that this new sensor benefits from six bands (purple, blue, green, red, red edge, NIR), meaning that it has two additional bands compared to Pleiades-1. Moreover, the purple band will be interesting to investigate for the sake of bathymetry extraction.  In a study that will be published in the near future, photogrammetric reconstruction from a Pleiades-Neo panchromatic triplet at the 0.30 m spatial resolution will be evaluated, given that this new sensor benefits from six bands (purple, blue, green, red, red edge, NIR), meaning that it has two additional bands compared to Pleiades-1. Moreover, the purple band will be interesting to investigate for the sake of bathymetry extraction.