1. a compressed sensing-based acoustic CT temperature field reconstruction method is characterized by comprising the following steps: the method comprises the following steps:

arranging ns sound wave emitting/receiving devices around a tested layer surface with a rectangular or square cross section, wherein ns is more than or equal to 8, forming M sound wave propagation paths effectively penetrating through a tested area, uniformly dividing the tested layer surface into N grids, and 2M is more than or equal to N and more than or equal to 1.5M; the position of the transmitter/receiver is to make the sound wave path formed between the transmitter/receiver pass through the tested layer surface as uniformly as possible and avoid the superposition of the sound wave path and the grid boundary;

(II) assuming that the propagation path of the sound wave between the two sound wave transmitting/receiving devices is a straight line segment connecting the two transmitting/receiving devices, calculating an mth propagation path, wherein the mth propagation path belongs to {1, 2.. multidot.M }, and the length a of the mth grid_mnN belongs to {1,2,. and N }, and a compressed sensing measurement matrix A of the M multiplied by N dimension acoustic CT system is formed by using the formula (1);

(III) calculating the N × N dimensional DFT basis matrix Ψ by using the formula (2)_DFTRow p and column q elements Ψ_DFT(p, q) calculating the N × N-dimensional DCT basis matrix Ψ using the formula (3)_DCTRow p and column q elements Ψ_DCT(p, q), wherein j in the formula (2) is an imaginary unit, p belongs to {1, 2.., N }, and q belongs to {1, 2.., N }; let U ═ A Ψ_DFTCalculating the measurement matrix A and the DFT base matrix psi by using the formula (4)_DFTCoherency ofIs mu (A Ψ)_DFT) Then let U ═ A Ψ_DCTCalculating the measurement matrix A and the DCT basis matrix psi by using the formula (4)_DCTThe coherence of (A) is expressed as mu (A Ψ)_DCT) (ii) a If μ (A Ψ)_DCT)＜μ(AΨ_DFT) Let the sparse basis matrix Ψ ═ Ψ_DCTIf μ (A Ψ)_DCT)≥μ(AΨ_DFT) Let the sparse basis Ψ ═ Ψ_DFT(ii) a U in formula (4)_pAnd u_qRespectively representing the p-th column and the q-th column of the M multiplied by N dimensional U matrix, wherein p belongs to {1, 2.,. N }, q belongs to {1, 2.,. N }, the superscript T is used for solving transposition, | | | represents solving absolute value, and | | represents zero₂2, calculating a norm;

sequentially turning on and off the transmitters in the sound wave transmitting/receiving devices in a measuring period to transmit sound waves in turn, and ensuring that only one transmitter transmits sound waves at most each time; when any transmitter transmits sound waves, all receivers receive the sound waves and convert the sound waves into electric signals; the electric signals are synchronously collected and enter a computer after passing through a signal conditioner and a data acquisition card; by adopting a time delay estimation algorithm, the propagation time of the sound wave along M effective propagation paths can be measured by the data, and the data are combined into an M-dimensional sound wave propagation time vector t;

and (V) reconstructing sparse representation theta of the acoustic slow vector f in the psi domain by using an improved OMP algorithm according to the acoustic propagation time vector t, the measurement matrix A and the sparse basis matrix psiMeans of weightThe obtained sparse representation estimated value is obtained, and the N-dimensional sound slow vector f is ═ f₁,f₂,…,f_N]^TWherein f is₁,f₂,…,f_NRespectively representing sound slow values at the central points of the N grids, wherein the sound slow is the reciprocal of the sound velocity;

(VI) estimation by sparse representation using equation (7)Computing voiced slow vector estimates

(VII) estimation from slowness of sound vector using equation (8)Calculating a temperature vector estimateAnd will vectorEach element in the CT-CT data acquisition system is endowed with a geometric center of each grid, and the temperature distribution of the whole measured layer is obtained by an interpolation operation method, so that the accurate and quick acoustic CT temperature field reconstruction based on compressed sensing of the measured layer is realized;

2. the compressed sensing-based acoustic CT temperature field reconstruction method according to claim 1, wherein: the specific operation of reconstructing theta by the improved OMP algorithm in the step (five) is as follows:

firstly, initializing: let residual vector r₀T, index vector Λ₀＝[]Supporting matrix V₀＝[]The iteration counter k is 1, the parameter P is 3, and the parameter α is 0.6, respectively]Representing an empty matrix;

selecting: computing the q-th column of the matrix U and the residual vector r_k-1Inner product g between_qQ ═ 1,2,. ·, N; definition ofSuppose that | g is satisfied in the matrix U_q|＞α|g_maxThe number of column vectors of is Q, if Q<P, then Q columns are all selected and added to V_k-1In (1), forming a new support matrix V_kOtherwise, only the most relevant (inner product absolute value maximum) P column is added to V_k-1In (1), forming a new support matrix V_k(ii) a At the same time, the column index of the selected column is added to Λ_k-1In (2), a new index vector Λ is formed_k(ii) a Finally, setting the column of the selected support matrix in the U as a zero vector;

updating residual errors: first, the formula (5) is used to calculate theta_kThen updating the residual error by using the formula (6);

r_k＝t-V_kθ_k (6)

judgmentIf yes, let α be P be 1, where σ is an estimate of the standard deviation of the noise signal in the acoustic wave propagation time vector, and σ can be taken as the sampling interval for the acoustic wave transceiver to output the electrical signal;

fifthly, whenWhen satisfied, or support matrix V_kNumber of middle column vectorsStopping iteration and jumping to the step (c) when the value reaches or exceeds M; if the two conditions are not met, making k equal to k +1, and returning to the step II to perform the next iteration;

sixthly, constructing an N-dimensional zero vectorThen order(Vector)Namely, the acoustic slow vector f reconstructed by the improved OMP algorithm is estimated in a sparse representation of the psi domain.

3. The compressed sensing-based acoustic CT temperature field reconstruction method according to claim 1, wherein: and (4) adopting a time delay estimation algorithm in the step (four) which is a cross-correlation time delay estimation method combined with wavelet noise suppression.

Background

The temperature of the medium is indirectly obtained by measuring the temperature by an acoustic method according to the propagation speed of sound waves in the medium. The speed of sound c (m/s) in a gaseous medium and the absolute temperature T (K) of the gaseous mediumThe relationship between the gas acoustic constant z determined by the gas composition can be expressed asFor air z is 20.05.

The acoustic CT temperature field measurement technology has the advantages of non-contact non-interference to the measured temperature field, wide measurement range, strong environment adaptability, convenient maintenance and the like. In order to obtain the temperature distribution of the measured layer by adopting an acoustic CT method, a plurality of sound wave transmitting/receiving devices are required to be arranged around the measured layer to form M sound wave propagation paths which effectively pass through a measured area, and the reconstructed layer is uniformly divided into N grids. The acoustic CT temperature field measurement technology measures the propagation time of sound waves on the sound wave propagation paths, and then a proper temperature field reconstruction algorithm is adopted to reconstruct a temperature value for each grid according to M sound wave propagation time measurement data and endow the temperature value to the geometric center of the grid. Then, a more detailed description of the temperature distribution on the measured layer is obtained by an interpolation method. The performance of an acoustic CT temperature field measurement system depends to a large extent on the performance of the temperature field reconstruction algorithm. Classical acoustic CT temperature field reconstruction algorithms, such as Least Squares (LSA), simultaneous iterative methods (SIRT), etc., all require N < M. Since the number M of available acoustic travel time measurements is usually small, only a coarse meshing of the measured layer is possible. The thicker the grid, the weaker the description of the temperature field, and the more information of the temperature field at the edge of the measurement area is lost, i.e., it is more difficult to provide complete and accurate temperature distribution information.

Disclosure of Invention

The invention aims to provide an acoustic CT temperature field reconstruction method based on compressed sensing, which can realize the rapid and accurate reconstruction of a two-dimensional temperature field according to a small amount of acoustic wave propagation time measurement data.

In order to achieve the purpose, the invention provides the following technical scheme: an acoustic CT temperature field reconstruction method based on compressed sensing comprises the following steps

Arranging ns sound wave emitting/receiving devices around a tested layer surface with a rectangular or square cross section, wherein ns is more than or equal to 8, forming M sound wave propagation paths effectively penetrating through a tested area, uniformly dividing the tested layer surface into N grids, and 2M is more than or equal to N and more than or equal to 1.5M; the position of the transmitter/receiver is to make the sound wave path formed between the transmitter/receiver pass through the tested layer surface as uniformly as possible and avoid the superposition of the sound wave path and the grid boundary;

(II) assuming that the propagation path of the sound wave between the two sound wave transmitting/receiving devices is a straight line segment connecting the two transmitting/receiving devices, calculating an mth propagation path, wherein the mth propagation path belongs to {1, 2.. multidot.M }, and the length a of the mth grid_mnN belongs to {1,2,. and N }, and a compressed sensing measurement matrix A of the M multiplied by N dimension acoustic CT system is formed by using the formula (1);

(III) calculating the N × N dimensional DFT basis matrix Ψ by using the formula (2)_DFTRow p and column q elements Ψ_DFT(p, q) calculating the N × N-dimensional DCT basis matrix Ψ using the formula (3)_DCTRow p and column q elements Ψ_DCT(p, q), wherein j in the formula (2) is an imaginary unit, p belongs to {1, 2.., N }, and q belongs to {1, 2.., N }; let U ═ A Ψ_DFTCalculating the measurement matrix A and the DFT base matrix psi by using the formula (4)_DFTThe coherence of (A) is expressed as mu (A Ψ)_DFT) Then let U ═ A Ψ_DCTCalculating the measurement matrix A and the DCT basis matrix psi by using the formula (4)_DCTThe coherence of (A) is expressed as mu (A Ψ)_DCT) (ii) a If μ (A Ψ)_DCT)＜μ(AΨ_DFT) Let the sparse basis matrix Ψ ═ Ψ_DCTIf μ (A Ψ)_DCT)≥μ(AΨ_DFT) Let the sparse basis Ψ ═ Ψ_DFT(ii) a U in formula (4)_pAnd u_qRespectively representing the p-th column and the q-th column of the M multiplied by N dimensional U matrix, wherein p belongs to {1, 2.,. N }, q belongs to {1, 2.,. N }, the superscript T is used for solving transposition, | | | represents solving absolute value, and | | represents zero₂2, calculating a norm;

and (IV) sequentially turning on and off the transmitters in the sound wave transmitting/receiving devices in a measuring period to transmit sound waves in turn, and ensuring that only one transmitter transmits sound waves at most each time. When any transmitter transmits sound waves, all receivers receive the sound waves and convert the sound waves into electric signals. The electric signals are synchronously collected and enter a computer after passing through a signal conditioner and a data acquisition card. The data can be obtained by adopting a time delay estimation algorithm; measuring the propagation time of sound waves along the M effective propagation paths, and combining the propagation time of the sound waves into an M-dimensional sound wave propagation time vector t;

and (V) reconstructing sparse representation theta of the acoustic slow vector f in the psi domain by using an improved OMP algorithm according to the acoustic propagation time vector t, the measurement matrix A and the sparse basis matrix psiRepresenting reconstructed sparse representation estimates with an N-dimensional acoustic slow vector f ═ f₁,f₂,…,f_N]^TWherein f is₁,f₂,…,f_NRespectively representing sound slow values at the central points of the N grids, wherein the sound slow is the reciprocal of the sound velocity; the specific operation of reconstructing θ using the modified OMP algorithm is briefly as follows:

firstly, initializing: let residual vector r₀T, index vector Λ₀＝[]Supporting matrix V₀＝[]The iteration counter k is 1, the parameter P is 3, and the parameter α is 0.6, respectively]Representing an empty matrix;

selecting: computing the q-th column of the matrix U and the residual vector r_k-1Inner product g between_qQ ═ 1,2,. ·, N; definition ofSuppose that | g is satisfied in the matrix U_q|＞α|g_maxThe number of column vectors of is Q, e.g.Fruit Q<P, then Q columns are all selected and added to V_k-1In (1), forming a new support matrix V_kOtherwise, only the most relevant (inner product absolute value maximum) P column is added to V_k-1In (1), forming a new support matrix V_k(ii) a At the same time, the column index of the selected column is added to Λ_k-1In (2), a new index vector Λ is formed_k(ii) a Finally, setting the column of the selected support matrix in the U as a zero vector;

updating residual errors: first, the formula (5) is used to calculate theta_kThen updating the residual error by using the formula (6);

r_k＝t-V_kθ_k (6)

judgmentIf yes, let α be P be 1, where σ is an estimate of the standard deviation of the noise signal in the acoustic wave propagation time vector, and σ can be taken as the sampling interval for the acoustic wave transceiver to output the electrical signal;

fifthly, whenWhen satisfied, or support matrix V_kStopping iteration and skipping to the step (c) when the number of the middle column vectors reaches or exceeds M; if the two conditions are not met, making k equal to k +1, and returning to the step II to perform the next iteration;

sixthly, constructing an N-dimensional zero vectorThen order(Vector)I.e. reconstructed by using modified OMP algorithmSparse representation estimation of the acoustic slow vector f in a psi domain;

(VI) estimation by sparse representation using equation (7)Computing voiced slow vector estimates

(VII) estimation from slowness of sound vector using equation (8)Calculating a temperature vector estimateAnd will vectorEach element in the CT-CT data acquisition system is endowed with a geometric center of each grid, and the temperature distribution of the whole measured layer is obtained by an interpolation operation method, so that the accurate and quick acoustic CT temperature field reconstruction based on compressed sensing of the measured layer is realized;

compared with the prior art, the invention has the beneficial effects that: the performance of an acoustic CT temperature field measurement system depends to a large extent on the performance of the temperature field reconstruction algorithm. Classical acoustic CT temperature field reconstruction algorithms, such as Least Squares (LSA), simultaneous iterative methods (SIRT), algebraic reconstruction methods (ART), etc., require the number N of meshes divided by the layer to be measured to be greater than the number M of available measurements of the propagation time of the acoustic wave. M is usually few, so that only a coarse grid division can be adopted for a measured layer, so that more temperature field information at the edge of a measurement area is lost, and complete and accurate temperature distribution information is difficult to provide. The invention discloses an acoustic CT temperature field reconstruction method based on compressed sensing, which allows N to be greater than M and can realize more comprehensive and more accurate measurement of a two-dimensional temperature field. The reconstruction method provided by the invention can obviously improve the reconstruction precision of the temperature field, the reconstructed temperature distribution is closer to the actual distribution, and the reconstruction method has the advantages of easy programming realization, high speed and the like.

Drawings

FIG. 1 is a grid division diagram of the acoustic CT temperature field reconstruction method based on compressed sensing according to the present invention.

Detailed Description

The technical solutions in the embodiments of the present invention will be clearly and completely described below with reference to the drawings in the embodiments of the present invention, and it is obvious that the described embodiments are only a part of the embodiments of the present invention, and not all of the embodiments. All other embodiments, which can be derived by a person skilled in the art from the embodiments given herein without making any creative effort, shall fall within the protection scope of the present invention.

The invention provides a technical scheme that: an acoustic CT temperature field reconstruction method based on compressed sensing comprises the following steps

Arranging ns sound wave emitting/receiving devices around a tested layer surface with a rectangular or square cross section, wherein ns is more than or equal to 8, forming M sound wave propagation paths effectively penetrating through a tested area, uniformly dividing the tested layer surface into N grids, and 2M is more than or equal to N and more than or equal to 1.5M; the transmitter/receivers are positioned so that the acoustic paths formed between the transmitter/receivers pass through the layer under test as uniformly as possible and so that the acoustic paths do not coincide with the boundaries of the grid.

For example: as shown in fig. 1, 12 sound wave transmitters/receivers are arranged around a square measured plane to form 54 effective sound wave paths, and the measured plane is divided into 81 grids, that is, ns-12, M-54, N-81, and fig. 1 shows the effective sound wave paths without considering the sound ray bending phenomenon, that is, straight paths.

(II) assuming that the propagation path of the sound wave between the two sound wave transmitting/receiving devices is a straight line segment connecting the two transmitting/receiving devices, calculatingThe mth propagation path, M ∈ {1, 2.., M }, the length a within the nth trellis_mnN belongs to {1,2,. and N }, and a compressed sensing measurement matrix A of the M multiplied by N dimension acoustic CT system is formed by using the formula (1);

for example: for 54 sound wave propagation paths and 81 grid divisions shown in fig. 1, a 54 × 81-dimensional compressed sensing measurement matrix a of the acoustic CT system can be obtained in this way.

(III) calculating the N × N dimensional DFT basis matrix Ψ by using the formula (2)_DFTRow p and column q elements Ψ_DFT(p, q) calculating the N × N-dimensional DCT basis matrix Ψ using the formula (3)_DCTRow p and column q elements Ψ_DCT(p, q), wherein j in the formula (2) is an imaginary unit, p belongs to {1, 2.., N }, and q belongs to {1, 2.., N }; let U ═ A Ψ_DFTCalculating the measurement matrix A and the DFT base matrix psi by using the formula (4)_DFTThe coherence of (A) is expressed as mu (A Ψ)_DFT) Then let U ═ A Ψ_DCTCalculating the measurement matrix A and the DCT basis matrix psi by using the formula (4)_DCTThe coherence of (A) is expressed as mu (A Ψ)_DCT) (ii) a If μ (A Ψ)_DCT)＜μ(AΨ_DFT) Let the sparse basis matrix Ψ ═ Ψ_DCTIf μ (A Ψ)_DCT)≥μ(AΨ_DFT) Let the sparse basis Ψ ═ Ψ_DFT(ii) a U in formula (4)_pAnd u_qRespectively representing the p-th column and the q-th column of the M multiplied by N dimensional U matrix, wherein p belongs to {1, 2.,. N }, q belongs to {1, 2.,. N }, the superscript T is used for solving transposition, | | | represents solving absolute value, and | | represents zero₂2, calculating a norm;

for example, when the layer under test is divided into 81 grids as shown in FIG. 1, a 81 × 81-dimensional DFT basis matrix Ψ can be obtained in this way_DFTAnd a 81 × 81-dimensional DCT basis matrix Ψ_DCT。

And (IV) sequentially turning on and off the transmitters in the sound wave transmitting/receiving devices in a measuring period to transmit sound waves in turn, and ensuring that only one transmitter transmits sound waves at most each time. When any transmitter transmits sound waves, all receivers receive the sound waves and convert the sound waves into electric signals. The electric signals are synchronously collected and enter a computer after passing through a signal conditioner and a data acquisition card. These data can be derived using some delay estimation algorithm, such as cross-correlation delay estimation combined with wavelet noise suppression. And measuring the propagation time of the sound wave on the M effective paths, and combining the propagation time of the sound wave on the M effective paths into an M-dimensional sound wave propagation time vector t.

For example: the number ns of the sound wave emitting/receiving devices is 12, the number M of the effective sound wave paths is 54, and the data acquisition unit consists of 2 PCI 6123 data acquisition cards. In a measuring period, firstly, the No. 1 sound wave emitting/receiving device sends out a sound wave signal, and the data acquisition device synchronously acquires output signals of all the sound wave receiving devices amplified by the signal conditioner at the sampling frequency of 500KHz and feeds the output signals into the computer. According to the sound wave data, the propagation time of the sound waves among sound wave emitting/receiving devices of No. 1-2, 1-3, 1-4, 1-5, 1-6, 1-7, 1-8, 1-9, 1-10, 1-11 and 1-12 is measured by a time delay estimation method based on cross correlation, then the sound wave emitting/receiving devices of No. 2, No. 3, No. 9, No. 10, No. 1-11 and No. 12 are sequentially enabled to emit a sound wave signal, the sound wave propagation time on the corresponding path is measured, and finally a 54-dimensional sound wave propagation time vector t can be formed.

And (V) reconstructing sparse representation theta of the acoustic slow vector f in the psi domain by using an improved OMP algorithm according to the acoustic propagation time vector t, the measurement matrix A and the sparse basis matrix psiRepresenting reconstructed sparse representation estimates with an N-dimensional acoustic slow vector f ═ f₁,f₂,…,f_N]^TWherein f is₁,f₂,…,f_NRespectively representing sound slow values at the central points of the N grids, wherein the sound slow is the reciprocal of sound velocity, and the specific steps of reconstructing theta by using the improved OMP algorithm are briefly described as follows:

firstly, initializing: let residual vector r₀T, index vector Λ₀＝[]Supporting matrix V₀＝[]The iteration counter k is 1, the parameter P is 3, and the parameter α is 0.6, respectively]Representing an empty matrix;

selecting: computing the q-th column of the matrix U and the residual vector r_k-1Inner product g between_qQ ═ 1,2,. ·, N; definition ofSuppose that | g is satisfied in the matrix U_q|＞α|g_maxThe number of column vectors of is Q, if Q<P, then Q columns are all selected and added to V_k-1In (1), forming a new support matrix V_kOtherwise, only the most relevant (inner product absolute value maximum) P column is added to V_k-1In (1), forming a new support matrix V_k(ii) a At the same time, the column index of the selected column is added to Λ_k-1In (2), a new index vector Λ is formed_k(ii) a Finally, setting the column of the selected support matrix in the U as a zero vector;

updating residual errors: first, the formula (5) is used to calculate theta_kThen updating the residual error by using the formula (6);

r_k＝t-V_kθ_k (6)

judgmentWhether or not, if so, thenLet α be P be 1, where σ is an estimate of the standard deviation of the noise signal in the acoustic propagation time vector, and σ can be taken as the sampling interval for the acoustic transceiver output electrical signal;

fifthly, whenWhen satisfied, or support matrix V_kStopping iteration and skipping to the step (c) when the number of the middle column vectors reaches or exceeds M; if the two conditions are not met, making k equal to k +1, and returning to the step II to perform the next iteration;

sixthly, constructing an N-dimensional zero vectorThen order(Vector)The sparse representation estimation of the acoustic slow vector f in the psi domain is reconstructed by using an improved OMP algorithm;

(VI) estimation by sparse representation using equation (7)Computing voiced slow vector estimates

(VII) estimation from slowness of sound vector using equation (8)Calculating a temperature vector estimateAnd will vectorEach element in the CT-CT data acquisition system is endowed with a geometric center of each grid, and the temperature distribution of the whole measured layer is obtained by an interpolation operation method, so that the accurate and quick acoustic CT temperature field reconstruction based on compressed sensing of the measured layer is realized;