A novel image-based approach for interactive characterization of rock fracture spacing in a tunnel face

Jiyo Chen, Yifeng Chen, Anthony G. Cohn, Hongwei Hung,*, Jinhong Mn,Lijun Wei

a Key Laboratory of Geotechnical and Underground Engineering of Ministry of Education, Department of Geotechnical Engineering, Tongji University, Shanghai,200092, China

b School of Computing, University of Leeds, Leeds, LS2 9JT, UK

c Department of Computer Science and Technology, Tongji University, Shanghai, 211985, China

d School of Civil Engineering, Shandong University, Jinan, 250061, China

e Luzhong Institute of Safety, Environmental Protection Engineering and Materials, and School of Mechanical and Electrical Engineering, Qingdao University of Science and Technology, Qingdao, 255000, China

Keywords:Fracture spacing Rock tunnel Image processing Fracture set grouping

ABSTRACT This paper presents a novel integrated method for interactive characterization of fracture spacing in rock tunnel sections. The main procedure includes four steps: (1) Automatic extraction of fracture traces, (2)digitization of trace maps, (3) disconnection and grouping of traces, and (4) interactive measurement of fracture set spacing, total spacing, and surface rock quality designation (S-RQD) value. To evaluate the performance of the proposed method, sample images were obtained by employing a photogrammetrybased scheme in tunnel faces. Experiments were then conducted to determine the optimal parameter values(i.e.distance threshold,angle threshold,and number of fracture trace grouping)for characterizing rock fracture spacing.By applying the identified optimal parameters involved in the model,the proposed method could lead to excellent qualitative results to a new tunnel face. To perform a quantitative analysis, three methods (i.e. field, straightening, and the proposed method) were employed in the same study and comparisons were made. The proposed method agrees well with the field measurement in terms of the maximum and average values of measured spacing distribution. Overall, the proposed method has reasonably good accuracy and interactive advantage for estimating the ultimate fracture spacing and S-RQD. It can be a possible extension of existing methods for fracture spacing characterization for two-dimensional (2D) rock tunnel faces.

1. Introduction

Rock fracture spacing, defined as the distance between two discontinuities, is a fundamental but important parameter to characterize the rock mass quality (Hudson and Priest,1979; Narr and Suppe, 1991; Wines and Lilly, 2002; Cai et al., 2021). This parameter is essentially used to ascertain the size of the basic blocks that make up the rock mass,since the stability of rock mass is strongly influenced by its constituent blocks(Priest and Hudson,1981; Moomivand et al., 2021). So far, fracture spacing has been widely applied to many rock mass classification schemes, such as rock mass index (RMi) (Palmstrøm, 1996), and rock mass rating(RMR)(Bieniawski,1988).The measurement of the exact spacing of rock joints in tunnel can help to estimate the overall fragmentation of the tunnel section, and thus to assess the stability and permeability before the next excavation phase (Xing et al., 2018). In this context, exploring an effective method to accurately measure and assess the fracture spacing of the excavation section is essential.

Prior to characterizing the fracture spacing, it is vital to collect samples and determine the specific approach of spacing measurement(Hudson and Priest,1979;Zhang and Einstein,1998;Hudson and Harrison, 2000; Azizi and Moomivand, 2021). The sampling methods are roughly categorized into contact and non-contact methods (Menegoni et al., 2019; Chen et al., 2020). The contact methods are usually carried out in the field with hand-held devices including the tapeline, laser range finder, and geological compass.They are simple and direct approaches, but are dangerous, timeconsuming, labor-intensive, and difficult to conduct in inaccessible areas(Franklin et al.,1988;Chen et al.,2021a,b;Zhang et al.,2021; Zhou et al., 2021). On the other hand, the non-contact methods, e.g. photogrammetry and three-dimensional (3D)remote sensor, are effective in characterizing the on-site rock fracture by acquiring image patterns or point clouds with considerable resolution(Riquelme et al.,2017;Fazio et al.,2019).The noncontact methods are gradually employed in engineering practice due to their high efficiency and safety(Wang et al.,2019;Wei et al.,2021).

The measurements of rock fracture spacing can be broadly divided into two types (Priest, 2012), i.e. fracture set and total spacings. Fracture set spacing is defined as the average distance between two adjacent fractures in the normal direction of the same group of rock fractures. It is scanned along a line with designative orientation and location. As a typical fracture set spacing, the normal set spacing requires the scanline to be perpendicular to the average direction of the measured fracture set. Total spacing is defined as the sum of the distance between all adjacent rock fractures along the drawn scanline.It is also scanned along a line with specified orientation and location. A complete characterization of rock fracture spacing requires the accurate and objective representation of both forms. Meanwhile, numerous scholars have focused on the mean spacing parameter to reflect the overall fragmentation of the whole rock mass (Azizi and Moomivand,2021; Moomivand et al., 2021). This work is a challenge since the mean spacing is heavily influenced by the deviation of both large and small spacings.It can vary significantly with different scanline orientations. Proposed models have to be calibrated for the mean spacing characterization of the whole section (Moomivand and Vandyousefi, 2020). In this study, we have analyzed the data characterized along the scanline (e.g. spacing distribution,maximum, medium, and mean spacing values), while the estimation of the mean spacing of the whole section is not our focus.Furthermore,rock quality designation(RQD),as another parameter of rock mass fragmentation, is defined as the core recovery of the wellbore, which includes only solid core fragments longer than 10 cm along the centerline of the rock core(Deere,1988).However,due to the limited construction time and high cost in sampling, a core hole cannot be drilled every time to measure the RQD value.Instead,engineers use the scanline to measure the surface RQD(SRQD) value (Boadu and Long,1994; Palmstrom, 2005). The S-RQD value of the excavated section can be roughly calculated from the spacing statistics obtained from the total spacing parameter(Esfahani and Asghari,2013).Thus,S-RQD value is evaluated as an additional index in this study in addition to the fracture set and total spacing indicators.

Quantitative assessment of fracture spacing is conceptually straightforward,but in practice it is a challenging task(Wines and Lilly,2002;Leng et al.,2021;Pi et al.,2021).This is not only due to the difficulty in collecting satisfactory data, but also in finding an effective method to achieve the automatic measurement of samples to obtain the ultimate spacing (Riquelme et al., 2015; Buyer and Schubert, 2017). Of course, the prerequisite of reasonable evaluation of fracture spacing is the accurate extraction of the fracture trace information from two-dimensional(2D)image or 3D point cloud samples (Hadjigeorgiou et al., 2003; Buyer and Schubert, 2017). However, due to the difficulty in post-processing of complex 3D point cloud data, the current fracture spacing measurements are still performed by projecting point clouds onto 2D slices (Li et al., 2016; Guo et al., 2019; Chen et al., 2021c).Although the 3D fracture traces, which contains real spatial information, can be extracted from the disordered point clouds,numerous researchers are still using 2D data to determine the fracture spacings(Leng et al.,2021;Stavropoulou et al.,2021).The main challenges in measuring fracture spacing using 2D samples can be summarized as follows:

(1) It is difficult to ensure that the pixel-based discontinuous traces are accurately extracted from field samples and effectively converted into digitized coordinate information(Deb et al., 2008; Battulwar et al., 2021; Chen et al., 2021d).

(2) Measurement of fracture spacing is generally incomplete,i.e.most studies focus on the partial spacing and ignore the other spacing parameters(Li et al.,2016;Zhang et al.,2020a).

(3) The extracted fracture traces are often characterized by the connecting lines of the start and end points(i.e.straightening method), and ignoring the natural curve shape of the trace(Healy et al., 2017; Zhang et al., 2020b).

(4) It is difficult to achieve fine automatic grouping of bending fracture traces (Azarafza et al., 2021).

These four challenges have received extensive research attention.For the first challenge,the representative extraction methods include manual mapping-, morphology-, and deep learning-based methods (Narr, 1996; Azarafza et al., 2019, 2021; Chen et al.,2021d), while representative digitization methods include manual labeling-, Hough transform- and chain code-based methods (Watkins et al., 2015; Healy et al., 2017; Chen et al.,2021d). Among these methods, the accuracy of manual labeling is the best approach, but the enormous work required for labeling might discourage engineers(Bao et al.,2019).Thus,a deep learning model integrated with a chain code-based method is chosen to extract the fracture trace map and process the coordinate information. To solve the second challenge, all parameters need to be integrated for a comprehensive statistical analysis of trace spacing,especially for the fracture set spacing. Note that extracting the accurate fracture set spacing from the image-based data is undoubtedly challenging as it is difficult to calculate the optimized trace grouping information(Li et al.,2016).For the third challenge,although the straightening method is highly simplified and has been widely accepted in engineering practice(Ge et al.,2018;Kong et al.,2020),the simplifications in the fracture characterization has inevitably affected the final statistical results due to the complex shape of the various fracture sets(Chen et al.,2021d).Thus,in cases where the accuracy of measurement is required, this method may no longer be applicable. It is urgent to find a trace map that can reasonably preserve natural curve shape (Zhang et al., 2018). For the fourth challenge, since the continuous traces in the form of pixels may not always belong to the same group,it is challenging to break up the traces and group them meaningfully according to the actual trace curve (Han et al., 2018; Leng et al., 2021). Ideally, to enable unobstructed measurement and statistical analysis of fracture spacing, it would be useful if the scanline of spacing can be adjusted according to engineers experience (Ruf et al., 1998).Admittedly,these methods may inevitably have certain restrictions for the difficult process from extracting fracture traces to measuring trace spacing parameters.It is crucial to develop a new method for automatic trace extraction, trace node digitization,disconnection, grouping and measurement.

In this study, a novel method for interactive measurement of rock fracture spacing using 2D image samples is proposed. The spacing measurement is a multi-step comprehensive process,where each step is realized by different algorithms, including a convolutional neural network(CNN)based model called FraSegNet(Chen et al., 2021d) for automatic trace extraction, a chain codebased polyline fitting algorithm for trace map digitization, an angle threshold-based algorithm for trace disconnection, a Kmeans++based algorithm(Arthur and Vassilvitskii,2007)for trace grouping, and an interactive measurement for estimating set spacing,total spacing and S-RQD value.A case study was conducted with images of tunnel sections taken at a tunnel project which was under construction in Yunnan, China. The spacings and S-RQD value of the contact method and the straightening method were used for comprehensive comparisons. The measurement errors of the three methods were discussed in the final evaluation.

2. Material and methods

The method used in assessing rock fracture spacing from 2D images is shown in Fig.1. The workflow consists of four steps: (1)Automatic extraction of the rock fracture map by applying the deep learning model FraSegNet(see Section 2.1);(2)Digitization of rock fracture map using the data in pixel format (see Section 2.2); (3)Disconnection and grouping of rock fracture traces using a Kmeans++ clustering algorithm (see Section 2.3); and (4) Quantitative evaluation of set spacing, total spacing, and S-RQD value using an interactive measurement algorithm (see Section 2.4).

2.1. Automatic extraction of the skeletonized rock fracture

As shown in Fig. 2, a deep learning model named FraSegNet-VGG19-ASPP and a subsequent skeletonization algorithm were used to extract the rock fractures and the corresponding skeleton map.The FraSegNet was originally developed by Chen et al.(2021d)for rock fracture segmentation, which contains an image input portal, a modified VGG19-based encoding module (i.e. five convolutional layers,four pooling layers and an Atrous spatial pyramid pooling module (ASPP)) (Simonyan and Zisserman, 2014; Chen et al., 2018), a decoding module, and a prediction output module.Among them, a layer is a basic container that usually receives weighted input, transforms it with a set of mostly nonlinear functions and then passes these values as output to the next-layer.The stochastic gradient descent (SGD) algorithm was applied as an optimizer to the deep learning model (Bottou, 2012; Zhao et al.,2021a, b). After the processes of image acquisition, cropping,screening and labeling, a database of rock surface images was created and used as input data.The side-output feature maps of the modified encoding module with different scaled features were acquired completely as they are important for the final prediction.Then the obtained feature maps with various scales were forwarded to the decoding module. The final fracture map was processed through a guided filter (He et al., 2010). A modified skeletonization algorithm that can remove the twigs from the main fracture trunk was used to extract a skeleton fracture map (Saha et al., 2016). The extracted fracture skeleton map is an image in black and white, where white pixel points correspond to the skeleton lines representing the fractures and the black pixel points correspond to the non-fractured regions. A more detailed description of the above algorithm can be found in Chen et al. (2021d).

2.2. Rock fracture image digitization

The main purpose of digitizing the 2D image is to accurately quantify the fracture spacing,because the coordinates of the pixels are crucial for obtaining the relative geometric information. The process involves four steps: (1) elimination of single and forked nodes, (2) extraction of fracture endpoint, (3) acquisition of the chain code information of each fracture,and(4)polyline fitting and determining of the key nodes. Steps (2) and (3) are described in Section 2.2.2.

2.2.1. Elimination of single and forked nodes

Although the skeletonized fracture map has eliminated twigs,the single and forked nodes still exist. From Fig. 3, it can be seen that the single and forked points could significantly affect the statistics of the fracture information, thus increasing the difficulty in statistical analysis of the fracture geometry. Undoubtedly, these points need to be effectively eliminated, by a chain code-based search algorithm (Ling, 2007; Healy et al., 2017). For a certain image,we define that there are eight adjacent points around a given pixel, so that uniformly distributed points can comprehensively represent all possible directions(Fig.3).We traverse all pixels from left to right and top to bottom, and count the number (N) of the neighboring points of each selected pixel.Single and forked points could be identified as pixels with N=0 and N ≥3,respectively.By deleting these identified points,the first step of image digitization is completed.

2.2.2. Extraction of fracture endpoint

In the first step,the pixel points with N=1 are searched again,as they are the endpoints of each fracture. Then, the relative coordinates of each endpoint are recorded.Each fracture is traversed from one endpoint until it reaches the end of the fracture.Using the chain code-based algorithm, the endpoint coordinates and the pixel point directions are calculated.In the end,the coordinates of each skeletonized fracture pixel can be calculated.

2.3. Polyline fitting and key node determining

Fig.1. Workflow of the proposed method for rock fracture spacing measurement.

Fig. 2. Sketch of the deep learning-based fracture trace extraction process.

Fig. 3. The definition of chain code and the relative position of pixel nodes.

Although all the coordinate information of the fracture pixels can be obtained, it is time-consuming and pointless to use the entire coordinates for a task where the geometric information of the original map is of the primary concern. Thus, the polyline approximation method shown in Fig. 4 was developed. In this method, several key nodes (i.e. P1, P2, P3and P4in Fig. 4) are conditionally generated to represent the subordinate fracture, i.e.the connecting lines of the subsequent key nodes are used to represent the original fracture traces.In this method,the parameter D refers to the Euclidean distance between the straight line connecting the two consecutive key nodes and the skeletonized pixel points within the curve connected by these two key nodes.Dmaxis the manually set distance threshold to control the density of node selection.The selection of Dmaxwill be presented in Section 3.1.For a fracture trace that has not been traversed,the endpoint is the first key node (P1). The second key node (P2) is defined as the point where D=Dmaxfirst occurs.The iteration to find the next key node is repeated until the entire traces are traversed.Note that no matter how the relation between D and Dmaxchanges, the start and end points must be the key nodes as they control the coverage of the entire trace.The relative coordinates of the identified key nodes are recorded as(Xi,Yi).If there is no pixel point between the start and end of a trace where D ≥Dmax,it can be concluded that the original curve is smooth, which can be approximated to a straight line.

2.4. Rock fracture image digitization

In deep learning-based rock trace extraction,it is inevitable that some traces from different groups are wrongly combined together,because the difference between pixels in the appearance is negligibly small, if we do not consider their engineering meaning. This problem will affect the assessment of set spacing and subsequent statistics.Thus,although the polyline approximation method could greatly reduce the amount of data,it is important to disconnect and re-group each trace in a correct way. In this section, an algorithm based on angle thresholds is presented for automatic fracture disconnection.Then, a K-means++clustering algorithm is applied to achieving the optimal grouping of traces.

2.4.1. Fracture trace disconnection

For fracture traces with large tortuosity, inaccurate representation of the actual trace is ineluctable if the setting value of Dmaxis randomly increased without a limit during the polyline fitting phase. Thus, an angle threshold-based algorithm was proposed to detect the points with large curvature which can be used as the cutoff points to break the trace.The main principle of the algorithm is to find a breakpoint in the continuous trace so as to minimize the angle between the two lines connecting this point and the other two endpoints (Leng et al., 2021). This is to find out the furthest point of the straight line connected by the start and end points of the continuous trace.As shown in Fig.5,the parameter ∂ispecifies the minimum angle between these two straight lines. ∂maxis the manually adjustable angle threshold to control the position of the breakpoints.The selection of ∂maxwill be presented in Section 3.2.

To illustrate the algorithm, Fig. 5 shows the process to find the node with the minimum angle ∂1for a typical trace. In this algorithm, if ∂1≤∂max, the searched node is the first breakpoint B1;otherwise,the whole trace does not need to be disconnected.Next,the first breakpoint B1is taken as the new starting point of a trace and the algorithm will retrieve the next potential breakpoint for both sides of the remaining fracture. In the remaining trace segments, the point that forms a minimum angle ∂2(Fig. 5b) is searched. The values of ∂2and ∂maxare compared again. If ∂2≤∂max, then the searched node is the second breakpoint B2; otherwise, the remaining trace will not be disconnected. Furthermore,the subsequent breakpoints are positioned in sequence until the disconnection requirements cannot be satisfied (Fig. 5c). An example of disconnecting trace is shown in Fig. 5d, from which it can be seen that the whole trace is divided into four parts as represented in different colors. The traversal process does not finish until all fulfilled traces in the fracture map are detected and disconnected.Note that the nodes identified as breakpoints are not deleted, as the deletion would affect the statistical index of the trace length.Instead,these breakpoints are represented as sharing trace endpoints in the fracture set spacing statistics. In this approach, the curved trace can be effectively partitioned via the angle threshold-based algorithm.

2.4.2. Fracture trace grouping

An accurate trace grouping is the first step of rock spacing measurements, especially for the set spacing. In this paper, a clustering model based on the K-means++algorithm has been applied to automatic fracture clustering. Compared with the original Kmeans algorithm,the proposed model shows some improvements in selecting initial centroids.A centroid is a point which is assumed to be the center of the cluster. The principle is to control the distance between centroids as small as possible (Chen et al., 2021a).Assuming that miinitial cluster centroids (0 < mi< K) have been selected, the point farther away from the current cluster centroids has higher probability of being selected as the (mi+1)th cluster centroid.Among them,K represents the optimal number of cluster centroids. Note that all cluster centroids were selected randomly,including the first one(m1).The corresponding pseudocode for the K-means++ algorithm is demonstrated in Appendix, where the angle of the trace (T) is considered as a key input parameter for clustering, since it represents the expansion trend of the discontinuous plane.Each trace can achieve the reliable disconnection of a significant curved region by ensuring the accurate characterization of the original contour. Hence, the K-means++ algorithm is a relatively reliable approach for automatic fracture grouping after the operations of polyline extraction and trace disconnection.

The Silhouette validity index(SVI)(Rousseeuw,1987)was used to select an optimal number of clusters (K) for grouping fracture traces. It was first proposed to evaluate the cases where actual category information is unknown:

where tirefers to a trace angle within the Kith cluster, a(ti)denotes the average distance between tiand one other randomly selected trace angle in the same category, and b(ti) is the average distance between tiand the randomly selected trace angle in another categories.The total SVI value of all samples in a cluster is the average value of the corresponding SVI of each sample. The value of SVI ranges from -1 to 1. A value close to 0 indicates an overlapping cluster, and a negative value indicates an incorrect cluster assignment. Larger SVI values indicate better clustering performance,which means that the distance between samples in the same category is relatively small, while that in different categories is large.Thus,the cluster value where the maximum value of SVI(ti)is obtained can be selected as the optimal grouping of fracture traces.The automatic optimization of SVI value will be presented in Section 3.3.

Fig.5. The process of angle threshold-based fracture trace disconnection,including:(a)Determining the first breakpoint of the trace,(b)Determining the second breakpoint of the remaining trace,and(c)Determining the rest breakpoint in the same way.(d)shows the diagram of the final trace disconnection.Note that the plot only shows the trace traversal on one side to simplify the process. P1-Pn are the fracture points selected as key nodes of one fracture. B1-B3 are the breakpoints in the fracture polyline.

2.5. Interactive measurement of spacing parameters and S-RQD value

2.5.1. Fracture set spacing measurement

The most commonly used methods for measuring the fracture set spacing are line-, area- and volume-based measurements. In practice,the spacing of fractures is generally determined using the line-based method. Then, the intersection of a scanline and the trace is used to characterize the geometric features of spacing within the group, such as spacing distribution, average, median,and maximum spacings. Thus, an interactive line measurement(ILM)method is proposed,which allows the experienced experts to draw scanlines and to measure the geometric features of the trace.The intersection statistics can be automatically extracted when a user locates the two endpoints of the scanline in the clustered fracture map.The fracture set with the mean trace angle closest to the perpendicular direction of the scanline is selected as the statistical object. For example, in Fig. 6, the intersections of the scanline and fracture group 1(blue trace)are extracted.The length of the connecting lines (i.e. l1, l2, l3, l4, …, lm) between successive intersections is measured to obtain the ultimate characterization of fracture set spacing.

2.5.2. Measurement of total spacing and S-RQD

The ILM method (Fig. 7) is also adopted to determine the total spacing of the tunnel section. Unlike the fracture set spacing, the total spacing is calculated by extracting all points that intersect the scanline. Then, the spacing between the connected intersection points(i.e.l1,l2,l3,l4,…,ln)is calculated to obtain the total spacing.In addition,the S-RQD value can be measured as the percentage of cores within the length of 10 cm (Deere,1988). It should be noted that this S-RQD value is of great significance for engineers to understand the degree of fragmentation of the current rock section,even though it is not directly derived from the borehole. Thus, to obtain the S-RQD value from the total spacing,the proportion of the sum length of intersection segments that are longer than 10 cm out of that of all segments is calculated. The S-RQD value can be calculated as follows:

Fig. 6. Sketch of the definition of set spacing.

Fig. 7. Sketch of the definition of total spacing.

3. Parameter optimization for grouping fracture traces

To achieve accurate and efficient trace grouping, three parameters need to be optimized, i.e. the threshold distance (Dmax) for polyline fitting, the threshold angle (∂max) for trace disconnection,and the SVI for ultimate grouping.To select the optimal range of the parameters,five tunnel section images(Fig.8)from Yunnan,China,were randomly selected and analyzed. The tunnel section 1 (TS-1)is applied to showing the results of different Dmaxvalues and select an optimal Dmaxvalue. Tunnel sections 2 (TS-2) and 3 (TS-3) are used to test the applicability of the selected Dmaxvalue and compare the results of ∂maxvalues.Based on the selected Dmaxand∂maxvalues,the tunnel section 4(TS-4)is employed to optimize the number of clusters by selecting a maximum SVI value. Tunnel section 5(TS-5)is finally applied to conducting the whole processes of spacing measurements.

3.1. Parameter selection for polyline fitting

Fig. 8. Sketches of tunnel geology and five selected tunnel sections from Yunnan,China.

Polyline fitting is a vital step in digitizing image information,because the extracted coordinates of fracture pixels are directly related to the statistics of the geometric information.The essence of digitization is to extract key nodes with different densities by adjusting the Dmaxvalue to characterize the changes in trace morphology.When the density of the node layout is high,i.e.Dmaxvalue is small, the similarity of the geometric shape between the extracted and the original fracture map is high. In contrast, when the density is small,the similarity is low.When the node density is excessively increased, the trace can be characterized fairly well.However,computing a large number of nodes is costly. In order to select an appropriate Dmaxvalue, we should maintain a balance between the similarity of the fracture map and the operating cost.Taking TS-1 as an example, the extracted polyline-shaped fracture and the corresponding number of key nodes with various Dmaxvalues are shown in Fig. 9. It can be seen from Fig. 9 that the curvature similarity of partial traces starts to decrease when Dmaxis larger than 6 pixels. In particular, details showing bending curves cannot be fully represented when Dmaxis greater than 20 pixels.On the contrary,when Dmaxis less than 6 pixels,decrease in Dmaxvalue does not significantly affect the similarity,but may lead to a sharp increase in the number of key nodes. Thus, Dmaxis qualitatively optimized to 6 pixels for the tunnel scale in this work.

3.2. Selection of parameters for trace disconnection

After polyline extraction, the fracture morphology is characterized, and the nodes are coordinated well. The next step is to reasonably disconnect the fracture nodes for further grouping.Hence, the threshold angle ∂maxis defined to control the layout of breakpoints. Since the rationality of disconnection cannot be evaluated quantitatively,the disconnection effects of different ∂maxvalues are compared qualitatively to select the optimal ∂max. TS-2 and TS-3 are used to select the optimal parameter for trace disconnection.The main evaluation principle is to see whether the obvious turning nodes are accurately broken in most traces in the fracture map.Fig.10 shows the results with ∂maxvarying from 120°to 175°in the disconnection phase of fracture traces. The segmented polylines with different colors represent different broken segments.It can be seen that selecting a small value of ∂max(less than 120°)may ignore a large number of key inflection points,resulting in unreasonable grouping of the traces.However,using a large value of ∂max(more than 160°) may overemphasize the curvature of the fracture traces and result in too many meaningless breakpoints. Thus, it is recommended to select a threshold angle value between 120°and 160°that allows the proposed algorithm to disconnect the fracture trace lines with reasonable accuracy.To this end, the ∂maxvalue is set to 150°in this work to form the final disconnected trace lines.

3.3. Parameter optimization for fracture grouping

The use of optimal parameters is of great significance to achieve accurate trace grouping.The optimized parameters can be directly used to measure the trace spacing for other rock mass samples with similar size of surface. However, it should be noted that these parameters are highly affected by the ratio of pixels to the actual sample size,i.e. they may not be fully applicable to rock sites with different sizes, such as slopes and open pits.

Before using the K-means++ based trace grouping, the optimized parameters,i.e.Dmax=6 pixels,and ∂max=150°,are used to obtain the disconnection maps. Then, SVI is applied to automatically selecting the number of optimal fracture group (GN). For a complete analysis, the number of fracture trace clusters is tested sequentially in 2-5 groups.As previously explained in Section 2.3.2,a larger SVI value indicates better clustering performance, and the number of clusters with the highest SVI value is considered as the optimal grouping of fracture traces.Fig.11 shows the results of the grouped fracture maps and the corresponding SVI values.In Fig.11,each grouped fracture trace is shown in one color,and the specific number of grouping can be clearly seen.Taking TS-4 as an example,the SVI reaches a maximum value of 0.614 when GN=3,followed by GN = 2, 4 and 5 with SVI values of 0.603, 0.571 and 0.56,respectively. Although there is no obvious fluctuation in the SVI value for different groups,the rationality of the grouped traces can still be ensured in practice by maximum SVI value and the visual effects. As a result, the optimal number of fracture trace grouping for the 2D rock tunnel face in this work is 3. In summary, it is reasonable to optimize the final number of fracture groupings in accordance with the maximum SVI value.

4. Interactive measurements and comparisons

The main purpose of the ILM method is to allow users to selectively delineate the measurement location according to the project requirements, and then calculate the corresponding parameters, i.e. the fracture set and total spacings, and the S-RQD value. It should be noted that the units of these three parameters should be the length scale rather than pixel. Therefore, it is necessary to determine the relationship between the pixel value and the actual size of the site to convert pixel to length.Following Chen et al. (2021e), the actual size of a pixel in this study is calculated by placing a steel tapeline on the rock tunnel surface.Through the artificially added benchmark, conversion relationship between pixel and length can be easily obtained. In addition, the user can also calibrate the size with the calibration board to obtain a more accurate conversion ratio.

Fig.9. Results of polyline fitting with different Dmax values varying from 1 pixel to 50 pixels:(a)Corresponding Dmax values and the numbers of corresponding extracted nodes,(b)Polyline sketches of the extracted fracture traces, and (c) Node maps of the extracted fracture traces.

Fig.10. Diagrams of the angle threshold-based optimizations of disconnected rock traces:(a)Corresponding angle thresholds,(b)Optimizations of the TS-2,and(c)Optimizations of TS-3.

Fig.11. Sketches of the K-means++ based rock trace grouping optimization:(a) Maps of grouped fracture traces,(b)Overlapping graphs of trace grouping,and (c)Corresponding numbers of grouping and SVI values.

4.1. Interactive measurement of fracture set spacing

The mapping procedure for ILM method is as follows: (1)selecting a group of fracture traces, (2) locating the head and tail points of the scanline along the perpendicular direction of the group of traces, (3) drawing a line segment with the head and tail points as the endpoints, and (4) counting the intersection coordinates of this line segment and the selected trace grouping.Currently,there are two common methods in practice to determine fracture spacing: the field method and the computer vision-based straightening method. The former can provide accurate measurement results by directly counting the spacings between scanline and real traces.Thus,the data obtained from the field measurement can be used as benchmark to evaluate the performance of different methods. For the latter, the main simplification in the current statistics of spacing is to straighten the curve traces,i.e. to ignore the curvatures of the traces.

In order to comprehensively analyze the rationality of the proposed interactive measurement, the field measurement and straightening method were used for comparison.A case study(TS-5) was carried out using the optimized parameters. For the case study,a threshold distance(Dmax)of 6 pixels,and a threshold angle(∂max) of 150°were chosen, and the maximum SVI value of 0.691 was identified to determine the optimal fracture grouping of three.Note that the grouping results of other SVI values are not shown due to space constraints. All statistical results are converted from pixel to length. Fig. 12 illustrates the visualization maps and the corresponding scanline locations.To achieve a proper comparison,the consistency of scanline position is ensured for each comparison during the measurement. It can be seen that the main difference between the proposed and straightening methods is the visualization of the key nodes along the fracture traces.Based on the field measurement, the parameters of Error-1 and Error-2 (see Eqs. (3)and (4)) are applied to characterizing the measurement performance of different methods.

Fig.12. Sketches of set spacing measurement and the corresponding scanline locations using(a)field,(b)proposed and(c)straightening methods.In(b),Group 1 is indicated in red,Group 2 is indicated in blue, and Group 2 is indicated in green.

For quantitative analysis, Fig.13 gives the distributions and the measurement errors of the counted fracture set spacing for three comparative methods. It can be seen from Fig. 13 that the distributions of the three trace groupings in the proposed method are in good agreement with that in the field measurement. For the straightening method, although it performs relatively well in the distribution statistics for Groups 2 and 3,large discrepancy exists in Group 1. In general, the performance of the two non-contact methods in the group with complex traces is slightly worse than that with simple traces, especially for the straightening method.The possible reason is that the original bending mode of the trace was changed after straightening, resulting in an offset of intersections between the traces and scanline. This affects not only the final distribution, but also the mean errors of fracture set spacing.The detailed statistical values of the fracture set spacing of different groups are summarized in Tables 1-3.The median,mean and maximum fracture set spacings obtained by the proposed method are close to the values measured in the field.The results of the straightening method for the statistical indicators deviate greatly from the other two methods,except for the average spacing.In summary,the mean error of the proposed method is less than 0.1(i.e. error values of 0.079, 0.062 and 0.033 for the three trace groupings), which is better than that of the straightening method(i.e. error values of 0.124,0.105 and 0.072).

4.2. Interactive measurements of total spacing and S-RQD value

The mapping of interactive measurements of the total spacing and corresponding S-RQD value is generally similar to that of fracture set spacing. It is worth noting that the scanline for total spacing statistics can be delineated in any direction according to engineering requirements.It is not necessary to be perpendicular to a particular group of fracture traces. The next step is counting the intersection coordinates of the scanline with all traces,recording all adjacent distances,and calculating the final S-RQD values.

In order to make a comprehensive comparison, the field and straightening methods are also employed in the analysis. Fig. 14 shows the visualization maps and the scanline locations of the three comparison methods.Note that for the proposed method,the scanlines can be drawn interactively by experts. To simplify the presentation, we only compare the statistical information generated by a relatively complex scanline. It can be seen that there are many differences between the proposed and straightening methods in the trace morphology.The proposed method appeared to be closer to the field measurement.

The statistical distributions of the total spacing and the corresponding error are indicated in Fig. 15. It can be found that the consistency between the proposed method and the field results is significantly higher than that of the straightening method. The detailed results of total spacing and S-RQD values are listed in Table 4. From this, it can be seen that the median, maximum and mean total spacings for the proposed method are 12.5 cm,34.1 cm and 16.1 cm, respectively. The results are close to the field measurements,i.e.12.3 cm,39.1 cm and 15.4 cm.The proposed method with a mean measurement error of 0.055 is better than that in the straightening method of 0.123. The reason is that using the straightening method changes the original shape and nesting of the fracture trace, which is different from the actual characterization.Furthermore, the S-RQD value obtained by the interactive measurement is 85.1%,which is close to 84.4%of the field measurement and much better than 88.8% of the straightening measurement.Thus,it can be seen that the proposed method performs well in SRQD prediction of the rock mass section. In conclusion, the proposed method can accurately measure and evaluate the fracture spacing and S-RQD value of the studied rock mass.

4.3. Comparative discussion

In this study, the filed and straighten methods have been employed for comparison.The filed measurement has been applied to a benchmark model since it measures the spacings between scanline and real in situ traces,in consideration of the real shape of rock fractures.Nevertheless,it is still controversial due to its limited measurement scale and huge labor and time requirements. The original shapes of traces are complex in practice,making it difficult to describe accurately using the obtained images. Thus, straightening method was gradually used by engineers. However, it is proven that the straightening method affects the accuracy of the statistical results although it is simple. Through comprehensive comparison, it was found that the proposed interactive measurement can not only characterize the accurate shapes of rock traces,but also achieve the automatic grouping of traces, to accurately measure the trace parameters, i.e. fracture set spacing, total spacing,and RQD value.

Fig.13. Statistical results and corresponding errors of the fracture set spacing of the three methods along the specified scanline for(a)Group 1,(b)Group 2,and (c)Group 3.

Table 1 Statistical results and errors of the estimated set spacing from Group 1 using three methods.

5. Conclusions

In this paper, an interactive measurement is proposed for automatic calculation of set and total spacings of fracture traces based on rock tunnel face images. Optimization methods were proposed to select the values of three parameters.By testing on real tunnel surface images, it was shown that the proposed methodcould provide comparable accuracy compared with the results of field measurement,and outperforms the existing computer visionbased straightening method in estimating fracture set and total spacings,as well as the S-RQD values.The conclusions of the study are drawn as follows:

Table 2 Statistical results and errors of the estimated set spacing from Group 2 using three methods.

Table 3 Statistical results and errors of the estimated set spacing from Group 3 using three methods.

Fig.14. Sketches of total spacing measurement and the corresponding scanline locations using (a) field, (b) proposed and (c) straightening methods.

(1) The proposed polyline fitting algorithm can represent the fracture traces in pixel form with key nodes in coordinate form. The threshold distance parameter (Dmax) can effectively control the density of nodes. For the tunnel samples,Dmaxwas qualitatively optimized to 6 pixels to balance the fitting effect and operating cost.

(2) The threshold angle (∂max) based algorithm can reasonably disconnect the fracture traces by automatically searching the breakpoints. Setting the ∂maxto 150°was considered to be the closest to the actual trace disconnection in this project.

(3) The K-means++ algorithm was proved to realize automatic clustering of traces. Meanwhile, the SVI can reasonably optimize the number of fracture trace grouping.

(4) The field and straightening measurements were conducted for a comprehensive comparison. Using the optimal parameters, the proposed interactive measurement model was proved to have good prediction performance for both spacing parameters and S-RQD value.

The proposed method is an extension of the geological characterization of 2D rock tunnel faces,especially for estimating fracture trace spacing and S-RQD.It is noted that this method is also suitable for measuring the trace spacing of other rock images,such as open pits and slopes. It would also be interesting to examine the applicability of the proposed method to different tunnels.It should also be mentioned that as the parameters involved in this method and the corresponding statistical trace spacing information are highly affected by the ratio of pixels to the actual sample size. Thus,for a rock surface with different scales,it is suggested to re-evaluate the optimal model parameters using the procedure recommended in this paper.

Fig.15. The statistical results and corresponding errors of the total spacing of the three methods along the specified scanline.

Table 4 Statistical results of total spacing and S-RQD values using three methods.

Declaration of competing interest

The authors wish to confirm that there are no known conflicts of interest associated with this publication and there has been no significant financial support for this work that could have influenced its outcome.

Acknowledgments

This study is supported by Key Innovation Team Program of Innovation Talents Promotion Plan by Ministry of Science and Technology(MOST)of China(Grant No.2016RA4059),Science and Technology Project of Yunnan Provincial Transportation Department(Grant No.25 of 2018),and Shanghai Science and Technology Committee Program (Grant No.20dz1202200).

Appendix A. Supplementary data

Supplementary data to this article can be found online at https://doi.org/10.1016/j.jrmge.2021.10.012.