**Miguel Hernández-Cabronero 1,***∗***, Jordi Portell 2,3, Ian Blanes 1,4 and Joan Serra-Sagristà 1,4**


Received: 30 July 2020; Accepted: 10 September 2020; Published: 11 September 2020

**Abstract:** The capacity of the downlink channel is a major bottleneck for applications based on remote sensing hyperspectral imagery (HSI). Data compression is an essential tool to maximize the amount of HSI scenes that can be retrieved on the ground. At the same time, energy and hardware constraints of spaceborne devices impose limitations on the complexity of practical compression algorithms. To avoid any distortion in the analysis of the HSI data, only lossless compression is considered in this study. This work aims at finding the most advantageous compression–complexity trade-off within the state of the art in HSI compression. To do so, a novel comparison of the most competitive spectral decorrelation approaches combined with the best performing low-complexity compressors of the state is presented. Compression performance and execution time results are obtained for a set of 47 HSI scenes produced by 14 different sensors in real remote sensing missions. Assuming only a limited amount of energy is available, obtained data suggest that the FAPEC algorithm yields the best trade-off. When compared to the CCSDS 123.0-B-2 standard, FAPEC is 5.0 times faster and its compressed data rates are on average within 16% of the CCSDS standard. In scenarios where energy constraints can be relaxed, CCSDS 123.0-B-2 yields the best average compression results of all evaluated methods.

**Keywords:** multispectral; hyperspectral; CCSDS; FAPEC; data compression; transform

### **1. Introduction**

Hyperspectral imaging (HSI) is one of the main components of modern remote sensing [1]. In recent years, new public and private organizations have joined the remote sensing ecosystem [2]. As a result, the number of deployed HSI sensors has grown significantly, especially for Earth Observation tasks, and this trend is predicted to continue in the future [3]. Parallel to this increment of active sensors, substantial efforts are being devoted to improving the analysis of scenes obtained with HSI sensors [4]. Classification, unmixing, and segmentation are among the problems currently more actively researched [4–15]. Thanks to the increased availability of HSI data and the accuracy of their interpretation, many commercial and scientific applications are becoming feasible, and attainable for a larger base of users. Some of the most widespread applications include resource management, e.g., in agriculture, mining. and forestry; geological structure observation; and disaster monitoring and response coordination [2,16–18].

One of the main bottlenecks when retrieving HSI information from spaceborne sensors is the downlink channel capacity [19]. This bottleneck is emphasized by the large size of HSI data. Retrieved scenes typically contain from tens to hundreds of spectral bands. For comparison, traditional color images contain only three bands. The inclusion of more than three bands is needed for a more detailed observation of the electromagnetic spectrum, but they multiply the number of samples produced by the sensor. Furthermore, newly deployed sensors also tend to exhibit larger spatial resolutions, i.e., a higher number of samples per spectral band. As a result, the latest sensors—such as NASA's HyspIRI—can produce up to 5 Terabytes of data per day [20], i.e., significantly more than can be transmitted through existing space links [19].

In order to maximize the effective capacity of downlink channels, compression is routinely employed to reduce the volumes of data to be transmitted [19]. In particular, lossless compression allows data volume reduction, while at the same time enabling perfect reconstruction of the scenes. Thus, decompressed scenes will produce the exact same result as the original scenes when subject to any deterministic analysis algorithm. This includes most classification, unmixing, and segmentation algorithms [4].

This paper investigates two aspects of lossless compression specific to HSI remote sensing—spectral band redundancy—and the performance of low-complexity algorithms. Redundancy among spectral bands in HSI is typically higher than in traditional color images. Thus, to attain competitive coding efficiency, it is paramount to effectively exploit this redundancy. In turn, low complexity in the compression pipeline, including the spectral decorrelation stage, is needed due to the limitations in spendable energy and available hardware in spaceborne sensors [19,21]. These two aspects of compression must be jointly considered so that obtained conclusions are relevant in practical remote sensing HSI scenarios.

Hereafter, a selection of the best performing, reversible spectral decorrelation transforms is studied, including the integer wavelet transform (IWT [22]), the pairwise-orthogonal transform (POT [23]), and the regression wavelet analysis (RWA [24]) transform. These transforms are evaluated in combination with some of the most competitive low-complexity compression algorithms including the Fully Adaptive Prediction Error Coder (FAPEC [25]), the CCSDS 122.1-B-1 standard [26], and the well-known JPEG-LS [27] algorithm. For completeness, the aforementioned transform-compressor combinations are also compared to CCSDS 123.0-B-2 [28] and M-CALIC [29], which exploit redundancy across bands by means of prediction. This work extends upon our previous work [30] by considering a wider selection of spectral decorrelation transforms, compression algorithms, and test scenes. This provides a more complete analysis of the state-of-the-art and increases the significance of the obtained results Note that, even though hardware implementations have been described for some of the aforementioned algorithms (see, e.g., in [31–33]), this study focuses on software implementations. This approach is consistent with the current trend in Earth Observation, in which constellations of small satellites and cubesats (mainly from the private sector) are being deployed [2]. These constellations often employ commercial-off-the-shelf hardware [3], as opposed to the custom hardware implementations traditionally used in larger satellites launched by space agencies. By analyzing the efficiency of software implementations, results become more relevant to the foreseeable future of Earth Observation.

The rest of this paper is structured as follows. Section 2 contains a summary of the spectral decorrelation transforms and compression methods considered for this study. Section 3 presents and analyzes experimental results based on a comprehensive corpus of real hyperspectral scenes. Conclusions are drawn in Section 4.

#### **2. Methods and Materials**

State-of-the-art data compression algorithms often consist of two conceptually separated parts: data decorrelation and entropy coding. Data decorrelation exploits the fact that neighboring pixels—either within a band or across several of them—tend to have similar values. Decorrelated samples normally exhibit less variability, i.e., lower entropy rates. The entropy coding stage then compacts decorrelated

samples in a reversible way. The output compressed data size depends directly on the entropy of those samples, hence the importance of effective decorrelation.

This work focuses on compression pipelines in which decorrelation is divided into two consecutive stages: spectral and spatial decorrelation. Typically, scenes are first processed along the spectral axis, although some compressors invert this order or skip some of the stages. The resulting transformed bands are often processed separately, i.e., they are spatially decorrelated and entropy coded. The compressed data obtained for each band is then packed together to form the final compressed bitstream, which can also be output to a file. A diagram of this general design is provided in Figure 1. This approach allows combining and comparing different methods for each stage in order to get the best compression results. At the same time, it helps to obtain an overall low-complexity algorithm. As discussed in the introduction, both aspects are critical for HSI remote sensing.

The remainder of this section provides some background on the spectral decorrelation methods considered in this study in Section 2.1. The spatial decorrelation and entropy coding algorithms with which they are combined are described in Section 2.2. In order to obtain representative results, a comprehensive corpus of hyperspectral scenes has been selected. This corpus is presented in Section 2.3.

**Figure 1.** Main stages of transform-based hyperspectral scene compression.

#### *2.1. Spectral Decorrelation Methods*

Spectral decorrelation transforms are applied to co-located samples across different spectral bands. This process is repeated until all spatial positions are decorrelated. Three spectral transforms have been selected for this work due to their competitive decorrelation-complexity trade-offs: the integer wavelet transform (IWT), the pairwise-orthogonal transform (POT), and the regression wavelet analysis (RWA) transform. These are briefly described next in chronological order of publication.

The IWT [22] is a well-known discrete wavelet transform (DWT). It is a linear, multilevel filter that produces integer coefficients. This allows reversible transformation, which is required to attain lossless compression. In each level, the IWT filter separates the input signal into so-called low-frequency and high-frequency sub-bands. After that, the low frequencies of one level are the input of the next level. In this work, five levels of the LeGall 5/3 integer wavelet decomposition are employed except for the FAPEC compressor described below, which applies three levels of 9/7-M integer wavelet decomposition. One important property of the IWT is that it can be implemented very efficiently using so-called lifting steps [34]. This approach significantly reduces the number of arithmetic operations needed per input sample, thus enhancing throughput.

The POT [23] is a low-complexity, reversible approximation of the Karhunen-Loève Transform (KLT). The KLT provides perfect decorrelation, i.e., zero covariance between different spectral bands, which normally results in very high compression performance. However, the KLT is generally considered too computationally intensive for low-complexity scenarios. In turn, the POT offers much lower computational cost, while trying to reach the compression performance gains of the KLT. To achieve this goal, the POT applies a divide-and-conquer strategy that decorrelates band pairs in a multilevel structure. This is opposed to the KLT, which decorrelates all spectral bands simultaneously.

The RWA [24] is a reversible transform that operates in two steps. First, a reversible wavelet transform is applied along the spectral axis. The Haar wavelet is used for this purpose due to its low computational cost and reasonable decorrelating power. In the second step, pyramidal nonlinear estimators are used to predict coefficients in high-frequency sub-bands using coefficients from low-frequency sub-bands.

#### *2.2. Low-Complexity Compression Methods*

In this study, three compression algorithms are applied after each of the spectral decorrelation transforms described in Section 2.1: CCSDS 122.1-B-1, FAPEC, and JPEG-LS. Two additional compressors are considered in the study, CCSDS 123.0-B-2 and M-CALIC, which do not exactly follow the structure depicted in Figure 1. Instead, they apply their own spectral decorrelation methods, so they are not combined with any of the aforementioned spectral transforms. FAPEC also follows a slightly different approach, as described later in this section.

The CCSDS 122.1-B-1 [26] standard (hereafter CCSDS 122.1) was proposed by the Consultative Committee for Space Data Systems (CCSDS) as an extension to the CCSDS 122.0-B-2 standard [28] (All CCSDS standards are available free of charge at https://public.ccsds.org/Publications/default. aspx). The CCSDS agglutinates the most important space agencies in the world, as well as other important space industry actors, with the goal of creating state-of-the-art standards for space communications. The CCSDS 122 standard follows closely the structure shown in Figure 1. Both the IWT and the POT are available in the spectral preprocessing stage, as well as a no-transform option that allows combination with other transforms including the RWA. In the spatial decorrelation stage, CCSDS 122.1 employs a 5-level 2D 5/3 IWT, which is then followed by a bitplane encoder. Note that the advantages of bitplane encoding are most evident for lossy compression, although this coding regime is not explored in this study.

FAPEC [25] contains a highly efficient entropy coder that works in a block-wise manner. In each block, the probability density function of its samples is computed, and used to determine which of the available coding modes to use. High-entropy blocks are output without further processing to reduce computational complexity and avoid data expansion. Samples in medium-entropy blocks are encoded using variable-length codes. Several code tables are available depending on the statistics of each block. Low-entropy blocks, where most samples are in {−1, 0, +1}, are divided in groups of six samples and compressed using fixed-length codes. Finally, runs of 5 or more zeros are encoded using run-length coding. The block-wise code selection of FAPEC is highly robust to outliers and other noise sources in the data, while being suitable for the complexity restrictions of space hardware. This entropy coder is applied to the samples after decorrelation. When no spectral transform is selected, decorrelation consists in 3 levels of 2D 9/7-M IWT. When the POT or the RWA are selected, they are applied before the same 2D IWT. For the IASI instrument, a unit-delay predictor is used instead of the 2D 9/7-M IWT, due to the better performance observed empirically. Note that, in its current version, FAPEC does not directly support the application of only spatial 2D IWT to scenes with more than one band. Therefore, to apply the 2D IWT in the aforementioned cases, scenes are compressed as a single band of width *w* and height *h* · *z*, where *w*, *h* and *z* are the width, height and number of bands of the original scene, respectively. This reorganization is performed after the spectral transform, when one is applied. Finally, when the IWT is selected, two different approaches are used based on the number of spectral bands. If the scene contains 22 bands or fewer, three spatial 2D 9/7-M

IWT decomposition levels are applied on each image band, and then the unit-delay differentiator is applied along the spectral axis. If the scene contains 23 or more bands, the unit-delay differentiator is substituted by three levels of 9/7-M IWT, also along the spectral axis. Hence, regardless of the number of bands, FAPEC does not exactly follow the diagram shown in Figure 1 when the IWT is selected. Details and justification of this design decision are provided in [25,35]. Note that this approach is in contrast to CCSDS 122.1 with the IWT, which always applies spectral decorrelation before the spatial 2D IWT. In all cases, FAPEC is invoked with parameters -chunk 128M -bl 128 -mt 1 -ow -noatt -be, adding -dwt when the spectral IWT is selected. This sets the block size to 128 samples, and forces a single thread of execution. In addition, parameter -signed is used when appropriate.

The JPEG-LS [27] is a well-known compression standard that offers lossless compression, not to be confused with JPEG. JPEG-LS uses a simple predictor to spatially decorrelate input samples. This predictor uses only three of the nearest neighboring pixels, namely, those at the North, West, and Northwest positions with respect to the pixel being coded. Entropy coding is performed using Golomb-Rice codes, in combination with context modeling to further improve compression performance. Context modeling uses one additional neighbor, situated at the Northeast position relative to the current pixel. An additional run-length coding mode is also available for low entropy distributions. In this work, the reference implementation (Available at https://github.com/thorfdbg/ libjpeg) is employed, using parameters -ls 0 -m 0. As JPEG-LS does not support compression of more than three spectral bands, scenes are compressed as a single band of with *w* and height *h* · *z*, as described above for CCSDS 123.0-B-2.

The CCSDS 123.0-B-2 [28] standard is the latest compression algorithm published by the CCSDS. Unlike FAPEC, CCSDS 122.1, and JPEG-LS, CCSDS 123.0-B-2 combines spectral and spatial decorrelation in a single prediction stage (Note that FAPEC provides a prediction-based multiband algorithm that combines spatial and spectral decorrelation, as described in [25]. Notwithstanding, this paper focuses on transform-based compression, so this alternative mode is not considered). This prediction stage uses so-called local sums and local differences, computed from neighboring pixels from the same and previously coded bands. These are combined to obtain a prediction of the current pixel based on a weighted average. Weights employed for prediction are automatically updated based on the encoded values. CCSDS 123.0-B-2 allows three different entropy coders to compress the prediction errors of the previous stage: block-adaptive, sample-adaptive, and hybrid. The hybrid is the latest entropy coder published by the CCSDS, and is designed to yield competitive performance for a wide range of data entropies. This coder updates two variables—the accumulator and the counter—whenever a sample is coded. These variables provide an adaptive estimation of the data source's entropy. This estimation is used to determine what coding mode is to be used for the current sample. In addition to Golomb power-of-two coding, 16 variable-to-variable codes are defined for low-entropy scenarios. In this work, the default parameters described in [36] are used for CCSDS 123.0-B-2.

M-CALIC [29] is the last compressor considered in this work. Spectral decorrelation is performed by linear prediction using co-located samples in previous bands. Prediction errors are then compressed using an adaptive arithmetic coder with 256 defined contexts. These contexts implicitly perform spatial decorrelation by using neighboring samples at the same spectral band to determine the coding context to be used.

#### *2.3. Test Corpus*

In order to obtain a representative corpus of test scenes, 47 scenes have been selected from those available in the "CCSDS MHDC" image corpus, which has been assembled over the years by the *Multispectral Hyperspectral Data Compression Working Group* (MHDC WG) of the CCSDS to provide a suitable testing framework for the compression of satellite images (See https://cwe.ccsds.org/sls/ docs/SLS-DC/123.0-B-Info/TestData/README.txt for futher details.). These images were obtained by 14 different instruments deployed in real missions. Both hyperspectral and multispectral scenes are

used so that results are illustrative of currently deployed sensors. The corpus includes mainly raw (unprocessed) scenes and some calibrated scenes. Dynamic ranges vary from 8 to 15 bits per sample, although most scenes employ at least 12 bits per sample. Note that all scenes are stored in files using 2 bytes per sample, regardless of the actual dynamic range. A summary of the selected scenes and their properties is provided in Table 1. For illustrative purposes, crops of two sample scene bands are shown in Figure 2. Moreover, the zero-order entropies of each spectral band of three different scenes are plotted in Figure 3. All scenes are publicly available (Download links for the test scenes can be found at at http://cwe.ccsds.org/sls/docs/sls-dc/123.0-B-Info/TestData), except for those produced by the IASI and MSG instruments due to licensing restrictions.

**Figure 2.** Crops (512 × 512) of two scenes from the test corpus. (**a**) Band 0 of the *vgt\_1b* scene from the Vegetation sensor; (**b**) band 5 of the *agriculture* scene from the Landsat sensor (over Northern Mexico). Scene brightness has been adjusted for display purposes.

**Figure 3.** Zero-order entropies in bits per sample for each band. (**a**) *gran9* scene, AIRS instrument (Pacific Ocean, Daytime); (**b**) *mantar\_Rad\_rmnoise* scene, SFSI instrument (human-made targets); (**c**) *Geo\_Sample* scene, Hyperion instrument (over Goldfield/Cuprite, Nevada).


**Table 1.** Summary of employed test corpus scenes and some of their properties.

#### **3. Experimental Results**

Experimental results for the aforementioned spectral decorrelation methods and compression algorithms are presented in this section. The POT implementation of CCSDS 122.1 is also applied to FAPEC and JPEG-LS. The IWT implementation of CCSDS 122.1 is also applied to JPEG-LS, while FAPEC uses its own. All implementations are in C/C++, except for the RWA transform, which is available as Matlab code (The RWA can be downloaded from http://gici.uab.cat/GiciWebPage/downloads.php). All time measurements are obtained as the average over 10 repetitions on a dedicated Intel(R) Core(TM) i7-8700 CPU @ 3.20 GHz machine, using a single thread for the experiments. The user and system CPU times are considered here (as opposed to real time) to remove measurement dependencies on the data storage hardware.

#### *3.1. Spectral Transforms*

Two important properties of the spectral decorrelation transform described in Section 2.1 are analyzed next: transformed entropy and dynamic range expansion.

Entropy rates before and after spectral decorrelation provide valuable information about how much redundancy is removed. Entropy is defined as − ∑*<sup>i</sup> pi* · log2(*pi*), where *pi* is the probability of the *i*-th sample value. This probability is computed empirically based on each symbol's observed frequency. The distribution of zero-order entropy for the original and transformed scenes are shown in Figure 4a. The vertical axis of the figure indicates the fraction of the images whose entropy falls within the range indicated by the width of each vertical bar. Note that a relative frequency of 1 is equivalent

to 100% of the images. No attempt is made to measure spatial redundancy present within the scene bands. Note that some of the analyzed spectral transforms produce side information, namely, the POT and the RWA. For these, the size of the side information is added to the zero-order entropy. Results are expressed in bits per sample (bps), where each pixel of each spectral band is considered a different sample. As expected, all spectral transforms reduce entropy rates. After the IWT, average entropy for all 47 scenes is 1.1 bps lower with respect to the average entropy of the original scenes. In turn, the POT and the RWA transform reduce entropy rates by 2.3 bps and 2.6 bps, respectively.

**Figure 4.** Property histograms for original and spectrally decorrelated scenes. (**a**) Entropy including side information; (**b**) dynamic range. Solid dots represent average values and horizontal lines span two standard deviations.

The dynamic range of an scene is hereinafter defined as the maximum number of bits required to store any of its samples. Transformed scenes typically exhibit larger dynamic ranges. This expansion is of interest because larger samples may slow down compression, e.g., when a bitplane encoder is employed or when additional bytes are needed to represent each sample. The distribution of dynamic ranges for the original and transformed scenes are shown in Figure 4b. The vertical axis indicates relative frequency as defined for Figure 4a above, considering dynamic range instead of entropy. All transforms introduce relatively small dynamic range expansions for the test corpus. The IWT introduces the largest expansion of all tested transforms, i.e., 0.60 bps. The POT and RWA transform both increase it by 0.23 bps. It is worth noting that all transformed scenes have dynamic ranges of 16 bits or less. Thus, compressors applied after these transforms need not support sample sizes beyond 2 bytes.

#### *3.2. Lossless Compression Rates*

Compressed data rates for the spectral transforms and compressors described in Section 2 have been obtained for all scenes of the test corpus. Table 2 shows average results for each instrument expressed in bps. Average results for all scenes are also shown for ease of comparison. Note that these averages are obtained as the arithmetic mean of the bps for each scene, and this mean is not weighted based on the number of samples of each scene. Consistent with the previous section, bits per sample are obtained assuming that each pixel of each spectral band is a different sample. For ease of comparison, Figure 5 displays the compression efficiency of all tested compressors, relative to that of CCSDS 123.0-B-2, averaged for each instrument. This relative efficiency is defined as the compressed data size of CCSDS 123.0-B-2 divided by the compressed data size of each algorithm, including side information when applicable.

**Figure 5.** Average compression efficiency relative to CCSDS 123.0-B-2 for each instrument. Values larger than 1 indicate better performance than CCSDS 123.0-B-2. (**a**) FAPEC compressors; (**b**) CCSDS 122.1 compressors; (**c**) JPEG-LS compressors; (**d**) M-CALIC compressor.

The lowest average compressed data rates, i.e., those with smallest compressed files, are obtained with the CCSDS 123.0-B-2 standard. This can be explained by the fact that CCSDS 123.0-B-2 uses the most complex decorrelation approach and entropy coding mechanism of all tested compressors. Both decorrelation and entropy coding are adaptive at the sample level, which improves compressed data rates. As discussed in the next section, this is also at the cost of slower execution times. When combined with the most advantageous spectral decorrelation strategy (when applicable), compressors sorted by compressed data rates averaged for all test images are CCSDS 123.0-B-2 (4.84 bps), JPEG-LS (5.01 bps, IWT), M-CALIC (5.27 bps), CCSDS 122.1 (5.46 bps, POT), and FAPEC

(5.54 bps, POT). Thus, on average, the best configuration of all codecs are within 16% of CCSDS 123.0-B-2.

**Table 2.** Average lossless compression rate in bits per sample. Column *No* indicates compression is applied to the original scenes. Results include side information overhead when applicable. The best average results per instrument are highlighted in bold face.


For some of the instruments considered separately, CCSDS 123.0-B-2 yields larger differences when compared to the most advantageous transform for each compressor. For example, rates 53% better than FAPEC and CCSDS 122.1 are observed for M3 scenes, 41% better than FAPEC and 37% better than CCSDS 122.1 for CRISM data, and 25% better than FAPEC for the AIRS instrument. On the other hand, some combinations outperform CCSDS 123.0-B-2 for other instruments, e.g., JPEG-LS with the POT is 13% better for MODIS and with the RWA transform is 2% better for AVIRIS.

Regarding the optimum transforms for each compressor, FAPEC typically performs best with the POT or the RWA transform. On average, the POT provides a compression gain of 0.72 bps as compared to not using any transform, 0.07 bps better than using the spectral IWT. In turn, the RWA transform yields average bitrates 0.05 bps higher than using the POT transform, although the RWA provides the best results for 8 of the 14 tested instruments. CCSDS 122.1 achieves the best rates when using POT for two-thirds of the instruments, with RWA leading on the rest. Finally, JPEG-LS benefits specially from RWA on half of the tested instruments, with the IWT on less than a third, and with the POT on just two instruments.

Consistent with the entropy reduction results shown in Section 3.1, all tested transforms typically improve average compressed data rates. There are only a few exceptions on specific instruments and compressors, such as MSG (where the use of a spectral transform always leads to worse results) or MODIS with several combinations (CCSDS 122.1 with IWT or RWA, and JPEG-LS with IWT or RWA). Moreover, PLEIADES and SPOT5 exhibit a few exceptions. Results suggest that instruments with a small number of bands or a high spatial resolution are prone to get worse lossless compression rates when prepending some of these spectral decorrelation transforms.

In general, the best transform for spectral decorrelation is compressor-dependent. This is due to the different assumptions made in the design of the different algorithms on the distribution of the samples. Furthermore, each transform spreads the noise present in the least significant bits of the original samples differently, e.g., due to sensor shot noise. Notwithstanding these differences, average results for all transforms are within 0.16 bps and 0.10 bps for CCSDS 122.1 and JPEG-LS, respectively. Note that results for CCSDS 123.0-B-2 and M-CALIC after spectral decorrelation are not presented, because both compressors are designed assuming pixel distributions typical for the original domain.

**Table 3.** Average compression time in seconds. Column *No* indicates compression is applied to the original scenes. Results include spectral transformation times when applicable. Best results for each instrument are highlighted in bold font.


**Table 4.** Average decompression time in seconds. Column *No* indicates decompression is applied to the original scenes. Results include inverse spectral transformation times when applicable. Best results for each instrument are highlighted in bold font.


### *3.3. Throughput Results*

The throughput yielded by the compressors described in the previous section is analyzed next. Average compression times for each instrument and for all scenes is provided in Table 3. Execution time is measured as the total user and system time reported by the operating system (Ubuntu 18.04 LTS). This is to avoid dependencies with the employed storage medium, and to provide a better surrogate for energy consumption in the comparison. It should be highlighted that all tested compressors are software implementations. Some differences in the level of optimization of these implementations were observed, also between different compression modes of the same algorithm.

The fastest of all tested compressors is FAPEC. When executed without any spectral decorrelation transform, it is 3.30 times faster than the next fastest non-FAPEC compressor, i.e., CCSDS 123.0-B-2. When the IWT is applied both spatially and spectrally with FAPEC, it is 5.0 times faster than CCSDS 123.0-B-2. The faster execution speed with the IWT, as compared to not using any spectral decorrelation, is due to the better decorrelation power of the IWT, and its low computational complexity. After applying the IWT, the entropy coder has to deal with a more favorable distribution of input symbols. This results in the emission of fewer bits per sample, and contributes to requiring shorter execution times. Thanks to the low complexity of the IWT, the execution time improvements due to its decorrelation power often outweigh the cost of applying this transform.

When compared to the fastest configuration of CCSDS 122.1 and JPEG-LS, FAPEC with the IWT is 37.1 times and 7.0 times faster, respectively. CCSDS 122.1 achieves average execution times significantly slower than FAPEC, JPEG-LS, and CCSDS 123.0-B-2. This is explained by the use of a bitplane encoder which imposes several passes on each sample when implemented on software. In turn, M-CALIC is the slowest of all tested compressors, 9.1% slower than the slowest configuration of CCSDS 122.1. These results can be explained by the use of an arithmetic coder in the entropy coding stage of M-CALIC. It is also interesting to analyze average compression time as a function of the employed spectral transform. As can be observed, compression using the IWT is faster than compression using the POT or the RWA transform. This is due to the lower complexity of the IWT. For the CCSDS 122.1 and JPEG-LS compressors, coding with the IWT takes respectively 3.7% and 17.0% less than with the fastest of POT and RWA. Significantly larger differences are observed for FAPEC and the IWT. As mentioned before, this can be explained by the use of a more efficient implementation of the IWT, more tightly integrated in the compression pipeline. For most instruments, the RWA attains higher throughput than the POT. Notwithstanding, the throughput of RWA for the IASI instrument is significantly lower due to the very high number of bands, which makes the throughput of POT better on average. For completeness, average decompression times—including inverse spectral transformation, when applicable—are presented in Table 4. Decompression times are generally lower than compression times. An exception to this is CCSDS 123.0-B-2, whose entropy decoder processes data in inverse order [28]. Another exception is FAPEC with the spectral IWT. This is due to an asymmetry between the optimization degree of the compressor (meant to be executed onboard spacecraft hardware) and the decompressor (meant to be executed on the ground).

Additional analysis is devoted to jointly compare compression performance and execution times. Figure 6a provides average compression speed as a function of the obtained compressed data rates for a selection of spectral transforms and compressors. For each compressor, the spectral decorrelation method that yields the best average compressed data rate is included. FAPEC with the IWT is also included for completeness. Figure 6b shows average decompression speed as a function of the obtained compressed data rates for the same selection as Figure 6a. For medium and low compression throughput, i.e., below 12 × 106 samples/s, CCSDS 123.0-B-2 dominates over JPEG-LS with spectral IWT, CCSDS 122.1 with spectral POT, and M-CALIC and FAPEC with the POT. It can also be observed that FAPEC with the spectral IWT yields significantly higher compression speeds, reaching 60 × <sup>10</sup><sup>6</sup> samples/s on average. As discussed above, the additional speed is traded for about 16% higher compressed data rates. Decompression throughput follows a pattern similar to compression, with the difference that most compressors achieve higher speeds. As execution speed is directly related to energy consumption, and assuming energy is a limited resource in the system where compression is to be performed, results suggest that FAPEC yields the most desirable trade-off when energy and compression performance are jointly considered.

**Figure 6.** Average (**a**) compression and (**b**) decompression speed as a function of the compressed data rate for all scenes in the corpus. The spectral decorrelation transform applied to the data is shown between brackets, and *Original* indicates no transform. When applicable, speeds include the spectral decorrelation stage times, and compressed data rates include side information overhead.

#### **4. Conclusions**

Hyperspectral imaging (HSI) is a staple tool of remote sensing and Earth observation. Data compression can significantly increase the amount of valuable industrial and scientific data retrieved on the ground. In particular, scene compression based on spectral decorrelation transforms enables important bitrate reductions at a reasonable computational complexity. This work presents a thorough analysis of the IWT, the POT, and the RWA transform as spectral decorrelators, in combination with FAPEC, CCSDS 122.1, and JPEG-LS coding techniques. Results for CCSDS 123.0-B-2 and M-CALIC are also shown for completeness. Compression rates and speeds on a large variety of scenes have been analyzed here, considering only lossless regimes. This analysis was based on software implementations, which enable more flexibility than their hardware equivalents (when available), especially on small satellites and cubesats. To the best of our knowledge, such a thorough comparison is not previously available in the literature. Experiments suggest that the use of a spectral decorrelation transform improves compression rates for most HSI instruments. The most advantageous transform depends on each kind of instrument and on the compressor to be used. Both the POT and the RWA transform can provide significantly better compression rates than IWT in specific cases, yet at a higher computing cost. Results indicate that the best trade-off is attained by the IWT followed by FAPEC. On average, this combination produces compressed data rates within 16% of CCSDS 123.0-B-2, which typically yields the smallest compressed data volumes. At the same time, FAPEC with its IWT is 5.0 times faster than CCSDS 123.0-B-2.

**Author Contributions:** Conceptualization, J.S.-S., J.P., M.H.-C., and I.B.; methodology, J.P., J.S.-S., M.H.-C., and I.B.; software, J.P., M.H.-C., and I.B.; validation, J.P. and J.S.-S.; formal analysis, M.H.C.; investigation, J.P., J.S.-S., and M.H.-C.; resources, I.B.; data curation, M.H.-C. and I.B.; writing—original draft preparation, M.H.C.; writing—review and editing, J.P., J.S.-S., M.H.-C., and I.B.; visualization, M.H.-C.; supervision, J.S.-S. and J.P.; project administration, J.S.-S. and J.P.; funding acquisition, J.P., J.S.-S., M.H.-C., and I.B. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by the Spanish Ministry of Economy and Competitiveness and the European Regional Development Fund under grants RTI2018-095287-B-I00, TIN2015-71126-R, RTI2018-095076-B-C21 (MINECO/FEDER, UE); by the ESA Business Incubation Programme and Barcelona Activa; by the Catalan Government under grant 2017SGR-463; by the postdoctoral fellowship programme Beatriu de Pinós, reference 2018-BP-00008, funded by the Secretary of Universities and Research (Government of Catalonia); and by the Horizon 2020 programme of research and innovation of the European Union under the Marie Skłodowska-Curie grant agreement #801370.

**Conflicts of Interest:** The authors declare the following conflicts of interest. J.P. is CTO of DAPCOM, original developer of the FAPEC algorithm. M.H-C., I.B., and J.S. have participated in the standardization of several CCSDS standards, including CCSDS 122.1-B-1 and CCSDS 123.0-B-2. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript; or in the decision to publish the results.
