CLASSIFICATION OF WATER SURFACES USING AIRBORNE TOPOGRAPHIC LIDAR DATA

Accurate Digital Terrain Models (DTM) are inevitable inputs for mapping areas subject to natural hazards. Topographic airborne laser scanning has become an established technique to characterize the Earth surface: lidar provides 3D point clouds allowing a fine reconstruction of the topography. For flood hazard modeling, the key step before terrain modeling is the discrimination of land and water surfaces within the delivered point clouds. Therefore, instantaneous shoreline, river borders, inland waters can be extracted as a basis for more reliable DTM generation. This paper presents an automatic, efficient, and versatile workflow for land/water classification of airborne topographic lidar data. For that purpose, a classification framework based on Support Vector Machines (SVM) is designed. First, a restricted set of features, based only 3D lidar point coordinates and flightline information, is defined. Then, the SVM learning step is performed on small but well-targeted areas thanks to an automatic region growing strategy. Finally, label probabilities given by the SVM are merged during a probabilistic relaxation step in order to remove pixel-wise misclassification. Results show that survey of millions of points are labelled with high accuracy (>95% in most cases for coastal areas, and >89% for rivers) and that small natural and anthropic features of interest are still well classified though we work at low point densities (0.5-4 pts/m). Our approach is valid for coasts and rivers, and provides a strong basis for further discrimination of land-cover classes and coastal habitats.


INTRODUCTION
1.1 Motivation for seashore and river monitoring Climate change and global warming should lead to an increasing number of severe storms, more significant winter rains, and sea level rise in the forthcoming years.Coastal and river areas are particularly at risk, but their physical characteristics are barely described, especially their accurate topography, yet essential data for forecasting and management purposes.Another need indeed arises from statutory provisions: the European Water and Floods Framework Directives (2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007) have influenced the strategies for establishing prevention and protection policies by imposing repeated area-wide surveying of all kinds of inland waters, rivers, and high-staked catchment basins.The characterization and quantification of coastal and river habitats have been improving over the last decades due to synergistic remote sensing techniques, that are able to deliver high-resolution spatio-temporal by-products (Yang, 2008).In addition to providing some initial maps, remote sensing is also essential to monitor and analyze the evolution of the measured physical characteristics.Repetitive and up-to-date measurements are also crucial for areas undergoing most changes that are flooding, erosion, accretion or retreating such as beaches, cliffs or unstable slopes (Miller et al., 2008;Addo et al., 2008).For this purpose, the small-footprint airborne lidar technology appears to be attractive because it provides fine-scale Digital Terrain Models (DTM) over large coverage.It allows to survey hundreds of kilometers of shoreline and rivers with a high spatial resolution within a few days only.Its very high vertical accuracy (<0.15 m) has opened up new possibilities of tackling very precise and specific problems, that were impossible to deal with before (Hladik and Alber, 2012;Collin et al., 2012).In France, in addition to European Directives, the National Institute of Geographic and Forest Information (IGN) and the Marine Hydrographic and Oceanographic Service of the Defence Ministry have initiate a national program, started in 2005, that aims to create a three-dimensional model of the French coastline, called Litto3D R (Pastol, 2011).Lidar is currently used to produce a continuous land-sea representation of the coast, allowing for accurate manual shoreline mapping.Generating DTM in coastal and river areas first requires the classification of water areas.Consequently, the aim of this paper is to propose a workflow for water/land classification in topographic lidar datasets, that is adapted to various landscapes.

Related works on water detection from lidar data
We focus on airborne topographic lidar data, that operates on the near-infrared channel (NIR).There is an abundant literature on shoreline extraction from Digital Terrain Models but only few papers exists about direct water detection.Methods based on lidar point clouds are divided into two main approaches.First, segmentation of water surfaces can be based on the recognition of specific patterns (breaklines, tidal channels or gullies).For instance, the approach developed in (Brzank et al., 2005) is based on the detection of vertical high frequencies of the relief: planimetric lidar point distribution should be sufficiently regular and water areas should be located in trenches.Such assumption is efficient for many habitats but cannot be generalized.Conversely, Höfle et al. (2009) look for atypically large triangles in a Triangular Irregular Network, which are a strong hint about the presence of water.Indeed, the reflection properties of water surfaces for NIR lidar beams are characterized by either significant absorption or specular reflection resulting in numerous non-recorded laser echoes.Secondly, more general classifiers can be adopted.Since merely height information is insufficient, additional information is inserted, namely lidar intensity and point density.Brzank et al. (2008) used both features in a fuzzy logic classifier.The method performs well for water/mudflat discrimination for various lidar point densities but necessitates a preliminary object based analysis.The authors also noted that calibration/correction steps should be carried out to normalize the intensity feature.Even when such steps are performed, with a rule-based approach, Höfle et International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XL-1/W1, ISPRS Hannover Workshop 2013, 21 -24 May 2013, Hannover, Germany al. (2009) note a poor discrimination between asphalt and water surfaces.Their method requires the knowledge on the position of the sensor, GPS timestamps and scan angle to model lidar drop-outs.To cope with poor intensity values, the RGB channels of an orthoimage are inserted in a Mean-Shift classifier (Lee et al., 2012).The approach is restricted since the timing of optical data acquisition is critical and difficult to achieve over large scales.More advanced lidar geometrical features are proposed in (Schmidt et al., 2011), coupled with full-waveform attributes.
The authors adopted Conditional Random Fields in order to introduce contextual knowledge, therefore improve classification, and discriminate new classes (mudflat and mussel bed) (Schmidt et al., 2012).Nevertheless, the approach is not conceivable for large scale mapping since graphical models require significant training and inference steps.Consequently, despite relevant existing works, none of these approaches can be adopted to deal with both the large-scale issue and the automatic adaptiveness to various landscapes.A three-step strategy is adopted, and is explained in Section 2. Section 3 describes the various areas of interest that have been processed.Then, experimental results are presented and commented in Sections 4 and 5. Finally, conclusions are drawn in Section 6.

Overall strategy
As we target an approach suitable both for coasts and rivers, a supervised classifier is adopted.The core of the workflow is designed in a 2D-raster mode in order to provide a fast classification, with a tailored learning step.Furthermore, the approach is fully automatic, parameter-free, and versatile: no multi-echo, intensity, full-waveform or multi-spectral information is necessary, making our work easy to reproduce.We note that our system requires the setting of parameters at the different steps.The parameter values provided in the following are set for all our experiments, i.e., they are not tuned on a per-dataset basis.The processing chain can be divided into three steps.
1. Computation of features of interest.Six attributes are designed for discriminating water areas.The 3D lidar attributes are then interpolated in a 2D regular grid for large scale processing (Section 2.2).
2. Learning procedure.The knowledge of the land/water interface is used to define the training pixels the most likely to belong to each class.However, as this delineation may not be up-to-date or accurate enough, more representative training pixels are automatically selected using spatial reasoning (Section 2.3).
3. Land/water classification.A supervised Support Vector Machines (SVM) classifier is run on the pixels of the raster grid.Results were finally regularized using probabilistic relaxation, taking into account the probabilities given by the SVM classifier (Section 2.4).
SVM was selected because they are adapted to deal with highdimensional spaces and have shown considerable potential in the supervised classification of remotely sensed data, with very limited training amount.In addition, they outperform standard classification methods, as demonstrated by several general or applicationdriven comparative studies (Mountrakis et al., 2011), and in particular with airborne lidar data (Braun et al., 2011).

Discriminative lidar features and 2D interpolation
We have retained three families of features.They are based on (i) the height, (ii) the local point density, and (iii) the local shape of the 3D point neighborhood.To deal with any kind of lidar point cloud, we have consciously discarded the intensity/amplitude echobased features, i.e., information about the position of the 3D point within the current emitted lidar pulse.As a result, our method applies only to the 3D point coordinates.
2.2.1 Feature computation and 2D interpolation.The size of a 3D point neighborhood is the single parameter to be tuned in the feature computation stage.Given d the point cloud density and n the minimal number of points needed to calculate robust 3D descriptors, the radius r is defined by r = (n/πd) 1/2 .In practice, 10 points are sufficient.The mean density was between 1-4 pts/m 2 for areas with equivalent surfaces of land and water.Thus, r Furthermore, the interpolation of the 3D point cloud features on a regular 2D grid allows to better handle the large data volume and high dimensionality of the raw point clouds.The main areas of interest are the land-water interface, and, in particular, shallow waters.The larger the grid cell, the more inaccurate the classification.Therefore, the point cloud is resampled at a resolution of 1 m.This is sufficient for retrieving the land/water interface with a horizontal accuracy inferior to 2 m.The mean value of each feature is selected as final descriptor.

Feature description.
A first obvious attribute is the height of the lidar point with respect to the geoid, HG.Nevertheless, it does not always help to discriminate water points in inland areas.However, near-infrared pulse penetrates very little in water volumes, and therefore the 2D distribution of the lidar points can be used to detect such surfaces.For wide laser pulse scan angles (>15 • ), part of the emitted radiation returning from water varies significantly and may not be distinguished from the background noise.Therefore, the point density is much lower in the tails of the swath comparing to inland areas.Several 2D density-based features can be derived.First, the local mean point density (D) is computed for each point, independently from the lidar strip they belong to.Several neighborhoods at different sizes can be computed, enhancing more or less land areas but also water areas at the nadir of the aircraft.However, this feature is preferable for single lidar strip analysis.The neighborhood increase allows to deal with overlap sensitivity (overlapped lidar strips locally increases point density).A 5×5 m offers a suitable trade-off.Nevertheless, this leads to the loss of small structures.To deal with such issue, the majority density (Dm) is computed: only lidar points from the strip with the highest number of lidar hits are taken into account.In order to avoid the sensitivity to the strip overlap, we define: Additionally, the density ratio (Dr) is introduced to favor areas acquired with several strips but not equally sampled, which is not the case of water surfaces.Dr is computed as follows: where Dmax and Dmin are the point densities corresponding, for each cell of the grid, to the strips with the highest and lowest number of points, respectively.Such a value is close to 0 on water areas.Although noisy values can be found at water strip borders, it is particularly suited for our problem.Finally, the 3D distribution of lidar points can also be evaluated International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XL-1/W1, ISPRS Hannover Workshop 2013, 21 -24 May 2013, Hannover, Germany through the computation of eigenvalue features (Demantké et al., 2011).A covariance matrix of the 3D coordinates is computed in a cylindrical vertical neighborhood of fixed radius (Section 2.2).Such a matrix provides three eigenvalues λ1, λ2, and λ3 (in descending order).These attributes allow to discriminate water surfaces, which are both horizontal and planar, from land, even for pixels lying on the ground.Two eigenvalue-based features were computed.The first one was the smallest eigenvalue, λ3 (called volume afterwards).Indeed, λ3 0 for planar elements, λ3 > 0 for ground pixels (due to ground microreliefs and surface roughness) and λ3 >> 0 for real 3D Earth surfaces.The second feature is the scatter, defined as S = λ3/λ1 ∈ [0, 1].After all, in case of lidar data with multiple strips, the feature set is: {HG, Dm, Dr, λ3, S }.For single strips, we use instead: {HG, D, λ3, S }.

Learning step
Many authors have shown that the size and composition of the training samples have a substantial impact on the classification accuracy (Foody and Mathur, 2004).For non-parametric classifiers such as SVM, only samples lying on the edges of a given class distribution in data space contribute to the analysis.Here, the subset generation is carried out automatically using the scatter and the volume features.Extreme values are representative of land and water classes, respectively.High scatter values correspond to areas with significant vertical scattering such as vegetation and buildings.Lowest volume values indicate flattest areas i.e., water surfaces close to the flightline nadir.Thus, we retrieve pure pixels (hereafter called seeds) by thresholding volume and scatter cumulative distribution functions (computed with 500,000 pixels randomly taken).Gradient is computed on both curves.For the volume (scatter) curve, areas inferior (superior) to the highest gradient value are labelled as water and land seeds, respectively.However, seed distribution is not satisfactory.For water areas, few nadir points are selected whereas they are the most similar to land surfaces, while for land regions, significant spatial heterogeneity exist and heterogeneous training pixels cannot be selected with such method.Since they are not representative of both classes, they are only used as a basis for designing a more adapted training set.This step is based on the historical coastline (HCL) or a rough landwater interface provided by an end-user.Discrepancies may exist between the dataset and such manually designed border, meaning that it cannot be used directly for enrichment.Consequently, we generate a buffer zone centered on the HCL, to prevent the aggregation of unreliable pixels that would deteriorate the final classification accuracy.A percentage of seed points within such buffer zone is used.Starting with pixels belonging to the HCL, the dilatation of such region is stopped when 40% of the seed points in both classes are included in the buffer (the [35-50%] interval leads to similar results).Then, connected regions are retrieved with a standard binary region growing procedure, and these regions are labelled as water or land according to the majority of seed points existing in such regions.At last, 1% of the pixels of these regions are randomly selected as training pixels.

Land-Water classification and regularization
For Support Vector Machines classification, the standard Gaussian kernel is selected and the SVM hyperparameters are optimized with a simple grid search.Despite their ability to handle high-dimensional feature spaces, SVM are limited to local classification of pixels, which results in noisy outputs.Several solutions are possible to overcome such limitation (Schindler, 2012).The simplest solution, which is adopted in this paper, is the filtering approach.For each pixel, a new label is computed with respect to the weighted incoming labels of a given neighborhood.Such a solution offers the advantage of linear growth of the computational cost with the number of pixels.For that purpose, we have adopted the relaxation probabilistic framework (Gong and Howarth, 1989).We take into account both neighborhood information and the probabilities of belonging to each class, as provided by libSVM (Wu et al., 2004).This is an iterative algorithm in which the probability values for each pixel are updated to make them closer to the probabilities of their neighbours.Thus, we compute the membership energy for both classes and assign the label corresponding to the lowest value.Such energy E for label l and pixel i is computed as follows: Gσ( i − j ) is the weight of pixel j, where Gσ corresponds to the zero-mean Gaussian density function with variance σ 2 (σ = 1 in our experiments) and V (i) is the vicinity window (here a 5×5 window).Finally, M is called compatibility matrix since it measures the compatibility between pixel i with label l and pixel j with label lj.We have: M defines a priori correlations between the probabilities of neighboring points and corresponds to conditional probabilities verifying: 0 ≤ Mi,j[l, lj] ≤ 1, and l Mi,j[l, lj] = 1.The coefficients have been empirically selected so as to enforce spatial homogeneity.However, in order to preserve water/land boundaries, the coefficients out of the matrix diagonal are not equal to 0.

DATASETS
Numerous datasets (Table 1) were processed in order to study the behaviour of the method for a large variety of landscapes.
• Perpignan area corresponds to 50km of the French Southern Mediterranean coast.Seven relevant areas (estuary, harbour, marshes, etc.) are extracted from this area.
• Martinique corresponds to the Northern part of the island located in the Lesser Antilles in the Caribbean Sea.It was selected for the sharpness of the seashore and the high density of the vegetation.
• Brittany is located in the North-West of France.It has been selected for our study because the survey was carried out during a high-tide period.Consequently, different water surfaces may exist at the same location depending on the hour of acquisition.
• Rhone corresponds to a part of the Rhone river, representative of a river area.Conversely to seashore areas, the point density is nearly constant on the water surface.
• A length of 6 km of the Moselle river in the North-East of France is selected.In this region, the river is meandering due to the flat topography.The area is very challenging since many small islands and sandbanks are scattered in the river bed, many lakes and reservoirs closely located to the river exist.
International Archives of the Remote Sensing and Spatial Information Sciences, Volume XL-1/W1, ISPRS Hannover Workshop 2013, 21 -24 May 2013, Hannover, Germany • Miami covers a district of the city of Miami (Florida, USA), representative of a dense city center, with complex structures (bridges, cranes, buildings) located close to water areas. •

RESULTS
For the Rhone, Moselle, Miami and Louisiana areas, since no historical coastline is available, the coarse knowledge of the landwater interface is substituted by a rough 2D manual plotting.Table 2 and Figure 1 show that very good results are achieved for land-water classification for all areas of interest.Since we deal with balanced classes, the standard Overall Accuracy (OA) is used as classification accuracy measure.The ground truth has been manually obtained.
The OA ranges from 95 to 99.5% for 10 out of 13 areas of interest.This proves both the SVM generalization ability and the efficiency of the designed learning procedure.In particular, even when the training step focuses on seashore areas or main river beds, the SVM is able to generalize to inland waters, harbours, other rivers and channels.Even with small training sets, the wide range of behaviours of both land and water areas is captured.Various anthropic and natural land cover classes are all correctly labelled, similarly to water surfaces with varying lidar point distributions.The method can deal with land regions with complex topography but also with terrains with high frequencies corresponding to cities. Indeed, land pixels are perfectly classified despite varying slopes, land covers and uses and city densities.Secondly, results show the robustness of the spatial features for quantifying the 2D and 3D distribution variability, namely fluctuations (due to the changing strip overlap, with only parallel and/or orthogonal configurations etc.) and low point densities.This simply shows the relevance of the proposed features.We can note that conversely to most of water classification existing methods, our approach is not limited to natural areas where no anthropic items exist.
The relevance of the filtering strategy using probabilistic relaxation is also assessed.Salt-and-pepper noise is removed and sharp boundaries preserved.Pixels misclassified by the SVM and located in the middle of both regions are corrected, as well as part of those lying on class edges.Contrary to a simple Gaussian weight that would have blurred such edges, the addition of the high confidence of the SVM classifier allows to preserve such borders.The Overall Accuracy is only slightly improved (mean value: +0.05%→0.8%),revealing that the class with the local highest confidence is dilated.Since for the 1 m resolution grid Table 2: Overall Accuracy (OA) for areas introduced in Section 3. Accuracy has been assessed of the whole areas.The left column corresponds to regions selected within the Perpignan dataset whereas the right column corresponds to various areas disseminated all around the world.a 5×5 neighborhood and a spatial Gaussian weight of σ=1 were adopted, the local erosion of either land or water classes is inferior to 1 m.This may remove very thin objects but does not affect objects of interest such as rivers, breakwaters and bridges (see Figure 1).Due to the higher homogeneity of behaviours of the water points, water areas are better discriminated than land points for both regions.Misclassification appears only very locally, meaning that this would not corrupt the subsequent DTM computation.Finally, the processing time of the method has been evaluated for 10 distinct areas within the seven datasets.The average computing time for the full workflow is around 30 minutes for 1 km 2 and approximately 4 million points (processor: 2.83 GHz with 3.9 GB RAM).Half of the time is spent for 3D feature computation and probabilistic relaxation.

DISCUSSIONS
A closer look to the classification of the thirteen datasets (Figure 1) allows a more in-depth analysis of the proposed method.Dense city centers, harbours, estuaries, cliffs, complex rocky shorelines, sandy beaches are correctly classified.Figure 1 shows the approach versatility in spite of both the significant heterogeneity of the land object 3D signatures and the very different sizes of land objects to be discriminated within each area.Moreover, the classifier is not perturbed by strong variations of the surface reliefs and is able to deal with various vegetation canopy covers or building densities (Miami).Such ability is illustrated for Complex land area and Cliffs datasets, where very high OA are reached (99.85 and 98.84%, respectively), despite steep slopes.The high vertical scattering of the lidar measurement is well modeled with the scatter feature and is used in conjunction with 2D density-based attributes, as proposed, to discriminate forested areas with dense canopy covers (Martinique) from water areas (OA=97.9%).Furthermore, small features are correctly preserved: rocks (Cliffs and Rocky sharped coast), bridges (Harbour), canal locks (Moselle) or ships of 5-10 m length (Harbour, Rhone, or Brittany), and small water surfaces on land regions (Brittany).In Miami dataset, bridge pillars are correctly labelled as land even if the point cloud is locally very sparse (mean density around 0.3 pt/m 2 for an overall accuracy of 97.15%, see Figure 2).In particular for DTM generation and subsequent flood modeling, rocks, breakwaters and bridges must be well classified (>99% in all cases).Despite low point densities (i.e.below 4 pts/m 2 ), these items are correctly retrieved enhancing the relevance of the proposed classification International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XL-1/W1, ISPRS Hannover Workshop 2013, 21 -24 May 2013, Hannover, Germany , where hydrographic networks are accurately delineated.In particular, in Moselle, one can show that the canal, small lakes and river branches are detected as water areas despite very few lidar points.Therefore, our method proposes a strong alternative for coastal object detection (tidal channels, trenches, gullies), assuming they are filled with water.This is no longer necessary to develop on-purpose tracking methods.Additionally, smooth transitions between both classes are well detected.This is the case for sandy beaches (Breakwater), river banks (Rhone, see Figure 2), and even for more complex landscapes such as Louisiana.Overall accuracies (OA) may be slightly inferior for such regions (95-96%), which comes from misclassification of a strip of 1-2 pixel width located on each side of the land/sea interface.This explains the difference in accuracy between regions with almost similar landscapes.For instance, Harbour and Estuary areas exhibit close topographic behaviors but sharp anthropic elements within the Harbour dataset allows to better delineate the land-water interface (OA=99.2%> 98% for the Estuary dataset, see Figure 2).A visual assessment confirms that errors mainly correspond to punctual misclassification.Wrongly assigned isolated pixels are present in both classes and are located when similar behaviours exist.In water areas, this corresponds to 3D points close to the strip centers, for which higher and regular point densities are reported.This is due to the specular properties of water surfaces.
Conversely, for land regions, pixels are labelled as water for flat and very irregular sampled surfaces.This is the case for asphalt roads, as noticed in (Höfle et al., 2009).Two main problems can be noticed, corresponding to the areas with lowest quality measures (Brittany, Moselle, and Louisiana), where 85%<OA<95%.On the one hand, classification can lo-cally fail in wetlands areas (marshes and mangroves).In the case of shallow waters, the laser pulse backscattering varies both with incidence angle and water depth.In addition, the presence of many small linear anthropic structures leads to difficult 2D and 3D point distribution modeling (Marshes).Consequently, the confusion between classes is higher (OA=85%).A similar behaviour is noticed for the Moselle dataset.While the canal is correctly classified, the river exhibits a larger number of erroneous labels (Water accuracy: 80.2%) due to very shallow depths and sandbanks coupled with narrow river branches.The accuracy of the results also depends on the tide strength, on the depth of the catchment basins, and on the sediment level.For instance, the Brittany dataset illustrates problems occurring when an area is acquired by several lidar strips with high tidal conditions.Several parallel sea surfaces exist in the point cloud at the same 2D location, and the classifier labels such areas as land (top left of the area in Figure 2).Moreover, for low point densities, the strong overlap between laser strips may result in densities similar in land and water areas.For nearly featureless areas, if the terrain is flat or shallowly slopped (e.g., Louisiana dataset), misclassification occurs in water surfaces, where the lidar point density is doubled (Figure 1).Nevertheless, since the land-water classification is mainly performed for shoreline derivation, such errors do not result in a decrease in accuracy of the land-water interface.Furthermore, the choice of a 2D-based classification can lead to other kinds of misclassification.First, for pixels mixing water and land points, both labels are conceivable.The selection of a 1 m cell size offers a good trade-off between computing time and accuracy but does not prevent erroneous labelling.It can be noticed for cliffs where high vertical shifts exist, resulting in a small dilatation of the land surfaces (Figure 2, Cliffs).This is also the case for harbour areas where small objects such as jetties and boats are partly classified as water.This is illustrated in Figure 2.

CONCLUSIONS
This paper proposes a workflow for land-water discrimination, which is efficient in terms of both classification accuracy and computing time.It applies to various landscapes: it is suitable for small or large area mapping and it is able to deal with single strips or complex survey configurations.The adoption of a stateof-the-art supervised classifier, the Support Vector Machines, has allowed to benefit from two strong features: a fast classifier with high generalization capability, and a method that has only to be fed with a list of attributes.Furthermore, we only need a coarse delineation of water and land regions or an existing coastline or river delimitation.Such inputs do not need to be geometrically accurate, and the effectiveness of our approach shows the potential for both coastline database refinement and updating.
The different landscape classifications prove the method reliability for both inland water surfaces, rivers and seashore areas.Objects of high relevance such as bridges are conserved, even with low point densities, allowing subsequent surface modeling valuable for many applications such as flood simulation purposes.
The method was established for areas with sufficient water surface but turns out to be efficient for linear features such as rivers and channels.A specific filtering step coupled with top-down knowledge (e.g., river banks are most of the time parallel to each other) would definitively improve the results for the latter areas.

Figure 1 :
Figure 1: Results for the various areas of interest.Both orthoimages and classifications ar provided (I Land -I Water -I No data).The ten areas on the left correspond to 1 km×1 km tiles.The orange line shows the historical coastline.features.Moreover, surveys with relatively low spatial resolution are sufficient for water monitoring since oceans, rivers, and shallow channels are correctly classified.This is enhanced in Figure 1 (Estuary, Rhone, Moselle, Louisiana), where hydrographic networks are accurately delineated.In particular, in Moselle, one can show that the canal, small lakes and river branches are detected as water areas despite very few lidar points.Therefore, our method proposes a strong alternative for coastal object detection (tidal channels, trenches, gullies), assuming they are filled with water.This is no longer necessary to develop on-purpose tracking methods.Additionally, smooth transitions between both classes are well detected.This is the case for sandy beaches (Breakwater), river banks (Rhone, see Figure2), and even for more complex landscapes such as Louisiana.Overall accuracies (OA) may be slightly inferior for such regions (95-96%), which comes from misclassification of a strip of 1-2 pixel width located on each side of the land/sea interface.This explains the difference in accuracy between regions with almost similar landscapes.For instance, Harbour and Estuary areas exhibit close topographic behaviors but sharp anthropic elements within the Harbour dataset allows to better delineate the land-water interface (OA=99.2%> 98% for the Estuary dataset, see Figure2).A visual assessment confirms that errors mainly correspond to punctual misclassification.Wrongly assigned isolated pixels are present in both classes and are located when similar behaviours exist.In water areas, this corresponds to 3D points close to the strip centers, for which higher and regular point densities are reported.This is due to the specular properties of water surfaces.Conversely, for land regions, pixels are labelled as water for flat and very irregular sampled surfaces.This is the case for asphalt roads, as noticed in(Höfle et al., 2009).Two main problems can be noticed, corresponding to the areas with lowest quality measures (Brittany, Moselle, and Louisiana), where 85%<OA<95%.On the one hand, classification can lo-

Figure 2 :
Figure 2: Classification results displayed in 3D for various areas of interest (see text for more details).

Table 1 :
Louisiana corresponds to an area in a bayou of Louisiana, USA (Passe A Loutre State Wildlife National Park, East part of the Timbalier Bay).It is a complex natural wetland, composed of many water meandering channels.Shallow slopes and little topographic relief exist.Description of the areas of interest.Miami and Louisiana datasets have been provided by the NSF through the OpenTopography portal (NSF, 2012) and acquired with Leica ALS50 and ALS40 respectively.Other datasets have been acquired by the French Mapping Agency using an Optech 3100EA device.