Research Article  Open Access
Fracture Segmentation Method Based on Contour Evolution and Gradient Direction Consistency in Sequence of Coal Rock CT Images
Abstract
The coal rock exhibits obvious heterogeneity and anisotropy after a longterm geological evolution in nature. Therefore, the coal rock CT (computed tomography) image shows uneven grey scale, low contrast, and the fractures with weak boundaries. The accurate segmentation of the coal rock fracture networks is challenging. In this paper, a segmentation method of fractures based on contour evolution and gradient direction consistency is proposed to accurately segment the fracture networks in the sequence of coal rock CT images. According to the contour variation rule of the fractures in the discrete 3D (three dimensional) space formed by the sequence of CT images, the fracture contour evolution model (FCEM) is constructed and the preliminary segmentation results of fractures are obtained from FCEM. A 3D adaptive median filtering (3DAMF) and a 3D bilateral filtering (3DBF) are proposed. The high density miscellaneous point noises in the coal rock CT images are filtered by the 3DAMF. And the boundaries of fractures are enhanced by 3DBF. According to the similarity of the preliminary segmentation results of fractures and the real contours of fractures, the preliminary segmentation results of fractures are optimized based on the gradient direction consistency model (GDCM) proposed in this paper to obtain the accurate boundaries of fractures. The fracture segmentation method proposed in this paper can obtain accurate boundaries of fractures with weak boundaries, and the experimental results show that the segmentation efficiency for sequence is high and adaptability is strong.
1. Introduction
With the increase in the depth of the underground coal mining, the gas content in the coal rock increases accordingly. The disasters such as gas outbursts and gas explosions have caused huge economic losses and a large number of casualties. The fractures of the coal rock are the storage place and the migration channels of the gas. The coal mining will cause changes in the stress field of the coal rock, which in turn will affect the extension, expansion, and redistribution of the fractures in the coal rock. Among them, the fractures with weak boundaries often indicate the most likely change position of the coal rock under mining stress. The accurate segmentation of the fractures in the coal rock is conducive to obtaining the shape parameters of the fractures (such as length, area, azimuth, etc.) and then can accurately give a quantitative description of the fractures in the coal rock. This is of great significance for mastering the evolution rule of the fractures in the coal rock under mining stress, as well as the coal rock seepage calculation and the damage calculation.
After many stages of the complex natural tectonic movement, the coal rock exhibits obvious heterogeneity and anisotropy, and the complex fracture networks were formed inside. Detecting the internal structure of the coal rock by industrial CT (computed tomography) is a nondestructive testing technology. Compared with other methods, CT scanning can detect the structure of fractures in the coal rock at the microscopic scale (micron or nanometre) and has the advantages of intuitive imaging, high resolution, and no limitation on the structure of the object to be detected. At the same time, a sequence of CT images contains thousands of images. The coal rock CT images are characterized by uneven grey scale and low contrast and contain high density miscellaneous point noises. Moreover, the fractures with weak boundaries are widely present in the coal rock, which makes it difficult to accurately segment the fractures in the coal rock.
For many years, the researches on segmentation methods of the fractures mainly have focused on the 2D image processing techniques based on threshold, edge detection, region growth, and texture analysis [1–9]. H Oliveira et al. [1] used dynamic threshold method and entropy to segment fractures of the road image. A Landstrom et al. [3] extracted longitudinal fractures from the image of road surface using morphological methods and logistic regression models. L Liu. et al. [5] studied the segmentation of fractures in CT image based on wavelet transform and CV model. Y. Fujita and Y. Hamamoto [7] applied the eigenvalues and eigenvectors of the Hessian matrix to discriminate the shape to distinguish fractures. The above methods only utilize the 2D features of the fractures, and these methods are sensitive to noises, uneven grey scale, and low contrast. Due to the lack of more reliable information of fractures, the segmentation results of the fractures with weak boundaries have serious misjudgement and misdetection by these methods.
In recent years, the segmentation methods of fractures based on 3D features have been widely concerned [10–12]. Since these methods can utilize the 3D grey scale and 3D structural features, a higher recognition rate and a lower misjudgement rate are obtained. But the corresponding 3D processing technologies of the fractures are rarely reported. MR. Jahanshahi [10] used depth information to reconstruct the image and then identify the fractures. M Voorn et al. [11] studied the fractures of CT image based on Hessian matrix multiscale planar filtering.
At present, the researches on the segmentation of the fractures in the coal rock are rarely reported. In this paper, by studying the 3D grey scale and 3D structure features of fractures in sequence of the coal rock CT images, a fracture segmentation method based on contour evolution and gradient direction consistency is proposed. The fracture contour evolution model (FCEM) is proposed based on the 3D variation rule of fracture contours, and the preliminary segmentation results of fractures are obtained by FCEM. According to the characteristics of the coal rock CT images, the 3D adaptive median filtering (3DAMF) is used to reduce the image noises, and the 3D bilateral filtering (3DBF) is used to enhance the boundaries of fractures. At the same time, the gradient direction consistency model (GDCM) is proposed and is used to optimize the preliminary segmentation results of fractures, and finally the accurate segmentation results of fractures are obtained.
2. Fracture Contour Evolution Model
The fractures are generated from natural evolution, and their features (such as shape, length, area, and azimuth) are random. Therefore, the structural features of the fractures are rarely mentioned in the studies of fracture segmentation. Most of the existing methods are based on the grey scale feature of the fractures. A discrete 3D space is formed by the sequence of coal rock CT images, where represents the 2D coordinates of the pixel points in the single CT image and represents the serial number of CT image. Discriminating whether a target in the image is a fracture, one often makes a judgment by analyzing the same position in the adjacent images. That is to say, in the 3D space, the fracture contours in the adjacent CT images are very similar, and the smaller the interlayer spacing is, the higher the contour similarity is. This characteristic of the fracture structure also clearly exists in the fractures with weak boundaries. That is, the fractures in the th image can be obtained by the slight shape deformation of the fractures among the adjacent CT images. Therefore, this paper proposes FCEM for obtaining the fractures in the CT image.
images of consecutive serial numbers are selected and added to the set . The images in the set require as much as possible including the complete fracture networks. The fracture networks in the set are segmented and added to the set . In this paper, the evolution rule is based on the curve evolution theory [13]. This theory mainly studies the process of curve evolution between two curves with time by the normal. During the process of curve evolution, only the point motion in its normal direction is considered.
Firstly, is discretized into a set [14]. At the same time, the normal lines of each point in the set are available, and these normal lines are added to the set . As shown in Figure 1, each normal line in and the contour curve intersect into the set . In the same process, the discrete points set and the normal line set on the curve can be obtained.
As shown in Figure 1, each point evolves along the normal direction and reaches the fracture contour of the adjacent image. The motion time of the point is related to . The speed of each point in the set is different. is obtained mainly by simulating the evolution process of . The evolutionary process of the points can be seen as a uniform acceleration motion. Set the speed of the point to zero. The acceleration formula for the point is calculated as follows:where is the Euclidean distance between the point and the point . The formula of the speed for the point is as follows:
Similarly, the acceleration and the speed of the point are calculated. Therefore, the distances of the points in the set moving along the normal direction are calculated by formula (3).
Finally a new set is obtained. The points in the set are fitted to a new curve [15, 16]. The normal line set of the is , and the corresponding normal direction angle set is . is very similar to the exact boundary of the fractures to be segmented, so is taken as the preliminary segmentation result of the image .
3. 3D Enhancement Methods of Coal Rock CT Image
The FCEM only utilizes the contour variation characteristic of the fractures in 3D space. In this paper, the result obtained by FCEM is optimized according to the 3D grey scale feature of the fractures. At the same time, for the characteristics of the coal rock CT images with low contrast, uneven grey scale, and high density miscellaneous point noises, this paper proposed 3DAMF and 3DBF, which can reduce noises and enhance the boundaries of fractures.
The coal rock is a natural body, and the coal rock CT image shows an obvious low contrast, as shown in Figure 2.
(a)
(b)
The red line in Figure 2(a) is located on a coal rock fracture with obvious fracture characteristics. The green line is located on the coal rock background and spans two fractures with weak boundaries. The curves of red and green in Figure 2(b), respectively, show grey scale value distribution of each point on the same color line in Figure 2(a). It can be seen from Figure 2(b) that the grey scale distribution of the weak boundaries is basically consistent with the coal rock background. Therefore, this paper uses the local histogram statistical method to perform contrast stretching on the coal rock CT images.
The coal rock image is a heterogeneous body with high density miscellaneous point noises. According to the median filtering principle, when the size of the filtering window is small, the image details can be well protected, but the filtering effect of noises is not very good. Moreover, if the number of the noise points in the filtering window is larger than the number of the coal rock pixels in the entire window, the noises cannot be well filtered. In order to protect the fractures with weak boundaries, this paper proposes 3DAMF. This method considers all information of the pixels in the 3D neighbour space of each pixel. So the size of the filtering window becomes a cube (e.g., 333), and the size of the filtering window can be adaptively adjusted according to the local information. The filtering method can protect the information of weak boundaries while filtering high density miscellaneous point noises, and the principle is as follows:
To Determine whether the median value in the filtering window is a noise point, and if so, increasing the size of window to find a nonnoise median output.
To determine whether the point to be processed is a noise point. If not, do not process; otherwise replace the point value with the median value of step .
For the feature of uneven grey scale in the coal rock CT images, the bilateral filtering is used to enhance the images. The bilateral filtering is a nonlinear filter composed of a spatial domain filter and a range filter. This method not only considers the distance between pixels, but also considers the similarity of grey scale in the pixel range domain. This method has characteristics of being simple, noniterative, and local. In order to improve the effect of bilateral filtering, the 3D domain information of the pixel to be processed is involved in the process of calculation, and the size of template becomes a cube (such as 333). The 3D range filter of this method is expressed aswhere is the coordinate of the central position within the template window and is the coordinates of the other positions within the template window. is the distance standard deviation of the Gaussian function. The 3D spatial domain filter is expressed aswhere represents the pixel value of coordinate and is the grey standard deviation of the Gaussian function. Then the 3DBF is expressed as
After being processed by the 3DBF, is the pixel value of the central point in template window. The 3D enhancement result of Figure 2(a) is shown in Figure 3.
4. Gradient Direction Consistency Model
According to the 3D grey scale feature of the fractures, the points in the set continue to evolve into a new point set . The point in the set will always be evolved into the point in . According to the curve evolution theory, the point is located on the normal line of the point . Since the points and are, respectively, located on two adjacent curves with similar shapes, the gradient direction of the point and the normal direction of the point are substantially identical. A segmentation method of fractures based on GDCM is proposed. This method finds a closest point on the normal line of the point , and the gradient direction of the point is similar to the normal direction of the point . In the fractures with weak boundaries, although the gradient value is small, the characteristic of gradient direction consistency is unchanged. Therefore, this method is effective for fractures with weak boundaries.
There is a set on the normal line , where the point and the point are the same point, and is the size of step when the discrete points are taken on the normal line , as shown in Figure 4.
Find a point in that is as close as possible to the point , and the gradient direction of the point is close to the direction of the normal . In this paper, the gradient values and gradient directions of the points in the set are calculated using the Sobel operator, and the convolution kernels are
For the pixel information of the point domain to be detected, is convoluted with two convolution kernels, respectively, to obtain gradient values and in the horizontal and vertical directions as follows
The gradient value of this point is calculated by and :
The gradient direction of this point can also be calculated.
The mathematical expression of GDCM iswhere , , and are the weights of each influencing factor, and is the distance from the point to the point . When is the smallest, the point is the point that the point is finally reaching along the evolutionary direction. This point is added into the set . Finally, the points in the set are fitted to form a curve which is the final segmentation result of the fracture.
5. Experimental Results and Analysis
In order to verify the effectiveness and superiority of the proposed segmentation method, the XRay (ACTIS300320/225) industrial CT detection system of the State Key Laboratory of Coal Resources and Safe Mining of China University of Mining and Technology (Beijing) was used to carry out CT realtime scanning experiments on different stages of coal rock triaxial loading. The coal samples used in the experiment were taken from the Liupanshui mining area in Guizhou Province. In this experiment, CT images of the coal rock under different loading conditions were obtained, and the CT scan scale was 50 microns. A sequence of CT images with relatively developed fractures was selected, and this sequence contained 876 CT images of the size .
As shown in Figure 2(a), the coal rock image contains an unrelated peripheral cavity. The coal rock images were processed by the Gaussian blur method and binarization method, and then the search method of contours is used to determine the coal rock region. Before the segmentation of the fractures, the peripheral cavity is first removed, so as to avoid the influence of unrelated information of the coal rock CT image on the segmentation of fractures.
Figure 5 is a set of two CT images with the peripheral cavity removed and located at different locations in the coal rock column.
(a)
(b)
The sequence of the coal rock images was processed by the traditional image processing methods. The median filtering [17] and bilateral filtering [18] methods were used to denoise and enhance the images, and the adaptive threshold method [19] was used for binarization. At the same time, morphological methods were used to reduce the connection of noises and fractures. The processed results are shown in Figure 6.
(a)
(b)
The processed results of Figure 6 contain a large number of nonfracture noises, and the loss and discontinuity of the fractures in the results of segmentation are serious. Due to the use of the morphological methods, the thickness of the fractures obtained does not match the actual fractures.
In this paper, the method in [11] was used to segment the fractures in the coal rock CT images. This method based on Hessian matrix multiscale planar filtering is used to segment the fractures from a sequence of CT images, and the main parameters used are as shown in Table 1. The results are shown in Figure 7.

(a)
(b)
As described in [11], the multiscale Hessian fracture filtering has a better effect on the segmentation of fractures in homogeneous materials, and the noises have a great influence on the segmentation of fractures with weak boundaries. Therefore, the fractures with weak boundaries in the rectangular mark area of Figure 7 are not recognized or incompletely segmented. At the same time, since the intersection position of multiple fractures is not a planar feature, this method cannot segment the fractures in the circular mark area in Figure 7(b). The Gaussian scale and the grey scale threshold in this method have a large influence on the noises, and the less the noise that is retained, the more the weak boundaries are lost. So the appropriate values of these parameters are difficult.
The processing results obtained by the method of this paper are shown in Figure 8, and the main parameters used are as shown in Table 2.

(a)
(b)
In the experiment, firstly, five continuous CT images of the complete fracture networks were selected from the sequence, and the corresponding fracture networks were obtained by manual segmentation. The images on both sides of the five serial numbers were iteratively processed using the method proposed in this paper. Therefore all results of segmentation are based on five a priori fracture networks. By the principle of FCEM and GDCM, the segmentation results of fractures are only the extension or contraction of the a priori fracture networks, so the fractures segmented do not appear at other positions of the image. That is to say, a small amount of noises in the results of segmentation will only be presented by the burrs in boundaries of the fractures. At the same time, a fracture that suddenly appears in an image area other than the a priori fracture networks cannot be segmented. Therefore, the a prior fracture networks are required to include all the fractures in the sequence as much as possible, or the sequence is divided into multiple subsequences with similar fracture networks for processing. Because FCEM and GDCM have strong applicability to fractures with weak boundaries and 3DAMF and 3DBF have strong ability to denoise and enhance images, the accurate results of segmentation were obtained in Figure 8.
6. Conclusions
The segmentation of fractures in the coal rock CT images is very important for other related basic researches. In this paper, a segmentation method of fractures based on contour evolution and gradient direction consistency is proposed, which provides an effective method to segment fractures from the sequence of coal rock CT images. By studying the contour change rule of fractures in 3D space, FCEM is proposed, and the preliminary segmentation results of fractures are obtained. Aiming at the characteristics of uneven grey scale, low contrast, and high density miscellaneous point noises in the coal rock CT images, 3DAMF and 3DBF methods are proposed to improve image quality. Finally, the GDCM based on the grey scale feature of the fractures is proposed to optimize the preliminary segmentation results, and the accurate segmentation results of fractures are obtained. The proposed method has strong antinoise ability and exhibits strong segmentation ability for fractures with weak boundaries. For a fracture that may suddenly appear, the sequence needs to be divided into multiple subsequences with similar fracture networks for processing.
Data Availability
Please contact the author by email to obtain the relevant data. The email address is 1173019984@qq.com.
Conflicts of Interest
The authors declare that there are no conflicts of interest regarding the publication of this paper.
Acknowledgments
The authors thank the State Key Laboratory of Coal Resources and Safe Mining of China University of Mining and Technology (Beijing) for its help.
References
 H. Oliveira and P. L. Correia, “Automatic road crack segmentation using entropy and image dynamic thresholding,” in Proceedings of the 17th European Signal Processing Conference, EUSIPCO 2009, pp. 622–626, UK, August 2009. View at: Google Scholar
 J. Tang and Y. Gu, “Automatic crack detection and segmentation using a hybrid algorithm for road distress analysis,” in Proceedings of the 2013 IEEE International Conference on Systems, Man, and Cybernetics, SMC 2013, vol. 60, no. 4, pp. 3026–3030, UK, 2013. View at: Google Scholar
 A. Landstrom and M. J. Thurley, “Morphologybased crack detection for steel slabs,” IEEE Journal of Selected Topics in Signal Processing, vol. 6, no. 7, pp. 866–875, 2012. View at: Publisher Site  Google Scholar
 S. A. Anwar and M. Z. Abdullah, “Microcrack detection of multicrystalline solar cells featuring an improved anisotropic diffusion filter and image segmentation technique,” Eurasip Journal on Image and Video Processing, vol. 2014, no. 1, pp. 1–17, 2014. View at: Google Scholar
 L. Liu, L. Zeng, and B. Bi, “A unified method based on wavelet transform and CV model for crack segmentation of 3D industrial CT images,” in Proceedings of the 6th International Conference on Image and Graphics, ICIG 2011, pp. 12–16, 2011. View at: Google Scholar
 S. Li, Y. Cao, and H. Cai, “Automatic pavementcrack detection and segmentation based on steerable matched filtering and an active contour model,” Journal of Computing in Civil Engineering, vol. 31, no. 5, 2017. View at: Google Scholar
 Y. Fujita and Y. Hamamoto, “A robust automatic crack detection method from noisy concrete surfaces,” Machine Vision and Applications, vol. 22, no. 2, pp. 245–254, 2011. View at: Publisher Site  Google Scholar
 H. N. Nguyen, T. Y. Kam, and P. Y. Cheng, “Automatic crack detection from 2D images using a crack measurebased Bspline level set model,” Multidimensional Systems and Signal Processing, vol. 29, no. 1, pp. 1–32, 2018. View at: Publisher Site  Google Scholar  MathSciNet
 X. Fan, J. Wu, P. Shi, X. Zhang, and Y. Xie, “A novel automatic dam crack detection algorithm based on localglobal clustering,” Multimedia Tools and Applications, pp. 1–19, 2018. View at: Publisher Site  Google Scholar
 M. R. Jahanshahi and S. F. Masri, “Adaptive visionbased crack detection using 3D scene reconstruction for condition assessment of structures,” Automation in Construction, vol. 22, pp. 567–576, 2012. View at: Publisher Site  Google Scholar
 M. Voorn, U. Exner, and A. Rath, “Multiscale Hessian fracture filtering for the enhancement and segmentation of narrow fractures in 3D image data,” Computers & Geosciences, vol. 57, no. 7, pp. 44–53, 2013. View at: Publisher Site  Google Scholar
 B. Peng, X. Cai, L. I. Shaobo, and Y. Zhang, “Automatic crack detection algorithm based on 3D virtual pavement,” Journal of Chongqing Jiaotong University, 2018. View at: Google Scholar
 S. Osher and J. A. Sethian, “Fronts propagating with curvature dependent speed: algorithms based on HamiltonJacobi formulation,” Journal of Computer Physics, vol. 79, no. 1, pp. 12–49, 1988. View at: Publisher Site  Google Scholar  MathSciNet
 S. B. Andersson, “Discretization of a continuous curve,” IEEE Transactions on Robotics, vol. 24, no. 2, pp. 456–461, 2008. View at: Publisher Site  Google Scholar
 Y. Jiang and Y. Li, “The research of the approximate algorithm based on cubic Bspline curves,” in Proceedings of the 2012 Fifth International Conference on Information and Computing Science (ICIC), pp. 23–26, IEEE, Liverpool, UK, July 2012. View at: Publisher Site  Google Scholar
 C. Weihua and Z. Tie, “The study of cubic uniform rational Bspline interpolation algorithm,” Machinery Design & Manufacture, pp. 3–5, 2010. View at: Google Scholar
 T. Galba, K. Romić, and A. Baumgartner, “Edgepreserving partial variable median filtering for fast noise reduction in CT slices,” in Proceedings of the International Symposium ELMAR 2014, pp. 1–4, Croatia, September 2014. View at: Google Scholar
 D. Edwin and S. Hariharan, “Projection space denoising with bilateral filtering and CT noise modeling for dose reduction in CT,” Medical Physics, vol. 36, no. 11, pp. 4911–4919, 2009. View at: Publisher Site  Google Scholar
 D. Edwin and S. Hariharan, “Liver and tumour segmentation from abdominal CT images using adaptive threshold method,” International Journal of Biomedical Engineering and Technology, vol. 21, no. 2, pp. 190–204, 2016. View at: Publisher Site  Google Scholar
Copyright
Copyright © 2019 Zhiwei Li and Guoying Zhang. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.