Study on Lossless Compression of Multispectral Remote Sensing Images

Study on Lossless Compression of Multispectral Remote Sensing Images

Kai Xu,Hiroshi Okayama
Remote Sensing and Image Research Center
Chiba University I-33 Yayoi-Cho, Chiba 260 Japan

Koji Kajiwara
Institute of Industrial Science, University of Tokyo 7-22 Roppongi
Minato-ku, Tokyo 106, Japan

Abstract
This paper describes a partially recoverable lossless data compression scheme for multispectral remote sensing images of raster format. The method takes both the intrepid and interurban redundancy into account and deals with intraband redundancy by converting each band of image to its differential image. The integrand correlations are handled by pre-coding the combination image (pairs of differential values) of two differential images. Then a line-by-line two bands of original images. The coding system is designed to maintain a sufficient line index information which ensures the possibility of reconstruction any part of the original image with no need to decode the whole compressed data ached. Thus we can access the compressed image data just as we do it with the original one. the method was applied to SPOT HRV and LANDSAT TM multispectral images an overall compression ratio up to 3.3. and 2.3 respectively was achieved.

Introduction
Multispectral satellite images, such as SPOT, LANDSAT TM or MSS images, are rich resources for earth observations and environmental analysis. However, satellite image data are normally of enormous sizes. For example, a single SPOT HRV multispectral (3 bands) scene takes up to a total 30 Mbytes of storage facilities. The steady growth implies a significant increase of storage requirements for such sensors, in the near future. from the point of view of image interpreters or image analysts who work on earth observation data, a global earth observation image data-base is highly desired. Thus, an easily accessible data format of digital images with low storage requirements should contribute somewhat to an efficient data-base management scheme.

For single band images, although a rather high sampling rate and quantization is required to obtain subjectively high-quality pictures, statistical analysis of images indicates that the entropy of most pictures is substantially less than what these high sampling and quantization rates requie (Tzay Y. Young, 1986). This implies that a rather high redundancy exists among the neighboring pixels due to the scanning mechanism of most sensors. For the multispectral remotely sensed images, interband redundancy should also exist as well because of the similar reflectance received from the objects for neighboring bands. And also, from an intuitive point of view, the fact that the objects in all bands share the same boundaries should lead to a considerable integrand correlations. These make the compression of multispectral images desirable and possible for storage reduction.

Generally, image compression schemes can be divided in two categories. One is the lossy approach by which a high compression ratio can be achieved when allowing a certain limited distortion in the reconstructed image provide that it leaves no apparent difficulty four visual interpretations. For this case the transform coding or a hybrid of DPCM (Differential Pulse Code Modulation) with transform coding scheme has been proved as practical and is widely utilized (Tzay Y. Young, 1986). The alternation is usually considered when a large amount of data is to be stored for future retrieval and analysis. For the lossless image compressions a reversible DPCM or its adaptive extension followed by a variable length encoder is widely used.

Our consideration for a lossless image compression concentrates on two objectives. One is to store the remotely sensed images in a lower storage volume for future retrieval. The other is to make the compressed file as easily accessible as the uncompressed original ones. Specifically, the compressed data must be partially recoverable so that we can extract any line, any block of nay band of original images with no need to extend the whole compressed files ahead. Thus, when the decoder is incorporated into some software developed in our center the users can just treat the compressed images the same as uncompressed ones. This is accomplished by maintaining a sufficient extra information as indices. Experimental results show that the overall ( including all the index tables ) compression ratio for SPOT HRV ( full scene 5369x3000 pixels/band x 3bands, 8 bits/pixel) and LANDSAT TM ( screen size 512x400.

Pixels/band x 4bands, 8bits/pixel) up to 3.3. respectively is achievable.

System design

1 Differential image
The most commonly used for of predictive coding system is DPCM, which uses a linear predictor to predict the incoming signal based on a weighted sum of the adjacent elements. The DPCM generates a set of uncorrelated signals, referred to as differential signals, from the correlated data, and the processing could either be 2-dimensional or I- dimensional. Since our goal is to make the compressed file partially decodable, we here choose the I-dimensional DPCM to reduce the redundancy degree among the neighboring pixels on each scan line. Thus given the image size (N pixels?line x M lines ) the inverse processor can easily reproduce any original data segment on any line. The prediction of pixel at ith column and jth line can be expressed like (1),

The optimal performance of the prediction ( in Mean Square Error sense ) depends on the optimal choose of the prediction coefficients, the aks in (tzay Y. Yung, 1986 ). Some statistical studies showed (Robert E Rice, 1971 ) that only the ui-lj contributed significantly to the prediction of ^u_ij. Thus (1) can be simplified to as

which is simply a differencing operation. We view the new image consists of the u_0js and the prediction error e_ijs as the differential image. Our next step is to code this differential image instead of the original, since the differential image, generally, has more compact dynamic range than that of the original image. Here a truncating operation is necessary to be carried out to the e_ijs before the coding procedure because, in most cases, remote sensing image data are usually quantized to an 8 bits/pixel level so that the value of e_ijs may vary from - 255 to + 255, which requires 9 bits for every e_ij, and this is not convenient for implementation on computers. We truncate the e_ijs according to (3) where the T(e_ij) is the truncated value of e_ij, and the table h_ij holds the transformed value of e_ijs which are out of the range [-126,+126] and the end-of-line (EOL) marks with value 129. The T (e_ij)s are to be recorded in a file as the actual differential image for further coding. And the h_ijs must be stored as a reference table for decoding.

So the inverse operation will simply be the inverse of T(e_ij) with a table looking up into hij whenever the value of 253 or 254 is encountered in the coded file. A line - by-line inverse of the truncating operation can be accomplished by counting the EOLs in h_ij.

2 Combination image
By converting the individual band image to its differential image the intraband correlations are significantly reduced. As explained above there should exist some correlations between adjacent bands, termed as the interband correlations. It is quite natural to consider a similar processing of interband prediction for corresponding pixels in two

Differential images then coding the prediction error image as proposed by (Soon Ing Yann, 1991) or simply the band-difference image. However, results of simulations showed that this did not work satisfactorily; Experiments indicated that simply taking the band-difference did not work well for correlation reduction. There must exist an optimal prediction coefficient for two band images, but to calculate it needs some intensive statistical computations especially when the data are of large volume

To reduce the correlations between two images we here try a method of pre-coding the combination value of corresponding pixel intensity values in two differential images. If the data file containing the combination values is viewed as an image file we call it the combination image. To guarantee the reversibility the 2-dimensional histogram of the underlying differential images must be stored. We allocate a 16 bits memory for each combination value in which the upper byte holds the pixel value of band one differential image and the lower byte for the corresponding pixel in

Band two. This can be written as (4) , where c_ij is the combination value for pixels at ith line jth column and uk_ij stands for the pixel intensity value of band k differential image.

Clearly, when the u^k_ijs have 256 possible values (gray levels ) the c_ijs will have 65536 possible values and to code the c_ijs directly will actually achieve nothing. Experiments showed that when (4) was applied to multispectral remote sensing images after the pixel differencing operation described above the number of different c_ij values was much less than 65536. This lead us to think of assigning a particular number to each different combination value was much less representation. A simple way of doing this is to first sort the 2-dimensional histogram calculated from the underlying differential images according to the decreasing of occurrence of each pixel combination and then use the index number of sort number to represent different combinations. For example, let L be the number of different pixel combinations so each of them can be represented by a integer in (0,1,2,…,L-1) where 0 stands for the most frequently occurred combination and 1 for the second and so on. By far a converting table is constructed and now we can scene the two differential images simultaneously and convert each pixel combination according to the converting table to generate a new data file. The reverse processing is also guaranteed as long as the sorted 2-D histogram of length L is stored. Consider each element as a point in a plane with the X-Y coordinates be its lower byte and upper byte respectively then this new data file will be the points covering a rectangular area in the X-Y plane with in 0£Y£L/255 and 0£X£255. It means that a histogram of this data calculated one every byte will have its distribution cover a full histogram of underlying data, such a data with a wide spreading histogram will not ( empirically) guarantee a good achieve a data file with more compact histogram. For L different combinations, the X-Y coordinates of L points in a square of area bounded by , ( as described above since L is much less than 65536 so that must5 be considerably less than 255) and highly peaked occurrence concentrating on the first few integers. Also, it is not necessary to store this converting table because the inverse processor can reconstruct it from the stored 2D sorted histogram.

Examinations one some 2-D histograms of remote sensing images indicates that they are usually of the shape as shown in Fig.1. A break-point can usually be found at B which separates the curve in two parts. The first B most frequent combinations usually have a much compact shape with fast decreasing and generally cover a 60% to 80% of the area under the curve. It means that for two differential images 60% to 80% of the pixel combinations are of the first B combinations. Because we assign the (x,y) values or the code words symmetrically, that is the sequence of (x,y) would be (0,0), (0,1), ( 1,0) , (1,1), (0,2), (2,0),…(0. [ ÖL] + l}, ([ÖL]+l,0), (1,2),(2,1(,… , and so forth, it is clear that if we first assign the (x,y) coordinates for the points in the square area bounded by 0£X£[ÖB] + 1 and 0£Y£[ÖB]+1 and then the rest of points in the square of the area 0£X£[ÖL]+1 and 0£Y£[ÖL]+1 by the same rule , and the resultant data file must have a more higher peak on the first few integers, shown in Fig. 2. This is the actual processing in our system. And the processing is applied to every tow differential images generated by 2.1, and the break-point B is the index of the 2-D histogram array where its second derivative is minimized ( here the position of B is selected empirically based on the observations of many sorted histograms of combination values.

3 Line-by-line Huffman coding
By far both the intraband correlation and the interband correlation have been significantly reduced. The resultant pre-coded combination image has relatively low entropy than the original images. Because we want the final compressed file to be partially decodable a line-by-line Huffman coding is carried out to the data file results from 2.2, The line-by-line Huffman coding in our system uses only one Huffman tree, calculated from the byte-based histogram of the data file, for each line and records the number of output bytes of each coded line to a line index file. Thus the line-by-line decoding is simply implemented by counting the number of the output bytes from the decoder up to the number recorded in the line index file for each individual line.

Results, Conclusions and discussions
We applied this method to the LANDSAT TM multispectral remote sensing image which consists of 4 bands each with a size of 512x400 bytes and a sampling rate of 8 bits / pixel. First, the differential images of all 4 bands were generated and then every combination image of two of the differential images was tested with the coding scheme to see the possible compression ratio. The results are listed in Table 1 and Table 2. A same experiment was also carried out to the SPOT HRV images of 3 bands each had 5369 x 3000 bytes with sampling rate of 8 bits/pixel, and the results are shown in Table 3 and Table 4. And the partly decoding processing were verified. As shown in the Tables, for correlated bands even the requirement of partly decidability of the compressed file was satisfied by holding several index tables a overall compression ratio up to 2.3 and 3.3 was still achieved, which could be concluded that the compression scheme described above is effective for bands of neighboring wavelengths. For the individual band such as band 4 of LANDSAT TM which has the last correlation with other bands, the procedure of taking the combination image of this band with any other band will not achieve a satisfactory compression ratio, as shown in the Tables, and thus is not encouraged. For this case, we can simply code the differential image of this band separately by the line-by-line Huffman coding to preserve that the whole compressed files are partly decodable. In the event that there are only odd bands of images, such as the 3 bands of SPOT HRV images, a simple way is to code the combination image of band 1 and 2 differential images, which are most correlated, and leave the differential image of band 3 to be coded separately, an other possible suggestion is that we can first insert the third band into band 1 and band 2 to produce the reproduced-band 1 and reproduced-band 2 images, then code the combination image of the differential images generated from the reproduced - 3 will be inserted into the upper part of band 1 such that they are arranged with the lines of band 1 alternatively, which results the upper part of the reproduced-band 1 contains line 1 of band 1, line 1 of band 3, line 2 of band 1, line 3 of band 3, and so forth. And the even lines of band 3 will be inserted into band 2 in the same way. We have tried this method on the SPOT HRV images and the results as shown in Table 3 and Table 4 with the entry named Reproduced Band 1 and Reproduced Band 2, showed a compression ratio of 3.1 was gained.

The lossless compression approach for multispectral remote sensing image described above has shown, experimentally, its usefulness or effectiveness. Since it is impossible to try this approach on all the images, a mathematical outline of this scheme must be rationally derived ad a theoretical estimation of its effectiveness and limitation on certainly modeled images should also be carried out. These leave for the topics of our further studies.

References

Robert F. Rice James R. Plaunt, "Adaptive Variable-length Coding for Efficient Compression of Spacecraft Television Data", IEEE Trans. On Communication Technology Vol. COM-19, No. 6. December 1971.

Tzay Y. Young, King-Sun Fu, "Handbook of Pattern Recognition and Image Processing," Academic Press, Inc. 1986.

Soon Ing Yann, "Lossless Compression for Remotely Sensed Images", Proc. of THe 12th ACRS , Vol. II, 1991.

Table 1. Statistics of LANDSAT TM image ( 4 bands )and their differential images

Sensor Name:LANDSAT TM	File Suze(bytes)	Entropy (bits/pixel)	Size of Line Index Table for Differential Image (bytes)	Enimpy of Differential Image (bits/pixel)
BAND 1	204800	4.627166	400	4.147409
BAND 2	204800	4.604271	400	3.454869
BAND 3	204800	5.167487	400	4.129300
BAND 4	204800	6.136996	400	4.824956

Table 2. Statistics of combinations between differential image of LANDSAT TM imageries

Table 3. Statistics of SPOT HRV image ( 3 bands and 2 reproduced - bands ) and their differential images

Table 4. Statistics of Combinations between differential images of SPOT HRV imageries

Fig A sorted 2-D histogram shows the occurrences of different pixel combinations of two differtial images

Fig 2. The corresponding histogram of the coded combination image according to the sorted 2-D histogram in Fig. 1. Shows the occurrences of dirrerent bytes(integers form 0 to 255)