A schlieren image shock wave recognition method based on image correlation and pattern search

Abstract

The present application provides a shock wave recognition method based on image correlation and pattern search. The method includes: adjusting the contrast of the image and performing cross-correlation operation to obtain a gradient image with enhanced shock wave boundary features, which is used for recognizing the shock wave; performing morphological processing on the image by using erosion and dilation methods to obtain a binary image used for region division, performing region division on the image by using a connected component algorithm to obtain a region image used for recognizing the shock wave; using an automatic initialization node method on the automatically divided region image to obtain the starting and ending nodes of the recognition; and performing correlation calculation on the starting and ending nodes by using a Dijkstra algorithm to obtain the path of the shock wave in the image.

Description

Technical Field

[0001] This application relates to the field of image processing, specifically to a schlieren image shock wave recognition method based on image correlation and graphic search. Background Technology

[0002] Schlieren technology is a non-contact optical observation method for visualizing density changes in transparent media. When light passes through a medium with uneven density, it is deflected due to changes in refractive index. Through a sophisticated optical system, these tiny light rays are deflected and converted into a visible image.

[0003] Conventional schlieren techniques capture flow field density gradients through changes in light refraction, enabling qualitative or semi-quantitative analysis of shock wave structure. Focused schlieren techniques can address the quantitative density measurement needs in shock wave wind tunnel experiments, achieving high-resolution observation and higher sensitivity measurements of local flow fields. With the development of high-sensitivity camera systems, schlieren techniques have expanded to three-dimensional flow field display and large-scale shock wave dynamic capture, significantly improving spatiotemporal resolution. Image processing in schlieren shock wave system identification mainly includes rapid preprocessing of the original schlieren image, eliminating interference noise through Gaussian filtering and thresholding, and improving image clarity and contrast. Classical operators such as Canny, Roberts, Prewitt, Sobel Canny, and Laplacian Gaussian are used for edge detection of the schlieren image to determine the shock wave location.

[0004] However, since the topological structure of shock waves in experimental schlieren images changes across different time periods, no single edge detection method is suitable for all situations, as a trade-off must be made between eliminating interference and improving the precision of recognition. Achieving rapid processing of dynamic flow field schlieren images, high-precision extraction of shock wave features from schlieren images, and rapid reconstruction of flow field topology are crucial problems that urgently need to be solved, and are of great significance for the extraction, diagnosis, and analysis of flow field features in experiments. Summary of the Invention

[0005] Purpose of the invention: To address the problem that existing image processing techniques, such as edge detection, cannot be applied to shock wave topologies that change over time, and cannot simultaneously eliminate interference and accurately identify shock waves, this invention discloses a schlieren image shock wave recognition method.

[0006] Technical solution: This invention discloses a schlieren image shock wave recognition method based on image correlation and image search, comprising the following steps:

[0007] The schlieren image to be identified is acquired, the schlieren image is cropped to highlight the display area of ​​the shock wave, and the schlieren image is L2 normalized.

[0008] The angle of the shock wave in the schlieren image is measured, and the kernel is set according to the angle to perform cross-correlation operation on the schlieren image to generate a gradient image that enhances the shock wave structure at a specific angle.

[0009] Set a threshold, binarize the gradient image, and process the binarized gradient image using erosion and dilation methods to obtain a binary image for automatic partitioning.

[0010] Based on the connected component algorithm, the binary image is divided into different connected branches, and the minimum bounding rectangle of each branch is used as the recognition region.

[0011] Based on the automated starting node method, additional node columns are added to both sides of the region to be identified. Each pixel in the image is used as a node and the weight between nodes is calculated to generate an adjacency matrix.

[0012] Based on Dijkstra's algorithm, the shortest path of the adjacency matrix is ​​calculated to obtain the shortest path node list, and the coordinates of the shortest path nodes are stored in a set.

[0013] Linear fitting is performed on the set to obtain the shock wave location and topology.

[0014] Furthermore, the L2 normalization includes: calculating the mean value of all pixel values ​​in the image, subtracting the mean value from each pixel value in the image; performing L2 normalization on the image after subtracting the mean value, and calculating the L2 norm of the image; dividing all pixel values ​​in the image by its L2 norm, and normalizing all pixel values ​​in the image to between 0 and 1.

[0015] Furthermore, the kernel for cross-correlation calculation based on the shock wave angle is first measured to determine the shock wave angle range to be identified, and a small matrix with a diagonal angle within the shock wave angle range is set as the kernel for cross-correlation calculation.

[0016] Furthermore, the binarized gradient image is processed using erosion and dilation methods. The structuring element used is a small matrix with an angle similar to that of the shock wave. After erosion and dilation, regions in the binary image that are different from the direction of the structuring element are erased, while regions with the same direction are retained.

[0017] Furthermore, the connected component algorithm divides branches based on the connectivity of each pixel's 8-neighborhood. Branches exceeding the shock wave size are not processed, while branches within the shock wave size have their minimum bounding rectangle region extracted.

[0018] Furthermore, the weights of the additional node columns are all set to 10. -5 .

[0019] Furthermore, the weight calculation method for the remaining nodes, excluding the extra node columns, is expressed as follows:

[0020] w ij =2-(G i +G j +10 -5

[0021] Among them, w ij G is the distance weight between nodes i and j. i G is the gradient value at node i. j It is the gradient value at node j.

[0022] Furthermore, the adjacency matrix is ​​a matrix symmetric along the main diagonal, which stores the distance weights between nodes.

[0023] Furthermore, the Dijkstra algorithm takes the adjacency matrix as input and outputs the path and distance with the lowest weight between the starting and ending nodes.

[0024] Furthermore, a linear fit is performed on the path coordinates, and the sum of squared residuals and coefficient of determination are calculated. Paths with a coefficient of determination less than 0.9 are excluded, and the remaining paths are used to fit a straight line to obtain the shock wave position. The formulas for calculating the sum of squared residuals and coefficient of determination are as follows:

[0025]

[0026] Among them, y i , These are the actual values, fitted values, and average values ​​for each node. RSS and TSS are the residual sum of squares and the total sum of squares, respectively. R 2 It is the coefficient of determination.

[0027] Beneficial effects:

[0028] This invention discloses a shock wave recognition method for schlieren images. Based on image correlation and graph search techniques, it improves the accuracy and anti-interference ability of the recognition path by calculating the residual sum of squares and the coefficient of determination of the path and by linear fitting. Compared with edge detection methods using conventional operators, the method provided by this invention can obtain more shock wave topology information, has higher accuracy, and solves the problem of weak anti-interference ability of conventional edge detection algorithms.

[0029] This invention utilizes the connected component algorithm and Dijkstra's algorithm to automatically identify shock waves in images, significantly reducing the processing time required for image segmentation and feature extraction. Compared to other shock wave identification methods, the shock wave identification method disclosed in this invention only requires the use of a suitable cross-correlation kernel and does not require the selection of start and end nodes on the image, greatly reducing the cost of feature extraction. It has a wide range of applications and can be extended to feature extraction in images of other imaging modes. Attached Figure Description

[0030] Figure 1 This is a flowchart illustrating a schlieren image shock wave recognition method based on image correlation and graphic search provided by the present invention.

[0031] Figure 2 This is a schematic diagram illustrating the enhancement of shock wave structure features by a specific angle kernel provided by the present invention;

[0032] Figure 3 This is a schematic diagram of the corrosion and expansion process for removing interference provided by the present invention;

[0033] Figure 4 This is a schematic diagram of the connected branch division region and automatic region segmentation process provided by the present invention;

[0034] Figure 5 This is a schematic diagram illustrating the automatic segmentation and recognition effect of shock wave structure in one direction provided by the present invention. Detailed Implementation

[0035] The present invention will be further illustrated below with reference to the accompanying drawings and specific embodiments. It should be understood that the following specific embodiments are for illustrative purposes only and are not intended to limit the scope of the invention. After reading the present invention, any modifications of the present invention in various equivalent forms by those skilled in the art will fall within the scope defined by the appended claims.

[0036] The schlieren image shock wave recognition method based on image correlation and image search provided by this invention has the following process: Figure 1 As shown, cross-correlation operations are performed using a kernel based on shock wave angle to enhance shock wave features. A connected components algorithm is used to segment the region containing each shock wave, and an image segmentation method based on Dijkstra's algorithm is used to extract shock wave features. The process includes the following steps:

[0037] The schlieren image to be identified is acquired, the schlieren image is cropped to highlight the display area of ​​the shock wave, and the schlieren image is L2 normalized.

[0038] The angle of the shock wave in the schlieren image is measured, and the kernel is set according to the angle to perform cross-correlation operation on the schlieren image to generate a gradient image that enhances the shock wave structure at a specific angle.

[0039] Set a threshold, binarize the gradient image, and process the binarized gradient image using erosion and dilation methods to obtain a binary image for automatic partitioning.

[0040] Based on the connected components algorithm, the binary image is divided into different connected branches, and the minimum bounding rectangle of each branch is used as the recognition region;

[0041] Based on the automated starting node method, additional node columns are added to both sides of the region to be identified, the weights between nodes are calculated, and an adjacency matrix is ​​generated based on the weight calculation results.

[0042] Based on Dijkstra's algorithm, the shortest path of the adjacency matrix is ​​calculated to obtain the shortest path node list, and the coordinates of the shortest path nodes are stored in a set.

[0043] Linear fitting is performed on the set to obtain the shock wave location and topology.

[0044] The specific implementation steps of the method of the present invention will be described below with reference to specific embodiments.

[0045] Step 1: In this embodiment, the schlieren image to be identified comes from the inlet mode conversion experiment of the external parallel TBCC inlet model with a design transition Mach number of 3.5 conducted in the supersonic wind tunnel of Nanjing University of Aeronautics and Astronautics. Based on the main imaging area of ​​the schlieren, the image to be identified is cropped to highlight the display area where the shock wave is located in the image to be identified.

[0046] Calculate the mean value of all pixels in the image. Subtract this mean value from each pixel value. Perform L2 normalization on the subtracted mean value and calculate the L2 norm of the image. Divide all pixel values ​​in the image by its L2 norm to normalize all pixel values ​​to between 0 and 1, thereby enhancing image contrast.

[0047] Step 2: Perform cross-correlation calculations on the L2-normalized image based on the shock wave angle to generate a gradient image that enhances the shock wave structure at a specific angle. A schematic diagram of this process is shown below. Figure 2 As shown, the leftmost image is the schlieren image to be identified in this embodiment.

[0048] The cross-correlation operation is performed using MATLAB software. The code is: R1 = imfilter(normalizeA, k, 'replicate'); where R1 is the image variable name after cross-correlation, normalizeA is the image variable name after L2 normalization, k is the kernel variable name used for cross-correlation, and 'replicate' is the boundary processing method.

[0049] In images, shock waves with the same development direction have small angular differences, generally within 20°, and therefore can be enhanced by the same kernel. Measuring the angles of some shock waves in the schlieren image serves as the basis for setting the kernel. When setting the kernel, the angle of the diagonal in the small matrix corresponds to the angle of the shock wave in the image. In this example, a small kernel with a main diagonal of 45° is chosen, corresponding to a shock wave angle of approximately 20° in the image. Because the small matrix is ​​much smaller than the shock wave, using 45° can also enhance oblique shock waves of approximately 20°, and it has greater versatility for oblique enhancement, effective between 15° and 75°. Setting the values ​​above and below the diagonal in the small matrix to -1 and 1, or vice versa, enhances the characteristics of the shock wave boundary. This is reflected in the image as areas with brightness variations at the angle corresponding to the kernel centerline becoming brighter, while the remaining areas become darker. Figure 2 The centerline angle of the kernel enhancement feature shown in the rightmost image is approximately the same as the shock wave development angle. This yields a gradient image that enhances the shock wave features in that direction. Shock structures and other schlieren regions opposite to the kernel direction will appear darker. In this embodiment, a (-1,1) type kernel is used. Specifically, it's a small 5x5 matrix with the values ​​on the diagonal set to 0, values ​​above the diagonal set to -1, and values ​​below the diagonal set to 1. This is an approximately 45° small kernel used to enhance oblique gradient changes. Figure 2 The black and white squares represent the image format using the kernel.

[0050] Step 3: Set a threshold based on image characteristics, binarize the gradient image, and perform morphological processing on the binary image using erosion and dilation methods to obtain a binary image for automatic region division.

[0051] like Figure 3 The diagram illustrates the morphological processing of the schlieren image in this embodiment. An appropriate threshold is selected based on image characteristics to binarize the gradient image. After cross-correlation calculation, a binarization threshold is set based on the shock wave region with lower pixel values ​​after feature enhancement. The set threshold should remove more noise interference while preserving the shock wave structure to the greatest extent possible. Figure 3 Using 60 / 255 as the binarization threshold, pixel values ​​below this value are set to 0, achieving more accurate shock wave localization.

[0052] Morphological processing of the image is performed based on erosion and dilation methods. A small matrix with an angle similar to that of the shock wave is used as the structuring element. After erosion and dilation, regions in the binary image that are not in the same direction as the structuring element are erased, while regions that are in the same direction are retained, resulting in a binary image for automatic region division. Figure 3The shock angle of the schlieren image is between 20° and 40°, so a small [3×7] matrix is ​​used for erosion and dilation. The values ​​of the small matrix (1,1), (1,2), (2,3), (2,4), (2,5), (3,6) and (3,7) are set to 1, and the values ​​of the rest of the matrix are 0. This is a structuring element with an approximate 23° angle, but its size is much smaller than that of the shock structure, so it can be used to erode and dilate shock regions with a wider angle range.

[0053] Step 4: Region segmentation of the image. Binary images use only two pixel values, typically 0 and 1, allowing for quick determination of whether adjacent pixels belong to the same region. Based on the connected components algorithm, connected branches are divided according to the connectivity of each pixel's 8-neighborhood. In binary images, the length and height of shock wave branches are significantly smaller than the length and height of surface branches within a certain range. Therefore, branches exceeding the shock wave length and height range are not processed. For branches within the shock wave length and width range, their minimum bounding rectangle region is extracted as the region to be automatically segmented for identification. For example... Figure 4 The diagram shows a flowchart of the connected components algorithm for dividing regions and automatically segmenting regions, according to an embodiment of the present invention.

[0054] Step 5: Based on the automated starting node method, add additional node columns on both sides of the image of the region to be recognized, with the weight of each node column set to 10. -5 Each pixel in the region to be identified is treated as a node, and weights are calculated between nodes.

[0055] An adjacency matrix is ​​generated based on the weight calculation results to reflect the distance relationship between each node and other nodes in the image. Specifically, the algorithm reads the height *h* and width *w* of the image region to be processed, and initializes a [hw×hw] matrix. For example, a 450×1000 image can typically generate a 450000×450000 adjacency matrix. All values ​​in the matrix are set to *inf*, representing that initially all nodes are unreachable. The diagonal of the matrix is ​​then set to 0, representing that the distance between each node and itself is 0. Next, each pixel node is traversed, and the weight calculation results of its eight neighboring nodes are written into the adjacency matrix, replacing the original *inf* values, representing that each node is reachable from its eight neighbors. The symmetric structure of the adjacency matrix means that different nodes have the same distance weight. For example, (a,b) and (b,a) have the same value in the same adjacency matrix, meaning that the distance between node *a* and node *b* is the same as the distance between node *b* and node *a*.

[0056] The formula for calculating the weights between nodes is as follows:

[0057] w ij=2-(G i +G j +10 -5

[0058] Among them, w ij G is the distance weight between nodes i and j. i G is the gradient value at node i. j It is the gradient value at node j.

[0059] Step 6: Use Dijkstra's algorithm to calculate the shortest path in the adjacency matrix to obtain the shortest path node list. Delete all added extra nodes from the shortest path node list and store all path coordinates into a set.

[0060] Because the path weights between nodes in the extra node column are much lower than the path weights between the node column and the region to be identified, start and end nodes are automatically initialized on both sides of the region to be identified during path calculation. The path and distance with the lowest weight between the start and end nodes are then output as the result. All added extra nodes are removed from the shortest path to obtain the automatic segmentation result.

[0061] Step 7: Perform linear fitting on the path coordinates in the set, calculate the sum of squared residuals and coefficients of determination for the paths, exclude paths with excessive fitting deviation, and use the remaining paths to fit a straight line to obtain the shock wave position.

[0062] Linear fitting is performed on the path coordinates, and the sum of squared residuals and coefficient of determination are calculated. Paths with a coefficient of determination less than 0.9 are excluded, and the remaining paths are used to fit a straight line to obtain the shock wave position. The expressions for calculating the sum of squared residuals and coefficient of determination are as follows:

[0063]

[0064] Among them, y i , These are the actual values, fitted values, and average values ​​for each node. RSS and TSS are the residual sum of squares and the total sum of squares, respectively. R 2 It is the coefficient of determination. The closer the coefficient of determination is to 1, the better the model fits.

[0065] like Figure 5 The diagram shown is a schematic of the process in this embodiment. The bottom image is the recognition result obtained by processing the schlieren image to be recognized by the method of the present invention. The black line in the figure is the shock wave topology structure recognized by the method of the present invention.

[0066] In summary, this application achieves shock wave feature extraction from schlieren images by finding the shortest path between automatically initialized nodes for image segmentation. This significantly reduces the processing time required for image segmentation and feature extraction, and improves the accuracy of shock wave recognition. Compared to other edge detection methods based on different operators, this automatic region partitioning and segmentation method can achieve high-precision shock wave recognition in schlieren images, derive more shock wave topology information from the image, and has a simple, reliable structure with wide applicability. It provides technical support for flow field topology reconstruction and rapid feature diagnosis.

Claims

1. A method for shock wave recognition in schlieren images based on image correlation and pattern search, characterized in that, Includes the following steps: The schlieren image to be identified is acquired, the schlieren image is cropped to highlight the display area of ​​the shock wave, and the schlieren image is L2 normalized. Measuring the angle of the shock wave in a schlieren image, and setting a kernel based on this angle to perform cross-correlation operations on the schlieren image to generate a gradient image that enhances the shock wave structure at a specific angle; including: setting a kernel for cross-correlation operations based on the angle of the shock wave, first measuring the range of shock wave angles to be identified, and setting a small matrix whose diagonal angle is within the range of shock wave angles as the kernel for cross-correlation operations; A threshold is set, and the gradient image is binarized. The binarized gradient image is then processed using erosion and dilation methods to obtain a binary image for automatic partitioning. The structuring element used in the erosion and dilation process is a small matrix with an angle similar to that of the shock wave. After erosion and dilation, regions in the binary image that are different from the direction of the structuring element are erased, while regions with the same direction are retained. Based on the connected component algorithm, the binary image is divided into different connected branches, and the minimum bounding rectangle of each branch is taken as the region to be identified. Based on the automated starting node method, additional node columns are added on both sides of the region to be identified. Each pixel in the image is used as a node and the weight between nodes is calculated to generate an adjacency matrix. Based on Dijkstra's algorithm, the shortest path of the adjacency matrix is ​​calculated to obtain the shortest path node list, and the coordinates of the shortest path nodes are stored in a set. To obtain the shock wave location and topology by performing linear fitting on the set, the following steps are taken: performing linear fitting on the path coordinates, calculating the sum of squared residuals and coefficients of determination of the paths, excluding paths with coefficients of determination less than 0.9, and using the remaining paths to fit a straight line to obtain the shock wave location.

2. The shock wave identification method according to claim 1, characterized in that, The L2 normalization includes: calculating the mean value of all pixel values ​​in the image, subtracting the mean value from each pixel value in the image; performing L2 normalization on the image after subtracting the mean value, and calculating the L2 norm of the image; dividing all pixel values ​​in the image by its L2 norm, and normalizing all pixel values ​​in the image to between 0 and 1.

3. The shock wave identification method according to claim 2, characterized in that, The connected component algorithm divides the branches based on the connectivity of each pixel's 8-neighborhood. Branches exceeding the shock wave size are not processed, while branches within the shock wave size have their minimum bounding rectangle region extracted.

4. The shock wave identification method according to claim 3, characterized in that, The weights of the additional node columns are all set to 10. -5 .

5. The shock wave identification method according to claim 4, characterized in that, The weights between nodes other than the extra node columns are calculated using the following expression: Among them, w ij G is the distance weight between nodes i and j. i G is the gradient value at node i. j It is the gradient value at node j.

6. The shock wave identification method according to claim 5, characterized in that, The adjacency matrix is ​​a matrix symmetric along the main diagonal, which stores the distance weights between nodes.

7. The shock wave identification method according to claim 6, characterized in that, The Dijkstra algorithm takes the adjacency matrix as input and outputs the path and distance with the lowest weight between the starting and ending nodes.

8. The shock wave identification method according to claim 7, characterized in that, The formulas for calculating the sum of squared residuals and the coefficient of determination are as follows: in, These are the actual values, fitted values, and average values ​​for each node. RSS and TSS are the residual sum of squares and the total sum of squares, respectively. R 2 It is the coefficient of determination.

Images