Range Versus Surface Denoising of Terrestrial Laser Scanning Data for Rock Discontinuity Roughness Estimation

Surface roughness represents a major component of rock discontinuity shear strength. To achieve comprehensive, accurate, and efficient estimates of in situ discontinuity roughness, the traditional contact measuring methods are being replaced by advanced remote-sensing technologies. Terrestrial laser scanner (TLS) is well suited for measuring large inaccessible discontinuities; however, inherent TLS range noise strongly influences the surface details and roughness estimation. The aim of this research is to establish an optimal wavelet-denoising procedure for the TLS data acquired with different scanning configurations (range and incidence angle), and for rock discontinuities having different roughness characteristics and surface reflectivity. The conventional discrete wavelet transform and stationary wavelet transform in combination with four threshold selection methods are applied on TLS data in the direction of range measurements (range denoising) and in the direction perpendicular to the best-fit plane (surface denoising). The performance of the denoising procedures is assessed by comparing the range and surface-denoised TLS surfaces with reference surfaces acquired with the Advanced TOpometric Sensor. Comparative analyses of the roughness calculated according to the angular thresholding method (Grasselli, in Shear strength of rock joints based on quantified surface description, Ph.D. thesis. EPF Lausanne, Lausanne; Grasselli, Shear strength of rock joints based on quantified surface description, Ph.D. thesis, EPF Lausanne, Lausanne, 2001) indicate that all the denoising methods improve the roughness estimated from the TLS data appreciably; however, the level of improve-ment depends intrinsically on geometrical characteristics of the rock surface and scanning configuration. Range denoising has been found to provide more

Maximum apparent asperity angle A * Total potential contact area ratio A 0 Surface area defined by an apparent dip greater than 0° normalized with respect to the total area of the surface C Dimensionless empirical fitting parameter calculated via a non-linear least-squares regression G Grasselli parameter Noise level e Noise-level estimation T Threshold value T 0 Threshold chosen according to one of the threshold selection methods for a signal model including white noise ( = 1) cD 1 First-level detail coefficients β Analysis direction for which the Grasselli parameter is computed G ATOS Grasselli parameter of reference ATOS surface G TLS Grasselli parameter of TLS surface

Introduction
Surface roughness represents a significant component of discontinuity shear resistance and influences the overall deformation behaviour of jointed rock masses. Since surface roughness is scale-and direction-dependent (ISRM 1978), it is preferable to perform measurements in the anticipated shear direction and at the engineering scale of interest. Larger scale roughness features are referred to as "waviness" and represent surface irregularities with a wavelength greater than about 10 cm (Priest 1993). Smaller scale features are referred to as "unevenness" and cover finer features, which are superimposed on the waviness. In the laboratory, roughness is traditionally measured on small rock samples (up to 50 cm length) using several measuring methods such as a mechanical profilometer (ISRM 1978), shadow profilometer (Maerz et al. 1990), laser profilometer (e.g., Huang et al. 1992), or a structured light projection method using a stereo-topometric camera, such as the Advanced TOpometric Sensor-ATOS (e.g., Grasselli et al. 2002). In situ discontinuity roughness is routinely estimated visually (using a catalogue of joint profiles), or measured using the traditional contact-based methods, e.g., compass and disc-clinometer method (Fecker and Rengers 1971), and linear mechanical profiling using a straight edge or contour gauge (Barton and Choubey 1977;Milne et al. 2009). The traditional in situ measurements provide analogue, discrete data that are restricted to an accessible area. In contrast, remote-sensing technologies such as photogrammetry and terrestrial laser scanning (TLS) allow the in situ acquisition of a large and remote surface in a short period of time (e.g., Poropat 2009). These methods are able to characterize the 3D surface structure, which is represented as a dense and precise point cloud that can be used for surface roughness estimation at different scales and along any direction.
The traditional stereo-photogrammetry has been used for roughness measurements for nearly 30 years (Lee and Ahn 2004;Haneberg 2007;Baker et al. 2008). New developments in computer vision methods, in particular Structure from Motion (SfM) combined with dense image matching (Hartley and Zisserman 2003;Carrivick et al. 2016), allow 3D reconstruction of a detailed surface from an array of unconstrained images. For SfM and dense image matching, the surface is preferably photographed at close range and from many different positions and orientations, which is not always possible in constrained terrain. Furthermore, reference points (visible on the acquired photos) need to be established and measured to define the linear scale, and setting the points can be time-consuming and precarious in rugged alpine terrain.
Terrestrial laser scanner has been successfully applied to in situ roughness measurements (Fardin et al. 2004;Renard et al. 2006;Rahman et al. 2006;Tesfamariam 2007;Haneberg 2007;Candela et al. 2009;Khoshelham et al. 2011;Pollyea and Fairley 2011;Mills and Fotopoulos 2013). While reasonable results have been obtained in quantifying waviness (Fardin et al. 2004), finer details of unevenness have been hindered by data precision and laser footprint size. The data precision mainly depends on the inherent random range error (noise), which results in an overestimation of surface roughness (e.g., Kulatilake et al. 2006;Poropat 2009;Khoshelham et al. 2011). The challenge is to eliminate noise, but preserves surface details.
With TLS, noise has typically been removed by surface interpolation techniques such as averaging the TLS range measurements (Schulz et al. 2008), orthogonal least squares (Fardin et al. 2004;Pollyea and Fairley 2011), the robust interpolation method RANSAC (Grasselli et al. 2002), or Fast Radial Basis Function (Rahman et al. 2006;Tesfamariam 2007). The disadvantage of interpolation methods is that irregular rock surfaces are smoothed and some topographic details are lost. By reducing the spatial complexity of a 3D randomly scattered point cloud (i.e., gridding to a regular 2.5D mesh), a wider range of image processing algorithms can be applied. An overview of image denoising methods and further references can be found in (Buades et al. 2005;Smigiel et al. 2011;Zhang et al. 2014).
Wavelet-denoising methods, including the discrete wavelet transform (DWT) and the wavelet packet (WP), have been applied to profile measurements subjected to fractal-based roughness length analysis and have been found to provide better roughness estimates than if no denoising was performed (Khoshelham et al. 2011). An important advantage of DWT compared to the other denoising techniques is that the TLS noise level is correctly estimated from the raw data (Bitenc et al. 2015a). Prior research has focused on wavelet denoising in a direction perpendicular to the best-fit plane (hereafter referred to as surface denoising). However, denoising in the range direction (hereafter referred to as range denoising) is considered preferable, since noise mainly relates to range.
Despite prior research showing the applicability of TLS to quantifying roughness, thorough investigations concerning the influence of noise and the choice of optimal denoising procedure under conditions of variable range, scanning direction (incidence angle), surface roughness, and reflectivity have been lacking. This contribution evaluates TLS data acquired for four rock samples having different roughness and reflectivity characteristics, and at different scanning ranges and directions. Reference data were measured with ATOS and the roughness was calculated according to the angular threshold method (Grasselli 2001). The objectives of this research are to: • investigate the influence of intrinsic TLS data noise on roughness (noise effect); • evaluate TLS noise estimations derived from the data itself using DWT and stationary wavelet transform (SWT), in the range and surface roughness measuring directions; • analyse the performance of eight denoising methods including DWT and SWT, in combination with four threshold selection methods that are applied in range and surface roughness measuring directions.
The success of wavelet-denoising procedures is evaluated through a comparative analysis of reference ATOS and denoised TLS surfaces and their derived roughness values. The goal is to identify optimal TLS denoising procedures for different rock surfaces and scanning configurations.

Parameterizing Rock Discontinuity Roughness
Discontinuity surface roughness describes the local departures of the actual surface from planarity or any higher order reference surface. Roughness can have a prevailing influence on the shear strength, particularly in cases of low normal stress combined with unfilled discontinuities. However, the parameterization of roughness, to fully capture the influence of roughness on shear strength, remains a challenge; parameterization needs to consider that roughness is direction and scale dependent.
The empirical roughness parameter developed by Grasselli (2001), hereafter referred to as the Grasselli parameter (G), has been adopted in this research, since it quantifies the direction dependence of roughness, is applicable to rock surfaces, and is least sensitive to data noise (Bitenc et al. 2015c). The Grasselli parameter is based on an angular threshold concept and was initially developed to identify potential contact areas during direct shear testing of artificial tensile rock fractures. Studying the damaged areas after shear testing, it was found that only the parts of the surface that face the shear direction and are steeper than a threshold inclination provide shear resistance. Surface data are first transformed into the average discontinuity plane (or shear plane) coordinate system and modelled with Delaunay triangulation. For each triangle, an apparent dip * is calculated according to the specific analysis direction and the shear plane. The sum of triangulated areas that are steeper than a certain value of * (denoted as A * and referred to as the total potential contact area ratio) is plotted against * (Fig. 1).
Based on curve fitting and regression analysis of A * as a function of * , G is calculated as follows: where * max is the maximum apparent asperity angle of the surface in the shear (analysis) direction and C is the dimensionless empirical fitting parameter calculated via non-linear least-squares regression. To account for surface anisotropy, the following correction has been proposed (Tatone and Grasselli 2009): where A 0 is the surface area defined by an apparent dip greater than 0°, normalized with respect to the total area of the surface. In the normal case that A 0 equals 0.5, the roughness matrix reduces to the simplified version in Eq. (1). The expression in Eq. (2) is used in this research to estimate the surface roughness. (1)

Denoising Laser Range Data Using Wavelets
Terrestrial laser scanner principles and limitations are enumerated below, along with the wavelet-denoising procedure.

Terrestrial Laser Scanning
Terrestrial laser scanner has a ranging and scanning unit for performing measurements. Range can be measured based on roundtrip time-of-flight (TOF scanner) or phase shift (phase scanner). Phase scanners are more precise (on the order of 1 mm), but have a limited range of a few tens of metres. TOF scanners have a precision in the order of 3-5 mm and range of up to a few kilometres, and are commonly used for rock mass characterization. The scanning unit consists of a motor-driven head with rotating mirrors to spatially redirect the laser beams to cover large areas. With TLS, the position and intensity of thousands of points can be recorded within a few seconds, at an angular resolution typically smaller than 0.01°. The point position is originally recorded in polar coordinates, which are, in real time, recalculated and stored as Cartesian coordinates. Intensity is a relative unit less number and is influenced mainly by the reflectance of the target surface and the laser scanning configuration (Wagner 2005). In general, the ability for TLS to capture surface roughness details depends on: the effective spatial resolution of TLS points and the error (noise) associated with the measurements. The effective resolution depends on sampling interval and laser beam footprint size (Lichti and Jamtsho 2006), and causes roughness underestimation (smoothing effect), while the presence of intrinsic TLS noise results in roughness overestimation (noise effect).
TLS point cloud errors arise from the following sources (Soudarissanane et al. 2011): (1) range and scan angle errors; (2) the scanning configuration and the material reflectivity; and (3) environmental conditions including lighting, humidity, and temperature. With systematic errors properly calibrated and removed, the remaining random error represents the measurement noise. Wujanz et al. (2017) showed that the range noise can be jointly considered as a function of backscattered signal strength or intensity I and proposed a one-term power series model = a × I b , where the unknown parameters a and b are determined within least-squares regression.

Wavelet Denoising
The wavelet transform (WT) is an extension of the classical Fourier transform. The WT basis functions (mother wavelets) exist in space, frequency, and amplitude, are irregular and asymmetric, and have limited duration. This makes the WT particularly suitable for natural rock surfaces that typically contain finite, non-periodic, and/or non-stationary signals that are characterized by discontinuities and sharp peaks. Scaling and shifting are central concepts of the WT. Scaling refers to the process of stretching or shrinking the wavelet in space, thus changing the frequency. Shifting a wavelet refers to delaying or advancing the onset of the wavelet along the length of the signal. The input signal is transformed into the space-frequency domain, where the particular unwanted or unneeded space or frequency components can be removed. The discrete wavelet transform (DWT) was developed to enable unique signal reconstruction using orthogonal wavelets and has, therefore, been widely applied for denoising and compression . The DWT has been adopted in this research.
The main steps of wavelet denoising involve multi-level signal decomposition, thresholding of detail coefficients, and signal reconstruction (Fig. 2). The signal is first decomposed (filtered) into several levels (j = 1,…,N) of approximation coefficients, cA j , and detail coefficients, cD j . This is performed using a DWT, a mother wavelet (filter), and N decomposition levels. An example of two-level DWT is shown in Fig. 3. Approximation and detail components contain low-and high-frequency contents, respectively. Next, an appropriate threshold value is calculated and cD j values below the threshold are discarded, because they are assumed to represent only noise. The cD j values above the threshold can be reduced, depending on the chosen threshold mode. The thresholding is intended to suppress the noise and preserve the original surface. Finally, the surface is reconstructed using the original cA N of the last level N and the modified cD j from levels 1 to N (j = 1,…,N). In this research, two DWT methods are tested; the conventional decimated DWT (e.g. Daubechies 1992; and its undecimated version, referred to as the stationary wavelet transform, or SWT (Coifman and Donoho 1995). Advantages of SWT are that it is shift invariant, provides more precise information regarding frequency localization, enables direct correlation of the space scale to the original data, reduces the edge effect (Coifman and Donoho 1995;Buades et al. 2004), and shows superior performance for image denoising (Gyaourova et al. 2002;Starck et al. 2004). For denoising 3D rock surfaces, 2D DWT and 2D SWT are used, where the input signal is a 2.5D grid (an image).
The threshold value (T) is a key to successful wavelet denoising, and is, in general, calculated as follows: where is the standard deviation of the noise and T 0 is the threshold chosen according to one of the threshold selection methods. In the work presented herein, the local fixedform universal (Donoho and Johnstone 1995) and the global penalised (Birgé and Massart 1997) threshold selection methods are tested. The local fixed-form universal threshold is used most commonly if the signal-to-noise ratio is small. The global-penalised threshold is a variant of the fixed-form and includes an adjustable (sparsity) parameter. A higher sparsity parameter returns a higher threshold, resulting in more coefficients being eliminated and a sparser (smoother) signal representation. In this research, the default sparsity parameter values for penalised low, penalised medium, and penalised high thresholds were tested. The noise level is estimated as a robust standard deviation of first-level detail coefficients cD 1 (Donoho and Johnstone 1995) and is defined as e :

Experiments and Results
The experiments described below were conducted to evaluate the effect of TLS noise on rock surface roughness and to analyse TLS noise details. The experiments were also designed to investigate optimal wavelet-denoising procedures in relation to: (1) threshold selection method (fixed form or penalised); (2) wavelet transform (DWT or SWT); and (3) denoising direction (range or surface denoising). Comparative analyses of noisy and denoised TLS data, and ATOS data based on the Grasselli parameter values were performed for different scanning configurations, and rock reflectivity and roughness characteristics.

Data Acquisition
The four rock samples comprising the experimental data set are shown in Fig. 4, along with the reference direction for calculating the Grasselli parameter.
Each sample was fixed on the wooden board equipped with four 10 cm square TLS black-white targets on top of which a 7 mm circular ATOS black-white target was pasted (see Fig. 5, right). The board with samples and targets was simultaneously scanned with the terrestrial laser scanner Riegl VZ400 and the optical 3D coordinate measuring sensor GOM ATOS I (Fig. 5). The key technical specifications of these sensors are summarized in Table 1. The principles of terrestrial laser scanner and 3D optical scanner are, in more detail, explained in, e.g., (Pfeifer and Briese 2007) et al. 2006), respectively. The black-white targets included on the sample mounting board were needed to coregister the TLS and ATOS data in postprocessing. Highprecision TLS target centers were measured in the point cloud by applying an algorithm based on image matching (Kregar et al. 2013) and ATOS target centers were identified automatically with built-in software.
For each sample, the TLS scans were acquired at 10 m intervals in the range of 10-60 m. For each range position, scans were performed perpendicular and obliquely to the mean rock plane (exception of an oblique scan that was not performed at the 60 m range, because of the low point resolution), resulting in 11 data sets. For the oblique orientations, the mean sample plane was rotated about a vertical axis and the angle between the plane normal and the TLS line of sight was between 25° and 47°. Table 2 summarizes the scanning configurations and is consistent with terminology presented in data graphs.
The scanning resolution on the rock surface ranged from approximately 0.2-2 mm. For a filtered 3D point cloud, Cartesian (X S , Y S , Z S ) and polar (Φ, Θ, R) coordinates along with raw intensity values were exported in the scanner coordinate system. ATOS data were acquired at a range of approximately 0.7 m in an indoor laboratory environment. Changing the position of sensor head with respect to the rock sample, shadow zones were eliminated, and a full and detailed 3D surface model was acquired. 3D images were automatically merged into a single point cloud with the help of the circular black and white targets (diameter 7 mm) sticked on the rock surface (see Fig. 4, basalt and clay stone samples). The ATOS high point density was reduced to the TLS point density of approximately 1 point/mm 2 using the subsample tool in CloudCompare (2018) to eliminate the influence of data resolution on roughness comparisons.

Data Processing
The overall workflow involved in wavelet denoising of TLS data and performing a comparative analysis of Grasselli    Fig. 6 and enumerated below. The processing was performed in Matlab. Wavelet denoising in 2D is an image processing technique requiring a full matrix (an image). Therefore, a range image was constructed from randomly scattered TLS points (Φ, Θ, R) as follows: 1. In the ΦΘ-plane, a rectangle including the rock surface to be denoised plus a surrounding buffer zone was defined. The buffer zone is necessary to prevent edge artefacts from contaminating the results of rock surface denoising. 2. Within the rectangle, a regular grid was created with an angular spacing equal to 1 mm at the corresponding Fig. 6 Data processing workflow for wavelet denoising of TLS data and performing a comparative analysis of Grasselli roughness parameters scanning range. The grid point spacing of 1 mm was selected to be equal to or larger than the scanning interval for most of the TLS data sets. Creating a regular grid assured uniform point density among the TLS data sets. 3. For each grid point, the range value was interpolated using the Nearest-Neighbour (NN) method.
In 2D wavelet denoising, range image was decomposed with the two wavelet transform variants, DWT and SWT, applying the most general and widely used Daubechies wavelet db3. The number of decomposition levels was determined according to the range image size and was set to 3. The DWT and SWT detail coefficients were then thresholded with the local fixed form and three global-penalised thresholds, penalised high, medium, and low, respectively. All thresholds were applied in the hard mode, following the previous findings that hard thresholding is more suitable for rock surface roughness estimation than soft thresholding (Khoshelham et al. 2011;Bitenc et al. 2015a). This denoising process resulted in a new set of range-denoised (R-den) images.
For surface reconstruction, only the grid points, which have an original TLS point closer than the defined angular spacing, were used. These subsets of grid points are referred to as valid grid points. For surface denoising, and roughness estimation and comparison, polar coordinates of valid grid points were transformed to Cartesian coordinates. These transformed TLS point clouds and the ATOS point cloud were then co-registered in a common coordinate system. It was defined by the targets on the wooden board and it is referred to as the mean plane coordinate system.
A surface image was created analogous to range image construction, except a rectangular area was defined in the mean plane. Following the wavelet-denoising procedure for range image, the surface-denoised (Z-den) images are obtained. To preserve point distribution in the mean plane, the denoised Z-values were linearly interpolated from Z-den images, resulting in the Z-den point cloud.
Noisy, R-den, Z-den, and ATOS point clouds within an identical rectangular area, which was defined in the mean plane coordinate system and is referred to as the common area, were used for roughness comparisons. The Grasselli parameter (G) was calculated for 72 analysis directions (β i , i = 1,…,72) spanning 5° increments. The resulting Grasselli parameters are denoted as G n for the noisy data set, G r for R-den and G z for Z-den data sets, and G ATOS for the ATOS data set. The accuracy of roughness estimates for noisy and denoised TLS surfaces G TLS was judged by comparing G TLS to G ATOS , and calculating a mean relative difference across all the analysis directions. This measure is referred to as the Grasselli parameter estimation error and is expressed as follows:

Results and Discussion
Details regarding rock surface representations, TLS noise effects, noise estimation, and optimal denoising procedure for rock discontinuity roughness estimation are summarized below.

Rock Surface Representations
For visual impressions of the rock samples, Fig. 7 shows the triangulated ATOS surfaces (Fig. 7a-d) and representative examples of noisy (Fig. 7e-h) and range-denoised (Fig. 7i-l) TLS surfaces. Noisy surfaces represent TLS data acquired at a range of 30 m in the perpendicular direction, and the range-denoised surfaces were obtained using the SWT with penalised high threshold, R-den SWTph. The differences of ATOS surface and TLS noisy (dZ noisy ) and denoised surfaces (dZ R-den ), are shown in Fig. 7, m-p and r-u, respectively. Values of dZ noisy are larger in areas of fast changing morphology, due to the lower resolution of TLS data. High dZ noisy values appear randomly across the area of sample 0936 and are probably due to variations in surface reflectivity. The pattern of dZ R-den values corresponds to the pattern of dZ noisy . In addition, in case of sample 0727, higher dZ R-den values can be observed along surface discontinuities (ridges), which indicate that wavelet denoising smoothes sharp edges. Medians and robust standard deviations (robust STD) of height differences for the 11 scanning configurations are presented in Table 3. Negative median height differences indicate systematic co-registration errors, which result from uncertainties of TLS target centre estimation (Kregar et al. 2013), but do not influence surface roughness values. The robust STD of dZ noisy values is an indicator of TLS noise level, and attenuates after range denoising is performed.

TLS Noise Effect on the Grasselli Parameter
The Grasselli parameter was computed for reference ATOS data (G ATOS ) to obtain insight to the sample roughness characteristics, and to facilitate roughness comparisons. Figure 8 depicts the direction-dependent G ATOS for the four rock surfaces, and Table 4 summarizes the common area size, and median robust STD values of G ATOS for all analysis directions. Sample 1309 is the roughest, with the median value almost three times higher than the smoothest sample 0727. Sample 1309 also has the least variability with analysis direction as indicated by the low robust STD. The robust STD for sample 0936 shows the highest roughness anisotropy.

Noise Estimation
The noise ( e ) was estimated using Eq. (4), where first-level detail coefficients were obtained from DWT and SWT transforms of range and surface images. The dependence of e on recorded intensity is for all data sets and decomposition combinations, as depicted in Fig. 10, together with the fitted curve of the one-term power-series (Wujanz et al. 2017) and the root-mean-square error (RMSE) as a measure of goodness of fit. Figure 10a, b shows these results for range and surface images, respectively.
These data show that σ e is inversely proportional to the mean intensity, that σ e using SWT is slightly higher than for DWT, and the σ e values calculated from range images are less disperse and have smaller RMSE than for surface images.
These results suggest that the dependence of noise on scanning configuration and surface reflectivity can be jointly explained by the TLS measurements of intensity. Highly reflective limestone (sample 1309) and schist (sample 0936) entail higher intensity values and lower σ e compared to dark basalt (sample 0727) and clay stone (sample 0934), which absorb laser light causing lower intensity values and, consequently, higher σ e . However, the strong reflection of sample 0936 scanned at 10 m range in the perpendicular direction results in higher σ e . Figure 11 shows the σ e dependence on the scanning configuration estimated using SWT (DWT estimates are similar and, therefore, excluded). The zig-zag pattern of σ e calculated for surface images (Z-den) shows that oblique scans have lower σ e than perpendicular scans. The σ e calculated for range images (R-den) tends to increase with the range and incidence angle. Exceptions are the high σ e for sample 0936 (configuration 8), and the low σ e for sample 1309 (configuration 11). Possible reasons for the 0936 anomaly   The results of Figs. 10 and 11 suggest that noise estimated from wavelet components of range images is more reliable than for surface images.

Optimal Wavelet-Denoising Procedure
To evaluate how the threshold selection methods, wavelet transform methods and denoising directions influence wavelet-denoising results; the errors for denoised surfaces are shown in Fig. 12 and discussed below.
Threshold selection method Low thresholds preserve most of the coefficients (conservative thresholding) and result in a rougher denoised surface. High thresholds, on the opposite, remove more coefficients, thus return a smoother surface. The error depends on the actual surface roughness and scanning configuration. With the Grasselli parameter overestimated up to 350%, the penalised low threshold is inefficient in case of sample 0727. This threshold is better suited for rougher sample 1309, where the error varies approximately between 70% and − 50%. The local fixed form and global penalised high thresholds return similar errors, which are the smallest for all the scanning configurations in case of smoother samples 0727 and 0934. However, those thresholds underestimate the actual surface roughness of rougher samples 0936 and 1309, especially for longer ranges.
Possible reasons for underestimated Grasselli parameters for rougher samples include the removal of some surface details due to high thresholds, and increase in the laser beam footprint size with scanning range. A solution for the first cause is to apply more conservative thresholds like penalised medium or penalised low (Bitenc et al. 2015b).
DWT versus SWT Denoised surfaces using the SWT are smoother than for the DWT. The reason is that the SWT reconstruction process averages slightly shifted DWTs. Because the SWT error is mostly smaller compared to the DWT for short ranges and their error differences are relatively small for longer ranges, the SWT is a preferred wavelet transform for its advantages as stationary invariance and reduced edge effect.
Range versus surface denoising The errors for range and surface denoising results are shown in Fig. 13. The results show a considerable amount of noise removal when compared to the error of noisy TLS surfaces shown in Fig. 9. The Grasselli parameter of denoised surfaces for samples 0727, 0934, and 0936 is over-or underestimated by maximum 39%, and for the roughest sample 1309, the parameter is underestimated up to 57%. For most samples and scanning configurations, range-denoised surfaces appear smoother than surface denoised, as a consequence of larger σ e and, thus, higher thresholds.

Summary and Conclusions
In this paper, advanced wavelet-denoising methods were applied to remove TLS range noise, to capture surface morphology details for quantifying rock surface roughness. The TLS data were acquired at different ranges and incidence angles, for four rock samples having different roughness and reflectivity. Two wavelet transforms, DWT and SWT, were used in combination with four threshold selection methods. Denoising was executed in the direction of range measurements (range denoising) and in direction of surface roughness measurements perpendicular to the mean plane (surface denoising). By systematically comparing reference ATOS surfaces to noisy and denoised TLS surfaces, the influence of noise and noise estimations from the TLS data has been quantified, and the success of wavelet-denoising procedures has been demonstrated. The analyses have shown roughness over-estimation due to the TLS noise, especially for smoother natural rock surfaces. By applying wavelet-denoising procedures, TLS data were substantially improved and more reliable estimates of rock surface roughness were obtained. However, the success depends on the threshold value defined by the threshold selection method and noise estimation.
The results of this study suggest that the optimal threshold selection method should be chosen based on surface roughness properties. High thresholds, including fixed-form or penalised high, successfully eliminated the high noise effect for smoother surfaces. More conservative thresholds (removing less coefficients), including penalised low, have shown to be more appropriate for rougher surfaces, where the noise is mixed with surface details.
The TLS range noise is not precisely known a priori and depends on surface reflectivity and scanning geometry. In this research, the noise (σ e ) was estimated with the Median Absolute Deviation of the first-level detail coefficients obtained from the DWT and SWT of range and surface images. The results show that σ e values calculated from both wavelet transforms are similar and are, especially when using range images, highly correlated with the raw intensity values. The observed inversely proportional functional relationship coincides with the intensity-based stochastical model presented in the previous research (Wujanz et al. 2017). However, the σ e depends on the image construction method. The grid (pixel) size should be equal to or bigger than the TLS sampling interval, and the Nearest-Neighbour (NN) interpolation method is preferred to the mean. Our initial studies showed that adopting the mean value within a radius results in lower σ e as compared to NN.
The comparative analysis of denoised surfaces using DWT and SWT shows that the SWT surfaces are smoother. This finding is attributed to the SWT reconstruction algorithm, which involves the averaging of slightly shifted DWTs. Since the differences of surface rouhgness are small, the SWT is preferred due to its stationary invariance and ability to reduce edge effects.
Finally, the main finding of this research is that range denoising outperforms surface denoising, because it returns more reliable σ e for arbitrary scanning configurations. If the TLS polar coordinates are available, denoising should preferably be performed in the range direction.
The Grasselli parameter is highly sensitive to TLS noise. Fractal representations of surface features have also been found to be sensitive to noise (Khoshelham et al. 2011), as well as roughness values obtained through virtual compass and disc-clinometer measurements (Bitenc et al. 2015c). Further research is warranted to assess the relative sensitivity of alternate roughness parameters to TLS noise.
To establish the optimal scanning configuration for discerning the rock surface details, the combined influence of TLS noise and the laser beam footprint size are deserving of future research. This is because the influence of noise is attenuated by the smoothing effect of the footprint size, particularly for rough surfaces combined with increasing scanning range and incidence angle. Constructing a full 2.5D image from a scattered 3D point cloud remains a challenge in applying image denoising methods to TLS data. A possible solution is to utilize diffusion wavelets (Coifman and Maggioni 2006), which are an extension of the classical Fig. 12 Error (%) of range (a-d) and surface (e-h) denoised surfaces using DWT (solid line with circles) and SWT (dashed line with triangles) in combination with the four thresholds (color-coded lines), local fixed form (fl), and global penalised high (ph), medium (pm) and low (pl), versus scanning configuration for sample a, e 0727, b, f 0934, c, g 0936, and d, h 1309 ◂ wavelet transform and have been shown to function well for discrete structures such as point.