CN103136758B - Rapid edge detecting method based on orthogonal polynomial fitting - Google Patents

Abstract

The fast edge detection method based on orthogonal polynomial fitting belongs to the field of sub-pixel image edge detection. This two-dimensional edge detection method uses the spatial moment to calculate the projection direction of the edge and combines the projection formula to calculate the projection direction of the edge every time the window image is numerically processed. The one-dimensional edge detection is transformed into one-dimensional edge detection, which avoids the complex operation required for fitting the gray surface. The analytical formula for calculating the position of one-dimensional edge points is obtained by using orthogonal polynomial fitting, which effectively overcomes the problem that the fitting method cannot usually obtain an analytical solution. The point position and edge normal direction determine the position of the edge line in the window image. This method greatly improves the calculation speed and can adapt to the two-dimensional edge detection process with high real-time requirements. It also has the characteristics of fast calculation speed of the moment-based method or interpolation method and the insensitivity of the fitting method to image noise, which can greatly improve the edge detection process. Accuracy and speed of detection.

Description

A Fast Edge Detection Method Based on Orthogonal Polynomial Fitting

technical field

The invention belongs to the field of sub-pixel image edge detection, in particular to a fast edge detection method based on orthogonal polynomial fitting.

Background technique

The most basic feature of an image is the edge. The so-called edge refers to the collection of pixels whose surrounding pixel gray levels have step changes or roof changes. The edge detection methods used in the early days are mostly pixel-level, and the positioning accuracy of these methods is an integer pixel. , that is, it can be judged that the edge position is in a certain pixel, but it cannot be further subdivided in this pixel. However, with the continuous improvement of the measurement accuracy requirements of industrial inspection, the pixel-level accuracy can no longer meet the actual detection requirements. The pixel method came into being.

The existing sub-pixel-level edge detection methods can be summarized into three types: the method based on moments, the method based on least squares and the method based on interpolation. Insensitive, so the detection edge is relatively accurate, but due to the need for nonlinear fitting operations, the calculation speed is relatively slow, which cannot meet the processing and monitoring occasions of real-time detection. The moment-based method or the interpolation-based method both have faster calculation speeds, but these two detection methods are more sensitive to image noise. In the case of more image noise, the detection accuracy of the two will be greatly affected.

There are many existing technologies that are relatively close to the present invention, such as "A Subpixel Subdivision Algorithm Based on Polynomial Interpolation Improvement" by Li Qingli et al. [J]. Journal of Beijing University of Science and Technology, 2003, 25 (3): 280-283 , this method applies the direction template constructed on the basis of the classic Sobel operator to process the gray image to obtain a gradient image, and then performs sub-pixel subdivision calculations on the gradient image along the gradient direction of the target edge by polynomial interpolation method. Sub-pixel accurate localization of object edges. However, this method inherits the characteristic that the interpolation method is sensitive to image noise, and its edge detection accuracy will be seriously disturbed when the image is noisy.

Wang Jianmin et al. in "Research on Sub-pixel Subdivision Algorithm of Spatial Moment" [J] Optical Technology. 1999, No, 4, 3-6. A sub-pixel edge detection algorithm based on spatial moments is proposed, which can achieve sub-pixel edge detection quickly. The basic model of this method is a step-edge model. However, in the actual imaging process, affected by imaging noise, the edge model is often not strictly step-like, but fuzzy. Therefore, this method is easily disturbed by noise, which affects the edge positioning accuracy.

Another example, Chen Jing et al. in "Research on Fast Polynomial Fitting Sub-pixel Edge Detection Algorithm" [J]. Applied Optics 2011, 32(1):91-95.DOI:10.3969/j.issn.1002-2082.2011.01.020 A fast sub-pixel edge detection algorithm based on polynomial fitting is introduced in this paper. According to the gray distribution of the image, the algorithm uses a cubic polynomial to fit the edge of the image to achieve sub-pixel positioning, so as to achieve high positioning accuracy. Although this method also uses polynomial fitting to eliminate noise interference, the calculation formula it gives contains more complicated inner product operations, which will lead to a slow solution process. Therefore, the solution process of this detection method is slow and cannot adapt to the Real-time edge detection with high effectiveness requirements. In addition, two-dimensional edge detection is the main application method to realize sub-pixel edge localization, but this method is limited to the discussion of one-dimensional edge detection method, and does not give how to apply its one-dimensional detection method to two-dimensional edge detection method.

Sun Qiucheng et al. proposed a sub-pixel edge detection algorithm based on the Bessel edge fitting model in "A Sub-pixel Accurate Edge Detection Method", Journal of Beijing University of Technology. 2009.10.35(10), 1332-1337. By fitting and correcting the Bessel edge model, the sub-pixel position of the image edge is obtained. The fitting method has good robustness to noise and high edge location accuracy. However, since there is no analytical formula for calculating the edge, a fitting operation is required in each detection process, which will cause the detection speed of the fitting method to be slow, especially when the method is applied to two-dimensional edge detection, it is necessary to fit the three-dimensional grayscale For curved surfaces, the detection speed will be slower due to the complexity of the fitting function.

Contents of the invention

In order to solve the existing sub-pixel edge detection methods, the methods based on moments or interpolation are relatively susceptible to interference from image noise, and the existing fitting methods do not have an analytical formula for calculating edge points, which leads to In the process of two-dimensional edge detection, fitting operation is required, especially in two-dimensional edge detection, surface fitting operation is required, which makes the operation process complicated, the solution speed is slow, and it cannot adapt to the two-dimensional edge with high real-time requirements. To solve the technical problem of detection, the present invention provides a fast edge detection method based on orthogonal polynomial fitting.

The technical solution adopted by the present invention to solve the technical problems is as follows:

A fast edge detection method based on orthogonal polynomial fitting includes the following steps:

Step 1: Image integer pixel edge detection and positioning: Use the Canny edge detection operator to perform edge detection on the image, and obtain the coordinates and gray levels of all integer pixel squares where the image edge is located;

Step 2: Extract the window image: take the integer pixel square containing the edge described in step 1 as the center, extract a local square grayscale matrix A _m×m of the area where the edge of the image is located, that is, the window image, and the window image contains a small segment Edge line to be detected;

Step 3: Numerical processing of the window image, specifically including the following steps:

Step 3.1: Taking the center of the window image described in step 2 as the origin, establish a xoy rectangular coordinate system, wherein the horizontal axis is the x axis and the vertical axis is the y axis; the unit length of the coordinate system is 1 pixel;

Step 3.2: Set the coordinates of the center point of each pixel grid in the window image in the Cartesian coordinate system as described in step 3.1 to be (x, y), then use the gray value of the entire pixel grid as the gray value corresponding to the point value, and add the Z axis perpendicular to the plane rectangular coordinate system to the original plane rectangular coordinate system to form a spatial three-dimensional coordinate system, and the Z axis represents the gray value corresponding to the pixel point, thus obtaining the discrete value of the gray value changing with the coordinate point Gray distribution function I(x,y);

Step 3.3: There is a straight line passing through the origin of the coordinate system and perpendicular to the straight line of the edge to be detected in the plane of the rectangular coordinate system. This straight line is the normal line of the edge, and the angle between the normal line of the edge and the horizontal horizontal axis is defined as θ;

Step 4: Determine the direction of the edge normal, the specific steps are:

Step 4.1: Using the calculation formula of spatial moments, calculate the spatial moments M _pq of each order corresponding to the inscribed circle of the square window image extracted in step 3.1:

${\mathrm{<span\; class="notranslate">\; m</span>}}_{\mathrm{<span\; class="notranslate">\; pq</span>}}<span\; class="notranslate">\; =</span>\underset{{\mathrm{<span\; class="notranslate">\; x</span>}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; the\; y</span>}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}<span\; class="notranslate">\; \&le;</span>\mathrm{<span\; class="notranslate">\; r</span>}}{<span\; class="notranslate">\; \&Integral;</span><span\; class="notranslate">\; \&Integral;</span>}{\mathrm{<span\; class="notranslate">\; x</span>}}^{\mathrm{<span\; class="notranslate">\; p</span>}}{\mathrm{<span\; class="notranslate">\; the\; y</span>}}^{\mathrm{<span\; class="notranslate">\; q</span>}}\mathrm{<span\; class="notranslate">\; f</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; x</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; the\; y</span>}<span\; class="notranslate">\; )</span>\mathrm{<span\; class="notranslate">\; dxdy</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; =</span>\mathrm{<span\; class="notranslate">\; 0,1</span>}<span\; class="notranslate">\; ,</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; N</span>}<span\; class="notranslate">\; \&Center\; Dot;</span><span\; class="notranslate">\; \&Center\; Dot;</span><span\; class="notranslate">\; \&Center\; Dot;</span><span\; class="notranslate">\; \&Center\; Dot;</span><span\; class="notranslate">\; \&CenterDot;</span><span\; class="notranslate">\; \&Center\; Dot;</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; )</span>$

In formula (1), f(x, y) represents the continuous gray scale distribution function of image pixels, p and q represent the order of spatial moments, which are positive integers, and r represents the radius of the inscribed circle;

Step 4.2: Rotate the continuous grayscale distribution function f(x, y) described in step 4.1 around the Z axis described in step 3.2 by the angle θ as described in step 3.3, and a new continuous grayscale distribution will be obtained in the window image function, at this time, the value of the spatial moment of each order corresponding to the inscribed circle of the window image changes;

Step 4.3: According to the theory of spatial moment rotation invariance, the relationship between the spatial moment M' _pq corresponding to the inscribed circle after rotation and the spatial moment M _pq before rotation is obtained as follows:

${\mathrm{<span\; class="notranslate">\; m</span>}}_{\mathrm{<span\; class="notranslate">\; pq</span>}}^{<span\; class="notranslate">\; \&prime;</span>}<span\; class="notranslate">\; =</span>\underset{\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; =</span>\mathrm{<span\; class="notranslate">\; 0</span>}}{\overset{\mathrm{<span\; class="notranslate">\; p</span>}}{\mathrm{<span\; class="notranslate">\; \&Sigma;</span>}}}\underset{\mathrm{<span\; class="notranslate">\; the\; s</span>}<span\; class="notranslate">\; =</span>\mathrm{<span\; class="notranslate">\; 0</span>}}{\overset{\mathrm{<span\; class="notranslate">\; q</span>}}{\mathrm{<span\; class="notranslate">\; \&Sigma;</span>}}}\left(\begin{array}{c}\mathrm{<span\; class="notranslate">\; p</span>}\\ \mathrm{<span\; class="notranslate">\; i</span>}\end{array}\right)\left(\begin{array}{c}\mathrm{<span\; class="notranslate">\; q</span>}\\ \mathrm{<span\; class="notranslate">\; the\; s</span>}\end{array}\right){<span\; class="notranslate">\; (</span><span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; the\; s</span>}}{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; cos</span>}\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; +</span>\mathrm{<span\; class="notranslate">\; the\; s</span>}}<span\; class="notranslate">\; \&times;</span>{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; sin</span>}\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; +</span>\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; the\; s</span>}}{\mathrm{<span\; class="notranslate">\; m</span>}}_{\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; +</span>\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; the\; s</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; +</span>\mathrm{<span\; class="notranslate">\; the\; s</span>}}<span\; class="notranslate">\; \&CenterDot;</span><span\; class="notranslate">\; \&CenterDot;</span><span\; class="notranslate">\; \&CenterDot;</span><span\; class="notranslate">\; \&CenterDot;</span><span\; class="notranslate">\; \&CenterDot;</span><span\; class="notranslate">\; \&CenterDot;</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 2</span>}<span\; class="notranslate">\; )</span>$

i represents an integer from 0 to p, s represents an integer from 0 to q, when the continuous gray distribution function f(x, y) is rotated by the angle θ, the edge line to be detected in the window image and the coordinates described in step 3.1 The horizontal axis x is vertical, that is, the grayscale function corresponding to the inscribed circle of the window image after rotation is symmetrical about the horizontal axis x. According to the calculation formula of the space moment, the space moment M' ₀₁ corresponding to the inscribed circle of the window image after rotation is 0 = 0. When p=0, q=1, use formula (2) to get:

M' ₀₁ =-sinθ·M ₁₀ +cosθ·M ₀₁ =0...(3)

Among them, M ₁₀ and M ₀₁ represent the space moment corresponding to the inscribed circle before rotation;

Step 4.4: Use the spatial moment calculation template to replace formula (1), and substitute the discrete grayscale function I(x, y) of the pixels in the inscribed circle of the window image into the above discrete template formula to calculate formula (3) The values of spatial moments M ₀₁ and M ₁₀ in ;

Step 4.5: Obtain formula (4) according to formula (3), and substitute the values of spatial moments M ₀₁ and M ₁₀ obtained in step 4.4 into formula (4) to obtain the value of angle θ:

$\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; =</span>\mathrm{<span\; class="notranslate">\; arctan</span>}\frac{{\mathrm{<span\; class="notranslate">\; m</span>}}_{\mathrm{<span\; class="notranslate">\; 01</span>}}}{{\mathrm{<span\; class="notranslate">\; m</span>}}_{\mathrm{<span\; class="notranslate">\; 10</span>}}}<span\; class="notranslate">\; \&CenterDot;</span><span\; class="notranslate">\; \&CenterDot;</span><span\; class="notranslate">\; \&CenterDot;</span><span\; class="notranslate">\; \&Center\; Dot;</span><span\; class="notranslate">\; \&Center\; Dot;</span><span\; class="notranslate">\; \&Center\; Dot;</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 4</span>}<span\; class="notranslate">\; )</span>$

Step 4.6: Set the edge normal as the new horizontal horizontal axis X, and take the θ angle in step 4.5 as the projection angle, then the pixel coordinates (x, y) of all points described in step 3.2 in the pixel plane and the edge normal The point coordinate X on the horizontal axis has the following projection formula:

X=xcosθ+ysinθ...(5)

The X value in formula (5) represents the coordinate reading value of the projected point on the new horizontal axis X described in step 4.6 after the projected point in the window image is projected according to the projection angle θ;

Step 5: Based on one-dimensional edge point detection based on orthogonal polynomial fitting, determine the intersection position of the edge normal and the edge line to be detected, which includes the following steps:

Step 5.1: reconstitute the Z-axis described in step 3.2 and the horizontal horizontal axis X described in step 4.6 into a plane Cartesian coordinate system XOZ required for one-dimensional edge point detection, that is, the projection plane;

According to the projection formula (5), the center point of each pixel grid in the window image is projected onto the horizontal horizontal axis X described in step 4.6; in the plane Cartesian coordinate system XOZ, the discrete gray value of the projected point on the edge normal is obtained Gray distribution function (X _i , Z _i ), i=1,...,2n+1;

Step 5.2: Normalization of the projection point coordinates, if the number of projection points formed on the horizontal axis X as described in step 5.1 is 2n+1, set $L=\frac{x_{2\mathrm{no}+1}-x_1}{2},$ and order $x_i=\frac{x_i}{L}i=1,...,2\mathrm{no}+1,$ Then there are x _i ∈ [-1,1];

Step 5.3: Calculate the orthogonal polynomial system, and use the schemite orthogonalization method to calculate the orthogonal polynomial system about the weight function ρ(x)=1 on the [-1,1] interval;

Step 5.4: Take the first four order orthogonal basis functions in step 5.3

Using formula (7) as the fitting basis to carry out cubic orthogonal polynomial fitting, the equation of fitting polynomial is:

The coefficients {a ₁ , a ₂ , a ₃ , a ₄ } of the fitting polynomial in formula (8) are the fitting parameters to be found;

Step 5.5: Using formula (8), the normal equations corresponding to the equations for fitting polynomials can be deduced as follows:

In formula (9) represents a discrete vector Inner product with itself, Z=[Z ₁ ,Z ₂ ,……,Z _2n+1 ] ^T is a discrete gray value vector, Indicates that Z and inner product;

Step 5.6: Use equation (9) to obtain the calculation formula of the coefficient {a ₁ , a ₂ , a ₃ , a ₄ } of the fitting polynomial:

Substituting each known quantity into formula (10), the results of a ₀ , a ₁ , a ₂ and a ₃ are obtained respectively;

Step 5.7: Substitute the result of the polynomial fitting described in step 5.6 into equation (8) as a known quantity, and the polynomial fitting equation becomes:

Step 5.8: Compute the second order derivative of the cubic polynomial (11) to zero, namely:

Solve and obtain the coordinates of the point where the second order derivative of the above formula (12) is zero;

Step 5.9: Solve formula (12) and denormalize the calculation results to obtain the position coordinate R of the one-dimensional edge point on the horizontal axis X, namely:

Equation (13) is the analytical formula for calculating the position of a one-dimensional edge point, which is the intersection point of the edge normal and the edge line to be detected;

Step 6: Use formula (13) to determine the intersection position of the edge normal passing through the origin in the window image and the edge straight line to be detected, use step 4.5 to obtain the direction angle θ of the edge normal, and determine the edge straight line to be detected in the xoy coordinate system equation, so as to realize the sub-pixel positioning of the edge line to be detected in the window image;

Step 7: The operator reselects the whole pixel grid described in step 1, and repeats steps 2 to 6 to obtain the sub-pixel positioning of the edge line to be detected in the newly obtained window image;

Step 8: When the operator determines that the sub-pixel positioning of the desired edge position on the entire image has been completed, the detection method ends.

The beneficial effects of the present invention are: the fast edge detection method of the present invention can be directly applied to two-dimensional edge detection, the method uses the spatial moment to calculate the projection direction of the edge each time the window image is numerically processed, and combines the projection formula to The two-dimensional edge detection is transformed into one-dimensional edge detection, which avoids repeating the complex operation required for fitting the gray surface in each one-dimensional edge detection process. The analytical formula for calculating the position of one-dimensional edge points is obtained by using orthogonal polynomial fitting, which effectively overcomes the problem that the fitting method usually cannot obtain analytical solutions. The analytical formula is simple in form and contains less inner product calculations. The position of the edge point and the direction of the edge normal determine the position of the edge line in the window image. This method greatly improves the calculation speed and can adapt to the two-dimensional edge detection process with high real-time requirements. It also has the characteristics of fast calculation speed of the moment-based method or interpolation method and the insensitivity of the fitting method to image noise, which can greatly improve the edge detection process. Accuracy and speed of detection.

Description of drawings

Fig. 1 is the flow chart of the fast edge detection method based on orthogonal polynomial fitting of the present invention;

Fig. 2 is the flowchart of step 3 of the present invention;

Fig. 3 is the flowchart of step 4 of the present invention;

Fig. 4 is the flowchart of step 5 of the present invention;

Fig. 5 is a schematic diagram of the application of projection of pixels on the window image to the horizontal axis X where the edge normal is located according to the projection angle θ.

Detailed ways

The present invention will be described in further detail below in conjunction with the accompanying drawings.

As shown in Figure 1, the fast edge detection method based on orthogonal polynomial fitting of the present invention comprises the following steps:

Step 1: Image integer pixel edge detection and positioning: use the Canny edge detection operator to perform edge detection on the image, and obtain the coordinates and gray levels of all the integer pixel squares where the image edge is located, so as to locate the pixel coordinates of the edge to be detected in the image The general orientation and range in the plane.

Step 2: Extracting the window image: centering on the whole-pixel square containing the edge described in Step 1, extract a local square grayscale matrix A _m×m of the area where the edge of the image is located, that is, the window image 2 . At the microscopic level, the edge to be detected in the window image 2 can be approximately regarded as a small straight line. As shown in FIG. 5 , in this embodiment, an image area formed by 5×5 integer pixel squares is used as a window image 2 , and a square image containing a small section of edge line 1 to be detected is extracted from the image.

Step 3: Numerical processing of window image 2, specifically including the following sub-steps as shown in Figure 2:

Step 3.1: As shown in Figure 5, take the center of the window image 2 described in step 2 as the origin, establish a xoy rectangular coordinate system, where the horizontal axis is the x axis and the vertical axis is the y axis. The unit length of the coordinate system is 1 pixel, then it can be seen from step 1 that the edge to be detected is located in the coordinate system.

Step 3.2: As shown in Figure 5, set the coordinates of the center point of each pixel grid in the window image 2 in the Cartesian coordinate system as described in step 3.1 to be (x, y), then set the gray value of the entire pixel grid As the gray value corresponding to the point, and adding the Z axis to the original xoy plane Cartesian coordinate system, it represents the gray value corresponding to the pixel point, thus obtaining the discrete gray distribution function I(x ,y).

Step 3.3: As shown in Figure 5, there is a straight line passing through the origin o of the coordinate system and perpendicular to the edge to be detected in the plane of the xoy rectangular coordinate system. The angle between the axes is θ.

Step 4: Determine the direction of the edge normal 3, which specifically includes the following sub-steps as shown in Figure 3:

Step 4.1: Using the spatial moment calculation formula, the spatial moment M _pq of each order corresponding to the inscribed circle of the square window image 2 extracted in step 3.1 can be calculated:

${\mathrm{<span\; class="notranslate">\; m</span>}}_{\mathrm{<span\; class="notranslate">\; pq</span>}}<span\; class="notranslate">\; =</span>\underset{{\mathrm{<span\; class="notranslate">\; x</span>}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; the\; y</span>}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}<span\; class="notranslate">\; \&le;</span>\mathrm{<span\; class="notranslate">\; r</span>}}{<span\; class="notranslate">\; \&Integral;</span><span\; class="notranslate">\; \&Integral;</span>}{\mathrm{<span\; class="notranslate">\; x</span>}}^{\mathrm{<span\; class="notranslate">\; p</span>}}{\mathrm{<span\; class="notranslate">\; the\; y</span>}}^{\mathrm{<span\; class="notranslate">\; q</span>}}\mathrm{<span\; class="notranslate">\; f</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; x</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; the\; y</span>}<span\; class="notranslate">\; )</span>\mathrm{<span\; class="notranslate">\; dxdy</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; =</span>\mathrm{<span\; class="notranslate">\; 0,1</span>}<span\; class="notranslate">\; ,</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; N</span>}<span\; class="notranslate">\; \&Center\; Dot;</span><span\; class="notranslate">\; \&Center\; Dot;</span><span\; class="notranslate">\; \&Center\; Dot;</span><span\; class="notranslate">\; \&CenterDot;</span><span\; class="notranslate">\; \&Center\; Dot;</span><span\; class="notranslate">\; \&Center\; Dot;</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; )</span>$

In formula (1), f(x, y) represents the continuous distribution function of image pixel grayscale, positive integers p and q represent the order of spatial moments, and r represents the radius of the inscribed circle.

Step 4.2: Rotate the continuous grayscale distribution function f(x, y) described in step 4.1 around the Z axis described in step 3.2 by the angle θ as described in step 3.3, and a new continuous grayscale will be obtained in window image 2 Distribution function, at this time, the value of the space moment of each order corresponding to the inscribed circle of window image 2 changes.

Step 4.3: According to the theory of spatial moment rotation invariance, the relationship between the spatial moment M' _pq corresponding to the inscribed circle after rotation and the moment M _pq before rotation can be obtained as follows:

${\mathrm{<span\; class="notranslate">\; m</span>}}_{\mathrm{<span\; class="notranslate">\; pq</span>}}^{<span\; class="notranslate">\; \&prime;</span>}<span\; class="notranslate">\; =</span>\underset{\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; =</span>\mathrm{<span\; class="notranslate">\; 0</span>}}{\overset{\mathrm{<span\; class="notranslate">\; p</span>}}{\mathrm{<span\; class="notranslate">\; \&Sigma;</span>}}}\underset{\mathrm{<span\; class="notranslate">\; the\; s</span>}<span\; class="notranslate">\; =</span>\mathrm{<span\; class="notranslate">\; 0</span>}}{\overset{\mathrm{<span\; class="notranslate">\; q</span>}}{\mathrm{<span\; class="notranslate">\; \&Sigma;</span>}}}\left(\begin{array}{c}\mathrm{<span\; class="notranslate">\; p</span>}\\ \mathrm{<span\; class="notranslate">\; i</span>}\end{array}\right)\left(\begin{array}{c}\mathrm{<span\; class="notranslate">\; q</span>}\\ \mathrm{<span\; class="notranslate">\; the\; s</span>}\end{array}\right){<span\; class="notranslate">\; (</span><span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; the\; s</span>}}{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; cos</span>}\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; +</span>\mathrm{<span\; class="notranslate">\; the\; s</span>}}<span\; class="notranslate">\; \&times;</span>{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; sin</span>}\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; +</span>\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; the\; s</span>}}{\mathrm{<span\; class="notranslate">\; m</span>}}_{\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; +</span>\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; the\; s</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; +</span>\mathrm{<span\; class="notranslate">\; the\; s</span>}}<span\; class="notranslate">\; \&CenterDot;</span><span\; class="notranslate">\; \&Center\; Dot;</span><span\; class="notranslate">\; \&Center\; Dot;</span><span\; class="notranslate">\; \&Center\; Dot;</span><span\; class="notranslate">\; \&CenterDot;</span><span\; class="notranslate">\; \&Center\; Dot;</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 2</span>}<span\; class="notranslate">\; )</span>$

i represents an integer from 0 to p, s represents an integer from 0 to q, when the continuous gray distribution function f(x, y) is rotated by the angle θ as shown in Figure 5, the edge to be detected in the window image 2 The straight line 1 is perpendicular to the horizontal axis x of the coordinate system described in step 3.1, that is, the grayscale function corresponding to the inscribed circle of the window image 2 after rotation is symmetrical about the horizontal axis x. According to the calculation formula of the space moment, the inscribed circle of the window image 2 after rotation corresponds to The spatial moment M' ₀₁ =0, at this time p=0, q=1, using formula (2) can also be obtained:

M' ₀₁ =-sinθ·M ₁₀ +cosθ·M ₀₁ =0...(3)

Where M ₁₀ , (p=1, q=0) and M ₀₁ , (p=0, q=1) represent the spatial moment corresponding to the inscribed circle before rotation.

Step 4.4: Since the pixel grayscale value I(x,y) in the window image is discrete, that is, there is no continuous grayscale function f(x,y), the spatial moment calculation template is used to replace the formula (1) , and substitute the discrete gray value I(x,y) of the pixel in the inscribed circle of the window image into the above-mentioned space moment calculation template formula to calculate the values of space moments M ₀₁ and M ₁₀ in formula (3) . The spatial moment calculation template is provided by LyversEP, Mitchell OR. in their book Precisionedgecontrastandorientationestimation[J], IEEETrans.PatternAnal.Mach.Intell, 1988, 10(6):927-937.[1].

Step 4.5: Obtain formula (4) according to formula (3), and substitute the values of spatial moments M ₀₁ and M ₁₀ obtained in step 4.4 into formula (4) to obtain the value of angle θ:

$\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; =</span>\mathrm{<span\; class="notranslate">\; arctan</span>}\frac{{\mathrm{<span\; class="notranslate">\; m</span>}}_{\mathrm{<span\; class="notranslate">\; 01</span>}}}{{\mathrm{<span\; class="notranslate">\; m</span>}}_{\mathrm{<span\; class="notranslate">\; 10</span>}}}<span\; class="notranslate">\; \&Center\; Dot;</span><span\; class="notranslate">\; \&CenterDot;</span><span\; class="notranslate">\; \&CenterDot;</span><span\; class="notranslate">\; \&CenterDot;</span><span\; class="notranslate">\; \&Center\; Dot;</span><span\; class="notranslate">\; \&Center\; Dot;</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 4</span>}<span\; class="notranslate">\; )</span>$

After the θ angle is obtained, the position of the edge normal 3 in the xoy coordinate system as described in step 3.3 is determined and becomes a known quantity.

Step 4.6: As shown in Figure 5, set the edge normal 3 as the new horizontal axis X, and take the θ angle in step 4.5 as the projection angle, then the pixel coordinates of all points in the pixel plane described in step 3.2 (x , y) and the point coordinate X on the horizontal axis of the edge normal 3 have the following projection formula:

X=xcosθ+ysinθ...(5)

The X value in formula (5) represents the coordinate reading value of the corresponding projected point 4 on the new horizontal axis X described in step 4.6 after the projected point in the window image 2 is projected according to the projection angle θ. Since the continuous grayscale distribution function f(x, y) does not change in this projection direction, the grayscale value of the point on the new horizontal axis is the same as the grayscale value 5 of the corresponding projected point.

Through step four, all the pixels in the two-dimensional sampling window can be projected onto the edge normal of the origin of the window image matrix, and then the two-dimensional edge detection problem is transformed into one-dimensional edge detection problem, and then the following steps are used: The position of the edge line can be located and obtained.

Step 5: Carry out one-dimensional edge point detection based on orthogonal polynomial fitting to determine the intersection position of the edge normal 3 and the edge line 1 to be detected, which includes the following sub-steps as shown in Figure 4:

Step 5.1: As shown in Figure 5, reconstruct the Z axis described in step 3.2 and the horizontal horizontal axis X described in step 4.6 into a plane Cartesian coordinate system XOZ required for one-dimensional edge point detection, that is, the projection plane 6, the new XOZ coordinates The system is a one-dimensional plane coordinate system, and its XOZ origin coincides with the origin of the original xoy coordinate system. Therefore, the coordinates of projected point 4 in the new coordinates are (X, Z).

According to the projection formula (5), project the center point of each pixel grid in the window image 2 onto the horizontal horizontal axis X described in step 4.6, since the gray values corresponding to the pixels along the projection direction θ in the window image are the same , so the gray value of the projected point is the same as the projected point, so in the plane Cartesian coordinate system XOZ, the discrete distribution function (X _i ,Z _i ), i =1,...,2n+1.

Step 5.2: Normalization of the projection point coordinates. If the number of projection points 4 on the horizontal axis X as described in step 5.1 is 2n+1, according to the symmetry of the square matrix projection, the above projection points 4 are on the horizontal axis Both sides of X are symmetrically distributed. if set $L=\frac{x_{2\mathrm{no}+1}-x_1}{2},$ and order $x_i=\frac{x_i}{L}i=1,...,2\mathrm{no}+1,$ Then there is x _i ∈ [-1,1], so that the range of the projection point on the X axis is [-1,1].

Step 5.3: Calculate the orthogonal polynomial system, use the schemite orthogonalization method to calculate the orthogonal polynomial system about the weight function ρ(x)=1 on the [-1,1] interval, the polynomial system is as follows:

Step 5.4: Take the first four order orthogonal basis functions in step 5.3

Using equation (7) as the fitting base to carry out cubic orthogonal polynomial fitting, the equation of fitting polynomial is

The coefficients {a ₁ , a ₂ , a ₃ , a ₄ } of the fitting polynomial in formula (8) are the fitting parameters to be obtained.

Step 5.5: According to the least squares fitting principle and the orthogonality property of orthogonal polynomials, the normal equations corresponding to the equations of fitting polynomials can be deduced as follows:

In formula (9) represents a discrete vector Inner product with itself, Z=[Z ₁ ,Z ₂ ,……,Z _2n+1 ] ^T is a discrete gray value vector, Indicates that Z and inner product.

Step 5.6: Use the equation (9) to obtain the calculation formula of the coefficients {a ₁ , a ₂ , a ₃ , a ₄ } of the fitting polynomial:

Substituting each known quantity into formula (10), the results of a ₀ , a ₁ , a ₂ and a ₃ are obtained respectively.

Step 5.7: Substitute the result of the polynomial fitting described in step 5.6 into equation (8) as a known quantity, and the polynomial fitting equation becomes:

Step 5.8: Since equation (11) is a cubic polynomial, there is only a second-order derivative that is zero, and the second-order derivative of the cubic polynomial (11) is zero, namely:

Solve and obtain the coordinates of the above point where the second derivative of equation (12) is zero.

Step 5.9: Solve formula (12) and denormalize the calculation results to obtain the position coordinate R of the one-dimensional edge point on the horizontal axis X, namely:

Equation (13) is an analytical formula for calculating the position of a one-dimensional edge point, which is the intersection point of the edge normal 3 and the edge line 1 to be detected.

Through the method described in step five, we can use this formula to calculate the position of the edge point on the X-axis, that is, the intersection position of the edge line and the normal.

Step 6: Use formula (13) to determine the intersection position of the edge normal 3 passing through the origin in the window image 2 and the edge line 1 to be detected, use step 4.5 to obtain the direction angle θ of the edge normal 3, and determine the edge line 1 to be detected at The equation in the xoy coordinate system, thereby realizing the sub-pixel positioning of the edge straight line 1 to be detected in the window image 2;

Step 7: The operator reselects the integer pixel grid described in step 1, and repeats steps 2 to 6 to obtain the sub-pixel positioning of the edge line 1 to be detected in the newly obtained window image 2;

Step 8: When the operator determines that the sub-pixel positioning of the desired edge position on the entire image has been completed, the detection method ends.

The fast edge detection method based on orthogonal polynomial fitting of the present invention can be directly applied to two-dimensional edge detection. It uses the invariance of space moment rotation to obtain the projection direction θ of the edge, and obtains the projection formula, so that the window image Pixels and their corresponding gray values are projected to the plane with the edge normal as the horizontal axis and the gray value as the vertical axis through the projection formula, and then the two-dimensional edge detection problem is transformed into a one-dimensional edge detection problem, through The one-dimensional edge detection method based on orthogonal polynomial fitting obtains the coordinate position of the intersection point of the edge normal and the edge in the window image, and finally combines the detection to obtain the edge normal direction θ to determine the position of the edge line in the window image.

In the detection method of the present invention, during the numerical processing of each window image, the spatial moment theory is used to first calculate the projection direction θ of the edge, and the two-dimensional edge detection problem is converted into a one-dimensional edge detection problem in combination with the projection formula, avoiding the problem of The complex operations required to fit the grayscale surface are repeated in each one-dimensional edge detection process. According to the theory of orthogonal polynomial fitting, an analytical formula for calculating the position of one-dimensional edge points can be obtained, which effectively overcomes the problem that the fitting method usually cannot obtain an analytical solution. Using this analytical formula, only the gray value of the edge position needs to be brought in. The analytical solution of the position of the relevant edge points can be obtained very conveniently. Since the analytical formula is simple in form and contains less inner product calculations, the calculation speed is greatly improved, so it can adapt to the edge detection process of real-time working conditions with high requirements for effectiveness.

In addition, the detection method of the present invention well integrates two sub-pixel detection methods based on moment and fitting method, not only inherits the characteristics of fast calculation speed of the method based on moment or interpolation method, but also has the advantages of fitting method for image The advantage of noise insensitivity can better reduce the interference and influence of image noise on the result, thereby greatly improving the accuracy and speed of edge detection.

Claims (1)

1. based on the fast edge detection method of orthogonal polynomial fitting, it is characterized in that: the method comprises the steps: Step 1: Image integer pixel edge detection and positioning: use the Canny edge detection operator to perform edge detection on the image, and obtain the coordinates and gray levels of all integer pixel squares where the image edge is located; Step 2: Extract the window image: take the integer pixel square with the edge obtained in Step 1 as the center, extract a local square grayscale matrix A _m×m of the area where the edge of the image is located, that is, the window image (2), the window image (2 ) contains a small section of edge straight line (1) to be detected; Step 3: numerical processing of the window image (2), specifically including the following steps: Step 3.1: take the center of the window image (2) described in step 2 as the origin, establish a xoy rectangular coordinate system, wherein the horizontal axis is the x axis, and the vertical axis is the y axis; the unit length of the coordinate system is 1 pixel; Step 3.2: Set the coordinates of the center point of each pixel grid in the window image (2) in the Cartesian coordinate system as described in step 3.1 to be (x, y), then use the gray value of the entire pixel grid as the corresponding point gray value, and add the Z axis perpendicular to the plane rectangular coordinate system to the original plane rectangular coordinate system to form a spatial three-dimensional coordinate system. Varying discrete grayscale distribution function I(x,y); Step 3.3: There is a straight line passing through the origin of the coordinate system and perpendicular to the edge line (1) to be detected in the plane of the rectangular coordinate system. This line is the normal line of the edge, and the edge normal line (3) and the horizontal horizontal axis are defined The included angle is θ; Step 4: Determine the direction of the edge normal (3), the specific steps are: Step 4.1: Using the calculation formula of spatial moments, calculate the spatial moments M _pq of each order corresponding to the inscribed circle of the square window image (2) extracted in step 3.1: ${\mathrm{<span\; class="notranslate">\; m</span>}}_{\mathrm{<span\; class="notranslate">\; pq</span>}}<span\; class="notranslate">\; =</span>\underset{{\mathrm{<span\; class="notranslate">\; x</span>}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; the\; y</span>}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}<span\; class="notranslate">\; \&le;</span>\mathrm{<span\; class="notranslate">\; r</span>}}{<span\; class="notranslate">\; \&Integral;</span><span\; class="notranslate">\; \&Integral;</span>}{\mathrm{<span\; class="notranslate">\; x</span>}}^{\mathrm{<span\; class="notranslate">\; p</span>}}{\mathrm{<span\; class="notranslate">\; the\; y</span>}}^{\mathrm{<span\; class="notranslate">\; q</span>}}\mathrm{<span\; class="notranslate">\; f</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; x</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; the\; y</span>}<span\; class="notranslate">\; )</span>\mathrm{<span\; class="notranslate">\; dxdu</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; =</span>\mathrm{<span\; class="notranslate">\; 0,1</span>}<span\; class="notranslate">\; ,</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; N</span>}<span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; )</span>$ In formula (1), f(x, y) represents the continuous gray scale distribution function of image pixels, positive integers p and q represent the order of spatial moments, and r represents the radius of the inscribed circle; Step 4.2: Rotate the continuous grayscale distribution function f(x, y) described in step 4.1 around the Z axis described in step 3.2 by the θ angle described in step 3.3, and a new continuous grayscale distribution function f(x, y) will be obtained in the window image (2). Gray distribution function, at this time, the value of the space moment of each order corresponding to the inscribed circle of the window image (2) changes; Step 4.3: The relationship between the spatial moment M' _pq corresponding to the inscribed circle after rotation and the spatial moment M _pq before rotation is as follows: ${\mathrm{<span\; class="notranslate">\; m</span>}}_{\mathrm{<span\; class="notranslate">\; pq</span>}}^{<span\; class="notranslate">\; \&prime;</span>}<span\; class="notranslate">\; =</span>\underset{\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; =</span>\mathrm{<span\; class="notranslate">\; 0</span>}}{\overset{\mathrm{<span\; class="notranslate">\; p</span>}}{\mathrm{<span\; class="notranslate">\; \&Sigma;</span>}}}\underset{\mathrm{<span\; class="notranslate">\; the\; s</span>}<span\; class="notranslate">\; =</span>\mathrm{<span\; class="notranslate">\; 0</span>}}{\overset{\mathrm{<span\; class="notranslate">\; q</span>}}{\mathrm{<span\; class="notranslate">\; \&Sigma;</span>}}}\left(\begin{array}{c}\mathrm{<span\; class="notranslate">\; p</span>}\\ \mathrm{<span\; class="notranslate">\; i</span>}\end{array}\right)\left(\begin{array}{c}\mathrm{<span\; class="notranslate">\; q</span>}\\ \mathrm{<span\; class="notranslate">\; the\; s</span>}\end{array}\right){<span\; class="notranslate">\; (</span><span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; the\; s</span>}}{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; cos</span>}\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; +</span>\mathrm{<span\; class="notranslate">\; the\; s</span>}}{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; sin</span>}\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; +</span>\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; +</span>\mathrm{<span\; class="notranslate">\; the\; s</span>}}{\mathrm{<span\; class="notranslate">\; m</span>}}_{\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; +</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; the\; s</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; +</span>\mathrm{<span\; class="notranslate">\; the\; s</span>}}<span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 2</span>}<span\; class="notranslate">\; )</span>$ i represents an integer from 0 to p, s represents an integer from 0 to q, when the continuous gray distribution function f(x, y) is rotated by θ, the edge line (1) to be detected in the window image (2) It is perpendicular to the horizontal axis x of the coordinate system described in step 3.1, that is, the grayscale function corresponding to the inscribed circle of the rotated window image (2) is symmetrical about the horizontal axis x. According to the calculation formula of the space moment, the rotated window image (2) is inscribed The space moment M' ₀₁ ＝0 corresponding to the circle, at this time p＝0, q＝1, can be obtained by formula (2): M' ₀₁ ＝-sinθ·M ₁₀ +cosθ·M ₀₁ ＝0...(3) Among them, M ₁₀ and M ₀₁ represent the space moment corresponding to the inscribed circle before rotation; Step 4.4: Use the spatial moment calculation template to replace the formula (1), and substitute the discrete grayscale function I(x, y) of the pixels in the inscribed circle of the window image (2) into the above-mentioned spatial moment calculation template formula, with Calculate the value of spatial moments M ₀₁ and M ₁₀ in formula (3); Step 4.5: Obtain formula (4) according to formula (3), and substitute the values of spatial moments M ₀₁ and M ₁₀ obtained in step 4.4 into formula (4) to obtain the value of angle θ: $\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; =</span>\mathrm{<span\; class="notranslate">\; arctan</span>}\frac{{\mathrm{<span\; class="notranslate">\; m</span>}}_{\mathrm{<span\; class="notranslate">\; 01</span>}}}{{\mathrm{<span\; class="notranslate">\; m</span>}}_{\mathrm{<span\; class="notranslate">\; 10</span>}}}<span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 4</span>}<span\; class="notranslate">\; )</span>$ Step 4.6: Set the edge normal (3) as the new horizontal axis X, and take the angle θ in step 4.5 as the projection angle, then the center point of each pixel grid in the window image (2) described in step 3.2 is at The pixel coordinates (x, y) in the pixel plane and the point coordinate X on the horizontal axis where the edge normal (3) is located have the following projection formula: X＝xcosθ+ysinθ...(5) In the formula (5), after the X value represents the projected point in the window image (2) by the projection angle θ projection, the coordinate reading value of the corresponding projected point (4) on the new horizontal horizontal axis X described in step 4.6; Step 5: Based on the one-dimensional edge point detection of orthogonal polynomial fitting, determine the intersection position of the edge normal (3) and the edge line (1) to be detected, which includes the following steps: Step 5.1: reconstitute the Z-axis described in step 3.2 and the horizontal horizontal axis X described in step 4.6 into a plane Cartesian coordinate system XOZ required for one-dimensional edge point detection, i.e. the projection plane (6); According to the projection formula (5), the center point of each pixel grid in the window image (2) is projected onto the horizontal horizontal axis X described in step 4.6; in the plane Cartesian coordinate system XOZ, the projection on the edge normal (3) is obtained Discrete gray distribution function (X _i , Z _i ) of point gray value (5), i=1,...,2n+1; Step 5.2: Normalization of the projection point coordinates, if the number of projection points (4) formed on the horizontal axis X as described in step 5.1 is 2n+1, set and order i=1,...,2n+1, then x _i ∈[-1,1]; Step 5.3: Calculate the orthogonal polynomial system, and use the schemite orthogonalization method to calculate the orthogonal polynomial system about the weight function ρ(x)=1 on the [-1,1] interval; Step 5.4: Take the first four order orthogonal basis functions in step 5.3 Using formula (7) as the fitting base to carry out cubic orthogonal polynomial fitting, then the equation of fitting polynomial is: The coefficients {a ₁ , a ₂ , a ₃ , a ₄ } of the fitting polynomial in formula (8) are the fitting parameters to be found; Step 5.5: Utilize formula (8) to deduce the corresponding normal equation group of the equation of fitting polynomial as follows: In formula (9) represents a discrete vector Inner product of (i=0,1,2,3) and itself, Z=[Z ₁ ,Z ₂ ,…,Z _2n+1 ] ^T is a discrete gray value vector, Indicates that Z and inner product; Step 5.6: Use equation (9) to obtain the calculation formula of the coefficient {a ₁ , a ₂ , a ₃ , a ₄ } of the fitting polynomial: Substituting each known quantity into formula (10), the results of a ₀ , a ₁ , a ₂ and a ₃ are obtained respectively; Step 5.7: Substituting the result of the polynomial fitting described in step 5.6 into formula (8) as a known quantity, the polynomial fitting equation becomes: Step 5.8: Calculate the second order derivative of the cubic polynomial (11) as zero, namely: Solve and obtain the coordinates of the point that makes the second order derivative of formula (12) zero; Step 5.9: Solve formula (12) and denormalize the calculation results to obtain the position coordinate R of the one-dimensional edge point on the horizontal axis X, namely: Equation (13) is an analytical formula for calculating the position of a one-dimensional edge point, which is the intersection of the edge normal (3) and the edge line (1) to be detected; Step 6: Use formula (13) to determine the intersection position of the edge normal (3) passing through the origin in the window image (2) and the edge line (1) to be detected, and use step 4.5 to obtain the direction angle θ of the edge normal (3), Determine the equation of the edge straight line (1) to be detected in the xoy coordinate system, thereby realizing the sub-pixel positioning of the edge straight line (1) to be detected in the window image (2); Step 7: The operator reselects the whole pixel grid described in step 1, and repeats steps 2 to 6 to obtain the sub-pixel positioning of the edge line (1) to be detected in the newly obtained window image (2); Step 8: When the operator determines that the sub-pixel positioning of the desired edge position on the entire image has been completed, the detection method ends.