New calibration noise suppression techniques for the GLORIA limb imager

. The Gimballed Limb Observer for Radiance Imaging of the Atmosphere (GLORIA) presents new opportu-nities for the retrieval of trace gases in the upper troposphere and lower stratosphere. The radiometric calibration of the measured signal is achieved using in-ﬂight measurements of reference blackbody and upward-pointing “deep space” scenes. In this paper, we speciﬁcally to calibrate GLORIA data exploiting the instrument’s imaging capability. The algorithms discussed here make use of the spatial correlation of parameters across GLORIA’s detector pixels in order to mitigate the noise and in the calibration measurements. This is achieved by combining a priori and empirical knowledge about the instrument background with noise-mitigating compression methods, speciﬁcally low-pass ﬁltering and principal component analysis (PCA). In addition, a new software package for the processing of GLORIA data is introduced which allows us to generate calibrated spectra from raw measurements in a semi-automated data processing chain.


Introduction
The chemical composition of the upper troposphere and lower stratosphere (UTLS) and its dynamic changes are known to have a significant influence on surface climate (Solomon et al., 2010;Riese et al., 2014).Key processes in the region occur on much smaller spatial scales than can be rendered by current climate models (Gettelman et al., 2011).More finely resolved observations of UTLS composition and temperature are therefore necessary to improve the parameterisations by which these processes enter into models.Most remote sensing instruments have, until recently, been lacking a sufficiently fine spatial resolution for this purpose.The typically better resolved in situ observations, on the other hand, lack the necessary spatial coverage.In order to narrow this gap, the Gimballed Limb Observer for Radiance Imaging of the Atmosphere (GLORIA) was developed jointly by the German research centres Forschungszentrum Jülich GmbH (FZJ) and Karlsruher Institut für Technologie (KIT) (Riese et al., 2014;Friedl-Vallon et al., 2014).
GLORIA is the direct successor of two successful lines of infrared limb sounders.One of these is the Cryogenic Infrared Spectrometers and Telescopes for the Atmosphere (CRISTA), which was developed originally as a satellite instrument and saw two orbital deployments during Space Shuttle missions in 1994 and 1997 (Offermann et al., 1999;Riese et al., 1997Riese et al., , 2002;;Grossmann et al., 2002).Due to the success of these missions, an airborne adaptation, called CRISTA-NF (for "New Frontiers"), was developed later and deployed on the M-55 Geophysica research aircraft (Kullmann et al., 2004).

T. Guggenmoser et al.: Calibration noise suppression for GLORIA
The other line of GLORIA's heritage is the Michelson Interferometer for Passive Atmospheric Sounding (MIPAS).Apart from the 10-year deployment of MIPAS on board the European Space Agency's Environmental Satellite (Envisat) (Fischer et al., 2008), MIPAS has also been highly successful as a balloon instrument (Friedl-Vallon et al., 2004).The airborne MIPAS-STR (for "Stratospheric Aircraft") adaptation (Piesch et al., 1996), also developed for the M-55 Geophysica, can be seen as GLORIA's direct precursor, alongside CRISTA-NF.
CRISTA-NF and MIPAS-STR were both conceived for observations in the UTLS, but with different priorities.The former was primarily designed for the study of transport and mixing processes.CRISTA-NF data from the African Monsoon Multidisciplinary Analyses (AMMA) and Reconciliation of Essential Process Parameters for an Enhanced Predictability of Arctic Stratospheric Ozone Loss and its Climate Interactions (RECONCILE) measurement campaigns have been used to derive trace gas mixing ratios with unprecedented vertical resolution and along-track sampling (Ungermann et al., 2012(Ungermann et al., , 2013;;Kalicinsky et al., 2013;Weigel et al., 2012).MIPAS-STR, on the other hand, was designed for a chemistry-centred analysis, offering a higher spectral resolution in order to derive more elusive trace gases (Woiwode et al., 2012), at the cost of lower spatial sampling.
Compared with CRISTA-NF and MIPAS-STR, GLORIA incorporates a number of advancements that greatly enhance and diversify its UTLS sounding capabilities.As an imager, it is capable of recording complete vertical profiles in a single measurement, without the need for elevation scanning.This greatly improves the spatial sampling of measurements.Another main feature of GLORIA is its ability to operate in either of two measurement modes: a chemistry mode (CM) with high spectral sampling (0.0625 cm −1 ) and a dynamics mode (DM) with lower spectral sampling (0.625 cm −1 ) but higher spatial sampling.
The CM can be used for the retrieval of trace gases with faint or not easily separable spectral signatures, such as ethane (C 2 H 6 ), using the Karlsruhe Optimized and Precise Radiative Transfer Algorithm (KOPRA) (Stiller, 2000;Höpfner et al., 1998).The DM measurements are 10 times as quick.Together with GLORIA's ability to pan the line of sight between measurements, it is possible to collect DM data sets suited for 3-D tomographic trace gas retrieval algorithms.This has been demonstrated via simulation (Ungermann et al., 2010(Ungermann et al., , 2011) ) and applied to measurement data (Kaufmann et al., 2015).Measurements recorded in DM are all processed using the JUelich Rapid Spectral Simulation Code (Hoffmann et al., 2008) version 2 (JURASSIC2), together with the JUelich Tomographic Inversion Library (JU-TIL) (Ungermann et al., 2015).
In this paper, several advancements in the processing of GLORIA spectra are presented, particularly the derivation of the instrument calibration function.Previously, the calibration parameters were found by treating each pixel of the fo-cal plane array individually.New techniques have now been developed to generate the calibration parameters which make use of the spatial homogeneity of the individual images.With this new approach, the noise present in the gain and offset parameters and resulting artefacts in the calibrated radiances are significantly reduced.
These new algorithms were developed in tandem with a new data processing system which incorporates previously separate processing stages, resulting in a semi-automated processing chain.Once the calibration parameters are found, it is now possible to process GLORIA raw interferograms into calibrated radiance spectra without the need for intermediate products or human intervention.
This paper is laid out as follows.Section 2 will introduce the basics of GLORIA data processing and provide the background and terminology for the following parts.We also briefly present the newly developed data processing system and evaluate its runtime performance.In Sect.3, an empirical spatial characterisation of GLORIA's radiometric background is given, and a method is discussed which exploits this knowledge in order to mitigate certain calibration artefacts.A more general approach for noise suppression which is applicable to all GLORIA calibration measurements, based on principal component analysis (PCA), is presented in Sect. 4. We conclude with a summary of the advancements achieved, as well as an outlook on work which builds upon it and is currently in progress.

L0/L1 processing and radiometric calibration
The heart of the GLORIA instrument is a linear Michelsontype interferometer (Friedl-Vallon et al., 2014).Incident infrared light is divided by a beam splitter, resulting in two rays which take different paths through the instrument until recombined at the detector.The difference between the distances travelled by the rays, the optical path difference (OPD), is varied continuously by a moving mirror.
One of the instrument's key features is that the recombined beam is captured with a two-dimensional focal plane array detector which provides images with up to 256 px × 256 px.Typically, a region of interest (ROI) of 48 horizontal and 128 vertical pixels is read out with a sampling rate of 6.281 kHz.The recorded intensities for each pixel are called the (timesampled) interferograms.
A co-modulated length measurement interferometer (LMI) provides the necessary reference data to interpolate these time-sampled recordings onto a spatial axis, i.e. express the radiance as a function of the OPD.This is referred to as the L0 processing step (Kleinert et al., 2014).The resulting spatial interferograms are, up to the instrument function, the Fourier transform of the spectrum of recorded radiances (Beer, 1992).
The L1 processing step derives absolute spectral radiances from the L0 data.The first and most straightforward step is to apply the Fourier transform.Because GLORIA records two-sided interferograms, the resulting spectra are complexvalued.If the two interferogram sides were perfectly symmetric, the imaginary part would vanish.This is not the case for three principal reasons: firstly, noise will not be symmetric around the point of zero OPD.Secondly, the point of zero OPD is never exactly sampled by the instrument in practice, although this effect can be compensated for by an interferogram alignment correction at the L0 stage.Thirdly and most importantly, because GLORIA does not operate near absolute zero temperature, the interferograms contain emission and absorption signatures from various parts of the instrument, which enter the signal at different phase angles relative to each other (Kleinert and Trieschmann, 2007).
The radiometric calibration, i.e. the inversion of the instrument function, is performed assuming C-linear behaviour.If y(ν) is the spectral radiance at the spectral sample ν, the spectrum obtained from the measurements is with a and b being the complex radiometric gain and offset parameters.The advantage of this complex calibration scheme is that the non-atmospheric signal components are automatically removed, without having to determine their individual magnitudes and phase angles (Revercomb et al., 1988).
As known reference radiation sources, the GLORIA instrument carries two blackbodies with it, BBH and BBC (for "hot" and "cold", respectively).They emit spatially homogeneous infrared radiation at their respective stabilised temperatures (Monte et al., 2014;Olschewski et al., 2013).Using measurements from these blackbody sources, as well as upward-pointing deep space (DS) measurements, the parameters a and b for each pixel and spectral sample can be obtained by either of two methods.The first method, called BB-DS calibration, uses averaged measurements of one blackbody source (y BB ) and averaged deep space measurements (y DS ): where B(T BB , ν) is the Planck function evaluated at the blackbody temperature T BB .The method assumes that the DS measurements are, to a good approximation, free of atmospheric emissions, and can therefore be identified with the radiometric offset.This assumption does not strictly hold true for airborne observations due to the non-negligible emissions from the upper stratosphere, especially ozone, along the line of sight.Before the DS measurements are used, these emissions are first removed by an iterative procedure based on the KOPRA forward model (Kleinert et al., 2014).This method is successful in most spectral regions but insufficient in a narrow band between 1020 and 1070 cm −1 , where BB-DS calibration can consequently not be used without introducing spectral artefacts.
The second method, called BB-BB calibration, uses measurements from both blackbody sources: where for the determination of the offset b any of the two blackbodies can be used.The BB-BB method avoids the problem of residual trace gas emissions, but at the cost of having to extrapolate the radiometric offset, which the BB-DS method can measure directly.The quality of this extrapolation depends on the absolute temperature of the blackbody sources, as well as on the temperature difference between them.
Having determined the radiometric gain and offset parameters for every pixel and spectral sample, each spectrum can be calibrated by application of where α(ν) is the inverse gain (IG) and β(ν) is the negative calibrated offset (NCO).The NCO is readily interpreted physically as the negative of the instrument background radiance.Its real part is dominated by thermal emissions from the components surrounding the detector, while the imaginary part is dominated by emissions from the beam splitter.This will be of importance in Sects.3 and 4.
During the first GLORIA measurement campaigns, a typical calibration sequence consisted of 20 measurements of each blackbody, followed by 10 DS measurements.These sequences were performed every 30 to 45 min and took about 5 min each.Only measurements with matching interferometer direction are processed together.Typical temperatures for the blackbodies were about 240 vs. 256 K, but values varied throughout the flights, with the difference between the two ranging from 15 to 25 K.These temperatures are both higher and less far apart than expected, which presents problems with error magnification for the BB-BB calibration.For this and other reasons, including an imperfect nonlinearity correction, the three calibration sources can currently not all be used together in a three-point fit (Kleinert et al., 2014).

Integrated data processing chain
Due to GLORIA's imaging capabilities and high sampling rate, the instrument produces a large amount of data.Assuming a 48 px × 128 px configuration, the raw detector frames alone generate over 70 MiBs (1 MiB = 2 20 bytes) per second.Accordingly, measurements from a long-range platform such as the High Altitude and Long Range Research Aircraft (HALO) can reach 1-2 TiBs (1 TiB = 2 40 bytes) per flight.

T. Guggenmoser et al.: Calibration noise suppression for GLORIA
Over 20 TiBs of data were recorded over the course of the Transport and Composition in the Upper Troposphere and Lowermost Stratosphere (TACTS) and Earth System Model Validation (ESMVal) measurement campaigns in 2012.To cope with this amount of data and to obtain L1 products in an efficient manner, a new data processing system has been developed, and upon it an automated L0/L1 processing chain has been built (Guggenmoser, 2014).
The new software is implemented mainly in the Python programming language.The optimised L0 algorithms (Kleinert et al., 2014) as well as some performance-critical L1 operations have been integrated in compiled form.Shared-memory parallelisation is used to take advantage of multi-core CPUs.One of the computationally most expensive operations, the Fourier transform of the whole array of L0-resampled interferograms, is delegated to the highperformance FFTW (fastest Fourier transform in the West) library.FFTW analyses the dimensionality of its input data and selects the best of a collection of Fourier transform algorithms (Frigo and Johnson, 1998).
The new design allows the previously separate components to interface seamlessly, without mediation through file input/output or human intervention.This, in turn, enables the development of a L0/L1 processing chain which produces calibrated spectra from raw interferogram data.
The runtime for the complete processing of a single measurement, and the effect of parallelisation, is shown in Fig. 1.The runtime scales linearly with the interferogram length, which is why chemistry and dynamics mode measurements have very similar parallelisation speedups.Single-thread runtime is approximately a factor of 10 longer than the data acquisition time: CM measurements (12 s acquisition time) are processed in slightly over 2 min, while DM measurements (1.2 s acquisition time) take about 13 s.The time is reduced by half using four parallel threads to about 5 times the acquisition time.Further parallelisation yields diminishing returns, mainly because of finite memory bandwidth and the non-parallel code segments.For comparison: the processing time for the L1 stage alone (i.e.Fourier transform and radiometric calibration), using the previous L1 processor, surpassed the acquisition time by a factor of 60.
In practice, the processing is usually performed on a cluster, with multiple measurements being processed simultaneously.This provides another level of parallelisation which is subject to different constraints.The optimal balance between the number of processes and number of threads depends on several factors such as mass storage bandwidth and memory available per compute node.

Spatial structure of the instrument background
When using DS measurements in combination with one BB source (BB-DS calibration), the radiometric offset can be measured with high quality using Eqs.( 2) and (3).The NCO function derived in this manner varies smoothly over the detector, as is to be expected of a quantity that represents the physical background radiation.
By contrast, the NCO derived from the BBC and BBH sources (BB-BB calibration) has turned out to be much less smooth than expected.This can be explained by a combination of noise and an imperfect nonlinearity correction for some pixels (Kleinert et al., 2014).Being able to use the BB-BB calibration method is desirable because it circumvents the problem of residual trace gas emissions in the DS measurements, and because it can serve as a consistency test for the parameters derived from BB-DS.
We found that the spatial variation of the real part of the NCO over the detector can be reproduced by a simple pseudo-hyperbola of the form The parameters x c , y c , a, b, c, d ∈ R describe the position and shape of the pseudo-hyperbola, while x and y are the coordinates of the detector pixels.A series of higherdimensional functions was also tested but did not yield better results than Eq. ( 7).While imperfect in certain spectral ranges, this function lets us correct the more serious artefacts in the following manner.First, the radiometric gain and offset are determined according to Eqs. ( 4) and ( 5), and from these two the NCO function.Then, individually for each spectral sample, a least squares fit is performed to find the pseudo-hyperbola that best reproduces its real part.This real part is then replaced by the smooth pseudo-hyperbola.Using the new NCO function and the original calibration measurements, the radiometric gain function is updated for consistency.
A representative calibration sequence from the TACTS measurement campaign has been chosen to illustrate the effects of this method.From the BB and DS measurements, the calibration parameters IG and NCO have been determined, both via BB-BB calibration and via BB-DS calibration.Both sets of parameters have then been subjected to the pseudohyperbola fit correction as described above, yielding a total of four sets of calibration parameters to compare.
The comparison of the NCO real parts for the spectral sample at 791.875 cm −1 is shown in Fig. 2. The non-fitted BB-DS calibration serves as a reference.The symmetric, ringlike shape of the NCO function which motivated the pseudohyperbola is recognisable in the uncorrected BB-DS and BB-BB data.However, the BB-BB image is much more noisy than the BB-DS reference, and differs from it both by a homogeneous bias and by a large number of small-scale structures.These structures, which caused artefacts in the L2 processing stage, are removed by the pseudo-hyperbola fit.The bias and smoothly varying difference still remain, but this kind of calibration error can be more easily handled by the Figure 1.Runtime for the automated L0/L1 processor for different levels of parallelisation.In order to eliminate the influence of file input/output, the raw data were read from and the calibrated spectra written to a memory-based file system.Note that the L0 runtime as shown here encompasses the entire L0 processing stage, including the interferogram alignment correction.Earlier analyses of the operational L0 processor (Kleinert et al., 2014) showed the runtime for single execution.trace gas retrieval algorithms.Because the IG and NCO functions are interdependent, updating the BB-BB IG with the corrected NCO removes similar structures from it, as shown in Fig. 3.
The pseudo-hyperbola function cannot reproduce the NCO for all spectral samples, as shown in Fig. 4 for the one at 840 cm −1 .Here, the difference between the BB-DS NCO with and without the correction shows ring-like structures, which are likely related to emissions of the outer interferometer germanium window that peak around this wave number.However, the improvements to the BB-BB calibration over the rest of the spectrum outweigh the cost of these artefacts, which only persist for a limited number of samples.Ongoing characterisation of the window emissions may eventually allow the inclusion of window emissions as dedicated terms in the fit formula.
Following the L1 processing step, the spectra are currently averaged over one row to reduce stochastically independent noise for the L2 processing.For this application, the average over each row of the NCO is an interesting diagnostic for the pseudo-hyperbola method.Figure 5 shows this quantity.Comparing the magnitude of the removed structures (of the order of 10 −7 W cm 2 cm sr −1 in the centre of the detector) with the magnitude of recorded atmospheric spectra, the artefacts would introduce a spatially uncorrelated error of up to nearly 8 % at this wave number for the BB-BB calibration unless corrected by the pseudo-hyperbola.
The size of this correction is attributable to two reasons.As mentioned above, the BB-BB calibration currently suffers from error magnification due to unexpectedly large and similar BB temperatures.A more fundamental reason is that the NCO itself is relatively large because the GLORIA spectrometer operates at 200-220 K (Friedl-Vallon et al., 2014), and the radiative background can be comparable in magnitude to atmospheric emissions, depending on the observed scene.
Let us assume that all pixels had the same response function.We could then pretend that the individual spectra within a calibration image were a series of measurements from the same detector, each time exposed to a slightly varying scene.
To put this in more mathematical terms, we have reason to suspect that the observed calibration scenes do not actually contain as many meaningful degrees of freedom as there are samples recorded by the instrument.If each of the n pixels' spectra consists of k partially correlated complex samples, this means we could find a way to instead represent it using a smaller but uncorrelated set of K < k complex numbers.Such a transformation can be achieved by PCA (e.g.Joliffe, 2002).The PCA is a linear basis transformation yielding the original data in terms of a new orthonormal basis, composed of the k eigenvectors of its covariance matrix and sorted in order of decreasing contribution to the total variance.These vectors are the "principal components".If the original data are noise normalised, such that the random noise becomes uncorrelated and uniform (white), any direction will carry the same amount of noise, so that reconstructing the data from a truncated set of K components will ideally reduce the sample noise SD σ (ν) by a factor of However, in practice, the relation holds true only asymptotically for large n because the ability to separate noise from signal depends on the number of measurements available for the decomposition (e.g.Pyatykh et al., 2013).A more realistic estimate, on the basis of a simulated white noise data set, will be shown in Sect.4.3.PCA-based methods for the removal of noise from Fourier transform spectrometer measurements have been investigated in great detail (Antonelli et al., 2004), and have been applied with great success for the operational data compression and noise suppression for L1 products, e.g. for the Infrared Atmospheric Sounding Interferometer (IASI, August et al., 2012).GLORIA's imaging capabilities, in combination with the homogeneity of calibration scenes, enable the use of single-or averaged-image decomposition, as opposed to decomposition based on a training set of previous measurements.Note that PCA is usually formulated in a way that applies to real-valued data sets.In the following, we will use a straightforward generalisation to complex vector spaces.This has the implication that the covariance matrix of the data set is complex-valued, which makes it less readily interpreted.However, because the covariance matrix is self-adjoint, its main diagonal and its eigenvalues remain real numbers, so the components can still be ordered in the regular fashion.

Data preconditioning and application of PCA
We start with a matrix X ∈ C n×k containing either a BB or DS data set in Fourier transformed but not radiometrically calibrated form.The calibration data sets are obtained from averaging over all available single measurements of the same calibration source during one calibration sequence.The DS measurements have additionally been subjected to the trace gas emission removal technique (Kleinert et al., 2014).Again, n is the number of pixels, while k is the number of spectral samples considered.In this section, we have chosen the spectral range between 780 and 1450 cm −1 .The calibration measurements are processed at the dynamics mode sampling of 0.625 cm −1 ; therefore k = 1072.For technical reasons, the topmost line of pixels is not considered; therefore n = 6096.
The relation with the true radiance image X as seen by the detector can be expressed as where G ∈ C n×k is the instrument gain matrix, and the symbol was chosen to represent the Hadamard (entrywise) matrix product.In the case where X is a DS data set, X is dominated by the NCO, potentially with minor contributions from insufficiently removed residual trace gas signatures.If X is a BB data set, the true image X contains the NCO and the nearly homogeneous blackbody radiation.
When PCA was introduced in Sect.4.1, it was assumed that all pixels had the same response function, so that each pixel measured the same statistical variables.More practically speaking, the set of pixel spectra resembled a series of measurements taken by the same instrument.This statistical interpretation provides the notion of a covariance matrix, whose eigendecomposition PCA is based on.In order to be able to apply PCA for noise suppression, the data set X therefore needs to be brought into a form which closely resembles this.To this end, the pixel response has to be homogenised; that is, the Hadamard inverse of the matrix G needs to be applied to the measurements.This poses a fundamental problem, as the purpose of the calibration measurements is to obtain this matrix, and it is therefore naturally not known at this point.
We solve this problem by the following process.First, a smoothing of the blackbody and deep space measurements is performed based on singular value decomposition (SVD) without gain homogenisation, using 20 singular vectors for reconstruction.These spectra, X BB and X DS are now used to infer an approximation of the gain matrix: where B contains the theoretically calculated inverse Planck radiances corresponding to the blackbody temperature.To avoid confusion with the ordinary matrix inverse, a tilde denotes inversion with respect to the Hadamard product.For the purposes of homogenisation, multiplication with the inverse Planck radiances is not strictly necessary, but facilitates the physical interpretation of intermediate results.The function is a Kaiser window which acts as a smoothing function on each row of the matrix.Its purpose is to keep noise and signal variance out of the homogenisation matrix.
Using the entrywise inverse of the matrix G , the blackbody and deep space measurements can now be homogenised by gain calibration: To further reduce the variation that needs to be decomposed by the PCA, the NCO pseudo-hyperbola from Sect. 3 is fitted to the gain calibrated deep space measurement at each spectral sample and then subtracted.Afterwards, the remaining mean value of each spectral sample is computed and subtracted as well, yielding the new matrix X .
This matrix now fulfills almost all the conditions necessary for the PCA.What is still missing is the normalisation of the individual spectral samples.Several methods exist to achieve this, the most simple of which being to divide each spectral sample by its sample SD.We opted for a method that approximates the spectral noise profile.A new reconstruction of the blackbody matrix is performed, based on the original decomposition, this time using 400 principal components.The difference between the original spectrum and the reconstructed spectrum is then used to normalise the noise variance of the matrix X , under the assumption that the reconstruction residual can be used as a proxy for the noise profile: where N ∈ R k×k contains the SD of the calculated difference on the main diagonal.Note that this is not an approximation of the total noise, but only of the spectral distribution for the purpose of rescaling.To arrive at the full noise, a more thorough analysis is necessary (Antonelli et al., 2004).The N matrix is computed separately for DS and BB measurements.Finally, the PCA is performed on Y by way of a complex SVD: Here, U ∈ C n×n and V ∈ C k×k contain the left-and rightsingular vectors, respectively, and S ∈ C n×k is a rectangular diagonal matrix which contains the singular values.The squares of the singular values are eigenvalues of the covariance matrix Y * Y, the rows of V * contain the orthonormal eigenvectors and the rows of the product T = US contain the linear factors, i.e. the coordinates of each row of Y in the new basis defined by V. A truncated reconstruction X K of X using K covariance eigenvectors can now be obtained by first reconstructing Y using truncated matrices S K and V K , then inverting the normalisation steps performed before.
After the procedure is finished, it is repeated once, starting with the original data set X.The new decomposition (13) for the BB and DS measurements is now used instead of the initial SVD to estimate the gain matrix and pseudo-hyperbola more accurately.An iterative approach using successive repetitions of the procedure was considered, but did not yield improved results.

Reconstruction from principal components
With the PCA calculated, the DS and BB spectra can be reconstructed using a chosen number of principal components.The more components are used, the more exactly the original data are reconstructed, the signal as well as the noise.Several diagnostics can be used to help decide on a sensible number of components to include.
The first diagnostic is visual inspection of the reconstructed spectra.As an example, a DS data set and a BB data set from the ESMVal campaign have been subjected to the PCA method.The original data sets were obtained by aver-aging the available calibration measurements from a calibration sequence recorded on 23 September 2012.The DS measurements were recorded at full resolution, and then the iterative DS correction procedure was applied to remove trace gas emissions from the upper atmosphere.Afterwards, the measurements were reduced to DM resolution, i.e. a spectral sampling of 0.625 cm −1 , the same as the blackbody measurements.Prior to the processing, both spectra were rotated so that their phase angle in the complex plane is approximately that of a calibrated spectrum.This is not necessary for the procedure, but helps with the physical interpretation of the real and imaginary parts.
Figures 6 and 7 show the full reconstructions of the DS and BB, respectively, from 10, 20, 40, and 70 principal components.Both spectra are reproduced well on this scale.In order for the differences to become more apparent, both between the different reconstructions and between the reconstruction and the original, we will concentrate on the DS spectrum within limited spectral ranges.
The DS reconstructions for the range between 800 and 850 cm −1 are shown in Fig. 8.In this region, the reconstruction is quite well behaved, and even from only 10 principal components the underlying signal can be reproduced without any apparent systematic residuals.This is not the case everywhere, as shown in Fig. 9 for the spectral range between 1000 and 1100 cm −1 .As has already been mentioned in Sect.2, this region contains artefacts from an incomplete removal of stratospheric trace gas emissions.These artefacts are apparently harder to reconstruct than the instrument background, presumably because they are less spatially homogeneous.In any case, the contaminated spectral range cannot be used, so the quality of its reconstruction is not crucial.

Covariance eigenvalues
Another avenue for diagnostics is provided by inspection of the singular values s i , i.e. the main diagonal elements of the matrix S from Eq. ( 13), whose squares are the covariance matrix eigenvalues.The cumulative sum of the squared singular values, normalised by the total sum, is the relative total variance of each variable restored in the reconstruction: As was mentioned in Sect.4.1, the distribution of white noise among the PCA vectors depends on the dimensionality of the data set.To estimate the performance of our method, we have generated a matrix X syn ∈ C 6096×1072 .The real and imaginary parts of the matrix entries were drawn independently from a standard normal distribution.
Figure 10 shows the first 70 normalised covariance eigenvalues for the DS spectrum from Fig. 6, the BB data set from Fig. 7, and the white noise matrix X syn .The cumulative sum of the first 30 eigenvalues is shown in Fig. 11.
For the DS spectra, the first five principal components already account for more than half of the variation in the measurement.Using 20 vectors, this figure increases to 60 %.It should be noted here that this is the variance after the homogenisation procedure; that is, most of the variation over the detector has already been removed.The spectrum reconstructed here is a combination of the detector noise, a correction to the assumed pseudo-hyperbola shape in the real part and the imaginary part of the spectrum.
The BB data set can be reconstructed to 90 % of the total variance from as few as 14 principal components.Because the BB spectra contain the same background as the DS as well as the much brighter blackbody radiation, the total vari-ance is higher.However, the blackbody radiation, which accounts for most of the signal, is easily reconstructed due to its homogeneity.
The synthetic noise is reconstructed much more linearly than the BB and DS measurements, but the covariance eigenvalues still decrease slightly for successive components.This is the expected non-uniform noise distribution which motivated the simulation.The contribution of successive principal components starts at 0.19 % and decreases to 0.16 % for the 70th component.The last eigenvalue (not shown in the plots) is 3.2 × 10 −4 , an order of magnitude smaller than the first one.The mean eigenvalue is 1/1072 ≈ 9.3 × 10 −4 , the expected contribution of each single component if the noise variance was distributed evenly.

Comparison with low-pass filtering
Another useful method for noise mitigation of GLORIA calibration measurements is a simple low-pass filtering (LPF) of the spectra.Note that, in this context, LPF refers to low and high frequencies in the Fourier space associated with the spectrum, not with the recorded signal.This space can be identified with the interferogram space.Applying an LPF to the spectrum is therefore equivalent to a shortening of the interferogram in combination with a Fourier interpolation of the spectrum.
Both methods employ the same general principle in that they are compression methods: the spectra are decomposed and components are discarded according to their expected signal-to-noise ratio.However, the mechanism by which this is achieved is different.The BB and DS measurements contain the instrument background and calibration blackbody radiation, both of which are spectrally slowly varying signals without discernible emission or absorption lines.This motivates and justifies the use of the LPF, which discards the higher frequencies, i.e. the rapidly varying components, because they can be expected to contain little or no signal but contribute equally to the noise level.The PCA method decomposes the data into uncorrelated components and retains only the most significant ones, discarding those that contribute little to the total variance of each spectral sample.
Because of the qualitatively different effects the two methods have on the data, a quantitative comparison of PCA and LPF proved difficult.In the following, we will make a qualitative comparison between a PCA reconstruction from 20 principal components (PCA-20) and a low-pass filter which dampens all but the 512 lowermost frequencies (LPF-512).As Fig. 11 shows, the first 20 PCA vectors reconstruct the simulated noise variance to 3.6 %, implying a reduction of the noise SD by over 80 %.
Figure 12 shows a reconstruction of a spectrum subjected to the PCA-20 and LPF-512 methods, as well as a third one, where the signal was first reconstructed from the principal components and subsequently low-pass-filtered.For a more detailed inspection of the differences, Figs. 14 and 13 show the same spectra within a reduced wave-number range.
The LPF-512 method effectively reduces the interferogram length from 8000 to 1024 real samples; that is, it retains 512 of 4001 complex modes.Assuming a uniform noise distribution among the modes, this implies a noise level reduction of 2.8 : 1, while for the PCA-20 method we expect a ratio of 5 : 1 from the simulation in Sect.4.4.A more fundamental difference, regardless of the number of principal compo- nents, is that the noise removed by the PCA method does not have a high-frequency bias, which is why the smoothing effect is less pronounced.
This difference between the PCA and low-pass methods is further illustrated by Fig. 15, which shows the magnitude of the complex reconstruction residual in Fourier space, for the LPF-512, PCA-20 and combined methods.The reconstruction residual for the pure low-pass method is an order of magnitude smaller for the lowermost frequencies, corresponding to the low-pass window employed.The PCA method, on the other hand, distributes the reconstruction residual more uniformly across the spectrum.This is also true for the combined method as long as the number of principal components is small.For larger numbers, the PCA-reconstructed spectrum approaches the original spectrum and the combined reconstruction according to the LPF-512 method.
The benefit of the combined approach is that most of the noise variance is already subtracted before the low pass is applied, and therefore the undesired interpolation of the lowfrequency noise modes is avoided.This results in a smooth calibration spectrum without the smoothing artefacts in the shape of slow oscillations.

Conclusions and outlook
One of the GLORIA limb sounder's key features is its imaging capability.In this paper, we presented two new ways to exploit this feature for the suppression of noise and certain artefacts in the calibration measurements.
The first technique was the empirical spatial characterisation of GLORIA's radiometric background, which helps identify and alleviate issues stemming from noisy calibra-   12) in the spectral range between 1100 and 1200 cm −1 .The pure low-pass method retains higher amounts of low-frequency noise, which manifest as slow oscillations in the spectrum.The PCA method subtracts noise variance irrespective of frequency (see also Fig. 15).tion measurements and an imperfect nonlinearity correction.Fitting a pseudo-hyperbola to the gain-calibrated radiometric offset allows us to use both blackbodies to generate calibration data of decent quality.Without the fit, the amplification of these errors by extrapolation of the offset would cause artefacts in the calibrated radiances, which would severely impact the L2 processing. .Spectrum of the reconstruction residual, i.e. the difference between the processed and original spectra, for different noise suppression methods.Shown are only the low frequencies.The LPF subtracts high-frequency noise from the spectrum, which means that low-frequency noise modes are retained.The PCA method introduces no such bias.The combined PCA + LPF method introduces no noticeable bias as long as the number of principal components is kept small.Note that the abscissa is in spatial units; this is because it is the Fourier space associated with the wave-number space.
We could significantly reduce the noise level of BB and DS calibration spectra by exploiting the spectral correlation of GLORIA's calibration scenes.Via truncated PCA, we were able to reconstruct the original measurements while reducing the noise level significantly.Compared with low-pass smoothing, an alternative option for GLORIA spectra, the PCA reconstruction is less smooth, but subtracts high-and low-frequency noise variance equally.By contrast, low-pass filters retain the low-frequency noise components, resulting in slowly oscillating perturbations in the result.A combined approach, where a low-pass smoothing is performed on previously PCA-compressed spectra, yielded the best results with suppressed low-frequency noise.
We expect that the techniques presented in this paper can also be adapted to other instruments with a focal plane array and similar absolute radiance calibration.The utility of the PCA-based method for other instruments would likely depend on the number of pixels available for each measurement, as well as the homogeneity of the calibration scenes.
In parallel with the new noise suppression techniques, a new L0/L1 data processing chain was developed.This chain enables high-performance processing of raw data into L1 spectra once the calibration parameters have been provided.Determining the parameters is difficult to automate, and is therefore still done under human supervision.Work on the integration of this task into the processing chain is currently in progress.

Figure 2 .
Figure 2. Real part of the NCO for different calibration methods at 791.875 cm −1 .The leftmost image shows the BB-DS reference without pseudo-hyperbola fit.The following two are the values for the original BB-BB calibration and the corresponding pseudohyperbola, and the last panel shows an application of the pseudohyperbola to the BB-DS calibration.The images in the second row show the difference between the NCO values and those for the BB-DS reference.An offset value given above a plot indicates that all values are to be understood as with this value subtracted.

Figure 3 .Figure 4 .
Figure 3. Magnitude of the IG for different calibration methods at 791.875 cm −1 .The panels are laid out the same way as in Fig. 2. The variation from the lower left to the upper right corner is caused by a composition change in the detector material.An offset value given above a plot indicates that all values are to be understood as with this value subtracted.

Figure 5 .
Figure 5. Averaged real part of the NCO over the detector rows for the different calibration methods.Both the unmodified data (suffix nh) and the pseudo-hyperbolas (suffix h) are shown.A discrepancy remains between the BB-DS and BB-BB methods which has to be compensated for in the trace gas retrievals.

Figure 9 .
Figure9.Reconstructions of the same DS spectrum as in Figs.6 and 8, in the spectral range between 1000 and 1100 cm −1 .The obvious discrepancies between the reconstructions and the original spectrum are due to residual emissions from stratospheric trace gases.

Figure 10 .
Figure 10.Normalised covariance eigenvalues of a DS, BB, and synthetic white noise data set.The eigenvalues indicate the rate at which each successive principal component adds to the total variance.

Figure 11 .
Figure 11.Cumulative normalised covariance eigenvalues, i.e. the cumulative sum of the values shown in Fig. 10.These indicate the relative variance restored after reconstructing the image from a given number of principal components.

Figure 12 .
Figure 12.Comparison between the PCA, low-pass and combined noise suppression methods, applied to the same DS measurement as in Fig. 6.

Figure 13 .
Figure13.Comparison between the PCA, low-pass and combined noise suppression methods (cf.Fig.12) in the spectral range between 1100 and 1200 cm −1 .The pure low-pass method retains higher amounts of low-frequency noise, which manifest as slow oscillations in the spectrum.The PCA method subtracts noise variance irrespective of frequency (see also Fig.15).

Figure 14 .
Figure14.Comparison between the PCA, low-pass and combined noise suppression methods (cf.Fig.12) in a limited spectral range.Between 920 and 960 cm −1 , a visible difference appears between the PCA methods and the pure LPF method.This difference is still within the original spectrum's noise range.