A&A 386, 1143-1152 (2002)
DOI: 10.1051/0004-6361:20020225

Lossy compression of scientific spacecraft data using wavelets. Application to the CASSINI spacecraft data compression

L. Belmon¹ - H. Benoit-Cattin² - A. Baskurt³ - J.-L. Bougeret¹

1 - Space Research Department of Paris Observatory (DESPA), Meudon, France
2 - CREATIS, Research Unit associated to CNRS (# 5515) and affiliated to INSERM, Lyon, France
3 - LIGIM, Computer Graphics, Image and Modelling Laboratory, Université Claude Bernard, Lyon, France

Received 28 July 2000 / Accepted 20 December 2001

Abstract
This paper presents a lossy coding scheme designed to be board on the Cassini spacecraft which has been launched in 1997 to study Saturn planet. To deal with specific time-frequency data and numerous constraints, our coding algorithm is based on wavelet decomposition associated to an adaptive bit allocation procedure. The performances of our approach are validated by an adapted scientific evaluation.

Key words: techniques: image processing

1 Introduction

The need for lossy compression in spacecraft data transmission has been growing in the last ten years consequently with the increasing number of space exploration missions. Table 1 gives some of the major spacecraft missions and the boarded compression approaches (Belmon 1998). One can see in this table that the boarded coding schemes are mainly based on Discrete Cosine Transform (DCT) or on Differential Pulse Code Modulation (DPCM) which can be found quite simple techniques compared to the last up to date techniques. That is justified by the need of robustness and by the severe hardware constraints imposed by spacecraft missions.

  

Table 1: Overview of several scientific spacecraft missions and the corresponding lossy coding approach on board.

$$\begin{tabular}{lcccc} \hline\hline Space Missions & \multico... ...t$ & & & DCT \\ [-1mm] (1998) & & & & \\ \hline \end{tabular}$$

The Cassini mission (Gibbs 1996; Kellog et al. 2001) is a scientific spacecraft which has been launched on 15 October 1997, to study Saturn and Titan and in particular their ionised and magnetised environment and the related radio and plasma phenomena. It is made in collaboration with the European Space Agency (ESA), the National Aeronautics and Space Agency (NASA) and some European national agencies including the French Centre National d'Études Spatiales (CNES).

The Radio and Plasma Waves Science instrument (RPWS) is one of the major instrument (Kellogg et al. 2001) on board the Cassini Spacecraft. It has been built by a consortium of laboratories led by the university of Iowa, and including CNET/CRPE France, NASA Goddard Space Flight Center, space research department of Paris Observatory (DESPA). The DESPA takes in charge the high frequency radio receiver of the RPWS, named Kronos. The Kronos analysis covers range from 4 kHz up to 16 MHz with numerous modes, many of them producing 450 bps bitrates of data which are beyond the allowed limit of the project.

Lossless coding algorithms are implemented on the Kronos instrument (Belmon 1998; Rice 1991), for compression ratios less than 2:1. For other needs, such as telemetry problem, the Kronos Project required an answer concerning the lossy compression at a global bit-rate of 2 bpp. Our work on the coding algorithm began in 1994. The decision to develop, implement and assess a whole compression method based on the wavelet theory validated by the Antonini et al.'s works (1992) has been taken at this time. As a consequence, our approach is similar to the existing work at this period and quite different from the more recent works like the SPIHT (Said et al. 1996) and EBCOT (Taubman 1999) algorithms.

At the beginning of the work (1994), the acceptability of a lossy compression scheme in the Kronos context has to be proved. That is the reason why the qualitative and quantitative quality of the images obtained with our lossy compression algorithm were evaluated. Saturn's Electrostatic Discharges (SED) statistics were analysed in terms of details such as the event duration, the occurrence distribution and the energy distribution (see Sect. 4.2). The images were also visually evaluated by 7 experts in order to appreciate the pattern and texture recognition in the compressed images. These results allow us to conclude that the lossy compression with our coding scheme is acceptable for the Kronos Project.

This paper presents the coding algorithm, the assessment procedure and its results. In part 2, we present the Kronos data we deal with, and we introduce the link between the scientific need (in terms of time and frequency resolution) and the compression constraints. The part 3 concerns the compression scheme adapted to Kronos constraints. It is based on a wavelet decomposition associated with an adaptive bit allocation procedure. In part 4, our coder performances are evaluated for several types of data. The data set used for both learning and testing is issued from past experiments (such as WIND spacecraft, Bougeret et al. 1995) when available, or simulated. We compare our approach with the well-known JPEG (Pennebaker et al. 1993) image coding standard for both numerical tests and dedicated tests such as physical parameters extracted from decompressed data.

	Figure 1: Example of a dynamic spectrum with kilometric signals from Wind spacecraft (1995), in 64-256 kHz frequency range in Earth's neighborhood.
Open with DEXTER

Notations

DPU: Data Processing Unit;
bpp: bit per pixel;
bps: bit per second;
KWC: Kronos Wavelet Coder;
JPEG: Joint Photographic Picture Expert Group;
RPWS: Radio and Plasma Waves Science instrument;
SNR: Signal to Noise Ratio;
SKR: Saturn's Kilometric Radiation;
SED: Saturn's Electrostatic Discharges;
SSR: Solid State Recorder;
$\langle A, B \rangle$ holds for the hermitian scalar product of function A and B.

2 Scientific framework and associated constraints

2.1 Scientific framework

The main scientific objectives of RPWS observations are the following:

a)

The study of Saturn's lightning, through the detection and analysis of short-lived sporadic broadband radio emissions associated to lightning (so-called SED for Saturn's Electrostatic Discharges, actually lightning-associated electromagnetic radio emissions) (Zarka et al. 1983), and their search at Titan (a moon of Saturn);

b)

The study of the aurora kilometric radiation from Saturn's magnetic poles (so-called SKR);

c)

The accurate localisation of the sources along aurora field lines to constrain the radiation mechanisms and the sources of accelerated particles involved in these emissions. It will be based on the monitoring of natural plasma waves and plasma density in Saturn's environment during Saturn's tour 2004-2007 and the remote studies of Jovian radio emissions during the Jupiter flyby in December 2000.

The first two objectives are the key ones for the RPWS Kronos instrument. They require non negligible amounts of data, due to the sporadicity of the studied events and to the requirements for direction-finding measurements (7 quantities must be measured for deriving one k-vector).

Data concerning SED and SKR will be generally obtained as a function of time and frequency (i.e. spectral and temporal variations), under the form of "pseudo-images'' of the intensity and/or polarisation versus time and frequency. Figure 1 gives an example of such a dynamic spectrum.

The Kronos instrument has been designed to work in three major modes (Zarka et al. 1996) adapted to the different post scientific analysis. The first one, so-called survey mode, allows the tracking of many events but is not dedicated to any special analysis. The time and frequency resolution are medium in order to cover the largest part of the spectrum. Frequency resolution is typically $\delta f/f =$10% in the 3.5-319 kHz frequency range, $\delta f = 25$ KHz up to 16 MHz. Time resolution, which is determined by integration time for the whole spectral range, is about 10 s. Kronos survey mode bitrate is 509 bps. In this paper, the dynamic spectra in survey mode come from WIND spacecraft experiment (Bougeret et al. 1995).

The second mode is dedicated to the lightning detection. The study of these phenomena requires from 1000 up to 2000 bps bitrate. For these wide band phenomena, only time resolution is crucial. For monitoring, fine structure study, and direction finding analysis, short sweep duration analysis (down to 0.3 s) are performed in the 100-16 KHz range. Frequency resolution is about 100-200 kHz. SED appear on dynamic spectra as short segments parallel to the frequency axis, because they are short lived broadband emissions, as detected by a swept-frequency analysis. A simple detection scheme based on subtraction of a boxcar average (on 3 pixels) and application of a threshold was applied successfully on Voyager observations of such events (Zarka et al. 1983). It can be used to detect a large number of events in the original and lossy compressed/decompressed dynamic spectrum, and to compare the properties of their statistical distributions. We simulate dynamic spectra (see Fig. 2) containing SED events (as well as interference) on the basis of the knowledge of SED properties (distribution of intensities, duration, and occurrences, and SED spectrum) derived from Voyager 1 and 2 observations (Zarka et al. 1983).

	Figure 2: Simulated dynamic spectrum with SED. The blanck lines are constant electrostatic parasites (horizontal) or sudden parasite discharges.
Open with DEXTER

The third mode is adapted to SKR direction finding. Time and frequency resolution depend on scientific goals: polarisation study, source location. Needed bitrate for SKR DF mode is about 1000 bps. In order to retrieve the direction of arrival of the wave at the studied frequency, 7 measurements performed successively by 2 pairs of antennas (+X, Z) and (-X, Z), including auto and cross-correlation coefficients, can be combined. Within the exploration of the whole frequency range, the SKR source will be seen as a position slowly variable with frequency, as it is distributed at the gyrofrequency along aurora field lines. It is thus possible to generate simulated SKR data corresponding to a fixed distribution of source positions. Performing the direction-finding analysis before and after lossy compression gives a quantitative measure of the compression error, and allows the comparison of different coding schemes (cf. Sect. 4.2).

The compression of Kronos data is mainly required to maximise the duty cycle (from about 100 ms up to several seconds depending on chosen operating mode) of RPWS measurements such as the lightning monitoring or the following of the motions of instantaneous aurora radio sources. Lossy compression should also help to reduce the loss of capabilities of the instrument in case of a major failure. As an example, a failure on the high gain antenna of the Galileo spacecraft has reduced the data rate by a factor 1000.

	Figure 3: Connections between the high frequency receiver Kronos (including the Digital Signal Processor) and the central Data Processing Unit (DPU) of the RPWS. On board storage is assumed by the Solid State Recorder (SSR). The storage capacity as well as the amount of data may vary in time.
Open with DEXTER

2.2 Coding constraints

Spacecraft boarded coding scheme are generally strongly constrained. For the Kronos instrument, the constraints which take into account the specific data structure and the instrument context are the following:

• Dynamic bitrate allocation: due to radio events sporadicity, the information is non-uniformly spanned over the spectral bandwidths. Consequently, bit allocation will be different from one band to another;
• Small and size-varying pictures: because of memory limitation and in order to minimise the loss due to transmission errors, the size of the images should be up to 32 acquisition time (spectra), giving pictures up to 32 channels $\times$ 32 spectra. The number of channels (frequency resolution) changes according to analysis mode requested by the DPWS data processing unit (DPU) (see Fig. 3);
• High dynamic range data, depending on the event's nature. So the algorithm should deal with non-stationary signals;
• Fast algorithm: because of real-time context, the data compression step should not consume more resources than the analysis one;
• Robust scheme: each picture has to be self-information relevant, i.e., we limit error propagation within each block of data. Accordingly, we deal with a coarse statistical modelling;
• Algorithm adapted to the ADSP 2100 Digital Signal Processor architecture of Kronos.

	Figure 4: Kronos Wavelet Coder overview. The same statistical Laplacian model is used for allocation and coding.
Open with DEXTER

3 Compression scheme

Figure 4 gives an overview of our coding scheme. It can be separated in 4 successive steps: the 2D discrete wavelet transform (DWT) which provides a subband decomposition of the original image, the bit allocation procedure which fixes a coding quality adapted to each subband, the scalar quantization and the entropy coder.

3.1 Wavelet decomposition

Wavelet theory has been investigated for a few decades in many signal processing applications (Rioul et al. 1991): series expansions, multiresolution analysis, singularities detection and subband coding for speech and images. A short overview on wavelet theory is given in Annex. On the other hand, coding techniques using filter banks to get a subband decomposition has been proposed for speech (Crochiere et al. 1976) and then for image coding (Woods et al. 1986). Different works (Daubechies 1988; Mallat 1989) have clearly establish the link between wavelet transform and subband decomposition. Indeed, wavelet transform can be easily implemented using filter banks (Vetterli et al. 1992) matching several constraints: Regularity, orthogonality, continuity, frequency selectivity, compact-support. Several families of wavelets filter banks satisfying these conditions can be found in the literature (Daubechies 1992; Antonini et al. 1992; Rioul et al. 1994).

	Figure 5: Two-dimensional dyadic successive approximations of a picture. a) The wavelet analysis is depicted as a tree decomposition with use of a low pass filter (G) and a high pass filter (H) and decimation operators (factor of 2). b) The wavelet coefficients organized in 7 sub-bands with different resoultions and directions.
Open with DEXTER

Wavelet based coder leads to interesting performances compared to standardised compression algorithms (Antonini et al. 1992; Said et al. 1996; Taubman 1999) and because of the good time and frequency localisation of the wavelet function, the wavelet transform is well suited for the non-stationnary signals such Kronos ones. Consequently, we decide to base our coder on the wavelet decomposition.

The filter bank used to implement the wavelet decomposition appears not to be so influent on the coder performances (Rioul 1993; Benoit-Cattin 1995). To implement the wavelet transform, regarding the small size of Kronos pictures, and to avoid border distortions, we use Daubechies' (Daubechies 1992) biorthogonal symmetric filters (tap length of 5 for low-pass filter and 3 for high-pass filter) with the symmetric signal extension (Brislawn 1996).

We apply an iterative 2D-dyadic decomposition on each picture with an analysis depth of 2. Then, we obtain 7 subbands (Fig. 5), 1 approximation (LL2) and 6 different horizontal (LH1, LH2), vertical (HL1, HL2) and diagonal (HH1, HH2) detail subbands. Each subband is quantized and coded respecting the bit allocation procedure.

3.2 Bit allocation

The Laplacian probability density function (Eq. (1)) is used (Mallat 1989) to model the wavelet coefficients amplitude distribution in each subband.

$$% p(x) = \frac{1}{\sqrt{2}\sigma_x} {\rm e}^{-\sqrt{2}\vert x\vert/_{\sigma_x}}\cdot$$ (1)

Such a model is taken for both the bit allocation and the entropy coding of quantized coefficients. In each subband i = 0...6, we compute the sample variances $\sigma^2_0$, $\sigma^2_1$, ..., $\sigma^2_6$ from which we derive the Laplacian model. Then, and because the transform is orthogonal, the overall reconstruction error $\sigma^2_{\rm r}$ is given by (Akansu et al. 1992):

$$% \sigma^2_q = \displaystyle \sum^{L}_{k = 0} w_k \cdot \sigma^2_{qk}.$$ (2)

Where wk are the weighting coefficients which are proportional to the subband size (Fig. 5).

For a given bit budget B, one can find the optimal bit allocation by minimising Eq. (2) using the Lagrange multipliers method (Everett 1963):

$$% b_1\! \! B + 1/2 \cdot \log_2 \left( \sigma^2_i \cdot \left( \left( \si... ...od^{L - 1}_{k = 0} \left( \sigma^2_k \right)^{1/2^{k + 1}} \right)^{-1} \right)$$ =

$$\! \! B + 1/2 \cdot \log_2 \left( \sigma^2_i \cdot P^{-1} \right).$$ = (3)

The quantity P is the geometric mean of the subband variances. It is often used to compute the coding gain, which measures the performance of the quantization method (Soman et al. 1993).

This bit allocation strategy may not be the optimal one, i.e. it is well known that wavelets coefficients are normally distributed only in the detail subbands. Moreover, we assume the high rate hypothesis which does not hold for our bitrates. Direct methods (Benoit-Cattin 1995; Ramchandran et al. 1993) lead to better bit allocation but they could not respect the current real-time context.

3.3 Quantization and coding

A residual correlation still exists between the approximation signal and the details signals. However in our approach, we propose to use "separate quantizers'' in order to have simplicity (independent quantizers) and flexibility (the possibility of using any kind of quantizer). In this implementation, we use the same uniform scalar quantizer.

As we do not deal with big pictures, we get some large magnitude coefficients that do not fit the statistical model. For these coefficients, the quantization step is larger.

Then, we use an Huffman coder (Huffman 1952) based on the Laplacian statistical model to code the coefficients included in the dynamic range [$- \sigma$, +$\sigma$]. For coefficients out of this range, we use a code with a length proportional to their magnitude.

Finally, the compressed picture includes a header giving the 7 rounded subband variance values that we need for decoding, and the coded coefficients for each sub-band.

4 Experimental results

The compression performance is the ability of the coding scheme to reduce the amount of data with a minimum loss of information. What is crucial in this context is that the dynamic spectra information is not directly relevant, i.e. the useful scientific parameters are extracted from the pictures with post-processing computations.

However, we first study the numeric performance in term of SNR. Indeed, SNR is a the distortion measure used for the compression management since we cannot envisage in-flight scientific processing, and SNR is a standard image compression measure.

Secondly, we evaluate the consequence of compression on the scientific parameters extraction on SED data (peak detection) and SKR data (source spatial localisation). Tests with "human experts'' are also performed to assess the limits of recognition of patterns or structures in the survey dynamic spectra at various compression rates.

These results have been obtained with two types of data: Survey dynamic spectra from WIND spacecraft (1993) in earth neighbourhood and simulated data for the SED detection and the source localisation.

	Figure 6: SNR performances for Kronos Wavelet Coder (KWC) and JPEG coder on 3 types of dynamic spectra form WIND TNR spacecraft. a) band A, mostly low frequency bursts, b) band C with noisy plasma ray, c) band E containing high dynamic range SKR.
Open with DEXTER

 

 

Table 2: SNR results for Kronos Wavelet Coder (KWC) and JPEG applied to SED simulated dynamic spectra.

	KWC	JPEG
Bit rate (bpp)	SNR (dB)	SNR (dB)
1.0	45.7	42.8
1.5	48.9	45.0
2.0	51.6	46.7

 

 

Table 3: SED occurrence analysis for reconstructed dynamic spectra. The first column gives the number of detected pixels whereas the second one gives the number of grouped pixels (events). The third column gives the percentage of missed events compared to the noisy spectrum (ii). The fourth and the last ones give the percentage of invented events in comparison to the noisy spectrum (ii) and with the original (i) spectrum respectively.

	SED pixels		SED events		Missed SED		invented SED		invented SED
					pixels/(ii)		pixels/(ii)		pixels/(i)
Original spectrum (i)	7838		5700		-		-		-
Noisy spectrum (ii)	3137		2427		-		-		0%
KWC/JPEG	99%	97%	97%	92%	4.0%	14.4%	2.6%	11.7%	1.1%	10.0%
2.0 bpp
KWC/JPEG	96%	106%	91%	94%	9.9%	17.8%	6.0%	22.7%	4.1%	20.5%
1.5 bpp
KWC/JPEG	101%	121%	84%	98%	18.0%	27.0%	18.5%	39.6%	16.2%	37.8%
1.0 bpp

4.1 Distortion results

For the survey dynamic spectra, Kronos Wavelet Coder (KWC) achieves better SNR results than JPEG for the planned average bitrate (2 bpp), especially on signals with large dynamic range (Figs. 6a and c). For low bitrates, the two methods are equivalent in term of SNR. Note that on a noisy band with poor information content (Fig. 6b), JPEG and KWC achieve the same performance at medium bitrate (2 bpp).

For the simulated SED spectra, KWC is much better than JPEG (Table 2). Indeed, when dealing with very sporadic events (1-3 pixels width), JPEG is not well suited to short-time variations characterisation because of the poor spatial localisation of the Discrete Cosine Transform (DCT). Consequently, the SED energy affects a large number of DCT coefficients and then the global quantization scheme used for JPEG produces lot of artefacts named "ringing'' effects around edges. At the opposite, the wavelet transform used in KWC packs peak energy in a very few coefficients and reduces this kind of distortion.

4.2 Scientific analysis

a) SED detection and characterisation

We use a simple SED detection based on a gradient algorithm (Zarka 1985). Then from the retained events, we derive statistics on the events duration, the occurrence distribution and energy distribution of the SED. We apply this analysis on 3 categories of spectrum:

i)

clean dynamic spectrum which contains only SED (uniform background);

ii)

noisy dynamic spectrum which is the clean spectrum with some added noise and signals simulating perturbations from electrostatic environment;

iii)

reconstructed noisy dynamic spectrum which is the lossy compressed-decompressed version of (ii), with the 2 coding methods, namely KWC and JPEG.

Table 3 shows that at the expected bitrate (2 bpp), KWC overcomes JPEG for the SED numbering and localisation. For example, 4.4% of missed SED pixels for KWC instead of 14.4% for JPEG, as well as 1.1% of invented SED pixels compared to 10% for JPEG. Furthermore, as the bitrate decreases, this difference increases.

At low bit rates, regarding only the SED events column, it seems that JPEG offers a better SED events numbering. But at the same bitrates, KWC gives less missed and invented SED pixels. In fact, JPEG induces false events numbering, which are detected by the missed and invented SED post-processing procedures which take into account the events localisation.

Finally, one notes that at high bitrates (2 bpp), wavelet compression performs a good denoising on the noisy spectrum. Indeed, only half of the invented pixels relative to (ii) does exist in the clean spectrum (i).

	Figure 7: Occurrence distribution defined as the population for a given delay (expressed in terms of number of samples) between two SED occurrences, for simulated SED spectrum versus JPEG and KWC coded spectra at a) 2.0 bpp and b) 1.0 bpp.
Open with DEXTER

Regarding the SED pixels detection at low bitrates, it appears that JPEG compression tends to add SED pixels (up to 121%). This could be explained by the ringing effect artefacts generated by JPEG around SED pixels. Such artefacts are attenuated (101%) when using wavelet coding, because of good spatial localisation of wavelet analysis.

Figure 7 shows the occurrence distribution in terms of number of samples between 2 events (related to the time between 2 SED), for a noisy simulated spectrum versus JPEG and KWC compressed spectra. It appears that occurrence distribution slope is better preserved with KWC than with JPEG. If we look more precisely at the distribution slops, we can see that most of the distortion is concentrated in the short-time events (delay between 2 events less than 30 pixels) as a consequence of the ringing effect artefact. At 1 bpp (Fig. 7b), JPEG compression introduces so many isolated edges that distribution profile is lost.

 

 

Table 4: Results of a visual test for dynamic spectra compressed-decompressed with JPEG and Kronos Wavelet Coder (KWC). Individual answers are binary (JPEG or KWC), i.e. the percentage of coder choice is computed from the 7 answers. Highlighted cells represent ambiguous results, where the method choice disagrees with SNR criteria.

	Kilometric signals				Small variations signal				Large variations signal
					(Plasma ray)				(Low Frequency radio bursts)
	JPEG		KWC		JPEG		KWC		JPEG		KWC
	SNR			SNR	SNR			SNR	SNR			SNR
Bit rate	(dB)			(dB)	(dB)			(dB)	(dB)			(dB)
1.0 bpp	40.7	29%	71%	40.6	46.3	0%	100%	46.6	41.2	15%	85%	40.4
1.5 bpp	43.1	43%	57%	42.6	48.4	50%	50%	48.2	43.6	43%	57%	43.4
2.0 bpp	45.1	43%	57%	45.2	50.2	25%	75%	52.6	45.2	29%	71%	45.1
2.5 bpp	46.6	57%	43%	48.0	51.9	25%	75%	54.1	46.7	71%	29%	47.3
3.0 bpp	47.9	0%	100%	50.1	53.5	50%	50%	56.6	48.0	29%	71%	49.5

b) Human experts evaluation of survey dynamic spectra

We also study the consequence of lossy compression on visual appreciation from experts, that is, patterns and textures recognition, for both JPEG and KWC methods. For three categories of radio signal and for 30 dynamic spectra in each category, 7 experts choose the best picture from 2 proposed compressed ones (blind test). Table 4 gives the percentage of chosen images compressed with JPEG and with KWC for different bit rates.

We notice that SNR criteria doesn't fit so well average choice, especially for bit rates greater than 2 bpp. Globally, this table shows that KWC psychovisually overcomes JPEG on every kind of spectra. More precisely, as the bitrate decreases, experts choose massively KWC.

At low bit rate, fine structures in kilometric signal (at 1 bpp) and low frequency bursts (at 1 bpp) are affected by JPEG blocking effect as shown in Fig. 8b.

At high bit rate, excepted for kilometric signals, it seems that experts cannot choose between JPEG and KWC (See Fig. 9). Indeed for such bitrates, the image quality is so high that expert could note make real differences between the two compressed images and as a consequence, their choice tend to be randomised.

c) Direction Finding (DF)

From autocorrelation and inter-correlation spectra obtained from 2 pairs of antennae, one can know precisely the source localisation (defined by azimuth and longitude) and Stockes parameters (power and polarisation of the radio source) (Zarka 1985). The desired precision is above 1 degree for both angles. Dispersion of intercorrelation data is very large, so these spectra are much harder to compress than autocorrelation spectra. Consequently, high compression ratio can not be achieved for such data. Thus we address the correlation between angle measurement quality and lossy coding using KWC which offers the best SNR performances. Table 5 presents the error variance of angle recovery results due to lossy coding at high bitrate. Autocorrelation and intercorrelation spectra are compressed at 2 bpp and 4 bpp respectively. One can see that KWC affects angle computation, even for good SNR performances. It shows how much angle computation method (singular value decomposition) is sensitive to input data precision.

	Figure 8: Example of survey dynamic spectra used for human expert evaluation. Original spectra a) and its coded version at low bit rate (1 bpp) using JPEG b) and KWC c) respectively. This spectrum is a small part (32 $\times$ 50 pixels) of a WIND dynamic spectrum (band E).
Open with DEXTER

	Figure 9: Example of survey dynamic spectra used for human expert evaluation. Original spectra a) and its coded version at high bit rate (2 bpp) using JPEG b) and KWC c) respectively.
Open with DEXTER

 

 

Table 5: Error variances on azimuth ($\theta$) and longitude ($\phi$) computations observed on 100 different spectra after high bit-rate KWC compression.

Data	Average bitrate	SNR	Longitude variance	Azimuth variance
	(bpp)	(dB)	$\sigma_{\Delta \phi}$	$\sigma_{\Delta \theta}$
Autocorrelation	2	42
			8.9$^\circ$	13.7$^\circ$
Intercorrelation	4	47

5 Conclusion

In this paper, we present a lossy coding algorithm designed for the coding of Kronos instrument data on board of the Cassini Spacecraft. It is based on wavelet decomposition associated to an adaptive bit allocation procedure and scalar quantization. We compare our coder with the JPEG standard using distortion criteria, automatic specific scientific data analysis as well as human expert analysis. Our wavelet coder suits very well the sporadicity of the expected radio signals and overcomes JPEG standard in most of the cases.

The results obtained show that on board lossy compression of scientific data is possible under limitations on specific data such as the dynamic spectra used for source localisation.

6 Appendix: Wavelet transform basis

We briefly review the basis of wavelet analysis. Theory of wavelets and further developments may be found in (Mallat 1989; Daubechies 1990).

A wavelet function is an oscillating function of the form:

$$% \Psi_{a,b} (t) = a^{-1/2} \cdot \Psi \left( \frac{t-b}{a} \right),~ a\;>\;0\cdot$$ (4)

The function $\Psi$ called the mother wavelet generates a family of functions $\Psi_{a,b}$ with a-dilatations and b-translations. Then the associated wavelet transform of a signal x(t), is:

$$% CWT (f)(a,b) = \frac{1}{\vert a\vert}\int\limits^{+\infty}\... ...(t) \cdot \Psi^\ast \left( \frac{t-b}{a} \right) {\rm d}t\cdot$$ (5)

It is also called continuous wavelet transform (CWT) because it derives a continuous decomposition of x(t) on an infinity of wavelets $\Psi_{a,b}$. In the discrete case, with writing a = 2^j and $b = k\cdot2^j$, it as been shown (Meyer 1986) that one can find functions $\Psi$ such that $\Psi_{j,k}$ is an orthonormal basis of L²(R). Then the discrete wavelet transform (DWT) of a continuous signal x(t) is:

$$\big\langle f,\Psi_{j,k} \big\rangle \Psi_{j,k}$$ DWT (f)(j,k) =

$$\int\limits^{+\infty}\limits_{-\infty}f(t) \cdot \Psi^\ast_{j,k}(t) \cdot {\rm d}t \cdot \Psi_{j,k}.$$ = (6)

Where quantities $c_{j,k} = \big\langle f,\Psi_{j,k} \big\rangle$ are the wavelet coefficients.

From the point of view of signal processing, the notion of wavelet analysis has many connections with the multiresolution analysis (Mallat 1989). Such an analysis is characterized by a subset of enclosed subspaces V_j, $j \in Z$, each subspace V_j being a basis of functions supported in $[-\pi/2^j,\; +\pi/2^j]$. There exists a function $\varphi(t)$ of V₀ such that the family of functions $\varphi(t-k)$ is a basis of V₀.

	Figure 10: Ideal dyadic partitioning of frequency space in approximation subspaces $^{V_j \subset V_{j-1}}$. At resolution j, the subspace Vj includes the average information while subspace Wj includes the detail information.
Open with DEXTER

Then, it has been shown (Meyer 1986) that there is a family of functions $\varphi_{j,k}$ that forms an orthonormal basis of V_j. That is, with dilating and translating a single function $\varphi$ called scaling function, one can define a set of orthonormal bases for the approximations subspaces V_j. The associated wavelet is then defined as:

$$% \Psi(t) = \sum^{+\infty}_{n = -\infty} (-1)^n \cdot \alpha_{n+1} \cdot \varphi (2t + n).$$ (7)

Now, the operator

$$% A_j \cdot X = \sum^{+\infty}_{k = -\infty} \big\langle x(t)... ...varphi_{j,k} (t) \big\rangle, j \in Z \;{\rm fixed},\; k \in Z$$ (8)

performs the projection of x(t) on subspace V_j spanned by the functions $\varphi_{j,k} \cdot A_j \cdot X$ is the average information about x(t) at resolution 2^-j, or the approximation of x(t) at resolution 2^j+1.

Consequently, the operator

$$% D_j \cdot X = \sum^{+\infty}_{k = -\infty} \big\langle x(t),~ \Psi_{j,k} (t) \big\rangle, j \in Z \;{\rm fixed},\; k \in Z$$ (9)

is the projection on the subspace W_j spanned by the functions $\Psi_{j,k} \cdot D_i \cdot X$ is the detailed information of x(t) at resolution 2j (Fig. 10).

One performs a m-DWT with these operators by retaining only the detail signal at different resolution (i=0...m), and the average signal at the maximum resolution $2^{\rm m}$. For discrete-time signals, this can be implemented by using a dyadic decomposition tree and a pair of appropriated half-band lowpass and highpass filters (Mallat 1989).

For image analysis, 2D-DWT can be achieved with computing the detail and the average projections called subbands for the raws (horizontal analysis), then for the columns (vertical analysis) at each resolution (Fig. 5) shows the 2-dimensionnal (2D) wavelet analysis implementation, as used in our image compression scheme. Vertical and horizontal analysis are successively performed for each scale j. Then at each analysis scale j four different signal are derived from the input picture:

• diagonal detail: $HH = \big\{ D_j \cdot \big\{ D_j \cdot X \big\}_{{\rm columns}} \big\}_{{\rm raws}}$
  
• horizontal detail: $LH = \big\{ D_j \cdot \big\{ A_j \cdot X \big\}_{{\rm columns}} \big\}_{{\rm raws}}$
  
• vertical detail: $HL = \big\{ A_j \cdot \big\{ D_j \cdot X \big\}_{{\rm columns}} \big\}_{{\rm raws}}$
  
• average: $LL = \big\{ A_j \cdot \big\{ A_j \cdot X \big\}_{{\rm columns}} \big\}_{{\rm raws}}\cdot$

Such an image dyadic decomposition has been implemented in KWC in order to get 7 subbands (see Fig. 5): one average subband at scale 2, three detail subbands at scale 2, and three detail sub-bands at scale 1.

References

 
• Akansu, A. N., & Haddad, R. A. 1992, Multiresolution signal decomposition (Academic Press Inc., New York) In the text  
• Antonini, M., Barlaud, M., Mathieu, P., & Daubechies, I. 1992, IEEE Trans. Image Proc., 1, 205 In the text  
• Baskurt, A., Odet, C., & Goutte, R. 1994, Signal Proc.: Image Comm., 6, 479  
• Belmon, L. 1998, Ph.D. Thesis, University of Paris XI Orsay In the text  
• Benoit-Cattin, H. 1995, Ph.D. Thesis, INSA Lyon, 258 In the text  
• Bougeret, J. L., et al. 1995, Space Sci. Rev., 71, 231 In the text NASA ADS  
• Brislawn, C. M. 1996, Appl. Comput. Harm. Anal., 3, 337 In the text  
• Crochiere, R. E., Webber, S. A., & Flanagan, J. L. 1976, Bell Syst. Tech. J., 55, 1069 In the text  
• Daubechies, I. 1988, Comm. Pure Appl. Math., 41, 909 In the text  
• Daubechies, I. 1990, IEEE Trans. Inf. Theory, 36, 961 In the text NASA ADS  
• Daubechies, I. 1992, Ten lectures on wavelets (Society for Industrial and Applied Mathematics, Philadelphia, Pennsylvania), 357 In the text  
• Everett, H. 1963, Oper. Res., 11, 399 In the text  
• Gurnett, D. A., et al. 2000, Space Sci. Rev., to be published  
• Huffman, D. A. 1952, IRE, 40, 1098 In the text  
• Mallat, S. G. 1989, IEEE Trans. Pattern Anal. Mach. Intell., 11, 674 In the text NASA ADS  
• Meyer, Y. 1986, Sémin. Bourbaki, 662 In the text  
• Pennebaber, W. B., & Mitchell, J. L. 1993, JPEG still image compression standard (Van Nostrand Reinhold, New York), 638  
• Ramchandran, K., & Vetterli, M. 1993, IEEE Trans. Image Proc., 2, 160 In the text  
• Rice, R. F. 1991, Jet Propulsion Lab. Publ., 3 In the text  
• Rioul, O., & Vetterli, M. 1991, IEEE Sign. Proc. Mag., 14 In the text  
• Rioul, O. 1993, in Proc. of IEEE ICASSP-93, Mineapolis, 550 In the text  
• Rioul, O., & Duhamel, P. 1994, IEEE Trans. Circ. Syst.-II, 41, 550 In the text  
• Said, A., & Pearlman, W. A. 1996, IEEE Trans. Circ. Syst. Video Techn., 6, 243 In the text  
• Soman, A. K., & Vaidyanathan, P. P. 1993, IEEE Trans. Sign. Proc., 41 In the text  
• Taubman, D. 1999, in Proc. IEEE, ICIP99 In the text  
• Vetterli, M., & Herley, C. 1992, IEEE Trans. Sign. Proc., 40, 2207 In the text  
• Woods, J. W., & O'Neil, S. D. 1986, IEEE Trans. Acous. Speech Sig. Proc., ASSP-4, 1278 In the text  
• Zarka, P., & Pedersen, B. M. 1983, J. Geophys. Res., 88, 9007 In the text NASA ADS  
• Zarka, P. 1985, Icarus, 61, 508 In the text NASA ADS  
• Zarka, P., & Manning, R. 1996, Obs. Paris Meudon: Internal Report In the text