Euclid: Estimation of the impact of correlated readout noise for flux measurements with the Euclid NISP instrument

The Euclid satellite, to be launched by ESA in 2022, will be a major instrument for cosmology for the next decades. \Euclid\ is composed of two instruments: the Visible (VIS) instrument and the Near Infrared Spectromete and Photometer (NISP). In this work we estimate the implications of correlated readout noise in the NISP detectors for the final in-flight flux measurements. Considering the multiple accumulated (MACC) readout mode, for which the UTR (Up The Ramp) exposure frames are averaged in groups, we derive an analytical expression for the noise covariance matrix between groups in the presence of correlated noise. We also characterize the correlated readout noise properties in the NISP engineering grade detectors using long dark integrations. For this purpose, we assume a $(1/f)^{\, \alpha}$-like noise model and fit the model parameters to the data, obtaining typical values of $\sigma = 19.7^{+1.1}_{-0.8}$ e$^{-} \rm{Hz}^{-0.5}$, $f_{\rm{knee}} = (5.2^{+1.8}_{-1.3}) \times 10^{-3} \, \rm{Hz}$ and $\alpha = 1.24 ^{+0.26}_{-0.21}$. Furthermore, via realistic simulations and using a maximum likelihood flux estimator we derive the bias between the input flux and the recovered one. We find that using our analytical expression for the covariance matrix of the correlated readout noise we diminish this bias by up to a factor of four with respect to the white noise approximation for the covariance matrix. Finally, we conclude that the final bias on the in-flight NISP flux measurements should still be negligible even in the white noise approximation, which is taken as a baseline for the Euclid\on-board processing


Introduction
Euclid 1 is a Medium Class satellite mission to be launched by ESA in 2022. The Euclid satellite is mainly devoted to cosmology and intends to unveil the nature of the Dark Energy and the Dark Matter (Euclid Collaboration: Blanchard et al. 2020). These two components dominate the content of the Universe today (Planck Collaboration et al. 2020) and are responsible for the accelerated expansion of the Universe and large scale structure formation, respectively.
Euclid will perform a survey of 15 000 deg 2 of the extragalactic sky with two instruments in the visible and near-infrared (NIR) domains (see Laureijs et al. 2012, for details). The visible instrument, VIS, operates in the visible regime providing high quality images in a wide-band to carry out precise weak lensing This paper is published on behalf of the Euclid Consortium. e-mail: jimenez@lpsc.in2p3.fr 1 https://www.euclid-ec.org/ galaxy shear measurements. The NIR instrument, NISP (NIR spectrometer and photometer), provides photometric measurements in three NIR bands (Y, J, and H) and slitless spectroscopy in a wide NIR band using blue (920-1250 nm) and red (1250-1850 nm) grisms. The NISP is designed to provide photometric redshifts of about 2 billion galaxies, as well as spectroscopic redshifts for more than 50 million galaxies. In terms of performance the NISP should be able to observe very faint IR distant galaxies obtaining a signal-to-noise ratio of at least 5 down to magnitudes of 24.0 in about 100 seconds per filter. The NISP focal plane (Maciaszek et al. 2016) holds 16 NIR sensitive H2RG (HgCdTe astronomical wide area infrared imager or HAWAII) detectors supplied by Teledyne (Beletic et al. 2008) and selected by NASA. Each one consists of an array of 2048 × 2048 pixels with 2040 × 2040 photosensitive pixels, the remaining ones, called reference pixels, being used for tracking biases and temperature variations over long exposures (Moseley et al. 2010;Rauscher et al. 2012). Generally, near-infrared detectors, as those of the NISP, perform non-destructive exposures of hundreds of seconds in an up the ramp (UTR) mode (left panel of Fig. 1). For each exposure, non-destructive frames containing the full array data (i.e., all the pixels) are read (without instrumental reset) at about 1 Hz and sent to the acquisition program. The set of frames, for a given exposure, forms a ramp and the total flux for the exposure can be computed for each pixel through the slope of the ramp via a simple linear fit (Rauscher et al. 2012;Kubik et al. 2016). The final uncertainties in the measured flux are given by the uncertainties in each of the frames, which are due both to the readout noise and the photon noise. In order to reduce the noise in the estimation of the flux the multiple accumulated (MACC, right panel of Fig. 1) sampling mode is applied by averaging the UTR frames into groups that are then sent to the acquisition system (Secroun et al. 2016).
In the case of the NISP instrument the typical exposure time varies from about 87 s in photometric mode to 547 s in spectrometric mode, and the frame readout frequency is about 0.692 Hz. In flight the MACC mode is used but as the amount of data is too large for the Euclid telemetry the individual groups are not transferred to Earth. Indeed, only the flux (ramp fitted slope) for each pixel is available for further data processing. Due to on-board CPU limitations the slopes are obtained from an analytical solution of the linear fit. As a consequence to obtain an accurate flux estimate it is necessary to have an accurate description of the photon and readout noise. In the current Euclid baseline both the readout and photon noise are described by a white noise approximation (Fowler & Gatley 1990;Kubik et al. 2015a). The main properties of the readout noise assuming white noise are characterized during ground calibration and used in-flight. However, it has been found that individual H2RGs may present some level of correlated noise in the form of (1/ f ) α -like noise (Smadja et al. 2010;Kubik et al. 2015b). Such noise correlation might bias Euclid in-flight flux estimates.
In this paper we use ground calibration data to characterize the readout noise of the NISP detectors in terms of correlated (1/ f ) α -like noise and study possible bias in the flux measurements induced by departures from the white noise approximation in the linear fit. The paper is structured as follows. In Sect. 2 analytical expressions for the noise covariance matrix of the MACC groups and for the maximum likelihood solution of the linear fit in the case of correlated readout noise are derived. Section 3 describes the characterization of the readout noise for the cali-bration data of the NISP detectors in terms of correlated noise. In Sect. 4, we estimate the expected flux bias in the case of realistic correlated readout noise. Finally, we conclude in Sect. 5.

Flux estimation
Here we concentrate on the MACC readout mode, for which the acquired data consist of n groups with m frames each, and n − 1 drops (non acquired frames) 2 with d frames each: MACC(n, m, d). During flight EUCLID uses MACC(4,16,4) and MACC(15,16,11) for photometry and spectroscopy, respectively. Following Kubik et al. (2015a) we estimate the total flux from the group differences, ∆G k = G k+1 − G k .
For each group k, we define S (k) i as the total signal in a frame, f 0 as the signal accumulated between the last reset and the first readout of the pixel, and f (k) i as the signal between the frames i (k) and (i + 1) (k) . Furthremore, we notate ρ (k) i as the readout noise of the frame i (k) , and, we represent the non-transferred signal between the last frame of a group and the first frame of the next one, which is accumulated in the drops, by D (k) . Thus, accounting for signal and readout noise contributions, G k , the averaged measured signal for group k, is given by and then 2 Drops are necessary to reduce the processing time on-board.
Using these group differences, we can then derive a Maximum Likelihood estimator for the total flux, g. We use the following Gaussian approximation for the likelihood function where ∆G = {∆G k , k = 1, n − 1} is a vector gathering all group differences. D is the covariance matrix for the group differences, and |D| is its determinant.

Correlated noise covariance matrix
We discuss here the computation of the group difference noise covariance matrix in the case of correlated readout noise. We refer to Kubik et al. (2015a) for the white noise case. We will account for both photon and readout noises. With respect to photon noise, the flux integrated over a frame is Poisson distributed and stochastically independent between frames. This applies both to fluxes of frames within a group and within a drop. We can then write: The correlated readout noise is assumed to be Gaussian distributed with a (1/ f ) α -like spectrum in Fourier domain as in Kubik et al. (2015b). In the time domain this is equivalent to a Gaussian distributed noise described by a correlation function 3 of the form: C [|δ t |], where δ t is the time interval between two given frames. Thus, we can write where t frame is the frame integration time. Using Eqs. (4) and (5) the group noise covariance matrix (see Appendix A for details) reads for the diagonal terms and for the off-diagonal ones. By constrast to the white noise case described by Kubik et al. (2015a) we find in the latter expression a contribution from the readout noise. 3 The correlation function is the inverse Fourier transform of the spectrum.
We also derive the group difference covariance matrix. The diagonal terms are given by and the off-diagonal terms are for k < l. By contrast to the white noise readout noise case presented in equations (7) and (8) of Kubik et al. (2015a), we observe that the contribution of the readout noise to the group difference covariance is not constant in the diagonal terms and it adds extra correlation in the off-diagonal ones.

Readout noise measurements
The readout noise of infrared detectors can be characterized from long exposure ramps in dark conditions. Here, we use dark test data obtained during the Euclid NISP detector characterization performed at the CPPM laboratory. We focus on one of the NISP engineering grade H2RG detectors, which was cooled down to a nominal operating temperature of 85 K. The testing facility was designed to achieve best possible dark conditions and special care was taken to achieve expected in-flight readout noise.
For proper dark measurement, long integration UTR ramps were acquired, typically, ramps of 8000 frames with a total exposure time of 3.21 hours corresponding to a frame exposure time of t frame = 1.445 s. For each frame and for each of the 2040 × 2040 photosensitive pixels we use the reference pixels to remove correlations in the readout noise induced by background variations (Moseley et al. 2010;Rauscher et al. 2012). After reference pixel corrections we find that correlations between pixels are reduced as expected. Furthermore, we compute the dark for each ramp using the Fowler & Gatley (1990) algorithm, for which the slope of the ramp (in this case the dark contribution) is computed from the difference of the average of blocks of frames at the end and at the beginning of the ramp. In our case we have considered blocks of 32 frames to reduce the uncertainties in the dark measurements. Every ramp is corrected for the dark by subtracting the median dark value of all of the pixels in a given detector. For the data used in this paper the median dark for all pixels in the array was about 6×10 −3 e − s −1 with a standard deviation of 2 × 10 −3 e − s −1 across pixels.
The left panel of Fig. 2 shows one of these ramps for one of the inner pixels in the array after correcting for the reference pixels (raw data, red line) and after subtraction of the dark contribution (dark corrected, violet line). We have used a conversion gain factor of f e = 2 e − ADU −1 . We can observe in the figure that the readout noise is not fully white. This can be better seen in the right panel of Fig. 2, where we show the power spectrum of the dark corrected data as a function of the time frequency in Hz.

Readout noise modelling and fitting
The readout noise power spectrum shows a (1/ f ) α -like spectrum with an excess of power at low frequencies. Similar patterns are found for all other pixels in the array. Thus, to characterize the correlated readout noise in the NISP detectors we assume that its power spectrum is given by with f the time frequency and σ, f knee and α the parameters of the model. The σ parameter gives us information about the flat part of the power spectrum, which corresponds to the white noise contribution. α and f knee inform us about the correlated noise contribution. Notice that for α = 0 this model converges to a flat power spectrum and thus to a white noise spectrum. For each pixel in the array we fit the readout noise power spectrum to this (1/ f ) α -like model. We use the python mpfit module, which gives the best-fit parameters and their uncertainties. Uncertainties on the data power spectrum are computed assuming Gaussian noise: σ P( f ) ∝ P( f ). In practice, the fit is performed in two steps. First, we estimate the uncertainties in the power spectrum from a smoothed version of the readout noise power spectrum (see red line in the right panel of Fig. 2) and compute the best-fit parameters for the (1/ f ) α -like model. Then, we use these first estimates of best-fit parameters to estimate the uncertainties in the power spectrum and perform a second fit to the readout noise power spectrum. The best-fit parameters obtained from this second fit are stored for further analysis. Using Monte Carlo simulations we have observed that this two-step procedure leads to non biased estimates of the best-fit parameters for the(1/ f ) α -like model. In Fig. 2 we show the best-fit (1/ f ) αlike model (black line) to the readout noise power spectrum (blue line) obtained from the second fit. The best fit-parameters and their uncertainties for this pixel are σ = 20.50 ± 0.23 e − Hz −0.5 , f knee = (5.5 ± 0.8) × 10 −3 Hz and α = 1.17 ± 0.15. We observe that the best-fit model is consistent with the data with χ 2 /N dof of 1.57 and N dof = 8000 − 3.

Readout noise properties
We present in the left column of Fig. 3 four maps representing the best-fit parameters and χ 2 /N dof values for all the 2040×2040 photosensitive pixels for one of the ramps of one of the tested detectors. The white dots in the maps correspond to either hot pixels or pixels for which we obtain a bad fit to the data. These pixels represent less than 0.1% of the total pixels and are uniformly distributed in the maps. We can observe in the maps vertical bands which are related to the 32 readout channels, for which we expect some correlations in the noise properties. We can also isolate some particular regions as the one in the f knee map for pixels around (2000,1400), which are also found when computing other characteristic quantities of the detectors as for instance the Correlated Double Sampling (CDS) noise. They mainly correspond to manufacturing defects.
The 1D distributions of the best-fit parameters and the χ 2 /N dof are shown in the right panels of Fig. 3 excluding hot and bad-fit pixels. We show in the figure four ramps of the same detector, for which we find consistent results. We observe that the distributions for the three parameters are skewed towards large values. We find that the median values for the best-fit parameters are σ = 19.7 +1.

Verification via simulations
In order to validate our analytical expressions for the group and group difference covariance matrices in the case of correlated readout noise we have performed Monte Carlo simulations. We have generated a large number of realizations of fake NISP readout noise using the (1/ f ) α -like model discussed in Sect. 3. The correlated readout noise simulations are obtained via three steps: 1) we produce realizations of Gaussian white noise in real space, 2) we take the Fourier transform of those and multiply each Fourier component by the square root of the value of the power spectrum model at the same frequency, and 3) we compute the inverse Fourier transform of the modified Fourier components of the readout noise simulation. From these simulations of readout noise we have constructed fake NISP ramps by adding a cumulative flux contribution as well as the corresponding photon noise assuming a Poisson distribution. We consider here the MACC(15, 16, 11) readout mode as chosen for the spectroscopic in-flight operations, for which the effect of the correlated noise is expected to be larger. As an example, we present in Fig. 4 the group difference covariance matrix as obtained from Eqs. (8) and (9) (left panel), and the relative difference with respect to the one obtained from Monte Carlo simulations (right panel). The estimates of the covariance matrix were computed for the values of σ, f knee and α found in Sect. 3 for the NISP detector data. The incident flux is set to 1 e − s −1 . We observe very good agreement between the two estimates. We have repeated this comparison for various values of the parameters σ , f knee and α, and for different input fluxes, and obtained the same results. We therefore validate our analytical expressions.

White and correlated readout noise covariance matrices
In this paper, we are interested in studying how using a white noise approximation in the case of a correlated readout noise can impact the on-board estimation of the total flux measured by the Euclid detectors. Therefore, it is interesting to compare the covariance matrix one would obtain for the same correlated input noise in the white and correlated readout noise approximations discussed in Sect. 2.2. The correlated readout noise is obtained via Monte Carlo simulations. We generate mock timelines using the (1/ f ) α -like model discussed above with the set of averaged best-fit parameters presented in Sect. 3. The covariance matrix for the correlated noise approximation is computed as described in Sect. 2.2. For the white noise approximation we start by computing the effective root mean square (rms) of the readout noise. In practice we deduce it from the CDS noise estimated from the simulated timelines of correlated readout noise and impose σ white = CDS/ √ 2. The covariance matrix in the white noise approximation is computed following Kubik et al. (2016). In terms of the signal contribution we have considered two cases: 1) dark conditions and 2) sky background on nominal Euclid flight operations. For the dark conditions we assume an incident flux of 10 −3 e − s −1 , while for the sky background we consider 2 e − s −1 , which is within th 10 % of the expected in-flight background in the NISP Y, J and H filters including zodiacal light and scattered light.
The main results of this analysis are presented in Fig. 5, where we represent the group difference covariance matrix for the white (left) and correlated (middle) readout noise approximations, and the relative difference between the two (right). In the top and bottom panels of the figure we show the covariance matrices for the dark conditions and nominal sky background cases, respectively. For all the covariance matrices displayed we observe strong correlation between adjacent group differences as one would expect for non-destructive exposures. For the dark conditions we observe important differences between the white and correlated readout noise approximations. This is due to the fact that the readout noise dominates with respect to the photon noise. We mainly find that the difference between the diagonal and adjacent terms is increased in the correlated readout noise approximation with respect to the white noise one. Furthermore, other off-diagonal terms are not strictly zero for the correlated approximation by contrast to the white noise one. These differences would increase in the case of either steeper power spectrum or a larger f knee frequency. For the nominal background conditions we expect the photon noise contribution to be dominant and therefore we find that the differences between the white and correlated readout noise approximations are smaller. As before, these differences would depend very much on the slope of the readout noise power spectrum and on its f knee frequency. Therefore, we think that an accurate characterisation of these parameters will be needed during ground calibration and in-flight commissioning and performance verification phases.

On-board flux estimation
As discussed above during Euclid flight operations the white noise approximation will be used to estimate the sky flux in each of the NISP array pixels and only the sky image will be transferred to Earth. In the presence of correlated noise we expect the estimate of the flux to be biased. To evaluate this bias we have constructed mock simulations of the sky emission by assuming a constant sky signal and adding realistic realisations of the readout noise in the NISP detectors. For the latter we use the (1/ f ) α -like model and best-fit parameters discussed in Sect. 3, and the simulation procedure described above. We explore the sky signal in the range 0.1-100 e − s −1 , i.e., from low dark values to bright objects. For each value of the sky signal we construct 10 000 mock ramps. For each of the simulated ramps we estimate the sky flux from the group differences using: 1) the white noise approxima- tion flux estimator (WNA, hereafter) developed by (Kubik et al. 2016), and 2) the maximum likelihood approach discussed in Eq. 3 assuming the noise group difference covariance matrix presented in Sect. 2.2, (CNA, hereafter). In practice, the maximum likelihood flux estimate for 1) is obtained using equation 11 in Kubik et al. (2016). For the CNA case the maximum likelihood flux estimate is obtained using a simple grid approach. For illustration, an example of the reconstructed likelihood function for a sky flux of 0.21 e − s −1 is presented as solid black curve in Fig. 6. The best-fit value and input value are indicated as vertical dashed red and dashed blue lines, respectively. We observe that even when using the expected group difference covariance matrix for the correlated readout noise there is still a small bias in the estimate of the sky flux.
In Fig. 7 we show the relative bias, F out −F in F in , in percent for both the CNA (red line and dots) and the WNA (green line and dots) cases as a function of the background flux F in . Uncertainties in the measured bias are given by the filled red (CNA) and green (WNA) areas as computed from the Monte Carlo simulations. For low background flux values (below 1 e − s −1 ) we find that the maximum bias for CNA is under 1% and about four times smaller than the WNA one. For fluxes above 1 e − s −1 the bias in both cases are equivalent and below 0.1 %. For sky observations we expect a background flux of about 1-2 e − s −1 , in the region where the bias is expected to be small. However, we have observed using the dispersion over the set of Monte Carlo simulations that the WNA systematically underestimates the uncertainties by a factor ranging from 2 to 5. However, for CNA the dispersion on the simulations and the measured uncertainties are consistent.

Summary and conclusion
The Euclid satellite mission will be a fundamental tool for cosmology and infrared astronomy thanks to its near infrared photometer and spectrometer instrument, NISP, which consists of 16 NIR sensitive H2RG detector arrays of 2048 × 2048 pixels each. The NISP is designed to observe very faint distant galaxies and therefore requires a low and well characterized readout noise.
In this paper we have studied the readout noise associated with the NISP detectors taking advantage of long exposures (few hours) performed during laboratory dark tests at the CPPM cryogenic facilities. We have found that the NISP readout noise is correlated and can be well characterized by a (1/ f ) α -like model with a typical knee frequency of f knee = 5.2 +1.8 −1.3 × 10 −3 Hz and a low frequency component with slope α = 1.24 +0.26 −0.21 . From this we conclude that the readout noise of the NISP detectors has non-negligible correlation within the typical in-flight NISP exposure time (575 s).
Infrared instruments, and in particular the NISP, acquire data using the MACC readout mode, which consists of a series of non-destructive exposures averaged into groups that form a ramp. The input flux in the detectors can then be obtained from the slope of the ramp using maximum likehood estimators, which generally assume white readout noise (Kubik et al. 2015a;Kubik et al. 2016). Here, we have extended these estimators to the case of correlated readout noise. Analytical expressions for the group and group difference covariance matrices are presented for the case of (1/ f ) α -like correlated readout noise. These have been validated via Monte Carlo simulations.
Finally, we have performed Monte Carlo simulations of the in-flight expected NISP detector signal and noise, including a realistic background signal and correlated readout noise as measured on the ground calibration tests. From these simulations we have been able to estimate the expected bias in the on-board flux estimates during in-flight operations, for which white readout noise is assumed. We find that for low background the flux bias can be up to four times larger than when accounting for the correlation in the readout noise. Nevertheless, this bias is negligible for typical sky background signals. Therefore, we expect no significant bias in the on-board fluxes measured by Euclid.
We detail here the computation of the group noise covariance matrix forcorrelated readout noise described by 5. We concentrate in the correlated readout noise terms. Other terms can be found in Kubik et al. (2015a).
The group noise covariance matrix is given by for the diagonal term, and for the off-diagonal elements. And the readout noise is For the diagonal terms, l = k , the readout noise term is and then For the off-diagonal terms, l k , we have and then

. Diagonal Terms
We can write the diagonal terms of the group differences covariance as (A.7) We compute now each term of the correlated readout noise contribution: Finally, the group difference diagonal covariance matrix is (A.11)

A.1.2. Off-diagonal terms
The off-diagonal terms are given by (A.12) We compute now each term of the correlated readout noise contribution: Finally, the group difference off-diagonal covariance matrix is