CN103501437A - Fractal and H.264-based hyper-spectral image compression method - Google Patents

Abstract

The invention provides a fractal and H.264-based hyper-spectral image compression method. The method comprises the following steps: firstly, converting a hyper-spectral image data cube into video of a YUV (Luma and Chroma) format and transmitting into a coder; performing intra-frame prediction on I frame of the hyper-spectral video by using a macro block flatness-based H.264 fast intra-frame prediction mode discrimination method to remove spatial correlation of a hyper-spectral image; performing block motion estimation/complementary fractal coding on P frame of the hyper-spectral video to remove spectral correlation of the hyper-spectral image; searching optimal matching blocks from a reference frame for various blocks of tree structures of the P frame coded macro block by an MSE (Mean Square Error) criteria, determining an iterative function system coefficient of each block, comparing the rate distortion costs of all block division modes, taking the block division mode with the least rate distortion rate as the final inter-frame coding mode, and recording a final fractal parameter; finally coding and writing residual data of the I frame and the P frame into a code stream by entropy coding CABAC after DCT (Discrete Cosine Transformation) and quantization, wherein the fractal parameter of the P frame is also subjected to CABAC entropy coding.

Description

A kind of based on method for compressing high spectrum image fractal and H.264

Technical field:

The invention belongs to the Compression of hyperspectral images field, for the spatial coherence existed in high spectrum image and Spectral correlation, propose a kind of based on method for compressing high spectrum image fractal and H.264, under the prerequisite that guarantees picture quality, greatly accelerate the compression speed of high spectrum image, improved compression ratio.

Background technology:

High spectral technique is frontier development and one of the focus paid close attention to of current remote sensing circle of 21 century remote sensing technology, it organically combines traditional two-dimensional imaging remote sensing technology and spectral technique, when by imaging system, obtaining the spatial information of measured object, by spectrometer system, the RADIATION DECOMPOSITION of measured object is become to the spectrum radiation of different wave length, can in a spectrum range, obtain an even hundreds of continuous carrier wave segment information of each pixel tens.Oneself is successfully applied high spectrum resolution remote sensing technique in geology, ecology, atmospheric research, a lot of fields of soil investigation by warp, demonstrates very large potentiality and vast potential for future development.

High spectrum image, on the basis of ground image two-dimensional signal, has increased third dimension spectral information.High-spectral data is the cube of a spectrum picture, and its spatial image dimension is described earth's surface two-dimensional space feature, and its spectrum dimension discloses the curve of spectrum feature of each pixel of image.In the spatial image dimension Existential Space redundancy of high spectrum image data, there is redundancy between spectrum in the spectrum dimension.Therefore the difference of Compression of hyperspectral images and normal image compression is to consider to remove spatial coherence, considers again to remove Spectral correlation.At present, method for compressing high spectrum image commonly used is divided into three major types: the method based on prediction, method and the method based on vector quantization based on conversion.During such as Chen Yu and Zhang Ye etc. the compression scheme based on the linear model optimum prediction has been proposed, utilize the bi-directional predicted thought of recurrence, by setting up the linear model between bands of a spectrum, shown that the optimum prediction under the signal to noise ratio is (referring to Chen Yushi, Zhang Ye, Zhang Jun is flat. the Compression of hyperspectral images based on the linear model optimum prediction [J]. and Nanjing Aero-Space University's journal, 2007,39 (3): 368-372.).Penna B counts and extracts the amount of calculation that strategy reduces KLT by employing, and decorrelation between this improved KLT is composed for high spectrum image, and reconstructed image quality is not subject to obvious impact (referring to Penna B, Tillo T, Magli E, et al.Transform coding techniques for lossy hyperspectral data compression[J] .IEEE Transactions on Geoscience and Remote Sensing, 2007,45 (5): 1408-1421.).Shen-En Qian has proposed a kind of method of vector quantity fast algorithm and has improved the efficiency of the generation of code book, at this algorithm without entirely searching for, reduced widely computing complexity (referring to Shen-En Qian.Hyperspectral data compression using a fast vector quantization algorithm[J] .IEEE Transactions on Geoscience and Remote Sensing.2004,42 (8): 1791-1798.).Except above-mentioned traditional coding method, fractal image also is applied to the compression of high spectrum image in recent years, as Xia Lili has designed a three-dimensional fractal Coding Compression Algorithm, removed the correlation between band image when removing the image space correlation, obtained preferably compression effectiveness (referring to Xia Lili. the research of the Hyperspectral image compression algorithm based on fractal theory [J]. computer and digital engineering, 2011,39 (9): 132-135.).

The present invention combines predictive coding, fractal image, dct transform coding, CABAC entropy coding, utilization reduces the spatial coherence of high spectrum image based on the H.264 quick intra-frame predictive encoding method of macro block flatness, utilize fractal image to reduce its Spectral correlation, and the coded residual data are carried out to dct transform, quantification, write code stream with CABAC entropy coding, fractal parameter also carries out CABAC entropy coding, has realized effective compression of high spectrum image.

Summary of the invention:

The present invention proposes a kind of based on method for compressing high spectrum image fractal and H.264.At first high spectrum image cubic data body is converted to the video of yuv format, sends into encoder.The I frame of high spectrum video, used the H.264 fast frame inner estimation mode method of discrimination based on the macro block flatness to carry out infra-frame prediction, to remove the spatial coherence of high spectrum image; The P frame of high spectrum video carries out block motion estimation/compensation fractal image, to remove the Spectral correlation of high spectrum image.At first the various tree piecemeals of P frame coded macroblocks find best matching blocks by the MSE criterion in reference frame, determine the iterated function system coefficient of each piece, then the rate distortion costs of more all dividing mode, using the piece dividing mode of rate distortion costs minimum as final interframe encoding mode, record final fractal parameter.Finally, the residual error data of I frame and P frame is through being encoded and write code stream by entropy coding CABAC after dct transform, quantification, and the fractal parameter of P frame also carries out CABAC entropy coding.

A kind of based on method for compressing high spectrum image fractal and H.264, it is characterized in that performing step is as follows:

Step 1: the video that the high spectrum image cube metadata is converted to yuv format;

Step 2: if (the first frame is necessary for the I frame to the I frame of high spectrum video, whether other frame can arrange is the I frame), the H.264 fast frame inner estimation mode method of discrimination of use based on the macro block flatness carries out Intra prediction mode selection, predicted macroblock, complete after the prediction of all the predictive frame that can obtain the I frame.Difference by primitive frame and predictive frame obtains the coding side residual frame, proceeds to the step 4 coded residual; If the P frame of high spectrum video, forward step 3 to;

Step 3: if the P frame of high spectrum video carries out fractal image to all macro blocks of present frame successively.In search window in reference frame, current macro is carried out to the piece coupling, the position of sub-block is as the initial search point of father's piece, and the size of father's piece is identical with the size of sub-block.Each macro block is carried out to the tree piecemeal, piecemeal can be divided into from big to small 16 * 16,16 * 8,8 * 16,8 * 8,8 * 8 down (inferior macroblock partition) can be divided into 8 * 4,4 * 8,4 * 4.At first the various tree piecemeals of coded macroblocks are found to best matching blocks by the MSE criterion, the iterated function system coefficient of determining each piece is the IFS coefficient; Then the rate distortion costs of more all dividing mode; Finally using the piece dividing mode of rate distortion costs minimum as final interframe encoding mode.Record final IFS coefficient, proceed to the reconstructed block that step 5 obtains this piece.If the reconstructed block that all macro blocks of present frame have all been encoded complete, all forms reconstructed image (being the reference frame of next frame), obtain the coding side residual image by the difference of original image and reconstructed image, forward the step 4 coded residual to.Described search window is the rectangular search zone in reference frame; Described IFS coefficient comprises the position skew of father's piece and sub-block, i.e. motion vector (x, y) and scale factor s, displacement factor o;

Step 4: the data of residual image are carried out Zig-Zag scanning on the one hand through the coefficient after DCT, quantification, then with entropy coding CABAC, encode and write code stream; Obtain the decoding end residual frame on the other hand after inverse quantization, anti-dct transform.Obtain reconstruction frames (being the reference frame of next frame) by predictive frame and decoding end residual frame sum.If the P frame also will carry out CABAC entropy coding to all IFS coefficients.Judge whether present frame is last frame, if last frame finishes coding; Otherwise, return to step 2 and continue to process the next frame image;

Step 5: by the IFS coefficient substitution decoding equation of preserving

r _i=s·d _i+o (1)

Calculate predicted value, by original block and the difference of prediction piece, obtain the coding side residual block, the coding side residual block obtains the decoding end residual block through dct transform, quantification, inverse quantization and anti-dct transform, then obtains reconstructed block by prediction piece and decoding end residual block sum.Proceed to next macro block of step 3 coding depth graphic sequence P frame.R in formula _ifor the pixel value of prediction piece, d _ipixel value for the corresponding father's piece of reference frame.

The video that in described step 1, the high spectrum image cube metadata is converted to yuv format comprises following two steps:

1) the high spectrum image data normalization is to [0,255] interval interior integer.The spectral reflectivity of high spectrum image data representation atural object, its bit depth is generally 16 bits/pixel, need to convert 8 bits/pixel to.Conversion method, at first the value of spectral reflectivity is normalized to [0,255] interval, then is converted into the integer without symbol 8 bits, and final numerical value is [0,255] interval integer;

2) using each wave band of high spectrum image data as a frame, the reflection coefficient after the normalization of each wave band is as the luminance elements (Y component) of this frame, and chromatic component Cb and Cr all are made as intermediate value 128.Circulation execution step 2) until all wave bands of high spectrum image be disposed.All frames form the video file of a yuv format.

H.264 fast frame inner estimation mode method of discrimination based on the macro block flatness in described step 2 comprises following five steps:

1) calculate current coding macro block two dimensional gray histogram, and record the maximum z of z axle in the two dimensional gray histogram _max.

The two dimensional gray histogram of image is based on the three-dimensional description figure of the Joint Distribution formation of image pixel gray scale and neighborhood of pixels gray average.If the size of image f (x, y) is M * N, gray scale is L.The neighborhood averaging gray level image that adopts k * k dot matrix smoothly to obtain by f (x, y) is g (x, y), is designated as

$g(x,y)=\frac{1}{k^2}\underset{m=-\frac{k}{2}}{\overset{\frac{k}{2}}{\mathrm{\Σ}}}\underset{n=-\frac{k}{2}}{\overset{\frac{k}{2}}{\mathrm{\Σ}}}f(x+m,y+n)---\left(2\right)$

In formula, 1≤x+m≤M, 1≤y+n≤N, k gets odd number.

The image of g (x, y) size is identical with f (x, y) with gray scale.Can form two tuples (i, j) by f (x, y) and g (x, y), each two tuple belongs to a point on two dimensional surface.The definition two-dimensional histogram is N (i, j), and presentation video f (x, y) grey scale pixel value is i, and the gray value of neighborhood of pixels average gray image is while being j, the number of times that two-dimensional points (i, j) occurs (i, j=0,1 ..., L-1).Z axle in the two dimensional gray histogram is N (i, j), represents the number of times that two-dimensional points (i, j) occurs;

2) capping threshold value Th _highwith lower threshold Th _low, Th _highand Th _lowbe the integer between [1,256];

3) if z _max>=Th _high, think that macro block is smooth, remove Intra4 * 4 patterns, select Intra16 * 16 patterns, and using the pattern of rate distortion expense minimum as the optimal frames inner estimation mode; Upgrade higher limit simultaneously

otherwise proceed to step 4);

4) if z _max≤ Th _low, think that the macro block details is abundant, remove Intra16 * 16 patterns, select Intra4 * 4 patterns, and using the pattern of rate distortion expense minimum as the optimal frames inner estimation mode; Upgrade lower limit simultaneously

otherwise proceed to step 5);

5) if Th _low<z _max<Th _high, think that macro block flatness feature is not remarkable, adopt H.264 standard intraframe prediction algorithm.

In described step 3, P frame macro block fractal image comprises following ten steps:

1) at first take 16 * 16 macro blocks carries out piece coupling to current macro as unit in the search window in reference frame.The position of sub-block is as the initial search point of father's piece, and the size of father's piece is identical with the size of sub-block.Travel through whole search window according to full way of search, find father's piece position of matching error MSE minimum.Recording this position iterated function system coefficient is IFS coefficient (comprising motion vector (x, y), scale factor s, displacement factor o), and calculates the rate distortion expense cost1 of this pattern;

2) fritter that is 2 16 * 8 by this 16 * 16 macroblock partitions, all according to the corresponding whole search window of full way of search traversal, find corresponding father's piece position of matching error MSE minimum to each fritter.Record the IFS coefficient of these two positions, and calculate the rate distortion expense sum cost2 of two pieces;

3) fritter that is 28 * 16 by this 16 * 16 macroblock partitions, all according to the corresponding whole search window of full way of search traversal, find corresponding father's piece position of matching error MSE minimum to each fritter.Record the IFS coefficient of these two positions, and calculate the rate distortion expense sum cost3 of two pieces;

4) fritter that is 48 * 8 by this 16 * 16 macroblock partitions, all according to the corresponding whole search window of full way of search traversal, find corresponding father's piece position of matching error MSE minimum to each fritter.Record the IFS coefficient of these 4 positions, and calculate respectively rate distortion expense cost4_1, cost4_2, cost4_3, cost4_4 and four rate distortion expense sum cost4 of these 4 pieces;

5) comparison step 1) to step 4) in rate distortion expense cost1, cost2, cost3, cost4, if minimum value is cost4, need further each 8 * 8 to be continued to divide, forward step 6 to); Otherwise the dividing mode of rate distortion expense minimum is the final P frame encoding mode of this macro block, retains corresponding IFS coefficient, finish;

6) at first be divided into 28 * 4 fritters by the 1st 8 * 8, each fritter, all according to the corresponding whole search window of full way of search traversal, is found to corresponding father's piece position of matching error MSE minimum.Record the IFS coefficient of these two positions, and calculate the rate distortion expense sum cost5 of these two pieces;

7) be divided into 24 * 8 fritters by the 1st 8 * 8, each fritter, all according to the corresponding whole search window of full way of search traversal, is found to corresponding father's piece position of matching error MSE minimum.Record the IFS coefficient of these two positions, and calculate the rate distortion expense sum cost6 of these two pieces;

8) be divided into 44 * 4 fritters by the 1st 8 * 8, each fritter, all according to the corresponding whole search window of full way of search traversal, is found to corresponding father's piece position of matching error MSE minimum.Record the IFS coefficient of these 4 positions, and calculate the rate distortion expense sum cost7 of these 4 pieces;

9) compare cost4_1, cost5, cost6, cost7, the piece dividing mode that minimum value wherein is corresponding, as these 8 * 8 final coding modes, retains corresponding IFS coefficient;

10) successively to each 8 * 8 of 16 * 16 macro blocks according to step 6) to step 9) and method select final coding mode, retain corresponding IFS coefficient.

Being calculated as follows of matching error MSE wherein:

$\mathrm{MSE}(R,D)=\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{[r_i-(s\&CenterDot;d_i+o)]}^2=\frac{1}{N}[\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{r}_i^2+s(s\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{d}_i^2-2\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{r}_i{d}_i+2o\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{d}_i^2)+o(N\&CenterDot;o-2\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{r}_i^2)]---\left(3\right)$

In formula, the calculating formula of s, o is respectively:

$s=\frac{[N\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{r}_i{d}_i-\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{r}_i\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{d}_i]}{[N\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{d}_i^2-{\left(\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{d}_i\right)}^2]}---\left(4\right)$

$o=\frac{1}{N}[\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{r}_i-s\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{d}_i]---\left(5\right)$

Wherein, r _ithe each point pixel value of current block, d _ithe each point pixel value of match block, the number of pixels that N is current block, s is scale factor, o is displacement factor.

Being calculated as follows of rate distortion expense:

cost=MSE(R,D)+λ _MODE·Rate(x,y,s,o,MODE|QP)(6)

Wherein, R represents current block, and D represents match block, and MODE means the macroblock partitions mode, and QP means quantization parameter, λ _mODEthe LaGrange parameter relevant with QP.

The advantage based on fractal and method for compressing high spectrum image H.264 proposed by the invention is:

(1) at first this method has carried out format conversion by the high spectrum image cube metadata, converts the video of yuv format to, the high spectrum image of so just can encoding as the coding ordinary video;

(2) this method adopts the H.264 fast frame intraprediction encoding based on the macro block flatness to reduce the spatial coherence existed in high spectrum image.Adopt the two dimensional gray histogram to characterize macro block flatness feature, and be provided with dual threshold, smooth macro block is removed to Intra4 * 4 patterns; The macro block abundant to details removed Intra16 * 16 patterns; The unconspicuous macro block of flatness feature is adopted to H.264 canonical algorithm.The dynamic iteration of threshold value is upgraded, and judgment threshold can be adjusted in real time according to the fluctuation of rebuilding the macro block flatness.Both guarantee the quality that high spectrum image decompresses, accelerated again compression speed;

(3) this method adopts fractal image to reduce the Spectral correlation existed in high spectrum image.Piece Matching supporting 16 * 16,16 * 8,8 * 16,8 * 8,8 * 4,4 * 8,4 * 4 is totally 7 kinds of piece dividing mode, at first every kind of piece dividing mode finds the match block of matching error MSE minimum, then the rate distortion costs of more all patterns, finally using the pattern of rate distortion costs minimum as optimum fractal image pattern.Improved the precision of fractal image.

The accompanying drawing explanation:

Fig. 1 is that the present invention is a kind of based on method for compressing high spectrum image flow chart fractal and H.264;

The 17th wave band grayscale sub-image that Fig. 2 is high spectrum image " DC Mall ";

Fig. 3 is that the present invention is a kind of based on the H.264 fast frame inner estimation mode method of discrimination flow chart based on the macro block flatness in fractal and method for compressing high spectrum image step 2 H.264

Fig. 4 is that the present invention is a kind of based on P frame macro block fractal image model selection flow chart in fractal and method for compressing high spectrum image step 3 H.264;

Fig. 5 is " DC Mall " high spectrum image the 17th wave band gray-scale map through the present invention is a kind of after rebuilding based on fractal and the compression of method for compressing high spectrum image H.264 again;

Fig. 6 (a) carries out the comparison diagram of the Y-PSNR of compressed encoding for the present invention is a kind of based on fractal and people's methods such as method for compressing high spectrum image H.264 and Lucana to " DCMall " high spectrum image;

Fig. 6 (b) is a kind of based on the comparison diagram fractal and bit number that people's methods such as method for compressing high spectrum image H.264 and Lucana are compressed " DCMall " high spectrum image for the present invention;

Fig. 6 (c) is a kind of based on the comparison diagram fractal and time that people's methods such as method for compressing high spectrum image H.264 and Lucana are compressed " DCMall " high spectrum image for the present invention.

Embodiment:

Below in conjunction with accompanying drawing, the inventive method is described in further detail.

The present invention proposes a kind of based on method for compressing high spectrum image fractal and H.264.At first high spectrum image cubic data body is converted to the video of yuv format, sends into encoder.The I frame of high spectrum video, used the H.264 fast frame inner estimation mode method of discrimination based on the macro block flatness to carry out infra-frame prediction, to remove the spatial coherence of high spectrum image; The P frame of high spectrum video carries out block motion estimation/compensation fractal image, to remove the Spectral correlation of high spectrum image.At first the various tree piecemeals of P frame coded macroblocks find best matching blocks by the MSE criterion in reference frame, determine the iterated function system coefficient of each piece, then the rate distortion costs of more all dividing mode, using the piece dividing mode of rate distortion costs minimum as final interframe encoding mode, record final fractal parameter.Finally, the residual error data of I frame and P frame is through being encoded and write code stream by entropy coding CABAC after dct transform, quantification, and the fractal parameter of P frame also carries out CABAC entropy coding.

As shown in Figure 1, a kind of based on method for compressing high spectrum image flow chart fractal and H.264.Shopping center, Washington DC (DC Mall) high spectrum image that spectral information technical applications center, Virginia provides of take is example, intercepts 256 * 1024 pixel subimages.Remove the wave band of 0.9 to the 1.4 μ m scope that its medium cloud coverage rate is larger, the final wave band number retained is 191.Accompanying drawing 2 is the 17th wave band grayscale sub-image of high spectrum image " DC Mall ".

Step 1: the video that the high spectrum image cube metadata is converted to yuv format.Concrete steps are as follows:

The 1st step, high spectrum image data normalization are to [0,255] interval interior integer.The spectral reflectivity of high spectrum image data representation atural object, its bit depth is generally 16 bits/pixel, need to convert 8 bits/pixel to.Conversion method, at first the value of spectral reflectivity is normalized to [0,255] interval, then is converted into the integer without symbol 8 bits, and final numerical value is [0,255] interval integer;

The 2nd step, using each wave band of high spectrum image data as a frame, the reflection coefficient after the normalization of each wave band is as the luminance elements (Y component) of this frame, chromatic component Cb and Cr all are made as intermediate value 128.Circulation is carried out the 2nd step until all wave bands of high spectrum image are disposed.All frames form the video file of a yuv format.

Step 2: judge whether the high spectrum frame of video of present encoding is that (the first frame is necessary for the I frame to the I frame, whether other frame can arrange is the I frame), if the I frame is used the H.264 fast frame inner estimation mode method of discrimination based on the macro block flatness to carry out Intra prediction mode selection, predicted macroblock.Accompanying drawing 3 is the H.264 fast frame inner estimation mode method of discrimination flow chart based on the macro block flatness.Specifically comprise following 5 steps:

The 1st step, calculate current coding macro block two dimensional gray histogram, and record the maximum z of z axle in the two dimensional gray histogram _max.

The two dimensional gray histogram of image is based on the three-dimensional description figure of the Joint Distribution formation of image pixel gray scale and neighborhood of pixels gray average.If the size of image f (x, y) is M * N, gray scale is L.The neighborhood averaging gray level image that adopts k * k dot matrix smoothly to obtain by f (x, y) is g (x, y), is designated as

$g(x,y)=\frac{1}{{\mathrm{<figure-callout\; id="2"\; label="k"\; filenames="HDA0000389629760000031.png"\; state="\{\{state\}\}">k</figure-callout>}}^2}\underset{m=-\frac{k}{2}}{\overset{\frac{\mathrm{<figure-callout\; id="2"\; label="k"\; filenames="HDA0000389629760000031.png"\; state="\{\{state\}\}">k</figure-callout>}}{2}}{\mathrm{\&Sigma;}}}\underset{n=-\frac{k}{2}}{\overset{\frac{k}{2}}{\mathrm{\&Sigma;}}}f(x+m,y+n)---\left(7\right)$

In formula, 1≤x+m≤M, 1≤ _y+ n≤N, k gets odd number.

The image of g (x, y) size is identical with f (x, y) with gray scale.Can form two tuples (i, j) by f (x, y) and g (x, y), each two tuple belongs to a point on two dimensional surface.The definition two-dimensional histogram is N (i, j), and presentation video f (x, y) grey scale pixel value is i, and the gray value of neighborhood of pixels average gray image is while being j, the number of times that two-dimensional points (i, j) occurs (i, j=0,1 ..., L-1).Z axle in the two dimensional gray histogram is N (i, j), represents the number of times that two-dimensional points (i, j) occurs;

The 2nd step, capping threshold value Th _high(empirical value is made as 175) and lower threshold Th _low(empirical value is made as 100), Th _highand Th _lowbe the integer between [1,256];

If the 3rd step z _max>=Th _high, think that macro block is smooth, remove Intra4 * 4 patterns, select Intra16 * 16 patterns, and using the pattern of rate distortion expense minimum as the optimal frames inner estimation mode; Upgrade higher limit simultaneously

otherwise proceed to the 4th step;

If the 4th step z _max≤ Th _low, think that the macro block details is abundant, remove Intra16 * 16 patterns, select Intra4 * 4 patterns, and using the pattern of rate distortion expense minimum as the optimal frames inner estimation mode; Upgrade lower limit simultaneously otherwise proceed to the 5th step;

If the 5th step Th _low<z _max<Th _high, think that macro block flatness feature is not remarkable, adopt H.264 standard intraframe prediction algorithm.

Complete after the prediction of all the predictive frame that can obtain the I frame.Difference by primitive frame and predictive frame obtains the coding side residual frame, proceeds to the step 4 coded residual; The P frame of high spectrum video, forward step 3 to if.

Step 3: the P frame of high spectrum video if, carry out fractal image to all macro blocks of present frame successively.In search window in reference frame, current macro is carried out to the piece coupling, the position of sub-block is as the initial search point of father's piece, and the size of father's piece is identical with the size of sub-block.Each macro block is carried out to the tree piecemeal, piecemeal can be divided into from big to small 16 * 16,16 * 8,8 * 16,8 * 8,8 * 8 down (inferior macroblock partition) can be divided into 8 * 4,4 * 8,4 * 4.At first the various tree piecemeals of coded macroblocks are found to best matching blocks by the MSE criterion, the iterated function system coefficient of determining each piece is the IFS coefficient; Then the rate distortion costs of more all dividing mode; Finally using the piece dividing mode of rate distortion costs minimum as final interframe encoding mode.Accompanying drawing 4 is P frame macro block fractal image model selection flow chart, specifically comprises following 10 steps:

The first step, at first take 16 * 16 macro blocks carries out piece coupling to current macro as unit in the search window in reference frame.The position of sub-block is as the initial search point of father's piece, and the size of father's piece is identical with the size of sub-block.Travel through whole search window according to full way of search, find father's piece position of matching error MSE minimum.Recording this position iterated function system coefficient is IFS coefficient (comprising motion vector (x, y), scale factor s, displacement factor o), and calculates the rate distortion expense cost1 of this pattern;

Second step, the fritter that is 2 16 * 8 by this 16 * 16 macroblock partitions, all according to the corresponding whole search window of full way of search traversal, find corresponding father's piece position of matching error MSE minimum to each fritter.Record the IFS coefficient of these two positions, and calculate the rate distortion expense sum cost2 of two pieces;

The 3rd step, the fritter that is 28 * 16 by this 16 * 16 macroblock partitions, all according to the corresponding whole search window of full way of search traversal, find corresponding father's piece position of matching error MSE minimum to each fritter.Record the IFS coefficient of these two positions, and calculate the rate distortion expense sum cost3 of two pieces;

The 4th step, the fritter that is 48 * 8 by this 16 * 16 macroblock partitions, all according to the corresponding whole search window of full way of search traversal, find corresponding father's piece position of matching error MSE minimum to each fritter.Record the IFS coefficient of these 4 positions, and calculate respectively rate distortion expense cost4_1, cost4_2, cost4_3, cost4_4 and four rate distortion expense sum cost4 of these 4 pieces;

Rate distortion expense cost1, cost2, cost3, cost4 in the 5th step, the comparison first step to the four steps, if minimum value is cost4, need further each 8 * 8 to be continued to divide, and forwards the 6th step to; Otherwise the dividing mode of rate distortion expense minimum is the final P frame encoding mode of this macro block, retains corresponding IFS coefficient, finish;

The 6th step, at first be divided into 28 * 4 fritters by the 1st 8 * 8, each fritter, all according to the corresponding whole search window of full way of search traversal, is found to corresponding father's piece position of matching error MSE minimum.Record the IFS coefficient of these two positions, and calculate the rate distortion expense sum cost5 of these two pieces;

The 7th step, be divided into 24 * 8 fritters by the 1st 8 * 8, each fritter, all according to the corresponding whole search window of full way of search traversal, is found to corresponding father's piece position of matching error MSE minimum.Record the IFS coefficient of these two positions, and calculate the rate distortion expense sum cost6 of these two pieces;

The 8th step, be divided into 44 * 4 fritters by the 1st 8 * 8, each fritter, all according to the corresponding whole search window of full way of search traversal, is found to corresponding father's piece position of matching error MSE minimum.Record the IFS coefficient of these 4 positions, and calculate the rate distortion expense sum cost7 of these 4 pieces;

The 9th step, comparison cost4_1, cost5, cost6, cost7, the piece dividing mode that minimum value wherein is corresponding, as these 8 * 8 final coding modes, retains corresponding IFS coefficient;

The tenth step, successively each 8 * 8 methods according to the 6th step to the nine steps of 16 * 16 macro blocks are selected to final coding mode, retain corresponding IFS coefficient, finish.

Being calculated as follows of matching error MSE wherein:

$\mathrm{MSE}(R,D)=\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{[r_i-(s\&CenterDot;d_i+o)]}^2=\frac{1}{N}[\underset{i=1}{\overset{\mathrm{<figure-callout\; id="2"\; label="N\; r\; i"\; filenames="HDA0000389629760000031.png"\; state="\{\{state\}\}">N</figure-callout>}}{\mathrm{\&Sigma;}}}{r}_i^{}{}_2^{}+s(s\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{d}_i^2-2\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{r}_i{d}_i+2o\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{d}_i^2)+o(N\&CenterDot;o-2\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{r}_i^2)]---\left(8\right)$

In formula, the calculating formula of s, o is respectively:

$s=\frac{[N\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{r}_i{d}_i-\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{r}_i\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{d}_i]}{[N\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{d}_i^2-{\left(\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{d}_i\right)}^2]}---\left(9\right)$

$o=\frac{1}{N}[\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{r}_i-s\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{d}_i]---\left(10\right)$

Wherein, r _ithe each point pixel value of current block, d _ithe each point pixel value of match block, the number of pixels that N is current block, s is scale factor, o is displacement factor.

Being calculated as follows of rate distortion expense:

cost=MSE(R,D)+λ _MODE·Rate(x,y,s,o,MODE|QP) (11)

Wherein, R represents current block, and D represents match block, and MODE means the macroblock partitions mode, and QP means quantization parameter, λ _mODEthe LaGrange parameter relevant with QP.

Obtain the final coding mode of P frame macro block by above-mentioned steps, record final IFS coefficient, proceed to the reconstructed block that step 5 obtains this piece.If the reconstructed block that all macro blocks of present frame have all been encoded complete, all forms reconstructed image (being the reference frame of next frame), obtain the coding side residual image by the difference of original image and reconstructed image, forward the step 4 coded residual to.Described search window is the rectangular search zone in reference frame; Described IFS coefficient comprises the position skew of father's piece and sub-block, i.e. motion vector (x, y) and scale factor s, displacement factor o.

Step 4: the data of residual image are carried out Zig-Zag scanning on the one hand through the coefficient after DCT, quantification, then with entropy coding CABAC, encode and write code stream; Obtain the decoding end residual frame on the other hand after inverse quantization, anti-dct transform.Obtain reconstruction frames (being the reference frame of next frame) by predictive frame and decoding end residual frame sum.If the P frame also will carry out CABAC entropy coding to all IFS coefficients.Judge whether present frame is last frame, if last frame finishes coding; Otherwise, return to step 2 and continue to process the next frame image.

Step 5: by the IFS coefficient substitution decoding equation of preserving

r _i=s·d _i+o (12)

Calculate predicted value, by original block and the difference of prediction piece, obtain the coding side residual block, the coding side residual block obtains the decoding end residual block through dct transform, quantification, inverse quantization and anti-dct transform, then obtains reconstructed block by prediction piece and decoding end residual block sum.Proceed to step 3 high next macro block of spectrum video P frame of encoding.R in formula _ifor the pixel value of prediction piece, d _ipixel value for the corresponding father's piece of reference frame.

This method is selected the implementation platform of Visual C++6.0 as described method, and CPU is Intel Core ^tM2Duo T8300,2.4GHz dominant frequency, memory size is 2G, 191 wave band cube metadatas of high spectrum image " DC Mall " have been carried out based on Compression of hyperspectral images experiment fractal and H.264, after converting the video file of yuv format to, it is 191 frames that totalframes equals total wave band number, every two field picture size is 256 * 1024 pixels, coded image group structure is IPPP ... I frame code period is 50, and the first frame is encoded to the I frame, afterwards every 1 I frame of 49 P frame codings, hunting zone is ± 7, and quantization parameter QP is 20.Accompanying drawing 5 is " DC Mall " high spectrum image the 17th wave band gray-scale map through the present invention is a kind of after rebuilding based on fractal and the compression of method for compressing high spectrum image H.264 again.

Adopt respectively people's methods such as Lucana (referring to Lucana Santos, Sebasti á n L ó pez, Gustavo M.Callic ó, et al.Performance evaluation of the is video coding standard for lossy hyperspectral image compression[J H.264/AVC] .IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, 2012, 5 (2): 451-461.) and the comparison diagram of the inventive method Y-PSNR that " DC Mall " high spectrum image is compressed as shown in accompanying drawing 6 (a), adopt respectively the comparison diagram of the bit number that people's methods such as Lucana and the inventive method compressed " DC Mall " high spectrum image as shown in accompanying drawing 6 (b), adopt respectively the comparison diagram of the time that people's methods such as Lucana and the inventive method compressed " DC Mall " high spectrum image as shown in accompanying drawing 6 (c).

Adopt respectively people's methods such as Lucana and the inventive method as shown in table 1 to 191 wave band compression performance mean value compare results of " DC Mall " high spectrum image.Wherein △ PSNR, △ bit rate, △ compression time are defined as follows:

△PSNR=PSNR _OURS-PSNR _Lucana (13)

△ bit rate=(bit rate _oURS-bit rate _lucana)/bit rate _lucana(14)

△ compression time=(compression time _oURS-compression time _lucana)/compression time _lucana(15)

Depth map sequence	△PSNR/dB	The △ bit rate	The △ compression time
				DC Mall	0.67	-90.50％	-90.00％

Table people's method such as 1Lucana and the inventive method are to 191 wave band compression performance mean value compares of " DC Mall " high spectrum image

From accompanying drawing 6 and table 1, can find out, the inventive method and traditional international video encoding standard H.264 standard testing model JM18.1 method are compared, Y-PSNR PSNR on average improves 0.67dB, the encoding code stream bit rate on average reduces by 17.73%, compression time decreased average 90.00%, show superior compression performance.This is because the inventive method on depth map sequence I frame coding, has adopted the improved H.264 quick intra-frame predictive encoding method based on the macro block flatness, removes the spatial coherence of high spectrum image; On depth map sequence P frame coding, adopted block motion estimation/compensation fractal coding algorithm, remove the Spectral correlation of high light image.

Claims (4)

1. one kind based on method for compressing high spectrum image fractal and H.264, it is characterized in that following steps:

Step 1: the video that the high spectrum image cube metadata is converted to yuv format;

Step 2: if (the first frame is necessary for the I frame to the I frame of high spectrum video, whether other frame can arrange is the I frame), the H.264 fast frame inner estimation mode method of discrimination of use based on the macro block flatness carries out Intra prediction mode selection, predicted macroblock, complete after the prediction of all the predictive frame that can obtain the I frame.Difference by primitive frame and predictive frame obtains the coding side residual frame, proceeds to the step 4 coded residual; If the P frame of high spectrum video, forward step 3 to;

Step 3: if the P frame of high spectrum video carries out fractal image to all macro blocks of present frame successively.In search window in reference frame, current macro is carried out to the piece coupling, the position of sub-block is as the initial search point of father's piece, and the size of father's piece is identical with the size of sub-block.Each macro block is carried out to the tree piecemeal, piecemeal can be divided into from big to small 16 * 16,16 * 8,8 * 16,8 * 8,8 * 8 down (inferior macroblock partition) can be divided into 8 * 4,4 * 8,4 * 4.At first the various tree piecemeals of coded macroblocks are found to best matching blocks by the MSE criterion, the iterated function system coefficient of determining each piece is the IFS coefficient; Then the rate distortion costs of more all dividing mode; Finally using the piece dividing mode of rate distortion costs minimum as final interframe encoding mode.Record final IFS coefficient, proceed to the reconstructed block that step 5 obtains this piece.If the reconstructed block that all macro blocks of present frame have all been encoded complete, all forms reconstructed image (being the reference frame of next frame), obtain the coding side residual image by the difference of original image and reconstructed image, forward the step 4 coded residual to.Described search window is the rectangular search zone in reference frame; Described IFS coefficient comprises the position skew of father's piece and sub-block, i.e. motion vector (x, y) and scale factor s, displacement factor o;

Step 4: the data of residual image are carried out Zig-Zag scanning on the one hand through the coefficient after DCT, quantification, then with entropy coding CABAC, encode and write code stream; Obtain the decoding end residual frame on the other hand after inverse quantization, anti-dct transform.Obtain reconstruction frames (being the reference frame of next frame) by predictive frame and decoding end residual frame sum.If the P frame also will carry out CABAC entropy coding to all IFS coefficients.Judge whether present frame is last frame, if last frame finishes coding; Otherwise, return to step 2 and continue to process the next frame image;

Step 5: by the IFS coefficient substitution decoding equation of preserving

r _i=s·d _i+o (1)

Calculate predicted value, by original block and the difference of prediction piece, obtain the coding side residual block, the coding side residual block obtains the decoding end residual block through dct transform, quantification, inverse quantization and anti-dct transform, then obtains reconstructed block by prediction piece and decoding end residual block sum.Proceed to next macro block of step 3 coding depth graphic sequence P frame.R in formula _ifor the pixel value of prediction piece, d _ipixel value for the corresponding father's piece of reference frame.

2. a kind of based on method for compressing high spectrum image fractal and H.264 according to claim 1, it is characterized in that: the video that in described step 1, the high spectrum image cube metadata is converted to yuv format comprises following two steps:

1) the high spectrum image data normalization is to [0,255] interval interior integer.The spectral reflectivity of high spectrum image data representation atural object, its bit depth is generally 16 bits/pixel, need to convert 8 bits/pixel to.Conversion method, at first the value of spectral reflectivity is normalized to [0,255] interval, then is converted into the integer without symbol 8 bits, and final numerical value is [0,255] interval integer;

2) using each wave band of high spectrum image data as a frame, the reflection coefficient after the normalization of each wave band is as the luminance elements (Y component) of this frame, and chromatic component Cb and Cr all are made as intermediate value 128.Circulation execution step 2) until all wave bands of high spectrum image be disposed.All frames form the video file of a yuv format.

3. a kind of based on method for compressing high spectrum image fractal and H.264 according to claim 1, it is characterized in that: the H.264 fast frame inner estimation mode method of discrimination based on the macro block flatness in described step 2 comprises following five steps:

1) calculate current coding macro block two dimensional gray histogram, and record the maximum z of z axle in the two dimensional gray histogram _max.

The two dimensional gray histogram of image is based on the three-dimensional description figure of the Joint Distribution formation of image pixel gray scale and neighborhood of pixels gray average.If the size of image f (x, y) is M * N, gray scale is L.The neighborhood averaging gray level image that adopts k * k dot matrix smoothly to obtain by f (x, y) is g (x, y), is designated as

$g(x,y)=\frac{1}{k^2}\underset{m=-\frac{k}{2}}{\overset{\frac{k}{2}}{\mathrm{\&Sigma;}}}\underset{n=-\frac{k}{2}}{\overset{\frac{k}{2}}{\mathrm{\&Sigma;}}}f(x+m,y+n)---\left(2\right)$

In formula, 1≤x+m≤M, 1≤y+n≤N, k gets odd number.

The image of g (x, y) size is identical with f (x, y) with gray scale.Can form two tuples (i, j) by f (x, y) and g (x, y), each two tuple belongs to a point on two dimensional surface.The definition two-dimensional histogram is N (i, j), and presentation video f (x, y) grey scale pixel value is i, and the gray value of neighborhood of pixels average gray image is while being j, the number of times that two-dimensional points (i, j) occurs (i, j=0,1 ..., L-1).Z axle in the two dimensional gray histogram is N (i, j), represents the number of times that two-dimensional points (i, j) occurs;

2) capping threshold value Th _highwith lower threshold Th _low, Th _highand Th _lowbe the integer between [1,256];

3) if z _max>=Th _high, think that macro block is smooth, remove Intra4 * 4 patterns, select Intra16 * 16 patterns, and using the pattern of rate distortion expense minimum as the optimal frames inner estimation mode; Upgrade higher limit simultaneously

otherwise proceed to step 4);

4) if z _max≤ Th _low, think that the macro block details is abundant, remove Intra16 * 16 patterns, select Intra4 * 4 patterns, and using the pattern of rate distortion expense minimum as the optimal frames inner estimation mode; Upgrade lower limit simultaneously

otherwise proceed to step 5);

5) if Th _low<z _max<Th _high, think that macro block flatness feature is not remarkable, adopt H.264 standard intraframe prediction algorithm.

4. a kind of based on method for compressing high spectrum image fractal and H.264 according to claim 1, it is characterized in that: in described step 3, P frame macro block fractal image comprises following ten steps:

1) at first take 16 * 16 macro blocks carries out piece coupling to current macro as unit in the search window in reference frame.The position of sub-block is as the initial search point of father's piece, and the size of father's piece is identical with the size of sub-block.Travel through whole search window according to full way of search, find father's piece position of matching error MSE minimum.Recording this position iterated function system coefficient is IFS coefficient (comprising motion vector (x, y), scale factor s, displacement factor o), and calculates the rate distortion expense cost1 of this pattern;

2) fritter that is 2 16 * 8 by this 16 * 16 macroblock partitions, all according to the corresponding whole search window of full way of search traversal, find corresponding father's piece position of matching error MSE minimum to each fritter.Record the IFS coefficient of these two positions, and calculate the rate distortion expense sum cost2 of two pieces;

3) fritter that is 28 * 16 by this 16 * 16 macroblock partitions, all according to the corresponding whole search window of full way of search traversal, find corresponding father's piece position of matching error MSE minimum to each fritter.Record the IFS coefficient of these two positions, and calculate the rate distortion expense sum cost3 of two pieces;

4) fritter that is 48 * 8 by this 16 * 16 macroblock partitions, all according to the corresponding whole search window of full way of search traversal, find corresponding father's piece position of matching error MSE minimum to each fritter.Record the IFS coefficient of these 4 positions, and calculate respectively rate distortion expense cost4_1, cost4_2, cost4_3, cost4_4 and four rate distortion expense sum cost4 of these 4 pieces;

5) comparison step 1) to step 4) in rate distortion expense cost1, cost2, cost3, cost4, if minimum value is cost4, need further each 8 * 8 to be continued to divide, forward step 6 to); Otherwise the dividing mode of rate distortion expense minimum is the final P frame encoding mode of this macro block, retains corresponding IFS coefficient, finish;

6) at first be divided into 28 * 4 fritters by the 1st 8 * 8, each fritter, all according to the corresponding whole search window of full way of search traversal, is found to corresponding father's piece position of matching error MSE minimum.Record the IFS coefficient of these two positions, and calculate the rate distortion expense sum cost5 of these two pieces;

7) be divided into 24 * 8 fritters by the 1st 8 * 8, each fritter, all according to the corresponding whole search window of full way of search traversal, is found to corresponding father's piece position of matching error MSE minimum.Record the IFS coefficient of these two positions, and calculate the rate distortion expense sum cost6 of these two pieces;

8) be divided into 44 * 4 fritters by the 1st 8 * 8, each fritter, all according to the corresponding whole search window of full way of search traversal, is found to corresponding father's piece position of matching error MSE minimum.Record the IFS coefficient of these 4 positions, and calculate the rate distortion expense sum cost7 of these 4 pieces;

9) compare cost4_1, cost5, cost6, cost7, the piece dividing mode that minimum value wherein is corresponding, as these 8 * 8 final coding modes, retains corresponding IFS coefficient;

10) successively to each 8 * 8 of 16 * 16 macro blocks according to step 6) to step 9) and method select final coding mode, retain corresponding IFS coefficient.

Being calculated as follows of matching error MSE wherein:

$\mathrm{MSE}(R,D)=\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{[r_i-(s\&CenterDot;d_i+o)]}^2=\frac{1}{N}[\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{r}_i^2+s(s\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{d}_i^2-2\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{r}_i{d}_i+2o\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{d}_i^2)+o(N\&CenterDot;o-2\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{r}_i^2)]---\left(3\right)$

In formula, the calculating formula of s, o is respectively:

$s=\frac{[N\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{r}_i{d}_i-\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{r}_i\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{d}_i]}{[N\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{d}_i^2-{\left(\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{d}_i\right)}^2]}---\left(4\right)$

$o=\frac{1}{N}[\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{r}_i-s\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{d}_i]---\left(5\right)$

Wherein, r _ithe each point pixel value of current block, d _ithe each point pixel value of match block, the number of pixels that N is current block, s is scale factor, o is displacement factor.

Being calculated as follows of rate distortion expense:

cost=MSE(R,D)+λ _MODE·Rate(x,y,s,o,MODE|QP) (6)

Wherein, R represents current block, and D represents match block, and MODE means the macroblock partitions mode, and QP means quantization parameter, λ _mODEthe LaGrange parameter relevant with QP.

Images