HALF-CAST MARK IDENTIFICATION AND DAMAGED FLATNESS EVALUATION AND CLASSIFICATION METHOD FOR BLASTHOLES IN TUNNEL BLASTING

Abstract

The present disclosure relates to a half-cast mark identification and damaged flatness evaluation and classification method for blastholes in tunnel blasting, including the following steps: S1-2: photographing first and second contrast images as well as a half-cast mark image after blasting; S3-6: performing denoising, gray-scale processing and binary processing on the above images, and identifying a boundary of a half-cast mark in each of the images; S7-9: determining a flatness damage variable, a quantitative relation among an area of a half-cast mark region, the damage variable and a fractal dimension, and a damage value of the half-cast mark image; S10-11: forming five-dimensional (5D) eigenvectors to obtain multi-dimensional digital information features of the images; and S12-13: selecting eigenvectors of 60 images as training data to input to a naive Bayes classifier (NBC), and taking eigenvectors of remaining 30 images as classification data to input the above well-trained NBC for classification.

Claims

1. A half-cast mark identification and damaged flatness evaluation and classification method for blastholes in tunnel blasting, comprising the following specific steps: S1: photographing standard contrast images, a size of each of the standard contrast images being comprehensively determined by an analysis scale, and the standard contrast images comprising a first contrast image of a surrounding rock surface with an ideal half-cast mark and a second contrast image that is unacceptable with excessive backbreak or overbreak; S2: acquiring a half-cast mark image after actual blasting to serve as a third analysis image, a size and a photographing environment of the analysis image being consistent with those of the above contrast images; S3: denoising the above three images by a two-dimensional (2D) Gaussian algorithm according to a characteristic that Gaussian noise in an acquired image obeys a normal distribution; S4: performing gray-scale processing on three denoised images according to a preset threshold of an image histogram to reduce an original data volume in each of the images; S5: performing binary processing on gray-scale processed images, and setting an optimal binary threshold by maximum entropy thresholding; S6: identifying a boundary and a related region of a half-hole mark in each of the images with a gradient vector flow (GVF)-Snake and active contour model (ACM) (GVF-Snake-ACM); S7: importing each of three half-cast mark identified images to ImageJ software, and determining a ratio of an area of a half-cast mark region to a total area of the image as a flatness damage variable; S8: normalizing an interval for an area ratio of a half-cast mark in each of the first contrast image and the second contrast image, and establishing a quantitative relation ω=(D).sup.β between a damage degree and a fractal dimension D, thereby determining a quantitative relation among the area ratio of the half-cast mark region, the damage variable and the fractal dimension; S9: determining a damage value of the third analysis image through linear interpolation calculation on an area ratio of a half-cast mark; S10: extracting, for each of a plurality of half-hole mark identified images with a gray-level co-occurrence matrix (GLCM) in 0°, 45°, 90°, and 135° directions or 180°, 225°, 270°, and 315° directions, five eigenvalues comprising an energy mean, an entropy mean, a contrast mean, a correlation mean, and a uniformity mean, thereby forming a five-dimensional (5D) eigenvector; S11: respectively photographing 90 blasting images for a tunnel excavation surface in three blasting plans, which specifically comprise conventional blasting, presplit blasting and smooth blasting and employ a same explosive charge during tunneling, to form a set C={c.sub.1, c.sub.2, c.sub.3}, analyzing target matrices corresponding to the conventional blasting, the presplit blasting and the smooth blasting, and obtaining 5D eigenvectors F={f.sub.1, f.sub.2, f.sub.3, f.sub.4, f.sub.5} in different blasting plans, wherein a feature attribute of each of the images is composed of a 5D eigenvector, comprising an energy mean, an entropy mean, a contrast mean, a correlation mean, and a uniformity mean, thereby obtaining a multi-dimensional digital information feature of the image; S12: randomly selecting 5D eigenvectors of 60 images from the image set in the three different blasting plans as training data to input to a naive Bayes classifier (NBC), and calculating, with five eigenvalues f.sub.i(i=1,2,3,4,5) of each of the images, a probability that a training image falls into a category c.sub.j(j=1,2,3) $P\left(c_j/F\right)=\frac{P\left(c_j\right)P\left(c_j/F\right)}{P\left(F\right)}=\frac{P\left(c_j\right)P\left(f_1,f_2,f_3,f_4,f_5/c_j\right)}{P\left(F\right)},$ wherein for an eigenvector F, a larger posterior probability P(c.sub.j/F) indicates a higher probability that the F falls into c.sub.j, and the category into which the F falls is calculated by: $c\left(F\right)=\arg \max P\left(c_j\right)\underset{j=1}{\overset{3}{.Math.}}P\left(f_i|c_j\right),$ and after P(c.sub.j) and P(c.sub.j/F) are calculated, a category with a maximum posterior probability is returned, thereby obtaining the category of the F; and S13: inputting eigenvectors of remaining 30 images in the image set in the different blasting plans to a well-trained NBC, and determining a state with a maximum probability in different given categories c.sub.j(j=1,2,3) based on a maximum a posteriori (MAP) principle of the NBC to take as a final result for evaluation and classification of the blasting plans, thereby implementing classification on blasting effects of the conventional blasting, the presplit blasting and the smooth blasting and on damaged flatness of surrounding rock surfaces.

2. The half-cast mark identification and damaged flatness evaluation classification method for blastholes in tunnel blasting according to claim 1, wherein in step S3, a Gaussian filter based model satisfies the following 2D Gaussian distribution: $G\left(x,y\right)=\frac{1}{2\pi {\sigma }^2}{e}^{-\frac{x^2+y^2}{2{\sigma }^2}},$ wherein x and y are transverse and longitudinal sizes of the image, and σ.sup.2 is a variance of the normal distribution.

3. The half-cast mark identification and damaged flatness evaluation classification method for blastholes in tunnel blasting according to claim 1, wherein in step S5, the maximum entropy thresholding specifically comprises: defining, according to a concept of a Shannon entropy, an entropy of an image having a gray scale in a range of 0, 1, 2, . . . , L−1 as $H=\underset{i=0}{\overset{L-1}{.Math.}}{\mathrm{iP}}_i,$ P.sub.i being a probability that a pixel having a value of i occupies the whole image, and dividing the image into a target O and a background B according to a set threshold t, probabilities of the pixel in the target and the background being: $\mathrm{region}O:\frac{P_i}{P_t},i=0,1,.Math.,t,\mathrm{wherein},P_t=\underset{i=0}{\overset{t}{.Math.}}{P}_i;\mathrm{and}$ $\mathrm{region}B:\frac{P_i}{1-P_t},i=t+1,.Math.,L-1.$

4. The half-hole mark identification and damaged flatness evaluation classification method for blastholes in tunnel blasting according to claim 3, wherein the target O and the background B each have the following entropy function: $H_o\left(t\right)=-\underset{i=0}{\overset{t}{.Math.}}\frac{p_i}{P_t}\ln \frac{p_i}{P_t}=\ln P_t+\frac{H_t}{P_t},\mathrm{wherein}{H}_t=-\underset{i=0}{\overset{t}{.Math.}}{p}_i\ln p_i;$ $H_B\left(t\right)=-\underset{i=t+1}{\overset{L-1}{.Math.}}\frac{p_i}{1-P_t}\ln \frac{p_i}{1-P_t}=\ln \left(1-P_t\right)+\frac{H-H_t}{1-P_t},\mathrm{wherein}H=-\underset{i=0}{\overset{L-1}{.Math.}}{p}_i\ln p_i;$ a total entropy of the image is: $H\left(t\right)=H_o\left(t\right)+H_B\left(t\right)=\ln P_t\left(1-P_t\right)+\frac{H_t}{P_t}+\frac{H-H_t}{1-P_t}$ the total entropy of the image is maximized through a maximum threshold T: T=arg max(H(t)); and the optimal binary threshold is determined according to a principle of maximizing the total entropy of the image.

5. The half-cast mark identification and damaged flatness evaluation and classification method for blastholes in tunnel blasting according to claim 1, wherein step S6 specifically comprises: S61: acquiring pixel points on a row and a column of an initial image I.sub.m×n (having a resolution of m×n) at a step length d=2, and performing downsampling for s times to obtain s+1 to-be-tested images (comprising the initial image I.sub.m×n), wherein a resolution in kth sampling is $I_{\frac{m}{2^k}\times \frac{n}{2^k}}.$ S62: preprocessing an image $I_{\frac{m}{2^s}\times \frac{n}{2^s}}$ in sth downsampling, performing edge detection with a Canny operator, defining an initial contour l.sub.s+1 on an edge image, a number of control points on the initial contour l.sub.s+1 being N.sub.s+1, and calculating a GVF external force field, wherein due to a diffusivity of the GVF external force field in a smooth region, the initial contour l.sub.s+1 approaches continuously to an edge of the image; and when 0.98N.sub.s+1 control points on the initial contour l.sub.s+1 reach a target contour (a termination criterion), it is considered that an initial contour curve reaches a target edge, and iteration is stopped to obtain a target contour curve l.sub.s at the resolution $\frac{m}{2^s}\times \frac{n}{2^s};$ S63: preprocessing an image $I_{\frac{m}{2^{s-1}}\times \frac{n}{2^{s-1}}}$ in (s−1)th downsampling with the Canny operator, taking the contour curve l.sub.s in the previous step as an initial contour curve, recalculating a GVF external force field on the image $I_{\frac{m}{2^{s-1}}\times \frac{n}{2^{s-1}}},$ and obtaining a target contour curve l.sub.s−1 according to the termination criterion in step S62 under pushing of an external force; and S64: repeating step S63, to until a target contour curve l of the identified image I.sub.m×n is obtained.

6. The half-cast mark identification and damaged flatness evaluation and classification method for blastholes in tunnel blasting according to claim 1, wherein in step S7, the flatness damage variable is calculated by $\omega =\frac{A^{*}}{A}=\frac{n^{*}}{n},$ wherein A* is the area of an identified half-hole mark region, A is the total area of the image, n* is a number of pixel points in the half-hole mark region of the image, and n is a total number of pixel points in the image.

7. The half-cast mark identification and damaged flatness evaluation and classification method for blastholes in tunnel blasting according to claim 6, wherein in step S7, ω=0 indicates that the half-cast mark is clear in a damage negligible state; ω=1 indicates that the half-cast mark is diffused in a completely damaged state; and 0<ω<1 indicates that the surrounding rock surface is damaged in different degrees.

8. The half-cast mark identification and damaged flatness evaluation and classification method for blastholes in tunnel blasting according to claim 1, wherein in step S12, a training number is set as 1,000 times; and the NBC is trained repeatedly with a plurality of sets of image data, until a classification error in training is less than 10.sup.−3.

Description

BRIEF DESCRIPTION OF THE DRAWINGS

[0041] FIG. 1 is a flowchart of a half-cast mark identification and damaged flatness evaluation and classification method for blastholes in tunnel blasting according to the present disclosure;

[0042] FIG. 2 illustrates a layout of a blasthole of a typical tunnel blasting;

[0043] FIG. 3 illustrates a third analysis image for a half-cast mark in tunnel blasting;

[0044] FIG. 4 illustrates a first photographed contrast image of a surrounding rock surface with an ideal half-cast mark of a blasthole;

[0045] FIG. 5 illustrates a second photographed contrast image of an unacceptable surrounding rock surface with excessive backbreak or overbreak;

[0046] FIG. 6 illustrates a binary processing result of an analysis image;

[0047] FIG. 7 illustrates an identified result of a GVF-Snake;

[0048] FIG. 8 illustrates an identified result of an ACM model for a half-cast mark;

[0049] FIG. 9 illustrates calculation of a fractal box dimension of an analysis image;

[0050] FIG. 10 illustrates division of a feature extraction direction in a GLCM;

[0051] FIG. 11 illustrates an eigenvector value of an analysis image;

[0052] FIG. 12 is a flowchart of an NBC;

[0053] FIG. 13 illustrates an evaluation and classification result of an NBC; and

[0054] FIG. 14 illustrates a half-cast mark and flatness identified result for a group of blastholes.

DETAILED DESCRIPTION OF THE EMBODIMENTS

[0055] The technical solutions of the present disclosure are further described in detail below with reference to specific embodiments.

[0056] Referring to FIGS. 1-14, a half-cast mark identification and damaged flatness evaluation and classification method for blastholes in tunnel blasting includes the following specific steps:

[0057] S1: Photograph standard contrast images, a size of each of the standard contrast images being comprehensively determined by an analysis scale, and the standard contrast images including a first contrast image of a rock surface with an ideal half-hole mark and a second contrast image that is unacceptable with excessive backbreak or overbreak.

[0058] S2: Acquire a half-cast mark image after actual blasting to serve as a third analysis image, a size and a photographing environment of the analysis image being consistent with those of the above contrast images.

[0059] S3: Denoise the above three images by a 2D Gaussian algorithm according to a characteristic that Gaussian noise in an acquired image obeys a normal distribution. A Gaussian filter based model satisfies the following 2D Gaussian distribution:

[00014]$G\left(x,y\right)=\frac{1}{2{\mathrm{\pi \sigma }}^2}{e}^{-\frac{x^2+y^2}{2{\sigma }^2}},$

where x and y are transverse and longitudinal sizes of the image, and σ.sup.2 is a variance of the normal distribution.

[0060] S4: Perform gray-scale processing on three denoised images according to a preset threshold of an image histogram to reduce an original data volume in each of the images, thereby preventing image distortion to lose information.

[0061] S5: Perform binary processing on gray-scale processed images, and set an optimal binary threshold by maximum entropy thresholding. The maximum entropy thresholding specifically includes: Define, according to a concept of a Shannon entropy, an entropy of an image having a gray scale in a range of 0, 1, 2, . . . , L−1 as

[00015]$\underset{i=0}{\overset{L-1}{.Math.}}{\mathrm{iP}}_i,$

P.sub.i being a probability that a pixel having a value of i occupies the whole image, and divide the image into a target O and a background B according to a set threshold t, probabilities of the pixel in the target and the background being:

[00016]$\mathrm{region}O:\frac{P_i}{P_t},i=0,1,.Math.,t,\mathrm{where},P_t=\underset{i=0}{\overset{t}{.Math.}}{P}_i;\mathrm{and}\mathrm{region}B:\frac{P_i}{1-P_t},i=t+1,.Math.,L-1.$

[0062] In the present disclosure, the target O and the background B each have the following entropy function:

[00017]$H_o\left(t\right)=-\underset{i=0}{\overset{t}{.Math.}}\frac{p_i}{P_t}\ln \frac{p_i}{P_t}=\ln P_t+\frac{H_t}{P_t},\mathrm{where}{H}_t=-\underset{i=0}{\overset{t}{.Math.}}{p}_i\ln p_i;$$H_B\left(t\right)=-\underset{i=t+1}{\overset{L311}{.Math.}}\frac{p_i}{1-P_t}\ln \frac{p_i}{1-P_t}=\ln \left(1-P_t\right)+\frac{H-H_t}{1-P_t},\mathrm{where}H=-\underset{i=0}{\overset{L-1}{.Math.}}{p}_i\ln p_i;$

[0063] a total entropy of the image is:

[00018]$H\left(t\right)=H_o\left(t\right)+H_B\left(t\right)=\ln P_t\left(1-P_t\right)+\frac{H_t}{P_t}+\frac{H-H_t}{1-P_t}$

[0064] the total entropy of the image is maximized through a maximum threshold T:

[0065] T=arg max(H(t)); and

[0066] the optimal binary threshold is determined according to a principle of maximizing the total entropy of the image.

[0067] S6: Identify a boundary and a related region of a half-cast mark in each of the images with a GVF-Snake-ACM.

[0068] S61: Acquire pixel points on a row and a column of an initial image I.sub.m×n (having a resolution of m×n) at a step length d=2, and perform downsampling for s times to obtain s+1 to-be-tested images (including the initial image I.sub.m×n), where a resolution in kth sampling is

[00019]$I_{\frac{m}{2^k}\times \frac{n}{2^k}}.$

[0069] S62: Preprocess an image

[00020]$I_{\frac{m}{2^s}\times \frac{n}{2^s}}$

in sth downsampling, perform edge detection with a Canny operator, define an initial contour l.sub.s+1 on an edge image, a number of control points on the initial contour l.sub.s+1 being N.sub.s+1, and calculate a GVF external force field, where due to a diffusivity of the GVF external force field in a smooth region, the initial contour l.sub.s+1 approaches continuously to an edge of the image; and when 0.98N.sub.s+1 control points on the initial contour l.sub.s+1 reach a target contour (a termination criterion), it is considered that an initial contour curve reaches a target edge, and iteration is stopped to obtain a target contour curve l.sub.s at the resolution

[00021]$\frac{m}{2^s}\times \frac{n}{2^s}.$

[0070] S63: Preprocess an image

[00022]$I_{\frac{m}{2^{s-1}}\times \frac{n}{2^{s-1}}}$

in (s−1)th downsampling with the Canny operator, take the contour curve l.sub.s in the previous step as an initial contour curve, recalculate a GVF external force field on the image

[00023]$I_{\frac{m}{2^{s-1}}\times \frac{n}{2^{s-1}}},$

and obtain a target contour curve l.sub.s−1 according to the termination criterion in Step S62 under pushing of an external force.

[0071] S64: Repeat Step S63, to until a target contour curve l of the identified image l.sub.m×n is obtained.

[0072] S7: Import each of three half-cast mark identified images to ImageJ software, and determine a ratio of an area of a half-cast mark region to a total area of the image as a flatness damage variable: The flatness damage variable is calculated by

[00024]$\omega =\frac{A^{*}}{A}=\frac{n^{*}}{n},$

where A* is the area of the identified half-cast mark region, A is the total area of the image, n* is a number of pixel points in the half-cast mark region of the image, and n is a total number of pixel points in the image. ω=0 indicates that the half-cast mark is clear in a damage negligible state; ω=1 indicates that the half-cast mark is diffused in a completely damaged state; and 0<ω<1 indicates that the surrounding rock surface is damaged in different degrees.

[0073] S8: Normalize an interval for an area ratio of a half-cast mark in each of the first contrast image and the second contrast image, and establish a quantitative relation ω=(D).sup.β between a damage degree and a fractal dimension D, thereby determining a quantitative relation among the area ratio of the half-cast mark region, the damage variable and the fractal dimension.

[0074] S9: Determine a damage value of the third analysis image through linear interpolation calculation on an area ratio of a half-cast mark.

[0075] S10: Extract, for each of a plurality of half-cast mark identified images with a GLCM in 0°, 45°, 90°, and 135° directions or 180°, 225°, 270°, and 315° directions, five eigenvalues including an energy mean, an entropy mean, a contrast mean, a correlation mean, and a uniformity mean, thereby forming a 5D eigenvector.

[0076] S11: Respectively photograph 90 blasting images for a tunnel excavation surrounding rock surface in three blasting plans, which specifically include conventional blasting, presplit blasting and smooth blasting and employ a same explosive charge during tunneling, to form a set C={c.sub.1, c.sub.2, c.sub.3}, analyze target matrices corresponding to the conventional blasting, the presplit blasting and the smooth blasting, and obtain 5D eigenvectors F={f.sub.1, f.sub.2, f.sub.3, f.sub.4, f.sub.5} in different blasting plans, where a feature attribute of each of the images is composed of a 5D eigenvector, including an energy mean, an entropy mean, a contrast mean, a correlation mean, and a uniformity mean, thereby obtaining a multi-dimensional digital information feature of the image.

[0077] S12: Randomly select 5D eigenvectors of 60 images from the image set in the three different blasting plans as training data to input to an NBC, and calculate, with five eigenvalues f.sub.i(i=1,2,3,4,5) of each of the images, a probability that a training image falls into a category c.sub.j(j=1,2,3)

[00025]$P\left(c_j/F\right)=\frac{P\left(c_j\right)P\left(c_j/F\right)}{P\left(F\right)}=\frac{P\left(c_j\right)P\left(f_1,f_2,f_3,f_4,f_5/c_j\right)}{P\left(F\right)},$

[0078] where for an eigenvector F, a larger posterior probability P(c.sub.j/F) indicates a higher probability that the F falls into c.sub.j, and the category into which the F falls is calculated by:

[00026]$c\left(F\right)=\arg \max P\left(c_j\right)\underset{j=1}{\overset{3}{.Math.}}P\left(f_i|c_j\right),$

and

[0079] after P(c.sub.j) and P(c.sub.j/F) are calculated, a category with a maximum posterior probability is returned, thereby obtaining the category of the F. A training number is set as 1,000 times; and the NBC is trained repeatedly with a plurality of sets of image data, until a classification error in training is less than 10.sup.−3.

[0080] S13: Input eigenvectors of remaining 30 images in the image set in the different blasting plans to a well-trained NBC, and determine a state with a maximum probability in different given categories c.sub.j(j=1,2,3) based on an MAP principle of the NBC to take as a final result for evaluation and classification of the blasting plans, thereby implementing classification on blasting effects of the conventional blasting, the presplit blasting and the smooth blasting and on damaged flatness of surrounding rock surfaces.

[0081] The ImageJ software has a version of 1.8.0. The GLCM is implemented by extracting the eigenvalues in 0°, 45°, 90°, and 135° directions.

[0082] According to FIG. 12, the workflow of the NBC can be divided into a preparation stage, a training stage of the NBC, and an application stage. The preparation stage is to determine an eigenvector attribute F and obtain a training sample. The training stage is to determine P(F) for each category, and calculate conditional probabilities of all divisions for each eigenvector attribute. The application stage is to calculate P(c.sub.j)P(c.sub.j/F) for each category, and take a maximum term of the P(c.sub.j)P(c.sub.j/F) as the category to which the eigenvector attribute belongs.

[0083] The foregoing are merely descriptions of preferred specific embodiments of the present disclosure, but the protection scope of the present disclosure is not limited thereto. Any equivalent replacement or modification made within a technical scope of the present disclosure by a person skilled in the art according to the technical solutions of the present disclosure and inventive concepts thereof shall fall within the protection scope of the disclosure.