Abstract
Quantitative analysis of the flow field is an effective method to study hydrodynamics. As a flow field measurement technology, the Background Oriented Schlieren (BOS) is widely used. However, it is difficult to measure the complex transparent flow field (flow field with large refractive index gradient) using the BOS experiment. In order to overcome this difficulty and improve the accuracy of the BOS experiment, this paper presents a hybrid adaptive wavelet-based optical flow algorithm for the BOS. The current algorithm is a combination of the traditional optical flow algorithm and the wavelet-based optical flow algorithm. By adding the initial value constraints, the adaptive scale constraints, and the adaptive regularization constraints, the algorithm can effectively overcome the above-mentioned difficulty and also improve its accuracy. To further illustrate the feasibility of the proposed method, this paper uses the simulation data, the data of the DNS datasets, and the data of the BOS experiment to verify the performance of the algorithm. The experiment of comparing the proposed algorithm with the existing ones is carried out on the DNS datasets and the data of the BOS experiment. Finally, the proposed method is verified by a practical BOS experiment. The results show that the proposed algorithm can effectively improve the measurement accuracy of displacements.
1. Introduction
Background Oriented Schlieren (BOS) is widely used to measure the transparent flow field owing to the simplicity of its experimental device and its high measurement accuracy [1, 2]. BOS is an optical flow diagnostic technique used to measure the refractive gradients in a flow field by tracking the distortion of a target dot pattern. Due to the uneven refractive index distribution in the transparent flow field, the light deflects after passing through the flow field; therefore, an object viewed through the transparent flow field will appear distorted due to the refraction of light rays traversing the flow field [3]. BOS relies on this principle to measure the transparent flow field.
The accuracy of the displacement estimation algorithm directly determines the measurement accuracy of the BOS. The displacement field obtained by the traditional cross-correlation algorithm is not dense, and especially in the region of large refractive index gradient, the algorithm is not accurate. The traditional optical flow algorithm and the improved multiscale optical flow algorithm overcome the shortcomings of the cross-correlation algorithm, especially in the region of large refractive index gradient. However, the optical flow algorithm is limited in practical application, as it is based on the assumption of small displacement and constant brightness. In addition, compared with the cross-correlation algorithm, the optical flow algorithm is more complex and computationally expensive. Jagadesh-Ramaiah et al. propose the use of GPU acceleration to improve the computational efficiency of the algorithm [4]. Because the wavelet transform itself has the characteristics of data compression and multiscale analysis, the wavelet-based optical flow algorithm has attracted people’s attention [5–8]. This paper presents a hybrid adaptive wavelet-based optical flow algorithm for the Background Oriented Schlieren (BOS) based on previous studies [2, 5–7, 9, 10]. The wavelet-based optical flow algorithm can effectively reduce the calculation cost of the algorithm by using the data compression characteristics of the wavelet. At the same time, the algorithm can effectively deal with the problem of large-scale displacement in the BOS experimental data by using the multiscale analysis characteristics of the wavelet. The rest of this paper is organized as follows: Section 2 introduces the related works at home and abroad; Section 3 describes the proposed algorithm; Section 4 presents the simulation method of the BOS experimental device based on COMSOL and the description of the truncation error; Section 5 presents the experimental results and discussion; and Section 6 summarizes the research.
2. Related Works
The performance of the data processing algorithm directly affects the measurement accuracy of the BOS. In order to improve the measurement accuracy of the BOS, more and more researchers are working on optimizing the data processing algorithm. Initially, an algorithm based on cross-correlation was used to process the BOS data. However, this algorithm has some shortcomings, which can be tackled by PIVlab, a useful MATLAB toolbox, proposed by Thielicke et al. [11]. The MATLAB toolbox uses cross-correlation algorithm to calculate pixel displacement. However, the calculated displacement is minimal and imprecise. Especially in the region of high refractive index gradient, the calculation error of the algorithm is large. Therefore, the optical flow method is introduced into BOS data processing. The optical flow method does not rely on the size of the interrogation window. It estimates the displacement of each pixel according to the overall correlation between two adjacent images, rather than depending on the relationship between the pixels of a certain region or a certain pixel. Optical flow algorithm is based on global correlation estimation of the displacement field, so that it can obtain dense, high-resolution displacement fields [5], such as Horn-Schunck algorithm [12] and Lucas-Kanade algorithm [13]. However, the optical flow algorithm also has its own potential drawbacks. It is based on the assumption of small displacement. This also limits the use of the optical flow algorithm. In order to overcome the limitation of displacement, the theory of multiresolution analysis is introduced into the optical flow algorithm [9, 14]. Although the improved optical flow algorithm can overcome the limitation of displacement, the optical flow algorithm cannot break through the assumption of constant brightness. This makes it difficult to apply the optical flow algorithm to the actual experiment. It may exhibit difficulties with experimental noise or events, such as particles, entering or leaving the image plane that commonly occur in experiments [6]. In order to enhance the adaptability of the optical flow algorithm, it needs some kind of regularization constraints, which are adjusted based on experience in most cases, or according to the characteristics of the flow field. There are several optical flow methods in the literature that have been designed and optimized for divergence-free motions [15]. However, some optical flows use periodic boundary conditions that are not well suited for experimental data [16]. The multiscales optical flow methods based on coarse-to-fine blurred images sequentially also have drawbacks, such as large errors occurring in coarser estimates cannot be subsequently corrected in the finer estimates [17]. However, the wavelet-based optical flow method could cope with this difficulty. The wavelet-based optical flow method representation hierarchy is different from the typical image coarse-to-fine hierarchy. The hierarchy is built along motion resolution level, not image resolution level. The wavelet-based optical flow algorithm was first proposed by Wu et al. [18]. Based on this, Kadri-Harouna et al. [10] improved the wavelet-based optical flow algorithm and mainly explored the application of regularization constraints in the optical flow algorithm based on the wavelet. Schmidt et al. [6] improved the wavelet transform on the basis of the Dérian algorithm. The transform in the WOF-X is performed using biorthogonal 9–7 wavelets and a lifting transform, unlike in previous wavelet-based optical flow methods. However, when dealing with the complex flow field (the flow field with great change of refractive index gradient), the calculation error of the above-mentioned wavelet-based optical flow algorithm is still large [6].
3. Hybrid Adaptive Wavelet-Based Optical Flow Algorithm
3.1. Multiscale Optical Flow Representation Based on the Wavelet
and are images obtained by imaging the same background pattern under the same BOS experimental device; is the imaging result without the flow field, and is the image with the flow field. The optical flow algorithm is used to estimate the pixel displacement between the two images. Based on the conservation assumption of the image brightness, images and fulfill the following formula:
In order to calculate the displacement field , it is necessary to change formula (1) to the standard quadratic function, as shown in the following formula:
Here, represents the set of effective pixels of image .
By calculating the minimum value of formula (2), the pixel displacement is obtained. For a point in two-dimensional space, formula (2) is a quadratic equation of two-dimensional vector . However, the number of unknowns in the equation is more than that in the equation, which is an underdetermined equation. In order to better solve the above optimization problems, a convex regularization term needs to be added [19]. The objective function is given by the following formula:
The wavelet framework introduces the decomposition of a signal into several sets of details at various scales and a remaining coarse approximation. The projection of the two components, , , of the velocity field onto such multiresolution approximation spaces proves to be efficient within the context of optical flow estimation, enabling the handling of large displacements. Here, let’s recall the basic principle of the two-dimensional wavelet transform. For any function , it can be expressed linearly by , as shown in the following formula:The closed subspace sequence generated by the expansion and translation of the scale function is called the approximate space of ; the closed subspace sequence generated by the expansion and translation of the wavelet function is called the detail space of .
The projection of each component of the optical flow field on wavelet basis can be obtained from formula (4), as shown in the following formulas:
Here, , , , , and the subscripts , represent the horizontal and the vertical displacements of the wavelet basis function; , represent the resolution scale; and the superscripts 1, 2, and 3 represent the horizontal direction, the vertical direction, and the diagonal direction, respectively.
It is known from formulas (5) and (6) that the reconstruction of the optical flow field from the wavelet basis function is a linear transformation. Therefore, the relationship between the optical flow and the wavelet basis function can be expressed by the following formula:
Here, , and represent the wavelet transform coefficients.
Using formula (7), the solution of the optical flow problem can be converted into the wavelet coefficients of the optical flow field, and the objective function can be changed into the following formula:
Here, is the unknown wavelet coefficient used to estimate .
3.2. Implementation
Formula (12) gives the mathematical expression of the optical flow problem based on the wavelet transform, which transforms the optical flow problem into the wavelet coefficient of the optical flow. The L-BFGS algorithm is also selected to solve the above-mentioned optimization problem. This algorithm can effectively save memory space and improve the efficiency of the operation, and is suitable for solving large-scale optimization problems. L-BFGS algorithm needs to calculate the gradient of the objective function when solving the optimization problem. The gradient of the objective function to the wavelet coefficients can be expressed as
In order to solve , this paper first solves , which is expressed as
For the solution of , Kadri-Harouna et al. [10] made a detailed derivation and explanation. In this paper, the method of discrete approximation in [10] is used to approximate , and the expression of iswhere is a diagonal matrix and the diagonal element is .
In the existing multiscale optical flow method based on the wavelet, the optimal scale is set by prior information of the flow field or experience, or is adjusted according to the calculation results. This setting method has certain subjectivity, aiming at a certain flow field or a certain kind of flow field. In order to enhance the adaptability of the algorithm and reduce its dependence on prior information, an adaptive scale selection method is added to the algorithm. Firstly, calculate the initial displacement field . Secondly, calculate the variance of , the gradient of , and the variance of . Thirdly, the threshold judgment is used to determine the optimal scale. When , it indicates that the overall amplitude of changes greatly and the distribution of the flow field characteristics is complex, which requires the use of the largest scale as the optimal scale; when , it indicates that the overall amplitude of changes greatly, and that the distribution of the flow field characteristics is relatively simple. It is possible that is a gradual flow field, which requires a smaller scale as the optimal scale; when , it indicates that the overall amplitude of change is small, and the distribution of the flow field characteristics is complex, which requires as the optimal scale; when , it indicates that the overall amplitude of change is small and that the distribution of the flow field characteristics is simple, which requires the use of scale as the optimal scale. As shown in Figure 1, this paper analyzes the feasibility of adaptive scale selection using 1000 images in the DNS dataset. Since the image resolution is , . For example, the threshold is set to and .(i), . According to Figure 1, , .(ii), . According to Figure 1, , .(iii), . According to Figure 1, , .(iv), . According to Figure 1, ,.
(a)
(b)
(c)
(d)
PSNR and SSIM represent the similarity between the image at different scales and the original image. The larger the value of PSNR and SSIM, the higher the similarity between the image at this scale and the original image. Through the above analysis, the adaptive selection of the scale can be realized by setting an appropriate threshold.
In order to avoid over-smoothing of the velocity field and the impact on the accuracy of the velocity field reconstruction, we improve the regular term parameter , so that the objective function can adaptively adjust the smoothness of the velocity field according to the gradient variation of . The expression of is shown in the following formula:where is the variance of the gradient of the velocity field , is the variation of the gradient of the velocity field , and , controls the value range of , . According to formula (17), is controlled by the ratio of to . When is less than , it means that the distribution characteristics of the flow field tend to be concentrated, which indicates that the phenomenon of over-smoothing occurs, and the parameter of the regularization term should be reduced; when is greater than , it means that the distribution characteristics of the flow field tend to be dispersed, which indicates that there is an under constraint phenomenon, and the regularization parameter should be increased. By adjusting adaptively, the objective function can reach the convergence condition smoothly.
4. Simulation of BOS Experiment Based on the Ray Tracing Method
4.1. Simulation Description
In BOS, the density gradients in a flow are measured by the apparent shift of a dot pattern viewed through a variable density medium, where the displacement is evaluated using methods similar to PIV [1, 2]. In practical application, it is necessary to use simulation data to verify the accuracy of the algorithm. To assess the accuracy of the proposed displacement estimation algorithms in the paper, synthetic BOS images are required [3], which are simulated based on the ray tracing method [20, 21]. The layout of the BOS experimental setup modeled in the simulation experiment is shown in Figure 2, and the parameters describing the placement of the elements are summarized in Table 1.
4.2. Simulation Results
The spatial refractive index distribution function of the flow field used in the simulation experiment is shown in the following formula:
Here, is the center of flow. The pixel displacements can be easily obtained by ray tracing. The simulation results are shown in Figure 3.
(a)
(b)
(c)
(d)
In order to further illustrate the effectiveness of the simulation method in this paper, PIVlab [11] and LK algorithm [22] are used to calculate the displacement distribution of the simulation image. The results are shown in Figure 4.
(a)
(b)
(c)
5. Experimental Results and Analysis
In order to verify the effectiveness of the proposed algorithm, firstly, this paper uses the simulation data generated by the proposed simulation method to verify the effectiveness of the proposed algorithm. Secondly, the validity of the algorithm is verified by using the public DNS datasets [6, 10]. At the same time, it is compared with a variety of existing optical flow algorithms on the datasets. Finally, in order to better illustrate the effectiveness of the proposed algorithm, this paper adds a BOS experiment. Detailed analysis will be carried out below.
Experiment 1. Algorithm verification on simulation data.
The proposed method was first validated on the data generated by the simulation method in this paper. Using synthetic data does not address certain complexities of experimental data, including 3D effects, noise, and nonuniform illumination. However, these issues are intentionally avoided in the current manuscript to assess the accuracy of the method in an ideal environment where a ground truth velocity is known [6]. The results obtained by using the hybrid adaptive wavelet-based optical flow algorithm at different scales is shown in Figure 5, and Figure 6 shows the error distribution.
(a)
(b)
(c)
(d)
(a)
(b)
(c)
(d)
Experiment 2. Algorithm comparison on DNS datasets.
DNS datasets are open and recognized datasets for optical flow algorithm verification. In order to further illustrate the effectiveness of the algorithm, this paper uses DNS datasets to verify the proposed algorithm and makes a comparative experiment. Firstly, it is compared with a single-scale optical flow algorithm, such as the optical flow algorithm based on cross-correlation proposed by Thielicke and Eize [11]. Secondly, it is compared with the traditional multiscale optical flow algorithm, such as the multiscale LK algorithm proposed by Qi et al. [22] or the multiscale HS algorithm proposed by Meinhardt-Llopis et al. [23]. Finally, it is compared with the wavelet-based optical flow algorithms, such as the WOF algorithm proposed by Dérian et al. [7] and the WOF-X algorithm proposed by Schmidt et al. [6]. Figure 7 shows the real optical flow estimation results and the velocity field distribution. Figure 8 is the result of the optical flow algorithm based on cross-correlation. Figures 9 and 10 are the flow estimation results of two multiscale optical flow algorithms. Figures 11–13 are the flow estimation results of three wavelet-based optical flow algorithms.
Figures 11–13 show the calculation results of three wavelet-based optical flow algorithms. The coarse scale and the fine scale of these three optical flow algorithms are 3 and 6. The differences of the three algorithms are as follows: firstly, the WOF algorithm does not use the regular term constraint, the WOF-X algorithm uses the regular term constraint of fixed parameters, and the algorithm proposed in this paper uses the regular term constraint of adaptive parameters; secondly, the random values of the first two algorithms are used as the initial values of the velocity field, but the values calculated by the traditional optical flow algorithm in this paper are used as the initial values of the velocity field.
(a)
(b)
(a)
(b)
(c)
(a)
(b)
(c)
(a)
(b)
(c)
(a)
(b)
(c)
(a)
(b)
(c)
(a)
(b)
(c)
Experiment 3. BOS experiment.
In order to illustrate the effectiveness of the proposed algorithm, a practical BOS experiment is carried out. The experimental equipment is shown in Figure 14(a), including camera, hot air gun, and background pattern. Figure 14(b) is a reference image (background pattern). Figure 14(c) shows the experimental results (the image with the flow field).
Figure 15 is the result of the optical flow algorithm based on cross-correlation. Figures 16 and 17 are the flow estimation results of two multiscale optical flow algorithms. Figures 18–20 are the flow estimation results of three wavelet-based optical flow algorithms.
(a)
(b)
(c)
(a)
(b)
(c)
(a)
(b)
(c)
(a)
(b)
(c)
(a)
(b)
(c)
(a)
(b)
(c)
(a)
(b)
(c)
5.1. Analysis of Experimental Results
In this paper, root-mean-square endpoint error (RMSE) is used to quantitatively evaluate the experimental results [6], as shown in the following formula:
The lower the values of RMSE, the better the displacement filed. The evaluation results of Experiment 1, Experiment 2, and Experiment 3 are shown in Table 2, Table 3, and Table 4, respectively.
Comparing the RMSE values of different scales, as shown in Table 2, it is obvious that the larger the scale, the higher the accuracy of the hybrid adaptive wavelet-based optical flow algorithm. This also accords with the principle of the wavelet transform. The larger the scale, the more abundant are the details contained in the wavelet coefficients. To ensure the accuracy and improve the efficiency of the algorithm, it is necessary to choose the appropriate resolution scale. In order to illustrate the effect of the resolution scale on the accuracy of the algorithm, the results of four scales are presented in Table 2. The optimal resolution scale calculated by the adaptive scale selection method proposed in this paper is scale = 6. Table 2 shows that the accuracy of the algorithm in scale = 6 is similar to that in scale = 7, but the calculation amount of the algorithm is reduced by four times. However, the accuracy of the algorithm in scale = 6 is twice as high as that in scale = 5, and four times as high as that in scale = 4. This shows that the adaptive scale selection method is effective and feasible. This not only ensures the accuracy of the algorithm but also improves the efficiency of the algorithm. If the error allows, a lower resolution scale can be selected, such as scale = 5 or scale = 4.
Comparing the RMSE values of different algorithms, as shown in Table 2, the RMSE value of our method is the smallest. As can be seen from Figure 11(c), the error of the WOF-X algorithm in some regions with a large refractive index gradient is large, which may be caused by the inability of the fixed regularization constraint parameters to adapt to some regions. Compared with the WOF-X algorithm, the algorithm proposed in this paper has higher calculation accuracy and smaller calculation error in the whole image area, as shown in Figure 12(c). This shows that the adaptive regularization constraint method proposed in this paper is effective and feasible.
In order to observe the accuracy of each algorithm more intuitively, this paper extracts the 50th, 100th, 150th, and 200th lines of the estimated velocity field of each algorithm and compares them with the corresponding lines of the true velocity field, as shown in Figure 21. According to Table 5, the velocity field estimated by the HS algorithm, the WOF-X algorithm, and the algorithm proposed in this paper has a high coincidence with the true velocity field. Through the above analysis, it is easy to draw the following conclusion: the accuracy of our algorithm is higher than that of other algorithms mentioned in this paper.
(a)
(b)
(c)
(d)
In this paper, the effectiveness of the proposed algorithm is verified by the actual BOS experiment, and the motion-compensated image is reconstructed according to the calculated displacement field , as shown in Figure 16(b). The motion compensation image was compared with the experimental image (image with flow field， as shown in Figure 14(c)) by RMSE, PSNR (Peak Signal to Noise Ratio), and SSIM (structural similarity). PSNR is the most common and widely used objective image evaluation index. However, it is based on the error between corresponding pixels, that is, the image quality evaluation based on error sensitivity. The larger the PSNR value, the smaller the distortion. SSIM, an image quality evaluation index, evaluates image quality based on structural similarity. It measures the image similarity from three aspects: brightness, contrast, and structure. The higher the SSIM value, the higher the image similarity. Comparing the RMSE values, the PSNR values, and the SSIM values of different algorithms, as shown in Table 4, the RMES value of our method is the smallest and the PSNR and SSIM values of our method are the highest. According to the RMSE, PSNR and SSIM values in Table 4, the proposed method is superior to other comparison methods.
6. Conclusions
In order to improve the accuracy of the BOS device in transparent flow field measurement, this paper presents the hybrid adaptive wavelet-based optical flow algorithm for Background Oriented Schlieren (BOS) experiments. The algorithm can effectively estimate the displacement field. Firstly, this paper simulates the BOS experiment and uses the simulation data to verify the effectiveness of the proposed algorithm. In Experiment 1, simulation data are used to verify the proposed algorithm, and the influence of the resolution scale on the accuracy and efficiency of the proposed algorithm is discussed. The results of Experiment 1 show that the algorithm proposed in this paper is effective. Secondly, in order to further prove the effectiveness of the algorithm proposed in this paper, in Experiment 2, this paper uses DNS datasets to verify the algorithm proposed in this paper and compares it with the existing related algorithms. The results of Experiment 2 show that the accuracy of this algorithm is better than the other algorithms mentioned in this paper. Especially, in the region of large refractive index gradient, the accuracy of the algorithm proposed in this paper is higher than that of other comparison algorithms. It also shows that the proposed algorithm is effective and feasible. Lastly, the method proposed in this paper is verified by the actual BOS experiment. The results of Experiment 3 show that the accuracy of this algorithm is better than other algorithms mentioned in this paper.
Data Availability
The data used to support the findings of this study are available from the corresponding author upon request.
Conflicts of Interest
The authors declare that there are no conflicts of interest regarding the publication of this paper.
Acknowledgments
The research was supported by the Key Laboratory of Transient Impact Technology of China under project no. 614260603030817.