CN106657713A - Video motion amplification method - Google Patents

Abstract

A video motion amplification method, which belongs to the image processing method, solves the problem that the existing linear Euler video motion amplification method will significantly increase the noise when the amplification factor is increased, and the motion distortion of the video motion amplification method based on complex operable pyramids. . The present invention includes video image decomposition, video frame color space conversion, N-layer pyramid decomposition, time domain bandpass filtering, second spatial frequency band group amplification, obtaining fourth spatial frequency band group, reconstructing brightness matrix, video frame color space restoration , output video image step; the present invention solves the problem that the existing linear Euler video motion amplification method based on increasing the amplification factor can significantly increase the noise and the existing video motion amplification method based on complex operable pyramids. Problems such as motion distortion, It is suitable for amplifying tiny space-time movements.

Description

A video motion amplification method

technical field

The invention belongs to an image processing method, and in particular relates to a video motion amplification method, which is used for amplifying tiny spatiotemporal motion.

Background technique

The human eye has limited temporal and spatial sensitivity, and many real movements or color changes cannot be captured by the human eye, such as blood flowing through the skin, the skin color will change slightly, and the pulse caused by the periodic beating of the heart, however It is difficult for human eyes to observe the change of skin color and the beating of pulse. Using the method of video motion magnification, these motions that cannot cause visual responses can be amplified so that we can better observe various phenomena.

Mallat and Meyer proposed wavelet multiresolution analysis (Multiresolution Analysis, MRA) in the late 1980s. Wavelet multiresolution analysis is based on the scaling and translation of wavelet function and scale function, in the Hilbert space L ² (R) A certain subspace of a base is established, and then the base of the subspace is expanded to L ² (R) through scaling and translation transformation. In the wavelet multi-resolution signal decomposition, the scaling function φ _j,k and the wavelet function ψ _j,k together constitute the wavelet decomposition, where j is the scaling factor, k is the translation factor, the low frequency part of the signal is characterized by the scaling function, and the wavelet function Characterize the details of the signal. Comparative documents 1 and 2, published in July 1989, the United States, SGMallat, A theory for multiresolution signal decomposition: the wavelet representation, IEEE TPATTERN ANAL, 11 (1989) 674-693; published in September 1989, the United States, SGMallat, Multiresolution Approximations and Wavelet Orthonormal Bases of L ² (R), T AM MATH SOC, 315 (1989) 69-87. Wavelet decomposition and reconstruction of different scales can be performed through wavelet scaling, and wavelet decomposition and reconstruction of different scales constitute wavelet pyramid decomposition and reconstruction.

Existing video motion amplification methods are mainly divided into two types. The first is the linear Euler video motion amplification method, such as comparative document 3, published on August 5, 2012, the United States, H.Wu, M.Rubinstein, E. Shih, J. Guttag, Fr, D. Durand, W. Freeman, Eulerian video magnification for revealing subtle changes in the world, ACM Trans. Graph., 31(2012) 1-8. It uses the Laplacian pyramid algorithm for space-time decomposition and reconstruction. There are two defects in this method. The first defect is that when the spatial size of the video is small, it only supports a small magnification factor, that is, magnification The magnification is limited, and the second drawback is that increasing the amplification factor will significantly increase the noise. The second video motion magnification method is a video motion magnification method based on complex operable pyramids, such as Comparative Document 4, published in 2013, the United States, N.Wadhwa, M.Rubinstein, Fr, D.Durand, W.T.Freeman, Phase-based videomotion processing, ACM Trans.Graph., 32(2013) 1-10; although this method solves the two defects of the linear Euler video motion amplification method, it introduces a new defect, the first defect is the enlarged motion Distortion occurs, disorder occurs; the second defect is that it is relatively sensitive to vibrations of the video capture system, which means that small vibrations of the camera will be reflected in the resulting file.

Aiming at the fact that Eulerian video motion amplification is sensitive to noise and that phase-based operable pyramid algorithm motion amplification is disordered in the process of dealing with wave propagation, the applicant is committed to seeking a method that can solve the above problems. Because wavelet has advantages in noise processing, the applicant applies wavelet decomposition to video motion amplification. The Coiflet wavelet filter has the characteristic of linear phase, which can greatly reduce the signal distortion caused by the phase when the signal is reconstructed. Therefore, the wavelet used in this method is the Coiflet wavelet. Therefore, based on the previous work, the applicant In this paper, a new motion amplification method is proposed, which is called video motion amplification method based on Coiflet wavelet. This method first utilizes the advantages of multi-scale analysis of wavelet, uses the wavelet pyramid algorithm to spatially decompose the video frame, and decomposes the video frame into spatial frequency bands of different scales, then in order to extract specific motion, adopts a motion frequency adaptive The band-pass filter performs time-domain filtering on the spatial frequency band, and then amplifies the filtered spatial frequency band by a given factor α, and then returns it to the original signal, and finally uses the wavelet pyramid algorithm to reconstruct the video, that is, Small invisible fluctuations magnify into visible fluctuations.

For convenience of understanding the present invention, relevant concept is explained now:

RGB color space: the color space of a general video frame belongs to the RGB color space, and the RGB color space is obtained by changing the three color channels of red (Red), green (Green) and blue (Blue) Get a wide variety of colors. Usually, the RGB color space uses three unsigned integers (commonly used uint8, 8-bit unsigned integers) to represent a color. RGB each has 256 levels of brightness, expressed digitally from 0, 1, 2... until 255, the brightness increases sequentially, and up to 256 ³ colors can be formed through various combinations, so the RGB color space includes almost all colors that human vision can perceive, and is currently the most widely used One of the color systems.

YIQ color space: YIQ color space is usually adopted by the North American TV system and belongs to the NTSC (National Television Standards Committee) system. In the YIQ color space, the Y component represents the brightness information of the image, and the I and Q components carry color information. The I component represents the color change from orange to cyan, while the Q component represents the color change from purple to yellow-green. Compared with other color spaces, the YIQ color space has the advantage of being able to separate and extract the luminance components in the image, and the luminance components can be used to identify moving targets under complex backgrounds collected under natural conditions, so it is necessary to convert the video by RGB Color space conversion to YIQ color space;

Common wavelet functions include Haar, Daubechies, Coiflets, Symlets, etc.; taking the Coiflets wavelet function as an example, the abbreviation of the wavelet function is coif, which is divided into coif1, coif2, coif3, coif4 and coif5 according to the length of the support; the so-called The support length means that the wavelet high-pass filter or low-pass filter is a time-limited signal sequence, the signal value is non-zero in the support length interval, the signal value is zero outside the support interval, and the length of the non-zero interval is its support length;

Wavelet high-pass filter coefficients and wavelet low-pass filter coefficients: wavelet high-pass filter coefficients are one vector, wavelet low-pass filter coefficients are another vector, and wavelet high-pass filter coefficients corresponding to different wavelet functions can be obtained by consulting relevant literature and wavelet low-pass filter coefficients; such as Fan Yanbin et al., "Wavelet Theoretical Algorithms and Filter Banks", pp. 160-192, Science Press, Beijing, first edition in June 2011.

When wavelet decomposition algorithm or wavelet reconstruction algorithm is used, performing low-pass filtering or high-pass filtering on rows and columns is essentially to convolve the elements of the row or column with wavelet low-pass filter coefficients or wavelet high-pass filter coefficients, When the length of wavelet low-pass filter coefficients or wavelet high-pass filter coefficients is less than the number of elements in the row or column, you need to fill the missing digits with 0.

Computer file: Any computer file contains two parts, the first part is the header file part, and the second part is the data part. The header file part of the video image generally includes key content such as the width, height, storage space, format type, frame rate, data compression method and description of the video frame, and the data part of the video image stores the data associated with the header file . When a computer reads a video image, it first reads the header file to obtain various key information, especially the composition form of the data file, such as the width, height, frame rate and other information of the video frame. When processing video frames, it is the processing of data, such as color space transformation, wavelet pyramid decomposition, wavelet pyramid reconstruction and time domain filtering.

Contents of the invention

The invention provides a video motion amplification method, which belongs to the linear Euler video motion amplification method, which solves the problem that the increase of the amplification factor will significantly increase the noise existing in the existing linear Euler video motion amplification method and the video based on the complex operable pyramid. Problems such as motion distortion in the motion amplification method.

A video motion amplification method provided by the present invention includes a video image decomposition step, a video frame color space conversion step, an N-layer pyramid decomposition step, a time domain bandpass filtering step, a second spatial frequency band group amplification step, and a fourth Spatial frequency band group step, reconstruction brightness matrix step, video frame color space restoration step, output video image step, it is characterized in that:

(1) Video image decomposition steps:

Record the RGB color space video image of the moving subject moving slightly in space, and then decompose the RGB color space video image into RGB color space video frames frame by frame according to time, and read all the RGB color space video frames;

(2) video frame color space conversion steps:

Each RGB color space video frame is represented by three two-dimensional matrices of R matrix, G matrix, and B matrix. The R matrix represents the red color intensity of the pixel point, represents the green intensity of the pixel point, and the B matrix represents the blue intensity of the pixel point;

Perform matrix operations on the R matrix, G matrix, and B matrix by the following formula to obtain three two-dimensional matrices of the Y matrix, I matrix, and Q matrix of the YIQ color space video frame:

Y=0.299R+0.587G+0.114B,

I=0.596R-0.275G-0.321B,

Q=0.212R-0.523G+0.311B;

Among them, the Y matrix represents the brightness value of each pixel in the YIQ color space video frame, the I matrix represents the color intensity of each pixel in the YIQ color space video frame from orange to cyan, and the Q matrix represents the color intensity of each pixel in the YIQ color space video frame from purple to yellow-green the color intensity;

(3) N layer pyramid decomposition steps:

For the Y matrix in each YIQ color space video frame, use the wavelet decomposition algorithm to decompose it into four matrices: the first approximate coefficient matrix cA ₁ , the first horizontal detail matrix cH ₁ , the first vertical detail matrix cV ₁ and the first Diagonal detail matrix cD ₁ ; this process is the first layer of pyramid decomposition, cH ₁ , cV ₁ and cD ₁ constitute the first layer of the pyramid;

Then, the wavelet decomposition algorithm is used to decompose the first approximate coefficient matrix cA ₁ into the second approximate coefficient matrix cA ₂ , the second horizontal detail matrix cH ₂ , the second vertical detail matrix cV ₂ and the second diagonal detail matrix cD ₂ ; this process is the second layer of pyramid decomposition, cH ₂ , cV ₂ and cD ₂ constitute the second layer of the pyramid;

Carry out N times of decomposition in this way, and finally get the Nth approximate coefficient matrix cA _N , the Nth horizontal detail matrix cH _N , the Nth vertical detail matrix cV _N and the Nth diagonal detail matrix cD _N ; cA _N , cH _N , cV _N and cD _N constitute the Nth layer of the pyramid; N≥3, the larger the picture size, the larger the N value, otherwise the smaller the N value;

After performing N times of pyramid decomposition on all YIQ color space video frames, the nth horizontal detail matrix cH _n of each YIQ color space video frame forms a horizontal detail space frequency band with a scale of n, and each frame of YIQ color space video frame The nth vertical detail matrix cV _n forms a vertical detail space frequency band with a scale of n, and the nth diagonal detail matrix cD _n of each YIQ color space video frame forms a diagonal detail space frequency band with a scale of n; where n is called Be scale, represent different pyramid levels, n=1, 2, ..., N; The approximate coefficient space frequency band is formed by the Nth approximation coefficient matrix cA _N of each frame YIQ color space video frame; Above-mentioned each space frequency band forms the first space frequency band group;

(4) Steps of time-domain bandpass filtering:

The first low-pass IIR (Infinite Impulse Response, Infinite Impulse Response) digital filter is used to perform time-domain filtering on each spatial frequency band forming the first spatial frequency band group, and each spatial frequency band obtained after filtering constitutes the first temporary space frequency band group;

The second low-pass IIR digital filter is used to perform time-domain filtering on each spatial frequency band forming the first spatial frequency band group respectively, and each spatial frequency band obtained after filtering forms a second temporary spatial frequency band group;

Subtract each spatial frequency band in the first temporary spatial frequency band group from each corresponding spatial frequency band in the second temporary spatial frequency band group, and the resulting spatial frequency bands together form a second spatial frequency band group;

The cut-off frequencies of the first low-pass IIR digital filter and the second low-pass IIR digital filter are f ₁ and f ₂ respectively, 0<f ₂ <f ₁ <f _s /2, f _s is the recording frame rate of the video image ;

(5) second spatial frequency band group amplification step:

Calculate the maximum magnification α(n) _max of each layer of the pyramid separately, and compare whether α(n) _max ≤ α, if yes, use α(n) _max as the actual magnification of the corresponding layer, otherwise use α as the actual magnification of the corresponding layer Multiple; where α is the set magnification, 5<α<30;

Each spatial frequency band forming the second spatial frequency band group is multiplied by the actual magnification factor of the corresponding layer according to the number of layers of the pyramid where it is located, and the resulting spatial frequency bands jointly form the third spatial frequency band group;

(6) Obtain the fourth spatial frequency band group step:

performing matrix addition operation on each spatial frequency band constituting the third spatial frequency band group and corresponding spatial frequency bands in the first spatial frequency band group, and the obtained spatial frequency bands together form a fourth spatial frequency band group; The matrix addition operation is a matrix operation of two groups of corresponding scales and corresponding positions;

(7) Steps of reconstructing brightness matrix Y:

Adopt the wavelet reconstruction algorithm to reconstruct the Y4 matrix of the YIQ color space video frame with the ^fourth spatial frequency band group, the specific process is as follows:

Each layer of the different pyramid layers in the fourth spatial frequency band group is reconstructed layer by layer. Starting from the Nth layer of the pyramid, the wavelet reconstruction algorithm will be used to form the YIQ color space video frame of each frame in the fourth spatial frequency band group. The Nth approximate coefficient matrix cA _N ⁴ and the Nth horizontal detail matrix cH _N ⁴ , the Nth vertical detail matrix cV _N ⁴ and the Nth diagonal detail matrix cD _N ⁴ reconstruct the N-1th approximate coefficient matrix cA _N-1 ⁴ ;

Then carry out the reconstruction of the N-1th layer of the pyramid, and use the wavelet reconstruction algorithm for each frame of YIQ color space video frame to combine the N-1th approximate coefficient matrix cA _N-1 ⁴ and the N-1th The horizontal detail matrix cH _N-1 ⁴ , the N-1th vertical detail matrix cV _N-1 ⁴ , the N-1th diagonal detail matrix cD _N-1 ⁴ reconstruct the N-2th approximate coefficient matrix cA _N-2 ⁴ ;

Proceed in this way until the first reconstructed approximate coefficient matrix cA ₁ ⁴ is obtained; then the wavelet reconstruction algorithm is used to reconstruct the first layer of the pyramid, and for each frame of YIQ color space video frame, the wavelet reconstruction algorithm is used to transform the first The approximation coefficient matrix cA ₁ ⁴ and the first horizontal detail matrix cH ₁ ⁴ , the first vertical detail matrix cV ₁ ⁴ , and the first diagonal detail matrix cD ₁ ⁴ are reconstructed into a Y ⁴ matrix;

(8) video frame color space restoration steps:

Carry out matrix operations to the Y matrix reconstructed in step (7), the I matrix and Q matrix in step (2) by the following formula to obtain new R matrix, G matrix and B matrix:

R= ^1.000Y4 +0.956I+0.621Q,

G=1.000Y4-0.272I ^- 0.647Q,

B＝1.000Y4-1.106I ⁺ 1.703Q;

The new R matrix, G matrix and B matrix synthesize new RGB color space video frames;

Carry out above-mentioned conversion to each new YIQ color space video frame successively, just YIQ color space video frame is restored to RGB color space video frame;

(9) Output video image steps:

Construct a video image file, including the header file part and the data part, first set the storage location of the video image file, then construct the header file of the output video image according to the header file format of the input video image, and write it to the specified location on the hard disk, and then Designate a location on the hard disk for each frame of video frame data, thereby completing the construction of the RGB color space video image.

The video motion amplification method is characterized in that:

In described step (3), to Y matrix, adopting wavelet decomposition algorithm to decompose it into four matrices comprises the following substeps:

(3.1) carry out low-pass filtering to each row of Y matrix, then carry out column down-sampling, then carry out low-pass filtering to each column, carry out row down-sampling at last to obtain the first approximation coefficient cA ₁ matrix;

(3.2) Low-pass filtering is carried out to each row of the Y matrix, then column down-sampling is performed, then each column is subjected to high-pass filtering, and row down-sampling is finally performed to obtain the first horizontal detail c H ₁ matrix;

(3.3) High-pass filtering is performed on each row of the Y matrix, then column downsampling is performed, then each column is low-pass filtered, and row downsampling is finally performed to obtain the first vertical detail cV ₁ matrix;

(3.4) Carry out high-pass filtering to each row of Y matrix, then carry out column down-sampling, then carry out high-pass filtering to each column, finally carry out row down-sampling to obtain the first diagonal detail cD ₁ matrix;

Adopt wavelet decomposition algorithm to decompose the nth approximate coefficient matrix cA _n-1 , the sub-steps are the same as above, n=1, 2, ..., N;

The column downsampling is to keep the even columns and discard the odd columns; the row downsampling is to keep the even rows and discard the odd rows;

The low-pass filtering of the row is to convolve the row element with the wavelet low-pass filter coefficient to obtain a row of new element values to replace the original row;

Carrying out high-pass filtering to row is exactly to convolve this row element and wavelet high-pass filter coefficient, obtain a row of new element value, replace original row;

The low-pass filtering of the column is to convolve the column elements with wavelet low-pass filter coefficients to obtain a column of new element values to replace the original column;

Carrying out high-pass filtering to the column is to convolve the column elements with wavelet high-pass filter coefficients to obtain a column of new element values to replace the original column;

The coefficients of the wavelet low-pass filter and the wavelet high-pass filter can be obtained by consulting related documents.

The video motion amplification method is characterized in that:

In said step (4), adopting the first low-pass IIR digital filter to carry out time-domain filtering respectively to each spatial frequency band forming the first spatial frequency band group includes the following sub-steps:

(4.1) Calculate the filter coefficient r ₁ of the low-pass IIR filter:

r ₁ =2πf ₁ /f _s ,

In the formula, f ₁ is the cut-off frequency of the first low-pass IIR digital filter, and f _s is the recording frame rate of the video image;

(4.2) Calculate the output Y(m) of the first low-pass IIR digital filter:

Y(m)=(1-r ₁ )×Y(m-1)+r ₁ ×X(m),

In the formula, X(m) is the input of the filter, m is the sequence number of the video frame, m=1, 2, ..., K; wherein K is the total frame number of the video frame, when m is 1, X(1) is the Known, Y(0) is:

The second low-pass IIR digital filter is used to perform time-domain filtering on the spatial frequency bands constituting the first spatial frequency band group, the process is exactly the same as the above process, the only difference is that f ₁ is replaced by f ₂ .

The video motion amplification method is characterized in that:

In the described step (5), the maximum magnification α (n) _max of calculating each layer of the pyramid comprises the following sub-steps:

(5.1) Calculate the spatial wavelength λ(n) of each layer of the pyramid:

In the formula, W _n and H _n are the width and height of the nth layer of the pyramid respectively, and the unit is pixel;

(5.2) Calculate the displacement function δ(t):

In the formula, λ _c is the spatial critical wavelength, and the value is 16 to 20 pixels;

(5.3) Calculate the maximum magnification α(n) _max of each layer of the pyramid:

The video motion amplification method is characterized in that:

In the step (7), the wavelet reconstruction algorithm is used to form the Nth approximate coefficient matrix cA _N ⁴ and the Nth horizontal detail matrix cH _N ⁴ and the Nth level detail matrix cH N 4 of each frame YIQ color space video frame in the fourth spatial frequency band group. The N vertical detail matrix cV _N ⁴ and the Nth diagonal detail matrix cD _N ⁴ reconstruct the N-1th approximate coefficient matrix cA _N-1 ⁴ , including the following sub-steps:

(7.1) Perform row upsampling on the matrix cA _N ⁴ , then perform low-pass filtering on each column of the row upsampled matrix, then perform column upsampling on the filtered matrix, and then perform row upsampling on each row of the matrix obtained by column upsampling A _N1 matrix is obtained by low-pass filtering;

(7.2) Perform row upsampling on the matrix cH _N ⁴ , then perform high-pass filtering on each column of the matrix after row upsampling, then perform column upsampling on the filtered matrix, and then perform low-pass on each row of the matrix obtained by column upsampling Obtain A _N2 matrix through filtering;

(7.3) Perform row up-sampling on the cV _N ⁴ matrix, then perform low-pass filtering on each column of the row-up-sampled matrix, then perform column up-sampling on the filtered matrix, and then perform on-column up-sampling to obtain each row of the matrix A _N3 matrix is obtained by high-pass filtering;

(7.4) Perform row upsampling on the cD _N ⁴ matrix, then perform high-pass filtering on each column of the matrix after row upsampling, then perform column upsampling on the filtered matrix, and then perform high-pass on each row of the matrix obtained by column upsampling Filter to get A _N4 matrix;

(7.5) carry out matrix addition operation above-mentioned A _N1 matrix, A _N2 matrix, A _N3 matrix and A _N4 matrix obtain the N-1th approximation coefficient matrix cA _N-1 ⁴ ;

The so-called matrix row upsampling is to double the number of matrix rows, insert a row of zero vectors between any two rows of the original matrix, and finally fill the last row of the newly constructed matrix with zeros; the so-called matrix column upsampling is to multiply the matrix Double the number of columns, insert a column of zero vectors between any two columns of the original matrix, and finally fill the last column of the newly constructed matrix with zeros;

The low-pass filtering of the row is to convolve the row element with the wavelet low-pass filter coefficient to obtain a row of new element values to replace the original row;

Carrying out high-pass filtering to row is exactly to convolve this row element and wavelet high-pass filter coefficient, obtain a row of new element value, replace original row;

The low-pass filtering of the column is to convolve the column elements with wavelet low-pass filter coefficients to obtain a column of new element values to replace the original column;

Carrying out high-pass filtering to the column is to convolve the column elements with wavelet high-pass filter coefficients to obtain a column of new element values to replace the original column;

The coefficients of the wavelet low-pass filter and the wavelet high-pass filter can be obtained by consulting related documents.

The phase gradient of the video image is directly related to the motion in the video image. The basic idea of the existing phase-based complex operable pyramid motion amplification method is to amplify the phase gradient of the image to achieve the purpose of amplifying the motion in the image. Because the phase gradient of the image is distorted, the movement that causes the magnification is also distorted. The present invention does not involve the phase gradient of the image, so there is no motion distortion.

Existing linear Euler video motion amplification method has the problem that increasing the amplification factor can significantly increase noise, the present invention in step (4), the difference equation of the first low-pass IIR digital filter is: Y(m)=( 1-r ₁ )×Y(m-1)+r ₁ ×X(m), in actual operation, the existing linear Euler video motion amplification method performs the initial input of the filter X(1) and Y(0) In order to give up, time-domain filtering is performed directly from X(2) and Y(1), where Y(1)=X(1), the specific calculation process is:

Y(1)=X(1)

Y(2)＝(1-r ₁ )*Y(1)+r ₁ *X(2)

Y(3)＝(1-r ₁ )*Y(2)+r ₁ *X(3)

Y(m)＝(1-r ₁ )*Y(m-1)+r ₁ *X(m)

Y(K)＝(1-r ₁ )*Y(K-1)+r ₁ *X(K)

The concrete calculation process of step (4) of the present invention is:

Y(1)＝(1-r ₁ )*Y(0)+r ₁ *X(1)

Y(2)＝(1-r ₁ )*Y(1)+r ₁ *X(2)

Y(m)＝(1-r ₁ )*Y(m-1)+r ₁ *X(m)

Y(K)＝(1-r ₁ )*Y(K-1)+r ₁ *X(K)

By analyzing the noise of the output video of the existing linear Euler video motion amplification method, it is found that the noise is initially small, then increases rapidly, and then slowly decays. The test results using the present invention show that the noise in the output video is very stable and basically remains unchanged. In addition, from the perspective of the size of the noise, the present invention greatly reduces the noise in the output video.

To sum up, the present invention solves the problem that increasing the amplification factor will significantly increase the noise in the existing linear Euler video motion amplification method and the motion distortion problem in the existing video motion amplification method based on complex operable pyramids.

Description of drawings

Fig. 1 is a schematic flow sheet of the present invention;

Fig. 2 is original video frame sequence;

Fig. 3 is a schematic diagram of wavelet pyramid decomposition;

Figure 4 (A) is the original picture;

Fig. 4 (B) is a schematic diagram of a three-layer wavelet pyramid structure;

5 is a schematic diagram of a first spatial frequency band group;

Figure 6 is a schematic diagram of wavelet decomposition,

Fig. 7 (A) is the schematic diagram of the frequency response characteristic of the first low-pass IIR digital filter;

Fig. 7 (B) is the frequency response characteristic schematic diagram of the second low-pass IIR digital filter;

Fig. 7 (C) is the frequency response characteristic schematic diagram of the band-pass filter that uses the first, the second low-pass IIR digital filter two low-pass filters of different cut-off frequencies to construct;

Figure 8 is a fourth spatial frequency band group;

Fig. 9 is a schematic diagram of wavelet reconstruction;

Figure 10 shows the reconstructed video sequence.

detailed description

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

As shown in Figure 1, the present invention includes video image decomposition step, video frame color space conversion step, N layer pyramid decomposition step, time domain bandpass filtering step, second spatial frequency band group amplification step, obtains the 4th spatial frequency band group Step, the step of reconstructing the luminance matrix, the step of restoring the color space of the video frame, and the step of outputting the video image.

An embodiment of the present invention comprises the following steps:

(1) Video image decomposition steps:

Record the RGB color space video image of the moving subject moving slightly in space, and then decompose the RGB color space video image into RGB color space video frames frame by frame according to time, and read all the RGB color space video frames; as shown in Figure 2 As shown, the applicant used a SLR camera to record a video of two ginkgo trees moving slightly in the air, the frame rate f _s of the video is 25 frames per second, and the video size is 512×512 pixels.

(2) video frame color space conversion steps:

Each RGB color space video frame is represented by three two-dimensional matrices of R matrix, G matrix, and B matrix. The R matrix represents the red color intensity of the pixel point, represents the green intensity of the pixel point, and the B matrix represents the blue intensity of the pixel point;

Perform matrix operations on the R matrix, G matrix, and B matrix by the following formula to obtain three two-dimensional matrices of the Y matrix, I matrix, and Q matrix of the YIQ color space video frame:

Y=0.299R+0.587G+0.114B,

I=0.596R-0.275G-0.321B,

Q=0.212R-0.523G+0.311B;

Among them, the Y matrix represents the brightness value of each pixel in the YIQ color space video frame, the I matrix represents the color intensity of each pixel in the YIQ color space video frame from orange to cyan, and the Q matrix represents the color intensity of each pixel in the YIQ color space video frame from purple to yellow-green the color intensity;

(3) N layer pyramid decomposition steps:

For the Y matrix in each YIQ color space video frame, use the wavelet decomposition algorithm to decompose it into four matrices: the first approximate coefficient matrix cA ₁ , the first horizontal detail matrix cH ₁ , the first vertical detail matrix cV ₁ and the first Diagonal detail matrix cD ₁ ; this process is the first layer of pyramid decomposition, cH ₁ , cV ₁ and cD ₁ constitute the first layer of the pyramid;

Then, the wavelet decomposition algorithm is used to decompose the first approximate coefficient matrix cA ₁ into the second approximate coefficient matrix cA ₂ , the second horizontal detail matrix cH ₂ , the second vertical detail matrix cV ₂ and the second diagonal detail matrix cD ₂ ; this process is the second layer of pyramid decomposition, cH ₂ , cV ₂ and cD ₂ constitute the second layer of the pyramid;

Carry out N times of decomposition in this way, and finally get the Nth approximate coefficient matrix cA _N , the Nth horizontal detail matrix cH _N , the Nth vertical detail matrix cV _N and the Nth diagonal detail matrix cD _N ; cA _N , cH _N , cV _N and cD _N constitute the Nth layer of the pyramid; N≥3, the larger the picture size, the larger the N value, otherwise the smaller the N value;

After performing N times of pyramid decomposition on all YIQ color space video frames, the nth horizontal detail matrix cH _n of each YIQ color space video frame forms a horizontal detail space frequency band with a scale of n, and each frame of YIQ color space video frame The nth vertical detail matrix cV _n forms a vertical detail space frequency band with a scale of n, and the nth diagonal detail matrix cD _n of each YIQ color space video frame forms a diagonal detail space frequency band with a scale of n; where n is called Be scale, represent different pyramid levels, n=1, 2, ..., N; The approximate coefficient space frequency band is formed by the Nth approximation coefficient matrix cA _N of each frame YIQ color space video frame; Above-mentioned each space frequency band forms the first space frequency band group;

Figure 3 shows the storage form of wavelet pyramid decomposition image data using the wavelet pyramid algorithm. This figure is to perform wavelet pyramid decomposition on a picture with a size of 512×512 pixels. C represents the number of decomposed matrices, and S represents the number of each matrix. Size Information. Figure 4(B) shows the schematic diagram of the three-layer wavelet pyramid structure after wavelet decomposition. _Firstly , the wavelet decomposition algorithm is used to decompose the picture into cA ₁ , cH ₁ , _{cV 1} and cD ₁ . Decompose into the second layer pyramid cA ₂ , cH ₂ , cV ₂ , cD ₂ , cA ₁ is replaced by cA ₂ , cH ₂ , cV ₂ , cD ₂ , decompose cA ₂ into the third layer pyramid cA ₃ , cH ₃ , cV ₃ , cD ₃ , and cA ₂ are replaced by cA ₃ , cH ₃ , cV ₃ , and cD ₃ , thus forming a three-layer pyramid structure; Figure 5 shows the first spatial frequency band group, and the layout of this figure is similar to that of Figure 4 (B), for simplicity, only one layer of pyramid decomposition is performed, and the total number of layers of the pyramid is 1;

As shown in Figure 6, using the wavelet decomposition algorithm to decompose the Y matrix into four matrices includes the following sub-steps:

(3.1) carry out low-pass filtering to each row of Y matrix, then carry out column down-sampling, then carry out low-pass filtering to each column, carry out row down-sampling at last to obtain the first approximation coefficient cA ₁ matrix;

(3.2) Low-pass filtering is carried out to each row of the Y matrix, then column down-sampling is performed, then each column is subjected to high-pass filtering, and row down-sampling is finally performed to obtain the first horizontal detail c H ₁ matrix;

(3.3) High-pass filtering is performed on each row of the Y matrix, then column downsampling is performed, then each column is low-pass filtered, and row downsampling is finally performed to obtain the first vertical detail cV ₁ matrix;

(3.4) Carry out high-pass filtering to each row of Y matrix, then carry out column down-sampling, then carry out high-pass filtering to each column, finally carry out row down-sampling to obtain the first diagonal detail cD ₁ matrix;

Adopt wavelet decomposition algorithm to decompose the nth approximate coefficient matrix cA _n-1 , the sub-steps are the same as above, n=1, 2, ..., N;

The column downsampling is to keep the even columns and discard the odd columns; the row downsampling is to keep the even rows and discard the odd rows;

The low-pass filtering of the row is to convolve the row elements with wavelet low-pass filter coefficients to obtain a row of new element values to replace the original rows;

Carrying out high-pass filtering to row is exactly to convolve this row element and wavelet high-pass filter coefficient, obtain a row of new element value, replace original row;

The low-pass filtering of the column is to convolve the column elements with wavelet low-pass filter coefficients to obtain a column of new element values to replace the original column;

Carrying out high-pass filtering to the column is to convolve the column elements with wavelet high-pass filter coefficients to obtain a column of new element values to replace the original column;

In the present embodiment, the wavelet function that adopts is coif5, by consulting reference document 1 and reference document 2, obtain low-pass filter and high-pass filter, as listed in table 1;

Reference 1: China, Jianping Xuan, "Construction and Application of Higher-Order Coiflet Wavelet Systems", Journal of Vibration Engineering, (2008) 521-529; Reference 2: Rujie Yu, "Wave Propagation Problems Based on Zero-moment Scaling Functions Research on Solving Method", Master Thesis of Huazhong University of Science and Technology, 1-52, published in 2013;

(4) Steps of time-domain bandpass filtering:

The first low-pass IIR digital filter (Infinite Impulse Response, infinite impulse response digital filter) is used to perform time-domain filtering on each spatial frequency band forming the first spatial frequency band group, and each spatial frequency band obtained after filtering constitutes the second spatial frequency band. a temporary set of spatial frequency bands;

The second low-pass IIR digital filter is used to perform time-domain filtering on each spatial frequency band forming the first spatial frequency band group respectively, and each spatial frequency band obtained after filtering forms a second temporary spatial frequency band group;

Subtract each spatial frequency band in the first temporary spatial frequency band group from each corresponding spatial frequency band in the second temporary spatial frequency band group, and the resulting spatial frequency bands together form a second spatial frequency band group;

In the present embodiment, the cut-off frequencies f1 and f2 of the _first low-pass IIR digital filter and the _second low-pass IIR digital filter are 0.4Hz and 0.2Hz respectively, and Fig. 7 (A) shows the first low-pass IIR digital filter The frequency response characteristic of the filter, Fig. 7 (B) shows the frequency response characteristic of the second low-pass IIR digital filter, Fig. 7 (C) shows the band-pass filter that is made of the first and the second low-pass IIR digital filter frequency response characteristics.

(5) second spatial frequency band group amplification step:

Calculate the maximum magnification α(n) _max of each layer of the pyramid separately, and compare whether α(n) _max ≤ α, if yes, use α(n) _max as the actual magnification of the corresponding layer, otherwise use α as the actual magnification of the corresponding layer Multiple; where α is the set magnification, 5<α<30;

Each spatial frequency band forming the second spatial frequency band group is multiplied by the actual magnification factor of the corresponding layer according to the number of layers of the pyramid where it is located, and the resulting spatial frequency bands jointly form the third spatial frequency band group;

In the present embodiment, calculating the maximum magnification α(n) _max of each layer of the pyramid includes the following sub-steps:

(5.1) Calculate the spatial wavelength λ(n) of each layer of the pyramid:

In the formula, W _n and H _n are the width and height of the nth layer of the pyramid, and the unit is pixel; taking the first layer and the second layer as an example, the width W ₁ and height H ₁ of the first layer are 256 and 256 respectively, The width W ₂ and height H ₂ of the second layer are 128 and 128 respectively, refer to FIG. 3 .

(5.2) Calculate the displacement function δ(t):

In the formula, λ _c is the spatial critical wavelength, which takes a value of 16 pixels, and sets the magnification factor α to 15;

(5.3) Calculate the maximum magnification α(n) _max of each layer of the pyramid:

Taking the first layer and the second layer as an example, the specific process of calculating the maximum magnification is as follows: First, the width W ₁ and height H ₁ of the first layer are known to be 256 and 256 respectively, and the width W ₂ and height of the second layer are H ₂ are 128 and 128 respectively, and they are substituted into the calculation formula to calculate the spatial wavelength of each layer to obtain λ(1)=361.984 and λ(2)=180.992, and then according to the set spatial critical wavelength λ _c =16 and Set the magnification α=15 to calculate δ(t)=0.125, and finally substitute these three calculated values to obtain the maximum magnifications of the first layer and the second layer respectively: α(1) _max =119.661, α (2) _max = 59.331.

(6) Obtain the fourth spatial frequency band group step:

performing matrix addition operation on each spatial frequency band constituting the third spatial frequency band group and corresponding spatial frequency bands in the first spatial frequency band group, and the obtained spatial frequency bands together form a fourth spatial frequency band group; The matrix addition operation is a matrix operation of two groups of corresponding scales and corresponding positions;

As shown in FIG. 8 , this figure is the superimposed spatial frequency band. This superposition process is the core approximation part of linear Euler. In linear Euler video motion amplification, its core idea can be expressed as:

Assuming that the intensity I(x) of any pixel in a one-dimensional image changes with time t as I(x,t), and the initial value is written as I(x,0)=h(x), then at any time The intensity of the pixel point can be expressed as I(x,t)=h(x+δ(t)), when the motion of the object is small, it can be expanded by Taylor approximation as shown in formula 1:

In formula 1, h(x+δ(t)) is the original pixel intensity of the image, which represents the first spatial frequency band, In order to characterize the change in intensity of moving pixels of the object, and to obtain the motion of interest, a band-pass filter may be used to perform time-domain band-pass filtering on the first spatial frequency band. The second spatial frequency band is obtained by this time-domain processing, and the second spatial frequency band can be used Indicates that the intensity value of the second spatial frequency band is multiplied by an amplification factor α, so that the third spatial frequency band is obtained, using Represent the third spatial frequency band, and then add it to the original first spatial frequency band, the addition process is represented by the following formula:

The right side of formula 2 can also be approximated as:

By adding back, the fourth spatial frequency band is obtained. By comparing the left side of formula 1 with the right side of formula 3, it can be found that the motion representing the change with time is amplified by (1+α) times.

(7) Steps of reconstructing brightness matrix:

Adopt the wavelet reconstruction algorithm to reconstruct the Y4 matrix of the YIQ color space video frame with the ^fourth spatial frequency band group, the specific process is as follows:

Each layer of the different pyramid layers in the fourth spatial frequency band group is reconstructed layer by layer. Starting from the Nth layer of the pyramid, the wavelet reconstruction algorithm will be used to form the YIQ color space video frame of each frame in the fourth spatial frequency band group. The Nth approximate coefficient matrix cA _N ⁴ and the Nth horizontal detail matrix cH _N ⁴ , the Nth vertical detail matrix cV _N ⁴ and the Nth diagonal detail matrix cD _N ⁴ reconstruct the N-1th approximate coefficient matrix cA _N-1 ⁴ ;

Then carry out the reconstruction of the N-1th layer of the pyramid, and use the wavelet reconstruction algorithm for each frame of YIQ color space video frame to combine the N-1th approximate coefficient matrix cA _N-1 ⁴ and the N-1th The horizontal detail matrix cH _N-1 ⁴ , the N-1th vertical detail matrix cV _N-1 ⁴ , the N-1th diagonal detail matrix cD _N-1 ⁴ reconstruct the N-2th approximate coefficient matrix cA _N-2 ⁴ ;

Proceed in this way until the first reconstructed approximate coefficient matrix cA ₁ ⁴ is obtained; then the wavelet reconstruction algorithm is used to reconstruct the first layer of the pyramid, and for each frame of YIQ color space video frame, the wavelet reconstruction algorithm is used to transform the first The approximation coefficient matrix cA ₁ ⁴ and the first horizontal detail matrix cH ₁ ⁴ , the first vertical detail matrix cV ₁ ⁴ , and the first diagonal detail matrix cD ₁ ⁴ are reconstructed into a Y ⁴ matrix;

As shown in Figure 9, in this embodiment, the wavelet reconstruction algorithm is used to combine the Nth approximation coefficient matrix cA _N4 and the _Nth horizontal detail matrix ^cHN of each frame YIQ color space video frame in the fourth spatial frequency band group ^4. The Nth vertical detail matrix cV _N ⁴ and the Nth diagonal detail matrix cD _N ⁴ reconstruct the N-1th approximate coefficient matrix cA _N-1 ⁴ , including the following sub-steps:

(7.1) Perform row upsampling on the matrix cA _N ⁴ , then perform low-pass filtering on each column of the row upsampled matrix, then perform column upsampling on the filtered matrix, and then perform row upsampling on each row of the matrix obtained by column upsampling A _N1 matrix is obtained by low-pass filtering;

(7.2) Perform row upsampling on the matrix cH _N ⁴ , then perform high-pass filtering on each column of the matrix after row upsampling, then perform column upsampling on the filtered matrix, and then perform low-pass on each row of the matrix obtained by column upsampling Obtain A _N2 matrix through filtering;

(7.3) Perform row up-sampling on the cV _N ⁴ matrix, then perform low-pass filtering on each column of the row-up-sampled matrix, then perform column up-sampling on the filtered matrix, and then perform on-column up-sampling to obtain each row of the matrix A _N3 matrix is obtained by high-pass filtering;

(7.4) Perform row upsampling on the cD _N ⁴ matrix, then perform high-pass filtering on each column of the matrix after row upsampling, then perform column upsampling on the filtered matrix, and then perform high-pass on each row of the matrix obtained by column upsampling Filter to get A _N4 matrix;

(7.5) carry out matrix addition operation above-mentioned A _N1 matrix, A _N2 matrix, A _N3 matrix and A _N4 matrix obtain the N-1th approximation coefficient matrix cA _N-1 ⁴ ;

The so-called matrix row upsampling is to double the number of matrix rows, insert a row of zero vectors between any two rows of the original matrix, and finally fill the last row of the newly constructed matrix with zeros; the so-called matrix column upsampling is to multiply the matrix Double the number of columns, insert a column of zero vectors between any two columns of the original matrix, and finally fill the last column of the newly constructed matrix with zeros;

The low-pass filtering of the row is to convolve the row element with the wavelet low-pass filter coefficient to obtain a row of new element values to replace the original row;

Carrying out high-pass filtering to row is exactly to convolve this row element and wavelet high-pass filter coefficient, obtain a row of new element value, replace original row;

The low-pass filtering of the column is to convolve the column elements with wavelet low-pass filter coefficients to obtain a column of new element values to replace the original column;

Carrying out high-pass filtering to the column is to convolve the column elements with wavelet high-pass filter coefficients to obtain a column of new element values to replace the original column;

In the present embodiment, the wavelet function adopted is coif5, by consulting aforementioned reference document 1 and reference document 2, obtain wavelet low-pass filter coefficient and wavelet high-pass filter coefficient as listed in Table 1;

(8) video frame color space restoration steps:

Carry out matrix operation to the Y matrix in step ( ⁷ ) reconstruction, I matrix and Q matrix in step (2) by following formula, obtain new R matrix, G matrix and B matrix:

R= ^1.000Y4 +0.956I+0.621Q,

G=1.000Y4-0.272I ^- 0.647Q,

B＝1.000Y4-1.106I ⁺ 1.703Q;

The new R matrix, G matrix and B matrix synthesize new RGB color space video frames;

Carry out above-mentioned conversion to each new YIQ color space video frame successively, just YIQ color space video frame is restored to RGB color space video frame;

(9) Output video image steps:

Construct a video image file, including the header file part and the data part, first set the storage location of the video image file, then construct the header file of the output video image according to the header file format of the input video image, and write it to the specified location on the hard disk, and then Designate a location on the hard disk for each frame of video frame data, thereby completing the construction of the RGB color space video image. The constructed RGB color space video image is shown in Figure 10. Comparing it with Figure 2, it is obvious that the motion in the video is amplified.

Table 1 coif5 low-pass filter coefficients and high-pass filter coefficients

Claims (5)

1. A video motion amplification method, comprising a video image decomposition step, a video frame color space conversion step, an N-layer pyramid decomposition step, a time domain bandpass filtering step, a second spatial frequency band group amplification step, and obtaining the fourth spatial frequency band Group step, reconstruct brightness matrix step, video frame color space restoration step, output video image step, it is characterized in that: (1) Video image decomposition steps: Record the RGB color space video image of the moving subject moving slightly in space, and then decompose the RGB color space video image into RGB color space video frames frame by frame according to time, and read all the RGB color space video frames; (2) video frame color space conversion steps: Each RGB color space video frame is represented by three two-dimensional matrices of R matrix, G matrix, and B matrix. The R matrix represents the red color intensity of the pixel point, represents the green intensity of the pixel point, and the B matrix represents the blue intensity of the pixel point; Perform matrix operations on the R matrix, G matrix, and B matrix by the following formula to obtain three two-dimensional matrices of the Y matrix, I matrix, and Q matrix of the YIQ color space video frame: Y=0.299R+0.587G+0.114B, I=0.596R-0.275G-0.321B, Q=0.212R-0.523G+0.311B; Among them, the Y matrix represents the brightness value of each pixel in the YIQ color space video frame, the I matrix represents the color intensity of each pixel in the YIQ color space video frame from orange to cyan, and the Q matrix represents the color intensity of each pixel in the YIQ color space video frame from purple to yellow-green the color intensity; (3) N layer pyramid decomposition steps: For the Y matrix in each YIQ color space video frame, use the wavelet decomposition algorithm to decompose it into four matrices: the first approximate coefficient matrix cA ₁ , the first horizontal detail matrix cH ₁ , the first vertical detail matrix cV ₁ and the first Diagonal detail matrix cD ₁ ; this process is the first layer of pyramid decomposition, cH ₁ , cV ₁ and cD ₁ constitute the first layer of the pyramid; Then, the wavelet decomposition algorithm is used to decompose the first approximate coefficient matrix cA ₁ into the second approximate coefficient matrix cA ₂ , the second horizontal detail matrix cH ₂ , the second vertical detail matrix cV ₂ and the second diagonal detail matrix cD ₂ ; this process is the second layer of pyramid decomposition, cH ₂ , cV ₂ and cD ₂ constitute the second layer of the pyramid; Carry out N times of decomposition in this way, and finally get the Nth approximate coefficient matrix cA _N , the Nth horizontal detail matrix cH _N , the Nth vertical detail matrix cV _N and the Nth diagonal detail matrix cD _N ; cA _N , cH _N , cV _N and cD _N constitute the Nth layer of the pyramid; N≥3, the larger the picture size, the larger the N value, otherwise the smaller the N value; After performing N times of pyramid decomposition on all YIQ color space video frames, the nth horizontal detail matrix cH _n of each YIQ color space video frame forms a horizontal detail space frequency band with a scale of n, and each frame of YIQ color space video frame The nth vertical detail matrix cV _n forms a vertical detail space frequency band with a scale of n, and the nth diagonal detail matrix cD _n of each YIQ color space video frame forms a diagonal detail space frequency band with a scale of n; where n is called Be scale, represent different pyramid levels, n=1, 2, ..., N; The approximate coefficient space frequency band is formed by the Nth approximation coefficient matrix cA _N of each frame YIQ color space video frame; Above-mentioned each space frequency band forms the first space frequency band group; (4) Steps of time-domain bandpass filtering: Using the first low-pass IIR digital filter to perform time-domain filtering on the spatial frequency bands forming the first spatial frequency band group respectively, and each spatial frequency band obtained after filtering constitutes the first temporary spatial frequency band group; The second low-pass IIR digital filter is used to perform time-domain filtering on each spatial frequency band forming the first spatial frequency band group respectively, and each spatial frequency band obtained after filtering forms a second temporary spatial frequency band group; Subtract each spatial frequency band in the first temporary spatial frequency band group from each corresponding spatial frequency band in the second temporary spatial frequency band group, and the resulting spatial frequency bands together form a second spatial frequency band group; The cut-off frequencies of the first low-pass IIR digital filter and the second low-pass IIR digital filter are f ₁ and f ₂ respectively, 0<f ₂ <f ₁ <f _s /2, f _s is the recording frame rate of the video image ; (5) second spatial frequency band group amplification step: Calculate the maximum magnification α(n) _max of each layer of the pyramid separately, and compare whether α(n) _max ≤ α, if yes, use α(n) _max as the actual magnification of the corresponding layer, otherwise use α as the actual magnification of the corresponding layer Multiple; where α is the set magnification, 5<α<30; Each spatial frequency band forming the second spatial frequency band group is multiplied by the actual magnification factor of the corresponding layer according to the number of layers of the pyramid where it is located, and the resulting spatial frequency bands jointly form the third spatial frequency band group; (6) Obtain the fourth spatial frequency band group step: performing matrix addition operation on each spatial frequency band constituting the third spatial frequency band group and corresponding spatial frequency bands in the first spatial frequency band group, and the obtained spatial frequency bands together form a fourth spatial frequency band group; The matrix addition operation is a matrix operation of two groups of corresponding scales and corresponding positions; (7) Steps of reconstructing brightness matrix: Adopt the wavelet reconstruction algorithm to reconstruct the Y4 matrix of the YIQ color space video frame with the ^fourth spatial frequency band group, the specific process is as follows: Each layer of the different pyramid layers in the fourth spatial frequency band group is reconstructed layer by layer. Starting from the Nth layer of the pyramid, the wavelet reconstruction algorithm is used to form the YIQ color space video frame of each frame in the fourth spatial frequency band group The Nth approximate coefficient matrix cA _N ⁴ and the Nth horizontal detail matrix cH _N ⁴ , the Nth vertical detail matrix cV _N ⁴ and the Nth diagonal detail matrix cD _N ⁴ reconstruct the N-1th approximate coefficient matrix cA _N-1 ⁴ ; Then carry out the reconstruction of the N-1th layer of the pyramid, and use the wavelet reconstruction algorithm for each frame of YIQ color space video frame to combine the N-1th approximate coefficient matrix cA _N-1 ⁴ and the N-1th The horizontal detail matrix cH _N-1 ⁴ , the N-1th vertical detail matrix cV _N-1 ⁴ , the N-1th diagonal detail matrix cD _N-1 ⁴ reconstruct the N-2th approximate coefficient matrix cA _N-2 ⁴ ; Proceed in this way until the first reconstructed approximate coefficient matrix cA ₁ ⁴ is obtained; then the wavelet reconstruction algorithm is used to reconstruct the first layer of the pyramid, and for each frame of YIQ color space video frame, the wavelet reconstruction algorithm is used to transform the first The approximation coefficient matrix cA ₁ ⁴ and the first horizontal detail matrix cH ₁ ⁴ , the first vertical detail matrix cV ₁ ⁴ , and the first diagonal detail matrix cD ₁ ⁴ are reconstructed into a Y ⁴ matrix; (8) video frame color space restoration steps: Carry out matrix operation to the Y matrix in step ( ⁷ ) reconstruction, I matrix and Q matrix in step (2) by following formula, obtain new R matrix, G matrix and B matrix: R= ^1.000Y4 +0.956I+0.621Q, G=1.000Y4-0.272I ^- 0.647Q, B＝1.000Y4-1.106I ⁺ 1.703Q; The new R matrix, G matrix and B matrix synthesize new RGB color space video frames; Carry out above-mentioned conversion to each new YIQ color space video frame successively, just YIQ color space video frame is restored to RGB color space video frame; (9) Output video image steps: Construct a video image file, including the header file part and the data part, first set the storage location of the video image file, then construct the header file of the output video image according to the header file format of the input video image, and write it to the specified location on the hard disk, and then Designate a location on the hard disk for each frame of video frame data, thereby completing the construction of the RGB color space video image. 2. video motion amplification method as claimed in claim 1, is characterized in that: In described step (3), to Y matrix, adopting wavelet decomposition algorithm to decompose it into four matrices comprises the following substeps: (3.1) carry out low-pass filtering to each row of Y matrix, then carry out column down-sampling, then carry out low-pass filtering to each column, carry out row down-sampling at last to obtain the first approximation coefficient cA ₁ matrix; (3.2) Low-pass filtering is carried out to each row of the Y matrix, then column down-sampling is performed, then each column is subjected to high-pass filtering, and row down-sampling is finally performed to obtain the first horizontal detail c H ₁ matrix; (3.3) High-pass filtering is performed on each row of the Y matrix, then column downsampling is performed, then each column is low-pass filtered, and row downsampling is finally performed to obtain the first vertical detail cV ₁ matrix; (3.4) Carry out high-pass filtering to each row of Y matrix, then carry out column down-sampling, then carry out high-pass filtering to each column, finally carry out row down-sampling to obtain the first diagonal detail cD ₁ matrix; Adopt wavelet decomposition algorithm to decompose the nth approximate coefficient matrix cA _n-1 , the sub-steps are the same as above, n=1, 2, ..., N; The column downsampling is to keep the even columns and discard the odd columns; the row downsampling is to keep the even rows and discard the odd rows; The low-pass filtering of the row is to convolve the row element with the wavelet low-pass filter coefficient to obtain a row of new element values to replace the original row; Carrying out high-pass filtering to row is exactly to convolve this row element and wavelet high-pass filter coefficient, obtain a row of new element value, replace original row; The low-pass filtering of the column is to convolve the column elements with wavelet low-pass filter coefficients to obtain a column of new element values to replace the original column; Carrying out high-pass filtering to the column is to convolve the column elements with wavelet high-pass filter coefficients to obtain a column of new element values to replace the original column; The coefficients of the wavelet low-pass filter and the wavelet high-pass filter can be obtained by consulting related documents. 3. video motion amplification method as claimed in claim 1, is characterized in that: In said step (4), adopting the first low-pass IIR digital filter to carry out time-domain filtering respectively to each spatial frequency band forming the first spatial frequency band group includes the following sub-steps: (4.1) Calculate the filter coefficient r ₁ of the low-pass IIR filter: r ₁ =2πf ₁ /f _s , In the formula, f ₁ is the cut-off frequency of the first low-pass IIR digital filter, and f _s is the recording frame rate of the video image; (4.2) Calculate the output Y(m) of the first low-pass IIR digital filter: Y(m)=(1-r ₁ )×Y(m-1)+r ₁ ×X(m), In the formula, X(m) is the input of the filter, m is the sequence number of the video frame, m=1, 2, ..., K; K is the total frame number of the video frame, when m is 1, X(1) is known , Y(0) is: $\mathrm{<span\; class="notranslate">\; Y</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 0</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; =</span>\frac{\mathrm{<span\; class="notranslate">\; 1</span>}}{\mathrm{<span\; class="notranslate">\; K</span>}}\underset{\mathrm{<span\; class="notranslate">\; m</span>}<span\; class="notranslate">\; =</span>\mathrm{<span\; class="notranslate">\; 1</span>}}{\overset{\mathrm{<span\; class="notranslate">\; K</span>}}{<span\; class="notranslate">\; \&Sigma;</span>}}\mathrm{<span\; class="notranslate">\; x</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; m</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; ,</span>$ The second low-pass IIR digital filter is used to perform time-domain filtering on the spatial frequency bands constituting the first spatial frequency band group, the process is exactly the same as the above process, the only difference is that f ₁ is replaced by f ₂ . 4. video motion amplification method as claimed in claim 1, is characterized in that: In the described step (5), the maximum magnification α (n) _max of calculating each layer of the pyramid comprises the following sub-steps: (5.1) Calculate the spatial wavelength λ(n) of each layer of the pyramid: $\mathrm{<span\; class="notranslate">\; \&lambda;</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; no</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; =</span>\sqrt{{{\mathrm{<span\; class="notranslate">\; W</span>}}_{\mathrm{<span\; class="notranslate">\; no</span>}}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}<span\; class="notranslate">\; +</span>{{\mathrm{<span\; class="notranslate">\; h</span>}}_{\mathrm{<span\; class="notranslate">\; no</span>}}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; no</span>}<span\; class="notranslate">\; =</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; 2</span>}<span\; class="notranslate">\; ,</span><span\; class="notranslate">\; ...</span><span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; N</span>}<span\; class="notranslate">\; ;</span>$ In the formula, W _n and H _n are the width and height of the nth layer of the pyramid respectively, and the unit is pixel; (5.2) Calculate the displacement function δ(t): $\mathrm{<span\; class="notranslate">\; \&delta;</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; t</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; =</span>\frac{{\mathrm{<span\; class="notranslate">\; \&lambda;</span>}}_{\mathrm{<span\; class="notranslate">\; c</span>}}}{\mathrm{<span\; class="notranslate">\; 8</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; +</span>\mathrm{<span\; class="notranslate">\; \&alpha;</span>}<span\; class="notranslate">\; )</span>}<span\; class="notranslate">\; ,</span>$ In the formula, λ _c is the spatial critical wavelength, and the value is 16 to 20 pixels; (5.3) Calculate the maximum magnification α(n) _max of each layer of the pyramid: $\mathrm{<span\; class="notranslate">\; \&alpha;</span>}{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; no</span>}<span\; class="notranslate">\; )</span>}_{\mathrm{<span\; class="notranslate">\; m</span>}\mathrm{<span\; class="notranslate">\; a</span>}\mathrm{<span\; class="notranslate">\; x</span>}}<span\; class="notranslate">\; =</span>\frac{\mathrm{<span\; class="notranslate">\; \&lambda;</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; no</span>}<span\; class="notranslate">\; )</span>}{\mathrm{<span\; class="notranslate">\; 3</span>}<span\; class="notranslate">\; \&times;</span>\mathrm{<span\; class="notranslate">\; 8</span>}<span\; class="notranslate">\; \&times;</span>\mathrm{<span\; class="notranslate">\; \&delta;</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; t</span>}<span\; class="notranslate">\; )</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 1.</span>}$ 5. video motion amplification method as claimed in claim 1, is characterized in that: In the step (7), the wavelet reconstruction algorithm is used to form the Nth approximate coefficient matrix cA _N ⁴ and the Nth horizontal detail matrix cH _N ⁴ and the Nth level detail matrix cH N 4 of each frame YIQ color space video frame in the fourth spatial frequency band group. The N vertical detail matrix cV _N ⁴ and the Nth diagonal detail matrix cD _N ⁴ reconstruct the N-1th approximate coefficient matrix cA _N-1 ⁴ , including the following sub-steps: (7.1) Carry out row upsampling to matrix _cAN ⁴ , then perform low-pass filtering on each column of the matrix after row upsampling, then perform column upsampling on the filtered matrix, and then perform each row of the matrix obtained by column upsampling A _N1 matrix is obtained by low-pass filtering; (7.2) Perform row upsampling on the matrix cH _N ⁴ , then perform high-pass filtering on each column of the matrix after row upsampling, then perform column upsampling on the filtered matrix, and then perform low-pass on each row of the matrix obtained by column upsampling Obtain A _N2 matrix by filtering; (7.3) Perform row up-sampling on the cV _N ⁴ matrix, then perform low-pass filtering on each column of the row-up-sampled matrix, then perform column up-sampling on the filtered matrix, and then perform on-column up-sampling to obtain each row of the matrix A _N3 matrix is obtained by high-pass filtering; (7.4) Perform row upsampling on the cD _N ⁴ matrix, then perform high-pass filtering on each column of the matrix after row upsampling, then perform column upsampling on the filtered matrix, and then perform high-pass on each row of the matrix obtained by column upsampling Filter to get A _N4 matrix; (7.5) carry out matrix addition operation above-mentioned A _N1 matrix, A _N2 matrix, A _N3 matrix and A _N4 matrix obtain the N-1th approximation coefficient matrix cA _N-1 ⁴ ; The so-called matrix row upsampling is to double the number of rows of the matrix, insert a row of zero vectors between any two rows of the original matrix, and finally fill the last row of the newly constructed matrix with zeros; the so-called matrix column upsampling is to multiply the matrix Double the number of columns, insert a column of zero vectors between any two columns of the original matrix, and finally fill the last column of the newly constructed matrix with zeros; The low-pass filtering of the row is to convolve the row element with the wavelet low-pass filter coefficient to obtain a row of new element values to replace the original row; Carrying out high-pass filtering to row is exactly to convolve this row element and wavelet high-pass filter coefficient, obtain a row of new element value, replace original row; The low-pass filtering of the column is to convolve the column elements with wavelet low-pass filter coefficients to obtain a column of new element values to replace the original column; Carrying out high-pass filtering to the column is to convolve the column elements with wavelet high-pass filter coefficients to obtain a column of new element values to replace the original column; The coefficients of the wavelet low-pass filter and the wavelet high-pass filter can be obtained by consulting related documents.