1.

INTRODUCTION

Satellite images may be modelled as follows in the Fourier domain:

where P and I refer respectively to the landscape and image and ^ the Fourier transform and MTF to the Modulation Transfer Function, Fourier transform of the point spread function.

Image restoration aims at both reversing MTF effects and keeping noise at low level. While restoration techniques usually assume a white constant noise and do not take into account compression artefacts, the restoration technique described hereunder manages both compression noise and signal dependent instrumental noise.

Performances of this new restoration chain have been assessed both visually and quantitatively on space products, including pansharpened images and Digital Surface Models.

2.

INSTRUMENT NOISE RESTITUTION

Past and current French High-Resolution Earth observation satellites such as SPOT5 and PLEIADES-HR (PHR) used onboard compression algorithms to optimize the mission capacity in terms of number of imaged scenes. These algorithms imposed a fixed compression rate to ease mission planning and mass-memory handling. To fulfill the required image quality, the fixed targeted bit rate was therefore calibrated on complex scenes such as urban areas. As explained in [1], fixed bit rate compression leads to under-optimized system with inhomogeneous image quality where some parts of the image are insufficiently encoded and present poor image quality (see Figure 1). Obviously, it is possible to improve image quality by increasing the targeted bit rate, but this brings also a major drawback: reducing the system imaging capacity. Another solution is to perform adequate on-ground post-processing of the compressed image to reduce compression artefacts described above.

Figure 1 –

Example of compression artefacts in a PLEIADES image of the Grand Canal in Venice compressed at 2.86 bits/pix on panchromatic (spatial resolution 0.70m) and 3.33 bits/pix in multispectral bands (blue, green, red and near infrared bands, spatial resolution 2.80m).

Thus, the instrumental noise restitution proposed by CNES [2] is a novelty that provides visual quality improvement, while preserving the overall signal to noise ratio. The instrumental noise, modified by the effect of quantization, is restituted where needed all over the image and consequently, compression artefacts such as butterflies and smoothed areas strongly reduced (see Figure 2 vs Figure 1).

Figure 2 –

Effect of instrument noise restitution for PLEIADES image presented Figure 1.

This instrumental noise restitution technique proposed by CNES is extensively suited for transform-based compression algorithms such as Discrete Cosine or Discrete Wavelet Transforms (resp. DCT and DWT). It can be performed during the decoding process or after the decoding process, the former being more efficient because it prevents further direct and inverse transform stages.

The general idea is to compare each quantized transformed coefficient to the local expected instrumental noise level computed in the transform domain. If the coefficient is equivalent or greater than the local instrumental noise, it remains unchanged. If it is lower, it is replaced by the local instrumental noise level.

The cornerstone of this technique is the instrumental noise model which is well-known and monitored in CNES systems: the noise variance depends linearly on the local signal intensity (see below), the latter being approximated by using the DC coefficient of the 8x8 block in case of DCT and DWT, following CCSDS convention [1]. For a wavelet based compressor as used on PHR satellite, the DC coefficient of the block can be used for a coarse estimation of the local mean value. For a better estimation, the low-pass wavelet coefficient of each DWT decomposition level can be used.

The pseudo-code of the “thresholding” step of transformed coefficients is the following:

The two multiplication factors k1 and k2 are used to tune noise thresholding and restitution levels. They are equal to 1 for a general use.

If the instrumental noise restitution is applied after decompression, the thresholding step should still be performed between the direct and the inverse transforms.

Another major advantage of this technique is that efficient denoising algorithms can be used during on-ground processing, because the structured compression noise is superseded by instrumental noise.

3.

DENOISING

Most state-of-the-art denoising methods only consider images corrupted with constant additive white noise, whereas in reality acquisition devices introduce signal dependent noise [3]. Its standard deviation may be modelled as

where:

Parameters A and B are well known in satellite imagery; they are measured on ground before launch on uniform target and in flight with for example steady-mode acquisition principle [4] [5] or on uniform landscape depending on the satellite agility. Signal to noise ratio, defined as

is thus a monotonous increasing function of S(x,y). When increasing the resolution, it becomes more and more difficult to maintain high SNR values. Low radiance landscapes (shadows, winter acquisitions…) are zones where denoising is particularly needed.

3.1

Variance Stabilizing Transform (VST)

In order to be able to use denoising methods using constant additive noise assumption, a variance stabilizing transform may be applied to the noisy image.

Anscombe transform [6][7] is particularly suited to the linear noise variance model. It is defined by the following formula:

where S(x,y) stands for the noisy pixel value, (A,B) representing the noise model parameters expressed in the same units as S. For instance, using electron units, B=1 and A is the darkness noise standard deviation. Constant 3/8 coefficient may usually be neglected when compared to the darkness noise standard deviation (65e- for PHR panchromatic band).

4.

DECONVOLUTION

Once noise has been sufficiently removed, deblurring merely consists in multiplying the image Fourier transform by a function close to MTF inverse. High resolution systems have usually low MTF values at Nyquist frequency mainly because of constraints on telescope diameter but also in order to minimize aliasing artefacts.

MTF is thus close to zero at Nyquist and 1/MTF could take very large values, which would mean an oscillating deconvolution filter inducing ringing artefacts around transients. Such a filter would also dramatically increase the noise remaining after the denoising stage. Compromise has thus to be done between deblurring and artefacts. Wiener-Tikhonov technique [11] provides a well suited filter. It is defined by the following formula:

where MTF* stands for the conjugate complex of MTF.

The parameter s monitors the weight of the regularization term which penalizes the L2-norm of the gradient (Figure 3). High λ values means sharp deconvolution while low λ values lead to smoother images.

Figure 3

Wiener-Tikhonov filters according to regularization parameter

5.

RESTORATION PIPELINE

As shown in

Figure 4, the basic stages of the restoration pipeline starts with instrument noise restitution. Whenever willing to restore a compressed image, any compression artefacts have to be masked, otherwise the classical restoration process (denoising and/or deconvolution) may enhance the artefacts. Then, the denoising process is applied before deconvolution, so that the noise is not being colored by deconvolution if it was done reversely (as previously made on the PHR satellites, high resolution Earth observation satellites [3]). Some stages are optional, such as VST and VST-1 if the instrument noise is already independent on the signal intensity; and such as deconvolution if the MTF of the input image is sufficiently high. Typically, on the PHR images, whose MTF for the panchromatic PAN (resp. multispectral XS) band is around 0.15 (resp. 0.32), the deconvolution step will be only applied on the panchromatic band, whereas the multispectral bands skip it since those images are sharp enough.

Figure 4

The restoration pipeline.

This restoration pipeline has been patented [12] and tested successfully on PHR and SPOT6-7 images, currently in-flight Earth observation satellites. Note that this pipeline has only few parameters, some of them are perfectly known because they only depend on the acquisition device (MTF and noise model are calibrated). Some of them are unknown: here they only depend on the chosen denoising algorithm. For NL-Bayes, we tune 6 parameters (Table 1) on a bunch of images, typically not more than 3 or 4 test images because satellite images are sufficiently big [>106 pixels each]. These parameters correspond to the two-stage NL-Bayes algorithm (stage 1: oracle computation, stage 2: final denoising computation): the patch half size (w), the search window half size (k), the smoothing parameter (β), the number of kept patches (N1) of stage 1 and the maximum number of kept patches (N2) of stage 2, (τ) the minimum threshold to determine similar patches during the second stage.

Table 1

NL-Bayes parameters for PHR images. The PAN band has its own parameters. The 4 multispectral bands share the same parameters since they have similar image quality.

Tuning has to be done for each spectral band of each satellite. Once these parameters are set, they apply to all images (Figure 5) which make it perfectly suitable for an on-ground automatic chain, such as an operational ground segment. This is on-going work for the PHR ground segment [13], which will be hopefully equipped with this new restoration pipeline early 2018.

Figure 5

Restoration results on two extracts of a PAN PHR image: restored images with the current ground segment (left), and with the new restoration pipeline (right). One can see on uniform areas that noise has been well removed, compression artefacts have disappeared and edges are better preserved on the right.

6.

APPLICATIONS

6.1

Fusion

Most of high resolution satellite systems acquire one high resolution panchromatic spectral band with a large spectral bandwidth and several multispectral bands with a narrower spectral bandwidth and a coarser resolution. Typically, the resolution ratio between panchromatic and spectral bands is equal to four. This is mandatory to comply with download bit rate constraints and on-board storage capacity. However, pansharpening techniques implemented within the ground segment allows to compute high resolution multispectral images.

Up to now, only panchromatic image was restored but we found out that restoring the multispectral bands significantly improves the final pansharpened products (Figure 6). The improvement of PléiadesHR products is shown Figure 8.

Figure 6

Pansharpening without (left) and with (right) restoration

Figure 7

DSM results without restoration (left), and with the new restoration pipeline (right). One can see on uniform areas that altimetric noise has been well removed and new 3D points have been computed on the right. The ground truth is represented below and the scale is in meter.

Figure 8

PleiadesHR detail example - raw image (left) improved image (right)

7.

CONCLUSION

In this paper, we have presented a new restoration pipeline dedicated to compressed and noisy satellite images. This pipeline starts with an instrumental noise restitution step, superseding structured compression noise by instrumental noise. Then, efficient denoising techniques are applied, together with a Variance Stabilization Transform. Finally, a deconvolution step is performed to take into account the satellite’s MTF. An optimal parametrization of this pipeline is possible since the instrument response (noise and MTF) and compression algorithm are known. It provides better restoration than state of the art algorithms on satellite imagery. This improvement in image quality is shown on two kinds of application: pansharpening and 3D extraction.