Study of Texture Indicators Applied to Pavement Wear Analysis Based on 3D Image Technology

Pavement texture characteristics can reflect early performance decay, skid resistance, and other information. However, most statistical texture indicators cannot express this difference. This study adopts 3D image camera equipment to collect texture data from laboratory asphalt mixture specimens and actual pavement. A pre-processing method was carried out, including data standardisation, slope correction, missing value and outlier processing, and envelope processing. Then the texture data were calculated based on texture separation, texture power spectrum, grey level co-occurrence matrix, and fractal theory to acquire six leading texture indicators and eight extended indicators. The Pearson correlation coefficient was used to analyse the correlation of different texture indicators. The distinction vector based on the information entropy is calculated to analyse the distinction of the indicators. High correlations between ENE (energy) and ENT (entropy), ENT and D (Minkowski dimension) were found. The CON (contrast) has low correlations with HT (macro-texture power spectrum area), ENT and D. However, the differentiation of ENE and HT is more prominent, and the differentiation of the CON is smaller. ENE, ENT, CON and D indicators based on macro-texture and the corresponding original texture have strong linear correlations. However, the microtexture indicators are not linearly correlated with the corresponding original texture indicators. D, WT (micro-texture power spectrum area) and ENT exhibit high degrees of numerical concentration for the same road sections and may be more statistically helpful in distinguishing the characteristics of the pavement performance decay of the road sections.


Introduction
Many scholars have tried to establish the relationships between pavement texture and pavement performance, such as surface ageing [1,2], skid resistance [1,3,4] and pavement noise generation [5,6]. These pieces of research show that the influence of pavement texture on pavement performance is very prominent. In addition, researchers believe texture features can be used in the design of pavement performance rather than just evaluation. Chu and Fwa [4] pointed out that different surface designs may be required for pavement skid resistance for different vehicle speeds and geometric alignments. Ktari et al. [7] found that the interface texture parameters significantly affect shear strength. Pratico and Astolfi [8] verified that macrotexture (MTD) and pendulum test values (PTV) could be explained and predicted based on simple physical and geometric models. Zou et al. [9] proved that pavement microtexture has a more significant effect on SFC than macrotexture.
Fwa [10] demonstrated the conventional single-point representation of the skid resistance threshold and the two-parameter skid resistance models of Penn State and PIARC could not provide the needed information for wet weather driving risk assessment. Though evaluating asphalt concrete (AC) and gypsum asphalt (SMA) polished with real car tires, Xiao et al. [11] found the texture roughness (material area ratio and two-dimensional power spectral density) varied asynchronously at all observed scales (0.1 to 80 mm).
Above all, texture feathers may have more possibilities to reflect pavement performance decay. We need to find a way to evaluate texture indicators, not just correlated with anti-slip or noise. Spatial pattern analysis of the digital images (fractal dimension) distinguished between light, moderate and heavy degradation but failed to make a fine distinction between moderate levels of wear [12]. Chen et al. [13] developed a cost-effective and relatively accurate image-based texture analysis method (ITAM) based on digital image processing and spectral analysis techniques. Mean Texture Depth (MTD) calculated with an image-based procedure can be used instead of the MTD measured with the sand patch method, as the coefficient of determination is 0.99, and the standard error is 0.061 mm [14]. It was found that the parameter ω (the ratio between the volume of the pavement model and the current cutting depth) can be used in combination with traditional characterisation parameters to evaluate the macro-texture quality of pavements [15]. 3D-ITAM was validated as an effective method to characterise asphalt mixture surface macro-and microtexture properties by comparing the results from 3D-ITAM with those from SPM and HFT methods [16]. Medeiros Jr. et al. [17] further confirmed the accuracy of pavement macro-texture parameters from the digital image processing technique. Three-dimensional models constructed by stereo infrared scanning distinguished between all but the smoothest texture classes [12].
The development of 3D scanning and image technology has provided a more convenient way to obtain 3D texture data. Texture changes are more pronounced than SFC, MTD, and other test indicators. Extracting richer and more comprehensive texture feature information can help carry out road performance-related studies such as interface relationship, performance decay, and mix design of asphalt pavements based on the depth and multidimensional profiling of the apparent texture.
After data acquisition and indicators construction, we need to find a suitable method to evaluate the rationality of the indicators. Pearson correlation analysis has been used to find the relation between texture and other performance indicators [18] and exclude indicators with strong correlations among texture indicators to eliminate potential multicollinearity in model development [19]. We processed texture data by four methods: texture separation, texture power spectrum, grey level co-occurrence matrix, and fractal theory. The relevant feature indicators are compared and analysed to obtain less correlated, more distinguishable, and more suitable indicators for texture separation and wear study.

Data Acquisition
This study used the Gocator high-speed 3D laser contour sensor with a maximum field of view (FOV) of 365 mm and 1280 laser line contour points to collect texture data. In this test, the resolution in the x-direction was 0.171 mm, the resolution in the y-direction was 0.116 mm, the resolution in the z-direction was 0.013 mm, and the set motion speed in the y-direction was 50 mm/s. We selected three road sections with different wear levels for 50 m each and tested every 1 m. Furthermore, rutting plate specimens of AC-13 and SMA-13, the commonly used mix types for road surface layers, were prepared for texture indicator analysis. The pre-processing of texture data was completed by data standardisation (DS), slope correction (SC), missing value (MP) and outlier processing (OP), and envelope processing (EP), which is shown in Figure 1.

Data Standardisation (DS)
The samples need to be the same size to eliminate the influence of the data size on the conclusion. For the measured pavement planes, their sizes need to be unified as a picture with width multiplied by a size equal to 1123 × 1024 pixels and a texture information matrix with rows multiplied by columns equal to 1024 × 1123. The measured rutting plate and Marshall specimens were cut around to obtain an inner joint rectangle with width multiplied by 451 pixels × 359 pixels and a texture information matrix with rows multiplied by columns equal to 359 × 451.

Slope Correction (SC)
The texture data need to be corrected if the measured pavement has a certain slope concerning the horizontal surface. The road surface may have a horizontal and vertical slope, and the horizontal slope is generally between 1.0% and 2.0%. The maximum vertical slope does not exceed 8%, and this limit decreases as the design speed increases, so the vertical slope correction is generally performed. It is necessary to fit a plane to the texture data first, which can be achieved using the least-squares method or the more accurate Random Sample Consensus (RANSAC) algorithm.
The RANSAC algorithm randomly selects points for plane fitting, and the plane with the highest frequency of occurrence is the final fitted plane, so its accuracy is higher [20]. The RANSAC algorithm is used for the immediacy plane fitting for the texture data obtained from the above measurements. Then its plane equation can be obtained as follows.
= + + z Ax By C (1) where z: information about the elevation of a point of the measured pavement, x, y: information about the plane position of a point of the measured pavement, A, B, C: constants, coefficients of the plane equation.
Since the maximum longitudinal slope does not exceed 8%, an approximate calibration can be made by processing only the elevation data to achieve the slope correction. That is, the calibrated elevation is: where i z is the elevation of a point before the calibration, and ' i z is its elevation after the calibration. The pavement before and after slope correction can be obtained using MATLAB for programming calculation, as shown in Figure 2a

Data Standardisation (DS)
The samples need to be the same size to eliminate the influence of the data size on the conclusion. For the measured pavement planes, their sizes need to be unified as a picture with width multiplied by a size equal to 1123 × 1024 pixels and a texture information matrix with rows multiplied by columns equal to 1024 × 1123. The measured rutting plate and Marshall specimens were cut around to obtain an inner joint rectangle with width multiplied by 451 pixels × 359 pixels and a texture information matrix with rows multiplied by columns equal to 359 × 451.

Slope Correction (SC)
The texture data need to be corrected if the measured pavement has a certain slope concerning the horizontal surface. The road surface may have a horizontal and vertical slope, and the horizontal slope is generally between 1.0% and 2.0%. The maximum vertical slope does not exceed 8%, and this limit decreases as the design speed increases, so the vertical slope correction is generally performed. It is necessary to fit a plane to the texture data first, which can be achieved using the least-squares method or the more accurate Random Sample Consensus (RANSAC) algorithm.
The RANSAC algorithm randomly selects points for plane fitting, and the plane with the highest frequency of occurrence is the final fitted plane, so its accuracy is higher [20]. The RANSAC algorithm is used for the immediacy plane fitting for the texture data obtained from the above measurements. Then its plane equation can be obtained as follows.
where z: information about the elevation of a point of the measured pavement, x, y: information about the plane position of a point of the measured pavement, A, B, C: constants, coefficients of the plane equation.
Since the maximum longitudinal slope does not exceed 8%, an approximate calibration can be made by processing only the elevation data to achieve the slope correction. That is, the calibrated elevation is: where z i is the elevation of a point before the calibration, and z i is its elevation after the calibration.
The pavement before and after slope correction can be obtained using MATLAB for programming calculation, as shown in Figure 2a

Data Standardisation (DS)
The samples need to be the same size to eliminate the influence of the data size on the conclusion. For the measured pavement planes, their sizes need to be unified as a picture with width multiplied by a size equal to 1123 × 1024 pixels and a texture information matrix with rows multiplied by columns equal to 1024 × 1123. The measured rutting plate and Marshall specimens were cut around to obtain an inner joint rectangle with width multiplied by 451 pixels × 359 pixels and a texture information matrix with rows multiplied by columns equal to 359 × 451.

Slope Correction (SC)
The texture data need to be corrected if the measured pavement has a certain slope concerning the horizontal surface. The road surface may have a horizontal and vertical slope, and the horizontal slope is generally between 1.0% and 2.0%. The maximum vertical slope does not exceed 8%, and this limit decreases as the design speed increases, so the vertical slope correction is generally performed. It is necessary to fit a plane to the texture data first, which can be achieved using the least-squares method or the more accurate Random Sample Consensus (RANSAC) algorithm.
The RANSAC algorithm randomly selects points for plane fitting, and the plane with the highest frequency of occurrence is the final fitted plane, so its accuracy is higher [20]. The RANSAC algorithm is used for the immediacy plane fitting for the texture data obtained from the above measurements. Then its plane equation can be obtained as follows.
where z: information about the elevation of a point of the measured pavement, x, y: information about the plane position of a point of the measured pavement, A, B, C: constants, coefficients of the plane equation.
Since the maximum longitudinal slope does not exceed 8%, an approximate calibration can be made by processing only the elevation data to achieve the slope correction. That is, the calibrated elevation is: where i z is the elevation of a point before the calibration, and ' i z is its elevation after the calibration. The pavement before and after slope correction can be obtained using MATLAB for programming calculation, as shown in Figure 2a

Missing Value Processing (MP)
Some missing values may be generated for various reasons when measuring texture data. Suppose it is caused by the objective existence of gaps on the road surface and other reasons. The fixed value filling can be replaced with a value much lower than the lowest point of the texture data. If it is caused by the failure of the instrument or human operation error, the interpolation filling method can be used. In the case of more missing values at the edges, the adjacent rows and columns of data can be deleted to reduce the impact of large area interpolation on the overall texture analysis processing. MATLAB R2020a is used to detect the location of the missing values first, and then the third spline interpolation method is used to process them, and their local before-and-after comparisons are shown in Figure 3a,b.

Missing Value Processing (MP)
Some missing values may be generated for various reasons when measuring texture data. Suppose it is caused by the objective existence of gaps on the road surface and other reasons. The fixed value filling can be replaced with a value much lower than the lowest point of the texture data. If it is caused by the failure of the instrument or human operation error, the interpolation filling method can be used. In the case of more missing values at the edges, the adjacent rows and columns of data can be deleted to reduce the impact of large area interpolation on the overall texture analysis processing. MATLAB R2020a is used to detect the location of the missing values first, and then the third spline interpolation method is used to process them, and their local before-and-after comparisons are shown in Figure 3a  The MP fills a small number of data gaps caused by the inability of the instrument to measure part of the seam due to the reflection of the laser, making the texture data more complete and the subsequent analysis more reasonable and accurate.

Outlier Processing (OP)
Outliers are one of the possible situations that can also occur in the measured texture data and can be divided into two types, pseudo-anomalies and true anomalies. The former can be judged in the field based on their measurement coordinates and reflect the actual condition of the pavement and usually do not need to be dealt with, while the latter is caused by instrument failure or human operation and need to be changed. Otherwise, the final results may be affected. Suppose the number of true anomalies is too high. In that case, it is necessary to consider whether the experiment has made an error to decide whether to discard the data and re-measure it.
Outliers can be judged by the median absolute deviation (MAD) algorithm, distance, clustering and density. Among them, the absolute median deviation algorithm is a commonly used method. For a data set { }  The MP fills a small number of data gaps caused by the inability of the instrument to measure part of the seam due to the reflection of the laser, making the texture data more complete and the subsequent analysis more reasonable and accurate.

Outlier Processing (OP)
Outliers are one of the possible situations that can also occur in the measured texture data and can be divided into two types, pseudo-anomalies and true anomalies. The former can be judged in the field based on their measurement coordinates and reflect the actual condition of the pavement and usually do not need to be dealt with, while the latter is caused by instrument failure or human operation and need to be changed. Otherwise, the final results may be affected. Suppose the number of true anomalies is too high. In that case, it is necessary to consider whether the experiment has made an error to decide whether to discard the data and re-measure it.
Outliers can be judged by the median absolute deviation (MAD) algorithm, distance, clustering and density. Among them, the absolute median deviation algorithm is a commonly used method. For a data set {x 1 , x 2 , · · · , x n }, first, calculate its median x m , then calculate the total deviation of the data points from the median d i = |x i − x m |, i = 1, 2, · · · , n, then calculate the median of the absolute deviation, d m and finally obtain the distance of all data points from the centre l i = d i d m . Usually, d m can be used as a consistent estimate of the standard deviation σ. The relationship between the two is σ = k·d m , where k is a constant factor with different values depending on the data distribution, k is usually equal to 1.4286 for normal distribution. The mean value µ is calculated, and then the values of the µ + 3σ·k·d m and µ − 3σ·k·d m are calculated to determine the range of statistically significant standard data points. Data points outside the interval [µ − 3σ·k·d m , µ + 3σ·k·d m ] are considered outliers [21].
The absolute median difference algorithm is a distance value method that is robust against outlier data and amplifies the effect of outliers. A method similar to that used to fill in missing values can be used when correcting outliers. The MAD algorithm is used to determine its position for the data here, and the interpolation method is used to process its local outliers. Its local OP results before and after comparison are shown in Figure 4a,b.  3 ] m m μ σ k d μ σ k d are considered outliers [21].
The absolute median difference algorithm is a distance value method that is robust against outlier data and amplifies the effect of outliers. A method similar to that used to fill in missing values can be used when correcting outliers. The MAD algorithm is used to determine its position for the data here, and the interpolation method is used to process its local outliers. Its local OP results before and after comparison are shown in Figure 4a,b.

Envelope Processing (EP)
The tire is more likely to contact the protruding part of the pavement, while the depressed part does not. Therefore, the envelope of texture data is needed to better calculate the real contact surface between the tire and the road surface to establish a connection with the road performance indicators. There are three commonly used methods for envelope calculation. The first one is to filter the texture data with a certain length of Hilbert FIR filter after Discrete Fourier Transform (DFT) to obtain the envelope. The second one is to return the root-mean-square envelope by a sliding window of a certain length. The third one is to find the wave peaks, then take out some peaks separated by a certain number of waves and obtain the envelope by applying a spline interpolation method. A suitable envelope is obtained by changing some parameters in these algorithms and then crossvalidating them with the pavement evaluation parameters [22]. For the data here, each wave crest was collected and then interpolated with the spline interpolation method and combined with the pavement evaluation parameters to obtain the envelope. The result of the profile's envelope based on partial data is shown in Figure 5.

Envelope Processing (EP)
The tire is more likely to contact the protruding part of the pavement, while the depressed part does not. Therefore, the envelope of texture data is needed to better calculate the real contact surface between the tire and the road surface to establish a connection with the road performance indicators. There are three commonly used methods for envelope calculation. The first one is to filter the texture data with a certain length of Hilbert FIR filter after Discrete Fourier Transform (DFT) to obtain the envelope. The second one is to return the root-mean-square envelope by a sliding window of a certain length. The third one is to find the wave peaks, then take out some peaks separated by a certain number of waves and obtain the envelope by applying a spline interpolation method. A suitable envelope is obtained by changing some parameters in these algorithms and then cross-validating them with the pavement evaluation parameters [22]. For the data here, each wave crest was collected and then interpolated with the spline interpolation method and combined with the pavement evaluation parameters to obtain the envelope. The result of the profile's envelope based on partial data is shown in Figure 5.
The absolute median difference algorithm is a distance value method that is robu against outlier data and amplifies the effect of outliers. A method similar to that used fill in missing values can be used when correcting outliers. The MAD algorithm is used determine its position for the data here, and the interpolation method is used to proce its local outliers. Its local OP results before and after comparison are shown in Figure 4a

Envelope Processing (EP)
The tire is more likely to contact the protruding part of the pavement, while the d pressed part does not. Therefore, the envelope of texture data is needed to better calcula the real contact surface between the tire and the road surface to establish a connectio with the road performance indicators. There are three commonly used methods for env lope calculation. The first one is to filter the texture data with a certain length of Hilbe FIR filter after Discrete Fourier Transform (DFT) to obtain the envelope. The second o is to return the root-mean-square envelope by a sliding window of a certain length. T third one is to find the wave peaks, then take out some peaks separated by a certain num ber of waves and obtain the envelope by applying a spline interpolation method. A suit ble envelope is obtained by changing some parameters in these algorithms and then cros validating them with the pavement evaluation parameters [22]. For the data here, ea wave crest was collected and then interpolated with the spline interpolation method an combined with the pavement evaluation parameters to obtain the envelope. The result the profile's envelope based on partial data is shown in Figure 5.  The envelope-processed pavement texture section better characterises the tire-pavement contact, making the subsequent analysis more accurate and meaningful for wear analysis. By considering the distance between the texture profile and the reference surface as a smooth random function, it is possible to define the distance between two recurring constructions as the wavelength of the texture. Different texture scales can affect pavement performance indicators differently, so it is necessary to classify the textures. Pavement textures are usually divided into three categories by wavelength, including micro-texture (0~0.5 mm), macro-texture (0.5~50 mm), and mega-texture (50~500 mm) [6]. The Fourier analysis method can decompose the texture into sinusoidal components of different wavelengths and amplitudes, which occupy different proportions of the total texture data. Therefore, the texture data indicators of different wavelengths can be calculated to reflect the effects of different wavelength textures and their distribution characteristics in the total texture data.
For a set of pavement texture data, each column (i.e., its y-direction data) is subjected to a Fast Fourier Transform (FFT), which transforms it from the spatial domain to the frequency domain with frequency units of mm −1 . Then a bandpass filter is designed to filter the texture data. The wavelength of micro-texture is 0~0.5 mm, and the wavelength of 0.5 mm corresponds to a frequency of 2 mm −1 . While the resolution in the z-direction at the measurement time is 0.013 mm, the maximum frequency after FFT corresponding to this size is 38.46 mm −1 . Therefore, a bandpass filter of 2~38.46 mm −1 can be designed to filter out the micro-texture. The wavelength of macro-texture is 0.5~50 mm, and the frequency corresponding to the 50 mm wavelength is 0.02 mm −1 , so a bandpass filter of 0.02~2 mm −1 could be used for the macro-texture. However, since the components with frequencies below 0.02 mm −1 are small and negligible, the lower cutoff frequency of 0.02 mm −1 is difficult to achieve in bandpass filters from 0.02 to 2 mm −1 . A low-pass filter with a frequency of 2 mm −1 is used to approximate the bandpass filter of 0.02~2 mm −1 , which is more effective and avoids large distortion. In addition, according to Nyquist's sampling theorem, the selected sampling rate should be greater than twice the highest frequency of the original signal, which can avoid the occurrence of spectral aliasing and thus restore the original signal without distortion. Finally, the Fourier inversion can reduce the frequency domain's texture data to the spatial domain's texture data. The computational flow chart is shown in Figure 6.
The envelope-processed pavement texture section better characterises the tire-pavement contact, making the subsequent analysis more accurate and meaningful for wear analysis.

Texture Separation
By considering the distance between the texture profile and the reference surface as a smooth random function, it is possible to define the distance between two recurring constructions as the wavelength of the texture. Different texture scales can affect pavement performance indicators differently, so it is necessary to classify the textures. Pavement textures are usually divided into three categories by wavelength, including microtexture (0~0.5 mm), macro-texture (0.5~50 mm), and mega-texture (50~500 mm) [6]. The Fourier analysis method can decompose the texture into sinusoidal components of different wavelengths and amplitudes, which occupy different proportions of the total texture data. Therefore, the texture data indicators of different wavelengths can be calculated to reflect the effects of different wavelength textures and their distribution characteristics in the total texture data.
For a set of pavement texture data, each column (i.e., its y-direction data) is subjected to a Fast Fourier Transform (FFT), which transforms it from the spatial domain to the frequency domain with frequency units of mm −1 . Then a bandpass filter is designed to filter the texture data. The wavelength of micro-texture is 0~0.5 mm, and the wavelength of 0.5 mm corresponds to a frequency of 2 mm −1 . While the resolution in the z-direction at the measurement time is 0.013 mm, the maximum frequency after FFT corresponding to this size is 38.46 mm −1 . Therefore, a bandpass filter of 2~38.46 mm −1 can be designed to filter out the micro-texture. The wavelength of macro-texture is 0.5~50 mm, and the frequency corresponding to the 50 mm wavelength is 0.02 mm −1 , so a bandpass filter of 0.02~2 mm −1 could be used for the macro-texture. However, since the components with frequencies below 0.02 mm −1 are small and negligible, the lower cutoff frequency of 0.02 mm −1 is difficult to achieve in bandpass filters from 0.02 to 2 mm −1 . A low-pass filter with a frequency of 2 mm −1 is used to approximate the bandpass filter of 0.02~2 mm −1 , which is more effective and avoids large distortion. In addition, according to Nyquist's sampling theorem, the selected sampling rate should be greater than twice the highest frequency of the original signal, which can avoid the occurrence of spectral aliasing and thus restore the original signal without distortion. Finally, the Fourier inversion can reduce the frequency domain's texture data to the spatial domain's texture data. The computational flow chart is shown in Figure 6.

Texture Power Spectrum Density
The power spectral density can be used to assess the level of macro-texture and micro-texture of pavement. First, the macro-texture and micro-texture power spectra are calculated separately using the PWelch function in Matlab. Then, the PWelch function uses the input signal power spectral density (PSD) estimates found by Welch's overlapping segment averaging estimator. For the pavement, the window function is specified as the Hamming window, and the length of the Hamming window is 1024. The number of overlapping samples, defined as a positive integer less than the length of the window, is taken here as 512. The number of points of DFT is 256, and the length of the segments is 2048. The sampling rate is defined as the number of samples per unit, taken here as 100 mm −1 .

Pre-processing Fast Fourier Transform
Band-pass filtering Fourier inversion

Texture Power Spectrum Density
The power spectral density can be used to assess the level of macro-texture and micro-texture of pavement. First, the macro-texture and micro-texture power spectra are calculated separately using the PWelch function in Matlab. Then, the PWelch function uses the input signal power spectral density (PSD) estimates found by Welch's overlapping segment averaging estimator. For the pavement, the window function is specified as the Hamming window, and the length of the Hamming window is 1024. The number of overlapping samples, defined as a positive integer less than the length of the window, is taken here as 512. The number of points of DFT is 256, and the length of the segments is 2048. The sampling rate is defined as the number of samples per unit, taken here as 100 mm −1 . For the specimens, the Hamming window length is 256, the number of overlapping samples is 128, the number of DFT points is 512, and the sampling rate is 100 mm −1 . The area under the power spectral density curve for the frequencies in its range is used as the characteristic indicators of macro-texture and micro-texture for this pavement section (i.e., this column of data in the Y-direction). Finally, the data of all road sections are averaged to obtain this pavement's power spectral density indicators. This method avoids the influence of the difference in the number of texture data columns on the conclusion.

Grey Level Co-Occurrence Matrix
Grey level co-occurrence matrix (GLCM) is a processing method for texture data used in texture feature representation in high-precision manufacturing. It represents the joint probability density distribution of the simultaneous occurrence of two pixels separated by a certain distance in the whole image according to a certain translation direction. The texture features can be extracted by analysing the spatial relationship of the grey values of pixels separated by a certain distance in space and rewriting them with texture descriptors [23,24].
For the texture image M × M, take any point A(x, y), transform its horizontal and vertical coordinates to obtain A (x + p, y + q), and set the grey value of the point pair A, A as (g 1 , g 2 ). Different grayscale values (g 1 , g 2 ) can be obtained for different points. If the number of grayscale levels is G, then there are several possible grayscale values G 2 . The occurrences are counted, normalised to frequency, and arranged as a square matrix to obtain the grayscale co-occurrence matrix P for each grayscale value. The coordinate transformation value (p, q) needs to be chosen according to the texture characteristics. A smaller value (1, 0) representing a horizontal pixel pair, i.e., a 0 • scan, is chosen for finer textures. It is also reasonable to choose other values such as (1, 1), (0, 1), (−1, −1), whose scanning angle will also change to 45 • , 90 • , 135 • , etc. Haralick [25] proposed 14 statistics such as energy, entropy, contrast, uniformity, correlation, etc. Contrast, entropy, and energy can be chosen to describe the texture features. These indicators can be computed based on the original or separated texture.
Energy represents the thickness and uniformity of texture distribution, and its calculation method is shown in (3).
Entropy reflects the complexity of the distribution, and its calculation method is shown in (4).
The contrast reflects the sharpness of the image and the groove depth of the texture. Its calculation method is shown in (5).
where P(i, j) is the element of the i th row and j th column of the grayscale co-generation matrix.

Fractal Theory
Mandelbrot has defined fractals twice [26,27]. In 1986, Mandelbrot defined a fractal as a shape composed of parts similar to the whole. This definition emphasises the selfsimilarity of fractal sets [28]. Fractal theory studies its properties and applications. It has two principles, the principle of self-similarity and iterative generation, characterising fractals under the usual geometric transformations independent of scale. Self-similar shapes can be identical or similar statistically, which implies recursion. The value of a graph with no fractal structure is exactly the spatial dimension of this structure, but it may not be an integer for a fractal pattern with infinite details [29,30].
The texture of a certain section of a road surface also has self-similarity, and the fractal dimension can be used to evaluate the texture condition of the road surface. There are usually two types of fractal dimensions, the Hausdorff dimension and the Minkowski dimension, the latter also known as the box-counting dimension. The Hausdorff dimension is the most commonly used dimension, and its expression is shown in (6).
The equation is equivalent to (7).
Taking the natural logarithm and sorting gives (8).
where L is the number of times a certain shape expands along each independent direction, and K is the number of times the resulting fractal occupies a range of space relative to the original shape. The idea of the Minkowski dimension is to cover the set using disjoint boxes. Firstly, each point's height is represented by a grey value with a grey level G. Secondly, the set is divided into several rectangles with the size of R × R × R , where R = R· G M , ensuring the number of rectangles in each direction is equal in the 3D coordinates. In the position (i, j) of the grid R × R, find the maximum grey value u, the minimum grey value, and then the minimum number of squares covering the minimum to the maximum grey value is (9).
Then summing over the entire grid of R × R is (10).
The fractal dimension of the whole texture image is shown in (11) [31]. Similar to the previous indicators, this indicator can be calculated based on the original and separated textures.

Indicators Evaluation Method
Reasonable indicators should have weak linear correlations to eliminate the effect of subsequent model building of multicollinearity, which the Pearson correlation coefficient matrix can measure. Furthermore, they need to have a significant degree of distinction to amplify the differences between indicators of different textures, which can be calculated by distinction vectors based on information entropy.
The data sets are first dimensionless to obtain the dimensionless matrix of the texture indicator set X = x ij n × m (n is the set number, and m is the indicator number). The Pearson correlation coefficients of each two indicators are calculated by (12).
where s pq is the element of the correlation matrix S. The meanings of p and q are the p th and the q th columns of the dimensionless matrix X. The cov(p, q) is the covariance between p and q, and σ p , σ q are the standard deviations of p and q. Next, the entropy of the i th texture indicator is calculated by (13).
Then the differentiation of the i th texture indicator is where g i is the element of the distinction vector G.

Texture Separation
The macro-texture and micro-texture filtering results of one of the pavement sections are shown in Figures 7 and 8. As seen in Figure 7a, the macro-texture is obtained by filtering the texture using low-pass filter with a pass frequency of 2 mm −1 , which shows the overall undulating tre of the texture and ignores its subtle local variations. The overlapping parts of the pow spectrum curves before and after filtering are shown in Figure 7b. Similarly, as seen Figure 8a, the micro-texture is obtained by filtering the texture using a bandpass fil As seen in Figure 7a, the macro-texture is obtained by filtering the texture using a low-pass filter with a pass frequency of 2 mm −1 , which shows the overall undulating trend of the texture and ignores its subtle local variations. The overlapping parts of the power spectrum curves before and after filtering are shown in Figure 7b. Similarly, as seen in Figure 8a, the micro-texture is obtained by filtering the texture using a bandpass filter passing through the frequency interval from 2 to 38.46 mm −1 , exhibiting subtle local variations and ignoring the overall undulating trend of its texture. The overlapping parts of the power spectrum curves before and after filtering are shown in Figure 8b. Figures 7 and 8 show that the filter and sample rate settings do well in separating the macro-texture and micro-texture portions of the pavement texture. The 3D comparison figures before and after texture separation are shown in Figure 9, which demonstrate the effect of texture separation.
As can be seen in Figure 9b, the macro-texture after texture separation characterises the overall undulating trend of the original texture. In Figure 9c, the micro-texture shows only the fine features of the texture. Therefore, the original texture in Figure 9a is well separated.
passing through the frequency interval from 2 to 38.46 mm , exhibiting subtle local vari-ations and ignoring the overall undulating trend of its texture. The overlapping parts of the power spectrum curves before and after filtering are shown in Figure 8b. Figures 7 and 8 show that the filter and sample rate settings do well in separating the macro-texture and micro-texture portions of the pavement texture. The 3D comparison figures before and after texture separation are shown in Figure 9, which demonstrate the effect of texture separation. As can be seen in Figure 9b, the macro-texture after texture separation characterises the overall undulating trend of the original texture. In Figure 9c, the micro-texture shows only the fine features of the texture. Therefore, the original texture in Figure 9a is well separated.

Texture Power Spectrum Density
The area distributions under the power spectral density curve of macro-texture and micro-texture for different roads and specimens are shown in Figure 10.

Texture Power Spectrum Density
The area distributions under the power spectral density curve of macro-texture and micro-texture for different roads and specimens are shown in Figure 10. As shown in Figure 10, the indicators' distribution of roads 1 to 3 is more concentrated, while the same indicators' distribution of specimens is more discrete. For road sections, the HTs are more discrete than the WTs. In contrast, the specimens' WTs demonstrate higher discreteness than the HTs.
It seems that the micro-texture of the same road section tends to be consistent, which is more representative of the pavement wear state of the road section. The specimens' HTs may have certain representativeness, but the representativeness of specimens' WTs is insufficient.  As shown in Figure 10, the indicators' distribution of roads 1 to 3 is more concentrated, while the same indicators' distribution of specimens is more discrete. For road sections, the HTs are more discrete than the WTs. In contrast, the specimens' WTs demonstrate higher discreteness than the HTs.
It seems that the micro-texture of the same road section tends to be consistent, which is more representative of the pavement wear state of the road section. The specimens' HTs may have certain representativeness, but the representativeness of specimens' WTs is insufficient.

Grey Level Co-Occurrence Matrix (GLCM)
The distributions of the grey level co-occurrence matrix indicators for texture images of different pavement and specimen images are shown in Figure 11. Similar to the texture spectrum indicator, the GLCM indicators of roads 1 to 3 show similar concentration trends, while the specimens show dispersion. Among the three indicators, ENT seems to perform best.

Fractal Theory
The Minkowski dimension D for texture images of roads 1 to 3 and specimens were distributed in the range of 2.1 to 2.7, as shown in Figure 12.  Similar to the texture spectrum indicator, the GLCM indicators of roads 1 to 3 show similar concentration trends, while the specimens show dispersion. Among the three indicators, ENT seems to perform best.

Fractal Theory
The Minkowski dimension D for texture images of roads 1 to 3 and specimens were distributed in the range of 2.1 to 2.7, as shown in Figure 12. Figure 12 implies that D exhibits a high degree of numerical concentration for the same road section, which may be more statistically helpful in distinguishing the characteristics of the pavement performance decay. D can reflect the wear on different scales simultaneously, which is more representative. In addition, Ds have similar distribution characteristics to the previous indicators, indicating that texture indicators' general trends are consistent. dicators, ENT seems to perform best.

Fractal Theory
The Minkowski dimension D for texture images of roads 1 to 3 and specimens distributed in the range of 2.1 to 2.7, as shown in Figure 12.  Figure 12 implies that D exhibits a high degree of numerical concentration fo same road section, which may be more statistically helpful in distinguishin  Figure 12. Ds of different roads and specimens.

Pearson Correlation Coefficients Matrix
Pearson correlation coefficients of the six indicators (HT, WT, ENE, ENT, CON, and D) were calculated, and the results are plotted in Figure 13.

Pearson Correlation Coefficients Matrix
Pearson correlation coefficients of the six indicators (HT, WT, ENE, ENT, CON, and D) were calculated, and the results are plotted in Figure 13. As seen in Figure 13, ENE and ENT have a strong negative correlation, while D and ENT have a strong positive correlation. The texture information reflected by the two GLCM indicators (ENE and ENT) and the fractal indicator (D) has certain consistency, which means that these indicators have certain substitutability for each other. There are no obvious linear correlations between CON and other indicators, implicating that CON may be suitably used together with other indicators for multidimensional texture analysis.

Distinction Vector
According to the information entropy theory, the differentiation degrees of HT, WT, ENE, CON, ENT, and D are calculated to obtain the distinction vector G1. The result is shown in (15).
The distinction vector G1 indicates that the two indicators with the largest differentiation are HT and ENE, whose differentiation is more than twice as large as the differentiation of the other indicators. HT and ENE better reflect the gap between different textures and may be more suitable as texture indicators for wear analysis. The indicator with the smallest differentiation is CON, which may be not suitable for evaluating pavement wear, As seen in Figure 13, ENE and ENT have a strong negative correlation, while D and ENT have a strong positive correlation. The texture information reflected by the two GLCM indicators (ENE and ENT) and the fractal indicator (D) has certain consistency, which means that these indicators have certain substitutability for each other. There are no obvious linear correlations between CON and other indicators, implicating that CON may be suitably used together with other indicators for multidimensional texture analysis.

Distinction Vector
According to the information entropy theory, the differentiation degrees of HT, WT, ENE, CON, ENT, and D are calculated to obtain the distinction vector G 1 . The result is shown in (15). G 1 = (0.420 0.100 0.369 0.084 0.007 0.020) The distinction vector G 1 indicates that the two indicators with the largest differentiation are HT and ENE, whose differentiation is more than twice as large as the differentiation of the other indicators. HT and ENE better reflect the gap between different textures and may be more suitable as texture indicators for wear analysis. The indicator with the smallest differentiation is CON, which may be not suitable for evaluating pavement wear, as the values of the specific indicators evaluated will be closer. CON represents some unique texture statistical information of pavement texture (Section 4.2.1), which may be used together with other texture indicators to study other performances (such as noise, anti-sliding, etc.) related to pavement texture features.

The Effect of the Texture Separation on These Indicators
The grey level co-occurrence matrix and fractal theory can also be used to calculate indicators based on the filtered macro-texture and micro-texture data. The indicators before and after separation are linearly fitted, as shown in Figure 14. As seen in Figure 14, the macro-texture indicators HENE, HENT, HCON, and HD have solid linear correlations with the original texture indicators ENE, ENT, CON, and D, and the R 2 values of the linear fits are above 0.9. Except for WD, the micro-texture indicators are not linearly correlated with the original texture indicators. The macro-texture indicators seem to overlap more with the meanings expressed by the original texture indicators. The distinction vectors are calculated to analyse the effect of the texture separation on these indicators, as shown in (16). As seen in Figure 14, the macro-texture indicators HENE, HENT, HCON, and HD have solid linear correlations with the original texture indicators ENE, ENT, CON, and D, and the R 2 values of the linear fits are above 0.9. Except for WD, the micro-texture indicators are not linearly correlated with the original texture indicators. The macro-texture indicators seem to overlap more with the meanings expressed by the original texture indicators. The distinction vectors are calculated to analyse the effect of the texture separation on these indicators, as shown in (16). G 2~G5 denote the differentiation vectors of ENE, ENT, CON and D indicators based on original, macro-texture and micro-texture, respectively. The differentiation degrees of the micro-texture indicators are much smaller than those of the original and macro-texture. For micro-texture indicators, the ENT and D differentiation degrees are more significant than those of ENE and CON.
The differentiation degrees of the macro-texture are similar to those of the original textures in G 2~G4 but significantly higher in G 5 . As D has a greater differentiation than other indicators on the micro-scale, and its macro-texture has a higher differentiation than the overall one, its performance in wear characteristics seems to be more prominent.
Although not significant enough, ENT also shows the same trend characteristics as D. It is more reasonable and effective to perform texture separation for ENT and D applied to the wear analysis.

Conclusions
We scanned texture data from 3 road sections and 14 rutting plate specimens of the AC-13 mix. Based on the texture data achieved, six statistical indicators were calculated. The Pearson correlation coefficient was used to analyse the correlation of different texture indicators. The distinction vectors based on the information entropy were calculated to analyse the distinction of the indicators. The effect of texture separation on calculating the grey level co-occurrence matrix and fractal indicators was investigated. The following conclusions were obtained:

1.
D, WT and ENT exhibit high degrees of numerical concentration for the same road sections. The three indicators may be more statistically helpful in distinguishing the characteristics of pavement performance decay. D reflects different information from all the other statistical variables. The specimens do not possess all the actual pavements' statistical characteristics.

2.
ENE, ENT and D have a certain consistency, which means that these indicators have certain substitutability for each other. CON represents some unique texture statistical information.

3.
ENT and D have greater differentiation than other indicators on the micro-scale, and its macro-texture has a higher differentiation than the overall one. Their performance in wear characteristics seems to be more prominent.

4.
With high degrees of numerical concentration for the same road sections and differentiation degrees, D, PSD indicators, and ENT should be further studied based on more road sections or abraded specimens.
The fractal and texture-based energy indicators at different scales can provide better technical means for evaluating pavement wear and further be used to assess noise and skid resistance by multi-scale texture analysis. However, the practical engineering level comparison can hardly be carried out due to the coarse means of existing noise and skid resistance testing.
Further study can be carried out in conjunction with microscopic simulation modelling of noise or tire skid resistance performance to construct metrics suitable for early wear performance decay analysis of pavements.