Multispectral light field imaging system, and three-dimensional temperature field measurement method and measurement system

Abstract

Disclosed in the present invention are a multispectral light field imaging system, and a three-dimensional temperature field measurement method and measurement system. The multispectral light field imaging system comprises: a primary lens, a micro-lens array, a first relay lens, a dichroic mirror, a second relay lens, a third relay lens, a first image sensor, a first optical filter, a second image sensor and a second optical filter, wherein the micro-lens array is located at an image plane of the primary lens; the first relay lens, the second relay lens and the third relay lens together form a 1:1 relay lens group imaging system on an optical path of a spectral band I and an optical path of a spectral band II; an angle of 45° is formed between the dichroic mirror and an optical axis; and the first image sensor is located at a focal plane of the second relay lens, and the second image sensor is located at a focal plane of the third relay lens. The present invention realizes the collection of multispectral light field data while maintaining the amount of light field data of each spectrum unchanged, such that an instantaneous three-dimensional temperature field of a flow field can be measured.

Description

Multispectral light field imaging system, three-dimensional temperature field measurement method and measurement system Technical Field

[0001] The present invention relates to the technical field of fluid measurement, in particular to high-precision measurement of a fluid's three-dimensional temperature field. Background Art

[0002] High-precision measurement of the three-dimensional temperature field of high-temperature combustion gases in narrow turbine passages is crucial for optimizing the structure and improving the performance of aircraft engine turbine blades. However, the complex mechanical structure of these blades severely limits the application of traditional contact flow field measurement techniques. Among non-contact flow field temperature measurement technologies, laser-induced phosphorescence (LIP) temperature measurement is a non-contact optical temperature measurement technique based on the thermal quenching effect of phosphorescence. It offers advantages such as minimal optical window requirements and no interference with the flow field structure. It has great potential for measuring complex three-dimensional flow fields with limited optical windows. Currently, the more established LIP temperature measurement methods include the absolute intensity method, the lifetime decay method, and the intensity ratio method. The intensity ratio method measures temperature based on the relationship between the phosphorescence spectral intensity ratio and temperature. Compared to the absolute intensity and lifetime decay methods, it is less susceptible to factors such as uneven phosphorescence particle distribution and is therefore suitable for temperature measurement in dynamic flow fields. The position and spectral intensity ratio of phosphorescent tracer particles are used to invert the flow field temperature field. Therefore, multispectral imaging of the phosphorescent particle field and the calculation of the spectral intensity ratio are key components of flow field temperature measurement.

[0003] For multispectral imaging of phosphorescent particle fields, laser-induced phosphorescence flow field temperature measurement systems primarily employ traditional imaging modes, which can only acquire multispectral images of phosphorescent particles within the focal plane. Consequently, they can only achieve instantaneous measurement of the two-dimensional flow field temperature, but not the instantaneous three-dimensional flow field temperature. Within computational optical imaging, spectral light field imaging, by coupling spectral and light field imaging principles, can capture multispectral light field information of three-dimensional objects. Combined with tomographic reconstruction algorithms, it can reconstruct the three-dimensional spectral intensity ratio of the target object. Therefore, combining spectral light field imaging with laser-induced phosphorescence intensity ratio temperature measurement technology is expected to enable measurement of the instantaneous three-dimensional flow field temperature.

[0004] Currently, there are two main architectures for spectral light field camera systems: one that couples the aperture of the main lens of the light field camera with a filter array, and the other that uses an RGB camera as the image sensor. Spectral light field imaging systems coupled with a filter array have the following three issues: ① Because the aperture is evenly divided by the spectral filters in the filter array, the amount of light field data for each spectrum is inversely proportional to the number of spectra; ② Due to limitations in the filter array's process conditions, multispectral data aliasing occurs at the intersection of the spectral sub-aperture images of the spectral light field camera, resulting in invalid spectral light field data at the intersection of the sub-aperture images; ③ The greater the distance the imaging area deviates from the optical axis, the greater the aberration of the imaging system, and the size of the spectral sub-aperture image changes with the position of the imaging area. Spectral light field imaging systems coupled with a color RGB camera have the following issues: The Bayer filter coupled to the color RGB camera limits its applicability to collecting multispectral light field data for specific objects.

[0005] After searching the prior art, the Chinese invention patent number is CN201910463400.9, and the invention name is Multispectral temperature measurement method and system for high-temperature components based on light field camera. The invention method adopts a spectral light field imaging system of the light field camera main lens aperture coupled with a filter array. This system structure has the following four problems: ① While increasing the number of spectra, the amount of light field data of a single spectrum is reduced; ② Due to the aliasing of multispectral data at the junction of spectral sub-images, some image data is invalid; ③ Affected by aberrations, the size of the spectral sub-aperture image changes with the change of the imaging area position; ④ Since the multi-spectral light field information is imaged on the same image sensor, the spectral data processing method is complex and the data processing process is cumbersome.

[0006] To calculate the phosphorescence spectral intensity ratio, the two-dimensional temperature field measurement of the laser-induced phosphorescence flow field simply requires taking the pixel-by-pixel ratio of the multispectral image of the phosphorescent particles to calculate the phosphorescence spectral intensity ratio in the two-dimensional plane. To calculate the spectral intensity ratio of the three-dimensional phosphorescent particle field based on the spectral light field image solution, the spectral intensity distribution of the three-dimensional phosphorescent particle field must first be reconstructed using a tomographic reconstruction algorithm. The spectral intensity ratio of the three-dimensional phosphorescent particle field is then calculated by taking the ratio voxel-by-voxel.

[0007] A prior art search revealed a Chinese invention patent numbered CN201710562875.4, titled "Three-Dimensional Flow Field Measurement Method Based on Dual Light Field Cameras." This invention reduces the stretching effect of the reconstruction results by taking the intersection of the solutions. However, this invention requires an additional light field camera system, increasing the optical window requirements. Therefore, it is not suitable for complex three-dimensional flow field measurements where the optical window is limited.

[0008] Summary of the Invention

[0009] The first technical problem to be solved by the present invention is to provide a multispectral light field imaging system in which the amount of light field data of each spectrum remains unchanged.

[0010] The second technical problem to be solved by the present invention is to provide a three-dimensional temperature field measurement method and measurement system with small intensity reconstruction error and high reconstruction quality.

[0011] In order to solve the above technical problems, the technical solution adopted by the present invention is:

[0012] The present invention first provides a multispectral light field imaging system, comprising: a main lens, a microlens array, a first relay lens, a dichroic mirror, a second relay lens, a third relay lens, a first image sensor, a first filter, a second image sensor and a second filter; the microlens array is located at the image plane of the main lens; the distance between the first relay lens and the equivalent surface of the microlens array is f+F, f is the focal length of the microlens array, and F is the focal length of the first relay lens; the first relay lens, the second relay lens, and the third relay lens respectively form a 1:1 relay lens imaging system on the optical path of spectral band one and the optical path of spectral band two; the dichroic mirror is at a 45-degree angle to the optical axis; the first image sensor is located at the focal plane of the second relay lens, and the second image sensor is located at the focal plane of the third relay lens.

[0013] The present invention also provides a method for measuring a three-dimensional temperature field of a laser-induced phosphorescence flow field, comprising:

[0014] Acquire light field images captured by the first image sensor and the second image sensor respectively;

[0015] According to the acquired light field image, a three-dimensional spectral intensity distribution corresponding to the light field image is obtained;

[0016] According to the obtained three-dimensional spectral intensity distribution, the three-dimensional temperature distribution is obtained by using the intensity ratio method.

[0017] According to the acquired light field image, the three-dimensional spectral intensity distribution corresponding to the light field image is obtained, including:

[0018] According to the acquired light field image, the corresponding relationship between the target object position, the target object spectral intensity and the light field image is established:

[0019] Where λ represents the spectral wavelength, (m,n) is the pixel coordinate, (i,j,k) is the voxel coordinate, represents the grayscale value of the pixel (m,n) of the λ spectrum light field image, represents the λ spectrum intensity value of the (i, j, k) voxel, The ratio of the grayscale contribution of the (i, j, k) voxel to the (m, n) pixel to the λ spectrum intensity value of the (i, j, k) voxel is the weight coefficient;

[0020] The three-dimensional discrete voxels of the control volume are rearranged into voxel column vectors in the order of "row first, then column, and finally vertical", and any original discrete voxel (i, j, k) is numbered as (k-1)·I·J+(j-1)·I+i in the rearranged column vector; the two-dimensional discrete pixels of the light field image are rearranged into pixel column vectors in the order of "row first, then column", and any original pixel (m, n) is numbered as (n-1)·M+m in the rearranged pixel column vector; taking r = (n-1)·M+m, c = (k-1)·I·J+(j-1)·I+i, R = M·N, C = I·J·K, then the correspondence between the control volume voxel column vector and the light field image pixel column vector is expressed as:

[0021] The weight coefficient matrix is ​​A λ , the pixel column vector of the light field image is B λ , the column vector of the control volume voxel is X λ , then formula (2) is expressed as A λ X λ =B λ ;

[0022] Based on the characteristics of the tangent circle imaging area of ​​the light field camera, the weight coefficient matrix and A are realized by eliminating invalid pixels of the light field image and their corresponding weight coefficients. λ X λ =B λ The simplification of the equation group; the criterion for invalid pixels is: the pixels outside the tangent circle of the light field image are invalid pixels; by eliminating B λ Invalid pixels and A λ Invalid weight coefficients corresponding to invalid pixels in A λ X λ =B λ The simplified equations transform the pre-processed light field tomography reconstruction equation into A λ 2X λ =B λ 2;

[0023] The weight coefficient matrix is ​​calculated by ray tracing method. The calculation formula of ray tracing method is:

[0024] Among them, (x, y) is the coordinate of the light source, (x', y') is the coordinate of the pixel where the light enters, is the solid angle of the outgoing ray, is the solid angle of the incident light; f is the focal length of the microlens unit, F is the focal length of the relay lens 1, S1 is the object distance, S2 is the image distance, (S x ,S y) is the displacement of the center of the microlens array from the optical axis of the imaging system;

[0025] Tomographic reconstruction of phosphorescent particle three-dimensional spectral intensity:

[0026] The preliminary reconstruction result X is obtained by SART algorithm hierarchical reconstruction λ 1. The iterative calculation formula is:

[0027] Where t is the number of iterations, μ is the relaxation factor, p is the voxel number, b is the pixel gray value, and o is the pixel number. is the ratio of the grayscale contribution of the p voxel to the o pixel to the λ spectral intensity value of the voxel, represents the sum of the forward projections of the λ spectral intensity values ​​of all non-zero voxels within the line of sight of pixel o, represents the sum of the weight coefficients of non-zero voxels within the line of sight of pixel o to pixel o, represents the difference between the grayscale value of the pixel o and the sum of the forward projections of the λ spectral intensity values ​​of all non-zero voxels within the pixel's line of sight.

[0028] The maximum inter-class variance method was used to binarize the preliminary reconstruction result X1, and the voxel intensity values ​​below the binarization threshold were set to zero to eliminate the overlapping effect of the tracer particles.

[0029] Using the maximum between-class variance method, X λ 1 The method for binarization is:

[0030] Assuming the intensity threshold is τ, then X λ 1 The voxel spectral intensity will be divided into two categories: A and B, where A is greater than τ and B is less than τ; traverse X λ 1 Voxel spectral intensity distribution, when the variance of the two categories A and B is the largest, the intensity threshold τ is X λ A binarization threshold of 1.

[0031] The specific calculation method is: First, use the frequency statistics function to count X λ 1, assuming that the intensities are arranged in ascending order as E1, E2, E3…E n , the probability of each intensity is p1, p2, p3…p n ; Set the intensity threshold to E k (1≤k≤n), the probability of the spectral intensity of voxels of type A is p A (k), the average voxel spectral intensity is m A (k), the probability of the spectral intensity of the voxel of type B is p B (k), the average voxel spectral intensity is m B (k), X λ The average spectral intensity of 1 voxel is mG ; Then the variances of types A and B are expressed as:

[0032] Go through E1, E2, E one by one 3… E n , when the value of the expression shown in formula (5) is the largest, E k That is X λ A binarization threshold of 1;

[0033] X λ 1 Medium strength lower than E k The voxel spectral intensity is set to zero to obtain the X without particle overlap. λ 2;

[0034] According to the corresponding relationship between the voxel column vector number and the voxel three-dimensional coordinate, X λ 2 is restored to a three-dimensional voxel matrix X λ 3. Use Gaussian fitting to obtain X λ 3D connected voxel where the center of gravity is located.

[0035] X λ Any three-dimensional connected body X in 3 λ The calculation method for the voxel where the center of gravity of 3(i:i+Δi,j:j+Δj,k:k+Δk) is:

[0036] Slice the three-dimensional connected body along the plane perpendicular to the Z axis to obtain X λ 3(i:i+Δi,j:j+Δj,n):

[0037] Where, k≤n≤k+Δk;

[0038] Identify X using a two-dimensional Gaussian fitting function λ The centroid of 3(i:i+Δi,j:j+Δj,n):

[0039] Among them, i≤x≤i+Δi, j≤y≤j+Δj, E λ is the amplitude of the Gaussian distribution of the voxel spectral intensity, (x0, y0, n) is the voxel coordinate of the centroid of the two-dimensional slice, σ x , σ y are the half-peak width of Gaussian distribution of spectral intensity of voxel in two-dimensional slices;

[0040] A two-dimensional Gaussian fitting function is applied to each slice of the three-dimensional connected body to obtain the centroid of all slices; a one-dimensional Gaussian fitting is performed on the spectral intensity of all voxels where the centroid is located to determine the center of gravity of the three-dimensional connected body:

[0041] Among them, (x0, y0, z0) is a three-dimensional connected body Xλ 3(i:i+Δi,j:j+Δj,k:k+Δk) voxel coordinates of the center of gravity, σ z is the half-maximum width of the Gaussian distribution of the voxel spectral intensity;

[0042] According to the center of gravity of the three-dimensional connected body, calculate X λ All three-dimensional connected voxels where the center of gravity is located within 3, mark the three-dimensional coordinates of the voxel where the center of gravity is located, and eliminate X λ All voxels other than the voxel where the center of gravity is located in 3 get X λ 4. X λ 4 is reduced to a one-dimensional voxel column vector X λ 5. Meanwhile, A λ 2 will be simplified to A λ 3. The light field tomography reconstruction equation is transformed into A λ 3X λ 5=B λ 2;

[0043] To A λ 3X λ 5=B λ 2 Perform the second SART tomographic reconstruction to obtain X λ 6;

[0044] Finally, according to the correspondence between the three-dimensional coordinates of the centroid voxel and the voxel column vector number, X λ 6 is restored to the three-dimensional flow field control volume X λ 7.

[0045] The spectral light field imaging system proposed in this paper utilizes a cage-type light field camera system with a relay lens coupled to a dichroic mirror. This allows light field data from different spectral bands to be imaged on separate image sensors, achieving spectral data separation. This avoids the loss of single-spectral light field data, the aliasing of multi-spectral light field data, and the inconsistency of spectral subaperture image size due to aberrations as the imaging area changes position. Furthermore, when the wavelength of the spectrum to be captured changes, only the dichroic mirror and corresponding filter need to be replaced, resulting in high system flexibility and convenient data processing.

[0046] Compared with the prior art, the present invention has the following advantages:

[0047] (1) The amount of spectral light field data does not decrease as the number of spectra increases; existing spectral light field cameras achieve multi-spectral light field data acquisition by adding a filter array to the aperture of the light field camera's main lens. Since the aperture aperture is evenly divided by the spectral filters of the filter array, the amount of light field data for each spectrum is inversely proportional to the number of spectra. The present invention achieves multi-spectral light field data acquisition by adding a dichroic mirror inside the 1:1 relay lens group of the cage-type light field camera, and the amount of light field data for each spectrum remains unchanged.

[0048] (2) There is no multi-spectral light field data aliasing. Due to the process conditions of the filter array, the existing spectral light field camera has multi-spectral light field data aliasing at the intersection of each spectral sub-aperture image, which makes the spectral data at the intersection of the sub-aperture image invalid. However, in the present invention, the light field data of different spectra are imaged on different image sensors, and there is no multi-spectral light field data aliasing.

[0049] (3) There is no problem of inconsistent spectral sub-aperture image sizes due to the influence of aberrations. In existing spectral light field cameras, multi-spectral light field data are imaged on the same image sensor. The greater the distance the photographed area deviates from the optical axis, the greater the aberration of the imaging system, and the sizes of the sub-aperture images of each spectrum will become inconsistent. In the present invention, light field data of different spectra are imaged on different image sensors, and there is no problem of inconsistent sizes of the sub-aperture images of each spectrum.

[0050] (4) The three-dimensional spectral intensity reconstruction process is more convenient. In existing spectral light field cameras, multi-spectral light field data are imaged on the same image sensor, and the light field data of each spectrum needs to be obtained by extracting sub-aperture images. In the present invention, the light field data of different spectra are imaged on different image sensors. The light field data on each image sensor is the light field data of the same spectrum, and there is no need to extract sub-aperture images, which makes the reconstruction calculation more convenient.

[0051] (5) The reconstruction results are of higher quality; the position and spectral intensity ratio of the tracer particles are used for the inversion calculation of the flow field temperature field. Therefore, the accurate reconstruction of the three-dimensional spectral intensity ratio of the phosphorescent particle field is the key core link of the flow field temperature measurement. When the SART algorithm is used to reconstruct the three-dimensional particle field of light field imaging by tomography, the reconstruction results are severely stretched along the depth direction of the control volume, resulting in a large error in the intensity reconstruction results, which seriously affects the accuracy of the flow field measurement results. The present invention utilizes the characteristic that the particle intensity reconstruction results of the SART algorithm are approximately Gaussian distributed along the stretching direction. Based on the principle of Gaussian fitting to obtain the center of gravity of the three-dimensional connected body, a light field imaging three-dimensional particle field reconstruction method combining Gaussian fitting positioning technology and SART algorithm is proposed. The initial three-dimensional particle field is reconstructed using the SART algorithm, and then the overlap effect between the tracer particles in the particle field is eliminated by the maximum inter-class variance method. Then, the Gaussian fitting principle is used to locate the voxel where the center of gravity of the stretched connected body is located, and the intensity value of the voxel where the center of gravity is located is reconstructed again by the SART algorithm, thereby reducing the intensity reconstruction error and improving the reconstruction quality.

[0052] (6) Realizing the measurement of the instantaneous three-dimensional temperature field of the flow field; Existing laser-induced phosphorescence flow field temperature measurement systems mainly use traditional imaging modes and can only obtain multispectral images of phosphorescent particles in the focal plane. Therefore, they can only realize the instantaneous measurement of the two-dimensional temperature field of the flow field, and cannot obtain the instantaneous three-dimensional temperature field of the flow field. The present invention proposes a method for measuring the instantaneous three-dimensional temperature field of the flow field by combining spectral light field imaging technology with laser-induced phosphorescence intensity ratio temperature measurement technology. BRIEF DESCRIPTION OF THE DRAWINGS

[0053] Figure 1: Laser-induced phosphorescence flow field three-dimensional temperature field measurement method and system based on a spectroscopic spectral light field camera;

[0054] Figure 2 is a schematic diagram of the structure of a spectroscopic spectral light field imaging system;

[0055] Figure 3 Schematic diagram of the tangent circle imaging area of ​​a light field camera;

[0056] Figure 4 Flowchart of the 3D spectral intensity reconstruction method for spectral light field imaging combining Gaussian fitting positioning technology with the SART algorithm. DETAILED DESCRIPTION

[0057] Below in conjunction with accompanying drawing, the present invention is described in detail:

[0058] The technical solution of the present invention is shown in Figure 1. A 355nm laser is expanded into a bulk laser by a laser beam expander to excite phosphorescent particles in the flow field to be measured. Two image sensors of the spectroscopic laser-induced phosphorescence spectral light field imaging system simultaneously collect phosphorescent light field images of phosphorescent particles in two different bands in the same flow field area to be measured. The spectral light field imaging three-dimensional spectral intensity reconstruction algorithm combining the Gaussian fitting positioning technology proposed in the present invention with the SART algorithm is used to simultaneously solve the phosphorescent light field images of the two bands and reconstruct the spectral intensity distribution of the three-dimensional phosphorescent particle field. The intensity ratio method is used to calculate the spectral intensity ratio of the three-dimensional phosphorescent particle field. Based on the temperature calibration result of the phosphorescent spectral intensity ratio, the three-dimensional temperature field distribution of the flow field to be measured is solved.

[0059] Based on the above technical solution, the present invention can be divided into three parts, namely: (1) the structure of the spectroscopic spectral light field imaging system; (2) the spectral light field imaging three-dimensional spectral intensity reconstruction algorithm combining Gaussian fitting positioning technology with the SART algorithm; and (3) the reconstruction of the three-dimensional temperature field of the flow field based on the intensity ratio method.

[0060] Structure of spectroscopic light field imaging system

[0061] The structure of the spectroscopic spectral light field imaging system is shown in Figure 2, which includes a primary lens, a microlens array, a relay lens 1, a dichroic mirror, relay lens 2-1, relay lens 2-2, image sensor 1, filter 1, and image sensor 2, filter 2. The microlens array is located at the image plane of the primary lens. Assume that the focal length of the microlens unit is f, the focal length of relay lens 1 is F, the distance between relay lens 1 and the equivalent surface of the microlens array is f+F, relay lens 1, relay lens 2-1, and relay lens 2-2 form a 1:1 relay lens imaging system on the optical path of spectral band 1 and spectral band 2, respectively. The dichroic mirror is oriented at 45 degrees to the optical axis. Image sensors 1 and 2 are located at the focal planes of relay lenses 2-1 and 2-2, respectively.

[0062] Three-dimensional spectral intensity reconstruction algorithm for spectral light field imaging combining Gaussian fitting positioning technology with SART algorithm

[0063] The position and spectral intensity ratio of phosphorescent tracer particles are used to invert the flow field temperature field. Therefore, calculating the three-dimensional spectral intensity ratio of phosphorescent particles is a key step in flow field temperature measurement. The present invention first constructs a spectral light field tomography reconstruction model based on the relationship between the position of phosphorescent tracer particles, the spectral intensity of phosphorescent tracer particles, and the light field image. The weight coefficient matrix of the spectral light field tomography reconstruction model is calculated based on the ray tracing method. The spectral light field tomography reconstruction model is then simplified by combining the tangent circle imaging area characteristics of the light field camera. Finally, the stretching effect of the three-dimensional spectral intensity reconstruction results of phosphorescent particles is weakened based on Gaussian fitting positioning technology.

[0064] (1) Spectral light field tomography reconstruction model

[0065] As shown in Figure 2, the light emitted from the target object in different directions is converged by the main lens, imaged twice by the microlens array, expanded by the relay lens 1, and split by the dichroic mirror. The light in band λ1 is converged by the relay lens 2-1, filtered by the filter 1, and then enters different pixels of the image sensor 1, completing the photoelectric conversion. Similarly, the light in band λ2 is converged by the relay lens 2-2, filtered by the filter 2, and then enters different pixels of the image sensor 2, completing the photoelectric conversion. The pixel position information and grayscale information of light field image 1 and light field image 2 correspond to the direction information, wavelength information, and spectral intensity information of the light emitted by the target object, respectively. Based on this correspondence, a correspondence can be established between the target object position, the target object spectral intensity, and the light field image.

[0066] Where λ represents the spectral wavelength, (m,n) is the pixel coordinate, (i,j,k) is the voxel coordinate, represents the grayscale value of the pixel (m,n) of the λ spectrum light field image, represents the λ spectrum intensity value of the (i, j, k) voxel, The ratio of the grayscale contribution of the (i, j, k) voxel to the (m, n) pixel to the λ spectrum intensity value of the (i, j, k) voxel is the weight coefficient;

[0067] For ease of calculation, the three-dimensional discrete voxels of the control volume are rearranged into voxel column vectors in the order of "rows first, then columns, and finally verticals". Then, any original discrete voxel (i, j, k) is numbered as (k-1)·I·J+(j-1)·I+i in the rearranged column vector. At the same time, the two-dimensional discrete pixels of the light field image are rearranged into pixel column vectors in the order of "rows first, then columns". Then, any original pixel (m, n) is numbered as (n-1)·M+m in the rearranged pixel column vector. Taking r = (n-1)·M+m, c = (k-1)·I·J+(j-1)·I+i, R = M·N, C = I·J·K, the corresponding relationship between the control volume voxel column vector and the light field image pixel column vector can be expressed as,

[0068] The weight coefficient matrix is ​​A λ , the pixel column vector of the light field image is B λ , the column vector of the control volume voxel is X λ , then formula (2) can be expressed as A λ X λ =B λ . Through the light field image pixel grayscale information B λ Reconstructing the voxel positions and spectral intensity information of the control volume tracer particles is an inverse problem typically solved iteratively using the more stable and parallel SART algorithm. However, since the SART algorithm requires repeated use of the weight matrix during calculation, and the light field weight matrix occupies a large amount of memory, the reconstruction process is time-consuming.

[0069] In order to reduce the memory occupied by the weight coefficient matrix without repeatedly calculating it, this paper implements the weight coefficient matrix and A by eliminating invalid pixels of the light field image and their corresponding weight coefficients based on the characteristics of the tangent circle imaging area of ​​the light field camera (see Figure 3). λ X λ =B λ The simplification of the equations can reduce the time of tomographic calculation. Among them, the criterion for judging invalid pixels is: no matter whether the photographed object is in focus or not, its image formed by the light field camera will not exceed the tangent circle area shown in Figure 3, so the pixels outside the tangent circle of the light field image are invalid pixels. By eliminating B λ Invalid pixels and A λ The invalid weight coefficients corresponding to the invalid pixels can be realized by λ X λ =B λThe simplification of the equations reduces the amount of tomographic calculations and thus reduces the time required for tomographic reconstruction. The pre-processed light field tomographic reconstruction equation is transformed into A λ 2X λ =B λ 2.

[0070] (2) Obtaining the weight coefficient matrix

[0071] The weight coefficient matrix is ​​calculated by ray tracing method. The calculation formula of ray tracing method is:

[0072] Among them, (x, y) is the coordinate of the light source, (x', y') is the coordinate of the pixel where the light enters, is the solid angle of the outgoing ray, is the solid angle of the incident light; f is the focal length of the microlens unit, F is the focal length of the relay lens 1, S1 is the object distance, S2 is the image distance, (S x ,S y ) is the displacement of the center of the microlens array from the optical axis of the imaging system.

[0073] (4) Three-dimensional spectral intensity tomographic reconstruction of phosphorescent particles

[0074] Aiming at the stretching problem of SART algorithm in spectral light field imaging 3D particle field reconstruction, this paper optimizes the reconstruction result of SART algorithm by Gaussian fitting to obtain the center of gravity of 3D connected particles. λ Invalid voxels and A λ The invalid weight coefficients corresponding to the invalid voxels in 2 are used to achieve A λ 2X λ =B λ The further simplification of the 2 equations eliminates the stretching phenomenon and improves the reconstruction accuracy of the three-dimensional position and intensity of the tracer particles.

[0075] Due to X λ is an unknown quantity, so it is necessary to use the SART algorithm to reconstruct the initial reconstruction result, that is, X λ 1. Its iterative calculation formula is,

[0076] Where t is the number of iterations, μ is the relaxation factor, p is the voxel number, b is the pixel gray value, and o is the pixel number. is the ratio of the grayscale contribution of the p voxel to the o pixel to the λ spectral intensity value of the voxel, represents the sum of the forward projections of the λ spectral intensity values ​​of all non-zero voxels within the line of sight of pixel o, represents the sum of the weight coefficients of non-zero voxels within the line of sight of pixel o to pixel o, represents the difference between the grayscale value of the pixel o and the sum of the forward projections of the λ spectral intensity values ​​of all non-zero voxels within the pixel's line of sight.

[0077] Since the reconstruction results of the SART algorithm are severely stretched in the depth direction, the reconstruction results of the multi-particle field are prone to particle stretching overlap, that is, two or more tracer particles become the same three-dimensional connected body due to the stretching effect. Therefore, before performing Gaussian fitting to extract the center of gravity of the three-dimensional connected body, this paper uses the maximum between-class variance method (also known as the Otsu method, abbreviated as OTSU) to perform X λ 1. Binarization is performed and the voxel intensity values ​​below the binarization threshold are set to zero, thereby eliminating the overlapping effect of the tracer particles. λ 1 The method of binarization threshold is: Assuming the intensity threshold is τ, then X λ 1 The voxel spectral intensity will be divided into two categories: A (greater than τ) and B (less than τ). λ 1 Voxel spectral intensity distribution, when the variance of the two categories A and B is the largest, the intensity threshold τ is X λ The specific calculation method is: first use the frequency statistics function to count X λ 1, assuming that the intensities are arranged in ascending order as E1, E2, E 3… E n , the probability of each intensity is p1, p2, p 3… p n ; Set the intensity threshold to E k (1≤k≤n), the probability of the spectral intensity of voxels of type A is p A (k), the average voxel spectral intensity is m A (k), the probability of the spectral intensity of the voxel of type B is p B (k), the average voxel spectral intensity is m B (k), X λ The average spectral intensity of 1 voxel is m G ; Then the variances of types A and B can be expressed as,

[0078] Go through E1, E2, E one by one 3… E n , when the value of the expression shown in formula (5) is the largest, E k That is X λ A binarization threshold of 1.

[0079] X λ 1 Medium strength lower than E k The voxel intensity of zero can be obtained to obtain X without particle overlap. λ 2.

[0080] According to the corresponding relationship between the voxel column vector number and the voxel three-dimensional coordinate, X λ 2 is reduced to a three-dimensional matrix X λ 3. Use Gaussian fitting to obtain X λ 3D connected voxel where the center of gravity is located. λ Any three-dimensional connected body X in 3 λ The calculation method of the voxel where the center of gravity of 3(i:i+Δi,j:j+Δj,k:k+Δk) is as follows: slice the three-dimensional connected body along the plane perpendicular to the Z axis to obtain the X λ 3(i:i+Δi,j:j+Δj,n), that is,

[0081] Where k≤n≤k+Δk. Then use the two-dimensional Gaussian fitting function to identify X λ The centroid of 3(i:i+Δi,j:j+Δj,n) is expressed as,

[0082] Where i≤x≤i+Δi, j≤y≤j+Δj, E λ is the amplitude of the Gaussian distribution of the voxel spectral intensity, (x0, y0, n) is the voxel coordinate of the centroid of the two-dimensional slice, σ x , σ y are the half-peak widths of the Gaussian distribution of the spectral intensity of the two-dimensional slice voxels. A two-dimensional Gaussian fitting function is applied to each slice of the three-dimensional connected body to obtain the centroid of all slices; finally, a one-dimensional Gaussian fitting is performed on the spectral intensity of all voxels where the centroids are located to determine the center of gravity of the three-dimensional connected body. The calculation expression is,

[0083] Among them, (x0, y0, z0) is a three-dimensional connected body X λ 3(i:i+Δi,j:j+Δj,k:k+Δk) voxel coordinates of the center of gravity, σ z is the half-maximum width of the Gaussian distribution of the voxel spectral intensity.

[0084] According to the three-dimensional connected gravity center calculation method, calculate X λ All three-dimensional connected voxels where the center of gravity is located within 3, mark the three-dimensional coordinates of the voxel where the center of gravity is located, and eliminate X λ All voxels other than the voxel where the center of gravity is located in 3 get X λ 4. X λ 4 is reduced to a one-dimensional voxel column vector X λ 5. Meanwhile, A λ 2 will be simplified to A λ 3. Therefore, the light field tomography reconstruction equation is transformed into A λ 3X λ 5=B λ 2. Pair A λ 3Xλ 5=B λ 2 Perform the second SART tomographic reconstruction to obtain X λ 6. Finally, according to the correspondence between the three-dimensional coordinates of the centroid voxel and the voxel column vector number, X λ 6 is restored to the three-dimensional flow field control volume X λ 7.

[0085] In summary, the process of the spectral light field imaging three-dimensional spectral intensity reconstruction method combining Gaussian fitting positioning technology with SART algorithm is shown in Figure 4.

[0086] Reconstruction of three-dimensional temperature field of flow field based on intensity ratio method

[0087] Based on the three-dimensional spectral intensity reconstruction process, the light field images of the λ1 and λ2 phosphorescence spectra are solved respectively, and the three-dimensional spectral intensity X of the phosphorescence particle field λ1 and λ2 phosphorescence spectra are reconstructed. λ1 7 and X λ2 7. To X λ1 7 and X λ2 7 Take the ratio of each non-zero voxel and calculate the spectral intensity ratio of the three-dimensional phosphorescent particle field. The calculation formula is:

[0088] Where (i, j, k) is the three-dimensional coordinate of the voxel, R (i,j,k) is the phosphorescence spectrum intensity ratio of voxel (i, j, k);

[0089] The relationship between the phosphorescence spectrum intensity ratio and temperature is obtained through temperature calibration experiment: T (i,j,k) =f(R (i,j,k) ) (10)

[0090] Where, T (i,j,k) represents the temperature of the (i, j, k) voxel;

[0091] Based on formula (10), the spectral intensity ratio result R of the three-dimensional phosphorescent particle field is (i,j,k) Solve and calculate the instantaneous three-dimensional temperature field of the flow field.

Claims

1. A multispectral light field imaging system, characterized in that: include: A main lens, a microlens array, a first relay lens, a dichroic mirror, a second relay lens, a third relay lens, a first image sensor, a first filter, a second image sensor and a second filter; the microlens array is located at the image plane of the main lens; the distance between the first relay lens and the equivalent surface of the microlens array is f+F, f is the focal length of the microlens array, and F is the focal length of the first relay lens; the first relay lens and the second relay lens and the third relay lens respectively form a 1:1 relay lens imaging system on the optical path of spectral band one and the optical path of spectral band two; the dichroic mirror is 45° with the optical axis; the first image sensor is located at the focal plane of the second relay lens, and the second image sensor is located at the focal plane of the third relay lens.

2. A method for measuring three-dimensional temperature field of laser-induced phosphorescence flow field based on the multi-spectral light field imaging system of claim 1, characterized in that: include: Acquire light field images captured by the first image sensor and the second image sensor respectively; According to the acquired light field image, a three-dimensional spectral intensity distribution corresponding to the light field image is obtained; According to the obtained three-dimensional spectral intensity distribution, the three-dimensional temperature distribution is obtained by using the intensity ratio method.

3. The method for measuring three-dimensional temperature field of laser-induced phosphorescence flow field according to claim 2, characterized in that: According to the acquired light field image, the three-dimensional spectral intensity distribution corresponding to the light field image is obtained, including: According to the acquired light field image, a corresponding relationship equation group between the position of the target object, the spectral intensity of the target object and the light field image is established; Based on the characteristics of the tangent circle imaging area of ​​the light field camera, the corresponding relationship equation group established in the previous step is simplified by eliminating invalid pixels of the light field image and their corresponding weight coefficients; Initial reconstruction was performed using the first SART algorithm; The maximum inter-class variance method is used to binarize the results of the preliminary reconstruction, and the spectral intensity values ​​of the voxels below the binarization threshold are set to zero to obtain the voxel column vector without particle overlap. According to the correspondence between the voxel column vector number and the voxel three-dimensional coordinate, the voxel column vector is restored to a three-dimensional voxel matrix, and Gaussian fitting is used to obtain the voxel where the centroid of the three-dimensional connected body is located in the three-dimensional voxel matrix; According to the center of gravity of the three-dimensional connected body, the voxels where the center of gravity of all three-dimensional connected bodies is located in the three-dimensional voxel matrix are calculated, the three-dimensional coordinates of the voxel where the center of gravity is located are marked, and all voxels except the voxel where the center of gravity is located in the three-dimensional voxel matrix are removed to obtain a simplified three-dimensional voxel matrix, and the simplified three-dimensional voxel matrix is ​​restored to a simplified voxel column vector, and the weight coefficient matrix is ​​synchronized to a simplified weight coefficient matrix to obtain a simplified light field tomography reconstruction equation; Perform a second SART tomographic reconstruction on the simplified tomographic reconstruction equation to obtain the target voxel column vector; According to the correspondence between the three-dimensional coordinates of the centroid voxel and the voxel column vector number, the target voxel column vector is restored to the target 3D voxel matrix of the flow field control volume.

4. The method for measuring three-dimensional temperature field of laser-induced phosphorescence flow field according to claim 3, characterized in that: Using the maximum inter-class variance method, the method for binarizing the initial reconstruction results of the SART algorithm is as follows: Assuming that the intensity threshold is τ, the voxel spectral intensity of the preliminary reconstruction result will be divided into two categories, A and B, where A is greater than τ and B is less than τ; traversing the voxel spectral intensity distribution of the preliminary reconstruction result, when the variance of categories A and B is the largest, the intensity threshold τ is the binarization threshold of the preliminary reconstruction result.

5. The method for measuring three-dimensional temperature field of laser-induced phosphorescence flow field according to claim 3, characterized in that: The corresponding relationship equations between the target object position, the target object spectral intensity and the light field image are established as follows: In the formula, represents the gray value of the pixel (m,n) of the λ spectrum light field image, λ represents the spectral wavelength, and (m,n) is the pixel coordinate; represents the λ spectrum intensity value of voxel (i,j,k), where (i,j,k) is the voxel coordinate; It represents the ratio of the grayscale contribution value of the (i,j,k) voxel to the (m,n) pixel to the λ spectral intensity value of the (i,j,k) voxel, which is the weight coefficient; The three-dimensional discrete voxels of the control volume are rearranged into voxel column vectors in the order of "row first, then column, and finally vertical", and any original discrete voxel (i, j, k) is numbered as (k-1)·I·J+(j-1)·I+i in the rearranged column vector; the two-dimensional discrete pixels of the light field image are rearranged into pixel column vectors in the order of "row first, then column", and any original pixel (m, n) is numbered as (n-1)·M+m in the rearranged pixel column vector; taking r=(n-1)·M+m, c=(k-1)·I·J+(j-1)·I+i, R=M·N, C=I·J·K, the corresponding relationship between the control volume voxel column vector and the light field image pixel column vector is expressed as: The weight coefficient matrix is ​​A λ , the pixel column vector of the light field image is B λ , the column vector of the control volume voxels is X λ , then formula (2) is expressed as A λ X λ =B λ .

5. The laser-induced phosphorescence flow field three-dimensional temperature field measurement method according to claim 4 is characterized in that, based on the characteristics of the tangent circle imaging area of ​​the light field camera, the invalid pixels of the light field image and their corresponding weight coefficients are eliminated to simplify the corresponding relationship equation group established in the previous step, and the judgment criterion of invalid pixels is: The pixels outside the tangent circle of the light field image are invalid pixels; by eliminating the light field image pixel column vector B λ Invalid pixels and weight matrix A λ The invalid weight coefficients corresponding to the invalid pixels in A are realized λ X λ =B λ The simplification of the equation group, the light field tomography reconstruction equation after preprocessing is transformed into A λ 2X λ =B λ 2.

6. The method for measuring three-dimensional temperature field of laser-induced phosphorescence flow field according to claim 5, characterized in that: The weight coefficient matrix is ​​calculated by ray tracing method. The calculation formula of ray tracing method is: Among them, (x, y) is the coordinate of the source of the light, (x', y') is the coordinate of the pixel where the light enters, is the solid angle of the outgoing ray, is the solid angle of the incident light; f is the focal length of the microlens unit, F is the focal length of the relay lens 1, S1 is the object distance, S2 is the image distance, (S x ,S y ) is the displacement of the center of the microlens array away from the optical axis of the imaging system.

7. The method for measuring three-dimensional temperature field of laser-induced phosphorescence flow field according to claim 6, characterized in that: The initial reconstruction is performed by the first SART algorithm, including: The preliminary reconstruction result X is obtained by SART algorithm hierarchical reconstruction λ 1. The iterative calculation formula is: Where t is the number of iterations, μ is the relaxation factor, p is the voxel number, b is the pixel gray value, and o is the pixel number. is the ratio of the grayscale contribution of the p voxel to the o pixel to the λ spectral intensity value of the voxel, represents the sum of the forward projections of the λ spectral intensity values ​​of all non-zero voxels within the line of sight of pixel o, represents the sum of the weight coefficients of non-zero voxels within the line of sight of pixel o to pixel o, Represents the difference between the grayscale value of the pixel o and the sum of the forward projections of the spectral intensity values ​​of all non-zero voxels within the pixel's line of sight.

8. The method for measuring three-dimensional temperature field of laser-induced phosphorescence flow field according to claim 7, characterized in that: X λ Any three-dimensional connected body X in 3 λ The calculation method of the voxel where the center of gravity of 3(i:i+Δi,j:j+Δj,k:k+Δk) is: Slice the three-dimensional connected body along the plane perpendicular to the Z axis to obtain X λ 3(i:i+Δi,j:j+Δj,n): Among them, k≤n≤k+Δk; Identify X using a 2D Gaussian fit function λ The centroid of 3(i:i+Δi,j:j+Δj,n): Among them, i≤x≤i+Δi, j≤y≤j+Δj, E λ is the amplitude of the Gaussian distribution of the λ spectrum intensity of the voxel, (x0, y0, n) is the voxel coordinates of the centroid of the two-dimensional slice, σ x , σ y are the half-peak width of the Gaussian distribution of the λ spectrum intensity of the two-dimensional slice voxel, respectively; A two-dimensional Gaussian fitting function is used for each slice of the three-dimensional connected body to obtain the centroid of all slices; a one-dimensional Gaussian fitting is performed on the λ spectrum intensity of all voxels where the centroid is located to determine the center of gravity of the three-dimensional connected body: Among them, (x0, y0, z0) is a three-dimensional connected body X λ 3(i:i+Δi,j:j+Δj,k:k+Δk) voxel coordinates of the center of gravity, σ z is the half-peak width of the Gaussian distribution of the λ spectral intensity of the voxel.

9. The method for measuring three-dimensional temperature field of laser-induced phosphorescence flow field according to claim 7, characterized in that: According to the obtained three-dimensional spectral intensity distribution, the three-dimensional temperature distribution is obtained by using the intensity ratio method, including: The three-dimensional spectral intensity X of the λ1 and λ2 phosphorescence spectra λ1 7 and X λ2 7. Take the ratio of each non-zero voxel and calculate the spectral intensity ratio of the three-dimensional phosphorescent particle field: Where (i, j, k) is the three-dimensional coordinate of the voxel, R (i,j,k) is the phosphorescence spectral intensity ratio of voxel (i, j, k); The relationship between the phosphorescence spectrum intensity ratio and temperature is obtained through temperature calibration experiment: T (i,j,k) =f(R (i,j,k) ) (10) Where, T (i,j,k) represents the temperature of the (i, j, k) voxel; Based on formula (10), the spectral intensity ratio result R of the three-dimensional phosphorescent particle field is (i,j,k) Solve and calculate the instantaneous three-dimensional temperature field of the flow field.

10. A laser-induced phosphorescence flow field three-dimensional temperature field measurement system, comprising a pulsed laser, a laser beam expander, a spectral light field image acquisition system and an image processing system, characterized in that: The spectral light field image acquisition system is the multi-spectral light field imaging system described in claim 1; the image processing system calculates the instantaneous three-dimensional temperature field of the flow field according to the steps of the laser-induced phosphorescence flow field three-dimensional temperature field measurement method described in any one of claims 2-9 based on the image acquired by the spectral light field image acquisition system.

Images