CN107895063A - One kind compression EO-1 hyperion mask optimization method - Google Patents

Abstract

The invention discloses one kind to compress EO-1 hyperion mask optimization method.Step of the present invention is as follows：Step (1) utilizes geometric optics by the position of mask in systems, obtains initial projections matrix P₀, P₀It is the matrix of a p × n in real number field, i.e.,Initial cycle iteration count variable q=0 is set；Step (2) is by matrix P_qD enters ranks standardization, obtains effectively dictionaryStep (3) constructs Gram matrixesStep (4) shrinks Gram matrixes G_q, that is, update Gram matrixes：Step (5) is decomposed the Gram matrixes after contraction using SVDOrder be reduced to p；Step (6) solves Gram matrixesSquare root S_qSo thatStep (7) finds P_q+1So that errorMinimum, q=q+1.The property consistent with each other of mask after present invention optimization is significantly lower than the property consistent with each other of random mask, and the signal to noise ratio of reconstruction image is also significantly improved.

Description

A Method for Optimizing Compressed Hyperspectral Masks

technical field

The invention relates to hyperspectral imaging technology and the optimization of compressed hyperspectral masks, in particular to providing a compression hyperspectral mask optimization method when applying compressed hyperspectral imaging.

Background technique

Hyperspectral imaging is a major research hotspot today. How to achieve fast hyperspectral acquisition is always an urgent problem to be solved. In this paper, the method of placing a mask in front of the sensor to realize compressed hyperspectral imaging is further optimized and expanded. Combining the relevant knowledge of optimizing the projection matrix and the relevant content learned in the matrix class, including matrix differentiation, projection analysis, etc., the optimization of the mask for this practical problem is realized.

Hyperspectral imaging technology is a combination of traditional two-dimensional imaging technology and spectral technology. In the spectral range from ultraviolet to near-infrared, using imaging spectrometer, dozens or hundreds of spectral bands within the spectral coverage range can continuously detect the target object. Imaging, while obtaining the spatial feature imaging of the object, also obtains the spectral information of the measured object, and has important applications in a wide range of fields such as remote sensing, monitoring, and spectroscopy. The main task of hyperspectral imaging is to effectively acquire a three-dimensional data cube that varies from two-dimensional space and one-dimensional spectral space. The most commonly used method is to record one or several data points at a time by means of mechanical scanning or sequential scanning. The snapshot imaging method acquires complete 3D data through a single imaging, which has obvious advantages over the scanning method when capturing dynamic scenes or aerial photography. Snapshot imaging is achieved by multiplexing high-dimensional signals to a two-dimensional CCD sensor, thus sacrificing image resolution. Examples of snapshot imaging methods include four-dimensional imaging spectrometer (4DIS) and computed tomography imaging spectrometer (CTIS), among others. A recent Coded Aperture Snapshot Spectrometer (CASSI) encodes optical signals by applying computational compression reconstruction. This computational imaging approach is able to overcome the trade-off in the obtained spatial and spectral imaging resolutions. The CASSI system can be improved by acquiring multiple images corresponding to an encoded mask that is constantly flipped on the piezoelectric stage. A more flexible option is to utilize a digital micromirror device (DMD) continuous-shot spectral imaging system (DMD-SSI). The trend of optically encoding recorded data and reconstructing data through compressive calculations has been very obvious, however, all methods of compressive spectral imaging encode color based on the assumption of spatial uniformity, which inherently limits compressive reconstruction algorithms with compressive sparsity constraints. the quality of. Therefore, in our previous work, we proposed to place a mask code in front of the sensor to achieve spatially inconsistent modulation, and what needs to be done after obtaining spatially inconsistent modulation is to perform sparse reconstruction. How to make the collected information most conducive to sparse reconstruction is how to optimize the projection matrix. For compressed hyperspectral imaging, however, this problem has not been studied in detail. Therefore, on the basis of compressed hyperspectral imaging, we further use the relevant knowledge of matrix analysis, including matrix differentiation, projection analysis, and matrix decomposition knowledge, to propose a mask optimization algorithm to achieve mask optimization. to achieve the best reconstruction results. In order to verify the feasibility of the proposed optimization algorithm, we conducted some simulation experiments to prove and demonstrate the optimization effect of the algorithm.

Contents of the invention

The technical problem to be solved by the present invention is: in order to make the collected information most conducive to sparse reconstruction, it is necessary to determine an evaluation standard for judging the quality of the projection matrix in compressed hyperspectral imaging, and based on this, through the proposed mask Optimize the algorithm to achieve the purpose of improving the quality of sparse reconstruction.

In the present invention, the mask is a 100×1 vector, and a 3100×100 projection matrix P is constructed from it.

In order to achieve the above purpose, the first application object is the spectral camera. Compared with the ordinary camera, the spectral camera places a grating on the image plane to generate a multi-spectral camera, which is equivalent to mapping the spectral information of the scene to the image plane. Where light rays emerge from different angles, that is, the angular dimension of the light field. First of all, we need to determine a standard to measure the irrelevance of the matrix. Here we choose the most commonly used μ _t to measure.

The present invention solves the technical scheme that its technical problem adopts The goal of the present invention: minimize μ _t (PD)

The input parameters are as follows:

t: Correlation degree threshold

n:overcomplete dictionary

p: number of measurements

γ: Correlation degree threshold

Iter: number of iterations

The goal of the invention: to minimize μ _t (PD)

The input parameters are as follows:

t: Correlation degree threshold

n:overcomplete dictionary

p: number of measurements

γ: Correlation degree threshold

P: Projection matrix, consisting of 100×1 numbers between 0 and 1.

D: An over-complete dictionary obtained through dictionary learning, with a size of 3100×387.

Iter: number of iterations

The specific implementation steps are as follows:

Step (1) Use geometric optics to obtain the initial projection matrix P ₀ through the position of the mask in the system. P ₀ is a p×n matrix in the real number field, namely Set initial loop iteration count variable q = 0;

The derivation of the initial projection matrix _P0 is as follows:

First construct a 100×100 zero matrix;

Then the 100 numbers of the mask are used as the diagonal elements of the zero matrix, and the storage form of the elements above the diagonal is a symbol variable form. This 100×100 matrix is copied 31 times to form a 3100×100 matrix, namely the projection matrix P ₀ . Among them, the elements on the diagonal in the 31 100×100 matrices are the same, but the order of the elements is different, and there will be a translation in position. This 3100×100 matrix is the projection matrix P ₀ constructed by us.

Step (2) Perform column standardization on the matrix P _q D to obtain an effective dictionary

The standardization is to normalize each column of the matrix P _q D into a unit vector, so that the modulus value of each column vector is one, where P _q represents the projection matrix P obtained by the qth cycle;

Step (3) Construct the Gram matrix

i.e. effectively a dictionary The transpose of is multiplied by itself to construct the Gram matrix.

Step (4) Shrink the Gram matrix G _q , that is, update the Gram matrix using the formula (1):

Among them: g _ij represents the element of row i and column j in the Gram matrix; Indicates the new value transformed by formula (1);

Step (5) Use SVD decomposition to decompose the contracted Gram matrix The rank of is reduced to p, specifically:

Depend on It can be seen that since the diagonal elements of the Σ matrix obtained by the singular value algorithm are all non-negative, the column vectors of the V and U matrices may have different signs. Since the U and V vectors are unit orthogonal, the matrix in actual calculation is also There are only a very small number of semi-negative fixed elements, so the absolute value is directly used for calculation, that is, V·Σ· ^UT is used to obtain the square root of the Gram matrix.

Step (6) Solve the Gram matrix the square root of S _q such that in

Step (7) Find the projection matrix P _q+ 1 obtained in the q+1th cycle, so that the error Minimum, q=q+1.

Using the gradient descent method to find problem, the detailed calculation process is as follows:

7-1. The definition is computed:

i.e. minimize

tr(S _q ^T S _q -S _q ^T PD-D ^T P ^T S _q +D ^T P ^T PD) (3)

Using the properties of the matrix to find the trace, the original formula can be changed to

f(x)=tr(S _q ^T S _q -2D ^T P ^T S _q +D ^T P ^T PD) (4)

Further find the gradient to get

Using the step-differential formula, according to the properties of the projection matrix generated by the previous mask, the gradient value for each unknown mask variable can be obtained from the formula (5).

f(x) is a function about the projection matrix P, but the projection matrix P is a 3100×100 two-dimensional matrix constructed from a 100×1 mask, and there are only 100 variables in it, and the other elements are all 0. Elements that are 0 must be guaranteed to remain unchanged. Therefore, it is not possible to directly derive the projection matrix P. For this reason, in MATLAB, 100 of the variables are defined as symbolic variables, and the function f(x) is used to derive the 100 symbolic functions in MATLAB, and this method is used to bypass the direct derivative of the projection matrix P bring about problems.

Step (8) repeat step (2) to step (7) Iter times, and output the optimized projection matrix P _iter . Since the operation mode is fixed, the threshold t is a fixed value.

Beneficial effects of the present invention:

The method of the present invention only needs to optimize the mask to reduce the mutual consistency of the projection matrix, so as to realize a more sparse projection matrix and facilitate reconstruction. Simulation results show that compared with the commonly used random mask, the mutual consistency of the optimized mask is significantly lower than that of the random mask, and the signal-to-noise ratio of the reconstructed image is also significantly improved.

Description of drawings

Fig. 1 is a flowchart of the present invention;

Fig. 2 is the mutual consistency convergence curve of the projection matrix P of the present invention.

Detailed ways

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

As shown in Figure 1-2, a compressed hyperspectral mask optimization algorithm proposed by the present invention uses an algorithm to reduce the mutual consistency of the projection matrix to obtain a more sparse projection matrix. The method of the present invention comprises the following steps:

Step (1) Use geometric optics to derive the projection matrix through the position of the mask in the system Set the loop iteration count variable q=0, and repeat steps 2 to 7 for a total of Iter times. Since the operation mode is fixed, the threshold t is a fixed value.

Step (2) Perform column standardization on the matrix P _q D to obtain an effective dictionary The normalization method here is to normalize each column of the matrix P _q D into a unit vector, so that the modulus value of each column vector is one.

Step (3) Calculate the Gram matrix Here we use matrix-related knowledge to make the effective dictionary The transpose of is multiplied by itself to construct the Gram matrix.

Step (4) shrink the matrix, and use the following formula to update the Gram matrix:

The purpose of doing this here is to shrink the Gram matrix so that the matrix can converge more easily.

Step (5) Use SVD decomposition to decompose the contracted Gram matrix The rank of is reduced to p. Depend on It can be seen that since the diagonal elements of the Σ matrix obtained by the singular value algorithm are all non-negative, the column vectors of the V and U matrices may have different signs. Since the U and V vectors are unit orthogonal, the matrix in actual calculation is also There are only a very small number of semi-negative definite elements, so here we directly use the absolute value to calculate, that is, use V Σ ^UT to find the square root of the Gram matrix. The actual calculation results show that for the two norms caused by this problem The error in the sense is around 1.

Step (6) Solve the Gram matrix the square root of S _q such that in Here we set Doing so can guarantee the solution speed of the algorithm by sacrificing finite solution accuracy.

Step (7) finds P _q+1 such that the error Minimum, q=q+1. seeking We tried gradient descent separately for the problem. The detailed calculation process is as follows:

1. We know The definition is computed:

i.e. minimize

tr(S _q ^T S _q -S _q ^T PD-D ^T P ^T S _q +D ^T P ^T PD)

Using the properties of the matrix to find the trace, the original formula can be changed to

f(x)=tr(S _q ^T S _q -2D ^T P ^T S _q +D ^T P ^T PD)

Further find the gradient to get

Here, using the step-by-step differential formula, according to the properties of the projection matrix generated by the previous mask, the gradient value for each unknown mask variable can be obtained from the above matrix gradient.

2. Here f(x) is a function of the projection matrix P, but the projection matrix P is a 3100×100 two-dimensional matrix constructed from a 100×1 mask, and there are only 100 variables in it, and the other elements are 0, it must also ensure that the elements that are 0 remain unchanged. Therefore, it is not possible to directly derive the projection matrix P. For this reason, in MATLAB, we define 100 of the variables as symbolic variables, and let the function f(x) take derivatives of these 100 symbolic functions in MATLAB, and use this method to bypass the direct calculation of the projection matrix P lead problems.

3. Using the gradient descent method to obtain the process is as follows:

Step1: given the initial point;

Step2:

Step3: Straight line search, determine the step size t by backtracking straight line search method;

Step4: x=x+Δx;

Repeat step2-step4 until the stop criterion is met.

Step (8) outputs the optimized projection matrix P _iter .

Use the new 100×1 mask after the convergence of the gradient descent method to construct the projection matrix P as the new projection matrix P _iter , output it, and check whether its mutual consistency is converged.

Claims (2)

1. one kind compression EO-1 hyperion mask optimization method, it is characterised in that be implemented as follows：

Target:Minimize μ_t(PD)

The parameter of input is as follows:

t:Degree of correlation threshold value

n:Cross complete dictionary

p:Pendulous frequency

γ:Degree of correlation threshold value

P：Projection matrix, by the array between 100 × 10 to 1 into；

D：The excessively complete dictionary obtained by dictionary learning, size are 3100 × 387

Iter:Iterations

Wherein mask is the vector of one 100 × 1, and the projection matrix P of one 3100 × 100 is constructed by it；

It is as follows to implement step：

Step (1) utilizes geometric optics by the position of mask in systems, obtains initial projections matrix P₀, P₀It is in real number field One p × n matrix, i.e.,Initial cycle iteration count variable q=0 is set；

Step (2) is by matrix P_qD enters ranks standardization, obtains effectively dictionary

Described standardization is exactly by matrix P_qD each row normalization is unit vector so that the modulus value of its each column vector is One, wherein P_qRepresent to circulate obtained projection matrix P the q times；

Step (3) constructs Gram matrixes

That is effectively dictionaryTransposition and itself be multiplied construction Gram matrixes；

Step (4) shrinks Gram matrixes G_q, i.e., update Gram matrixes using formula (1)：

<mrow> <msub> <mover> <mi>g</mi> <mo>^</mo> </mover> <mrow> <mi>i</mi> <mi>j</mi> </mrow> </msub> <mo>=</mo> <mfenced open = "{" close = ""> <mtable> <mtr> <mtd> <mrow> <mi>&amp;gamma;</mi> <mo>&amp;CenterDot;</mo> <msub> <mi>g</mi> <mrow> <mi>i</mi> <mi>j</mi> </mrow> </msub> </mrow> </mtd> <mtd> <mrow> <mo>|</mo> <msub> <mi>g</mi> <mrow> <mi>i</mi> <mi>j</mi> </mrow> </msub> <mo>|</mo> <mo>&amp;GreaterEqual;</mo> <mi>t</mi> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <mi>&amp;gamma;</mi> <mi>t</mi> <mo>&amp;CenterDot;</mo> <mi>s</mi> <mi>i</mi> <mi>g</mi> <mi>n</mi> <mrow> <mo>(</mo> <msub> <mi>g</mi> <mrow> <mi>i</mi> <mi>j</mi> </mrow> </msub> <mo>)</mo> </mrow> </mrow> </mtd> <mtd> <mrow> <mi>t</mi> <mo>&gt;</mo> <mo>|</mo> <msub> <mi>g</mi> <mrow> <mi>i</mi> <mi>j</mi> </mrow> </msub> <mo>|</mo> <mo>&amp;GreaterEqual;</mo> <mi>&amp;gamma;</mi> <mi>t</mi> </mrow> </mtd> </mtr> <mtr> <mtd> <msub> <mi>g</mi> <mrow> <mi>i</mi> <mi>j</mi> </mrow> </msub> </mtd> <mtd> <mrow> <mi>&amp;gamma;</mi> <mi>t</mi> <mo>&gt;</mo> <mo>|</mo> <msub> <mi>g</mi> <mrow> <mi>i</mi> <mi>j</mi> </mrow> </msub> <mo>|</mo> </mrow> </mtd> </mtr> </mtable> </mfenced> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>1</mn> <mo>)</mo> </mrow> </mrow>

Wherein：g_ijRepresent the element that the i-th row jth arranges in Gram matrixes；Represent the new value after formula (1) conversion；

Step (5) is decomposed the Gram matrixes after contraction using SVDOrder be reduced to p, specifically：

ByUnderstand, because the Σ matrix diagonals members that unusual value-based algorithm obtains are non-negative, therefore V and U matrixes Column vector may differ sign, because U, V vector are that unit is orthogonal, it is actual calculate in matrix also only have only a few Negative semidefinite member, therefore directly calculated using absolute value, that is, utilize V Σ U^TTo ask for the square root of Gram matrixes；

Step (6) solves Gram matrixesSquare root S_qSo thatWherein

Step (7) finds the q+1 times and circulates obtained projection matrix P_q+1So that errorMinimum, q=q+1；

Asked for using gradient descent methodProblem, detailed calculating process are as follows：

Definition be calculate：

<mrow> <mo>|</mo> <mrow> <msub> <mi>S</mi> <mi>q</mi> </msub> <mo>-</mo> <mi>P</mi> <mi>D</mi> </mrow> <mo>|</mo> <msubsup> <mo>|</mo> <mi>F</mi> <mn>2</mn> </msubsup> <mo>=</mo> <mi>t</mi> <mi>r</mi> <mrow> <mo>(</mo> <msup> <mrow> <mo>(</mo> <mrow> <msub> <mi>S</mi> <mi>q</mi> </msub> <mo>-</mo> <mi>P</mi> <mi>D</mi> </mrow> <mo>)</mo> </mrow> <mi>T</mi> </msup> <mo>(</mo> <mrow> <msub> <mi>S</mi> <mi>q</mi> </msub> <mo>-</mo> <mi>P</mi> <mi>D</mi> </mrow> <mo>)</mo> <mo>)</mo> </mrow> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>2</mn> <mo>)</mo> </mrow> </mrow>

Minimize

<mrow> <mi>t</mi> <mi>r</mi> <mrow> <mo>(</mo> <msup> <msub> <mi>S</mi> <mi>q</mi> </msub> <mi>T</mi> </msup> <msub> <mi>S</mi> <mi>q</mi> </msub> <mo>-</mo> <msup> <msub> <mi>S</mi> <mi>q</mi> </msub> <mi>T</mi> </msup> <mi>P</mi> <mi>D</mi> <mo>)</mo> </mrow> <mo>-</mo> <msup> <mi>D</mi> <mi>T</mi> </msup> <msup> <mi>P</mi> <mi>T</mi> </msup> <msub> <mi>S</mi> <mi>q</mi> </msub> <mo>+</mo> <msup> <mi>D</mi> <mi>T</mi> </msup> <msup> <mi>P</mi> <mi>T</mi> </msup> <mi>P</mi> <mi>D</mi> <mo>)</mo> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>3</mn> <mo>)</mo> </mrow> </mrow>

Using Matrix Calculating trace property, former formula can be changed to

<mrow> <mi>f</mi> <mrow> <mo>(</mo> <mi>x</mi> <mo>)</mo> </mrow> <mo>=</mo> <mi>t</mi> <mi>r</mi> <mrow> <mo>(</mo> <msup> <msub> <mi>S</mi> <mi>q</mi> </msub> <mi>T</mi> </msup> <msub> <mi>S</mi> <mi>q</mi> </msub> <mo>-</mo> <mn>2</mn> <msup> <mi>D</mi> <mi>T</mi> </msup> <msup> <mi>P</mi> <mi>T</mi> </msup> <msub> <mi>S</mi> <mi>q</mi> </msub> <mo>+</mo> <msup> <mi>D</mi> <mi>T</mi> </msup> <msup> <mi>P</mi> <mi>T</mi> </msup> <mi>P</mi> <mi>D</mi> <mo>)</mo> </mrow> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>4</mn> <mo>)</mo> </mrow> </mrow>

Further gradient is asked to obtain

<mrow> <mfrac> <mrow> <mi>d</mi> <mi>f</mi> <mrow> <mo>(</mo> <mi>x</mi> <mo>)</mo> </mrow> </mrow> <mrow> <msup> <mi>dP</mi> <mi>T</mi> </msup> </mrow> </mfrac> <mo>=</mo> <mo>-</mo> <mn>2</mn> <msubsup> <mi>DS</mi> <mi>q</mi> <mi>T</mi> </msubsup> <mo>+</mo> <mn>2</mn> <msup> <mi>DD</mi> <mi>T</mi> </msup> <msup> <mi>P</mi> <mi>T</mi> </msup> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>5</mn> <mo>)</mo> </mrow> </mrow>

Using substep differential formulas, according to the property of the generation projection matrix of mask before, it becomes possible to be directed to from formula (5) In the Grad of each unknown mask variable；

Step (8) repeat step (2) arrives step (7) Iter times, the projection matrix P after output optimization_iter；Because mode of operation is Fixed, therefore threshold value t is the definite value of setting.

A kind of 2. compression EO-1 hyperion mask optimization method according to claim 1, it is characterised in that described initial projections Matrix P₀Derivation it is as follows：

The null matrix of one 100 × 100 is constructed first；

Then the diagonal entry using 100 numbers of mask as the null matrix, and surface element storage form is on diagonal Symbolic variable form；By this 100 × 100 reproduction matrix 31 times, 3100 × 100 matrix, i.e. projection matrix P are formed₀；Its In element is identical on diagonal in this 31 100 × 100 matrix, but element sequence is different, has the translation on position； This 3100 × 100 matrix is exactly the projection matrix P that we construct₀。