Optical Energy Variability Induced by Speckle: The Cases of MERLIN and CHARM-F IPDA Lidar

: In the context of the FrenchGerman space lidar mission MERLIN (MEthane Remote LIdar missioN) dedicated to the determination of the atmospheric methane content, an end-to-end mission simulator is being developed. In order to check whether the instrument design meets the performance requirements, simulations have to count all the sources of noise on the measurements like the optical energy variability induced by speckle. Speckle is due to interference as the lidar beam is quasi monochromatic. Speckle contribution to the error budget has to be estimated but also simulated. In this paper, the speckle theory is revisited and applied to MERLIN lidar and also to the DLR (Deutsches Zentrum für Luft und Raumfahrt) demonstrator lidar CHARM-F. Results show: on the signal path, speckle noise depends mainly on the size of the illuminated area on ground; on the solar ﬂux, speckle is fully negligible both because of the pixel size and the optical ﬁlter spectral width; on the energy monitoring path a decorrelation mechanism is needed to reduce speckle noise on averaged data. Speckle noises for MERLIN and CHARM-F can be simulated by Gaussian noises with only one random draw by shot separately for energy monitoring and signal paths.


Introduction
Atmospheric methane is a greenhouse gas responsible for about 20% of the additional radiative forcing due to human activities since the industrial revolution [1]. In order to increase knowledge on atmospheric methane burden, the MEthane Remote LIdar missioN (MERLIN) is being developed jointly by CNES (Centre National d'Études Spatiales) and DLR (Deutsches Zentrum für Luft und Raumfahrt) [2][3][4] for a scheduled launch in 2024. It plans to deliver as level 2 products: atmospheric methane total weighted columns (XCH 4 ) with their associated weighting functions describing the vertical sensitivity of the measurement to CH 4 variation along the atmospheric column. It targets to achieve on XCH 4 systematic (deterministic or correlated) errors less than 3 ppb and random (stochastic or uncorrelated) errors less than 22 ppb for 50 km average when the mean methane content value is 1780 ppb. That corresponds to a signal-to-noise ratio (SNR) value (the inverse of the relative random error) of 81. This challenging accuracy is needed to improve the estimates of methane emissions and sinks [5] and implies to track all the sources of variability, uncertainties and biases in the measurement.
Speckle appears when coherent light is scattered by heterogeneities at the scale of its wavelength like variations in the surface roughness or in the refractive index. Scattered waves propagate along different optical paths and interfere in any observation plan showing patterns with granular structure of alternately light and dark spots [30,31]. For IPDA lidar speckle is a serious contributor to the SNR as it induces energy fluctuations on the detector [32][33][34][35][36][37][38][39][40][41].
In Section 2 of this paper, speckle contribution to SNR on XCH 4 is addressed and with MERLIN and CHARM-F specifications the corresponding speckle characteristics are determined. In Section 3, for the two instruments, speckle noise level and shot noise level are compared for one shot and after averaging. And in Section 4, a way to simulate energy fluctuations on a lidar detector due to speckle with some hypotheses on the speckle pattern decorrelation in time is described.

How Does Speckle Contribute to Signal-to-Noise Ratio on XCH 4 Measurements?
Double-pulsed IPDA lidar emits coherent light pulses and collect the energy scattered back by the ground surface. The atmospheric methane total weighted column (XCH 4 ) is estimated from the difference in atmospheric transmission between two laser pulses: one (ON) at a wavelength selected to have a high methane absorption, the other (OFF) as reference at a wavelength with significantly less methane absorption. The two pulses are close enough in wavelength and emitted with a small enough time interval in order to consider atmospheric, ground and instrumental optical properties to be identical except for CH 4 absorption. Nevertheless, differences in H 2 O and CO 2 absorption between the two beams are accounted. In a typical case, these differences are in the order of a few 10 −4 to compare with the DAOD due to methane in the order of 0.4 to 0.6.
From the vertical DAOD due to methane (DAOD CH4 ), XCH 4 is obtained as an averaged value of methane dry-air volume mixing ratio with associated weighting function WF(p): where µ refers to the incident angle departure from the nadir, P surf is the surface pressure, and where σ CH4 is the absorption cross sections of a mole of methane depending on wavelength, pressure p and temperature T, g is the acceleration due to gravity, M dry is the average molar mass of dry-air, M H2O is the molar mass of water vapor and X H2O is water vapor dry-air volume mixing ratio. WF(p), DAOD H2O and DAOD CO2 are computed from data provided by meteorological centers and spectroscopic data found in GEISA database [42] with specific improvements [43,44]. Slant DAOD is determined as the ratio of the backscattered energy measurements for pulses ON and OFF (P on and P off ) normalized by emitted energy measurements for each pulse (E on and E off ) to deal with the fluctuations of the energy delivered by the laser source [40,41].
The detection chain acquires these energies together with the solar flux energy P sun which is also acquired on its own. The actual measurements add random noises: electronic noise due to electric fluctuations of intensity or voltage caused by electrons movements in the electronics, shot noise due to fluctuations of time occurrence for electrons release from incident photons caused by the corpuscular nature of light and speckle noise due to fluctuations of the optical energy on detector caused by interference reflecting the wave nature of light. In this study, the focus is put on the speckle, the major contribution on noise coming from the electronic chain is not analyzed. Comprehensive analyses of IPDA lidar noises are developed in [11,33,38].
Mean energy of scattered light by ground is space distributed according to the ground Bidirectional Reflectance Distribution Function [45]. However, for quasi monochromatic light, it fluctuates because of interference. If the light is not quasi monochromatic, interference exists for each wavelength but with a negligible effect on the spatial and temporal distribution of the total energy. For geophysical surfaces speckle pattern is fully developed. For quasi monochromatic light, relative energy fluctuations through a given aperture S during a given time T can thus be linked to the observation geometry and the wavelength through the coherence area S c and time τ c . Signal-to-noise ratio due to speckle (SNR sp ) is then given, with P the polarization index of the light, by (cf. Appendix A) In this analysis speckle due to scattering from aerosol and turbulence is not taken into account [35,37]. For an instrument, SNR sp can be compared with the shot noise SNR due to the random noise related to the photoelectric effect and to its amplification in the avalanche photo-diode used as a detector. This noise (SNR sn ) can be expressed as a variation of the incoming flux [11]: where N is the number of photons constituting the optical energy flux, η is the detector quantum efficiency, F is the noise factor of the avalanche photo-diode, it refers to the noise due to the electrons multiplication process (It is one plus the ratio between the variance and the square of the mean of the probability for a primary electron to give secondary electrons). From speckle and shot noise variability, assuming that the uncertainties are not correlated, slant DAOD uncertainty is given by: where, for each flux, the variance σ 2 Xxx is σ contribution to the other measurements and its own measurement. Neglecting uncertainties coming from the integrated weighting function WF or from DAOD H2O and DAOD CO2 estimates and without taking into account the full relation involving the electronic detection response function, the data sampling for ground transmission and the ground processing, the SNR for XCH 4 due to uncertainties of optical fluxes on the detector (SNR of ) is given by: Atmosphere 2019, 10, 540 4 of 13

MERLIN and CHARM-F Data Sets
For satellite or airborne IPDA lidar observation, the ground is counted as a plane rough surface illuminated by a coherent source and it is modelled as an extended source of chaotic light of intensity I G (∆x, ∆y, t) producing a fully-developed speckle pattern resulting from interference in the propagating space ( Figure 1). The field of scattered light is observed at distance z of the rough surface and the area of coherence at this location is determined using results from Appendix A. The instrument has relative speed to ground v. It emits wavelengths in (λon) and near (λoff) a methane multiplet. The emitted beam is polarized and its divergence is divbeam. The receiver telescope is an afocal design with huge magnification. It consists of two conical mirrors M1 and M2 plus an achromatised ocular lens which generates an image of the entrance pupil. The collecting mirror M1 is elliptical. To obtain a compact design, the secondary mirror M2 is close to M1, and M2 partly vignettes the entrance pupil. Up to the detector the focal length of the receiver chain is frec. A band pass filter is incorporated on the optical path to block the light outside a narrow window with size Lfilter. The incoming signal light scattered back by the earth is focused on the detector of diameter dAPD which collect also the calibration light from the energy monitoring path [12]. The use of the same detector for signal estimate and energy monitoring avoids variations in the optical-to-electrical response between the two, but the energy extracted from the emitted beam has to be reduced by several orders of magnitude to match the energy level of the lidar returns as the detector performances are for a limited energy range. For this purpose, integrating-spheres with fiber coupling are used [40]. The cut off frequency of the electronic detection is smaller than the sampling frequency νsample. Table 1 gives the values provided by manufacturers [6,14] for parameters which allow to compute speckle impact.  The instrument has relative speed to ground v. It emits wavelengths in (λ on ) and near (λ off ) a methane multiplet. The emitted beam is polarized and its divergence is div beam . The receiver telescope is an afocal design with huge magnification. It consists of two conical mirrors M1 and M2 plus an achromatised ocular lens which generates an image of the entrance pupil. The collecting mirror M1 is elliptical. To obtain a compact design, the secondary mirror M2 is close to M1, and M2 partly vignettes the entrance pupil. Up to the detector the focal length of the receiver chain is f rec . A band pass filter is incorporated on the optical path to block the light outside a narrow window with size L filter . The incoming signal light scattered back by the earth is focused on the detector of diameter d APD which collect also the calibration light from the energy monitoring path [12]. The use of the same detector for signal estimate and energy monitoring avoids variations in the optical-to-electrical response between the two, but the energy extracted from the emitted beam has to be reduced by several orders of magnitude to match the energy level of the lidar returns as the detector performances are for a limited energy range. For this purpose, integrating-spheres with fiber coupling are used [40]. The cut off frequency of the electronic detection is smaller than the sampling frequency ν sample . Table 1 gives the values provided by manufacturers [6,14] for parameters which allow to compute speckle impact.

Speckle Contributions to Signal-to-Noise Ratio on MERLIN or CHARM-F XCH 4 Measurements
Laser pulses energy are Gaussian distributed on the ground. The energy of scattered light is thus modeled as: where I 0 the mean intensity and σ R the spatial standard deviation. That gives by integration in the field of view (FOV) of the receiver with Equation (A16) the effective area of the footprint looked as a secondary source for laser returns: Solar energy is uniform in the FOV. The energy of scattered solar flux is thus modeled by a top hat distribution: I sun G(∆x,∆y) = I sun 0 inside the FOV, 0 elsewhere, so, (still using Equation (A16)) the effective area for sun light is given by: Because of the movement of the satellite (or of the aircraft) during the signal energy estimation there is time speckle renewal due to changes in the illuminating of the ground at the wavelength scale [36]. With d eff the characteristic size of the coherent area (the diameter for a circular area), the characteristic time of this renewal is d eff /v [32] which remains considerable even for space lidar (in the order of 0.1 ms). Moreover, as pulsed lidar are used, regardless of the laser beam characteristic time, there is a full coherence inside one pulse. There is no averaging of the interference pattern with increasing integration time as there is no more signal. Table 3 provides both for MERLIN and CHARM-F the effective surface for incoherent source on ground, the coherence area, the coherence time and the number of speckles. The number of temporal speckle for solar flux is estimated for a discretisation time δt dis settled between a tenth of the sampling time to a sampling time. Even if for simulation purposes the time signal is discretised at a tenth of the sampling time, for each discretisation interval there are several hundreds of temporal speckles for solar flux.
Subjective speckle for laser and solar fluxes are not to be counted as the detection optics is built in such a way that all the photons collected by the entrance pupil reach the detector. So, speckle only modulates the spatial distribution of the energy on the detector and not its integrated value.
Along the optic calibration path (see Figure 2), speckle induces fluctuations of the energy amount at each location where there is energy dilution, specifically the output of the integrated spheres and the entrance of the optical fibers.

Coherence surface for solar flux
Sc sun = 4/π (λ × frec/dAPD) 2  The number of temporal speckle for solar flux is estimated for a discretisation time δtdis settled between a tenth of the sampling time to a sampling time. Even if for simulation purposes the time signal is discretised at a tenth of the sampling time, for each discretisation interval there are several hundreds of temporal speckles for solar flux. Subjective speckle for laser and solar fluxes are not to be counted as the detection optics is built in such a way that all the photons collected by the entrance pupil reach the detector. So, speckle only modulates the spatial distribution of the energy on the detector and not its integrated value.
Along the optic calibration path (see Figure 2), speckle induces fluctuations of the energy amount at each location where there is energy dilution, specifically the output of the integrated spheres and the entrance of the optical fibers. For CHARM-F it is assumed that the detector does not truncate the fiber output [14], but for MERLIN in order to avoid alignment issues, the illuminating area is larger than the detector size. Therefore, subjective speckle also occurs. For MERLIN, AIRBUS has done a dedicated study on preliminary design for this calibration path and estimated SNREon/Eoff around 43 mainly from detector and fiber entrance face [personal communication]. For CHARM-F, the SNREon/Eoff could be estimated from [40] around 59 mainly from fiber entrance face. Table 4 summarizes SNR estimates obtained with Equation (4) and data from Table 3 for the different fluxes in the case of MERLIN and CHARM-F. Pon and Poff are polarized, Psun is not, Eon and Eoff neither because integrated-spheres depolarize [46]. The maximum number of backscattered photons for MERLIN has been estimated to less than 18000 and for CHARM-F it is about 3550 time more (because the distance between ground and the receiver is 8.5 km instead of 506.3 km). The relative random error (RRE) is the inverse of the SNR and can be expressed in percentage form. For CHARM-F it is assumed that the detector does not truncate the fiber output [14], but for MERLIN in order to avoid alignment issues, the illuminating area is larger than the detector size. Therefore, subjective speckle also occurs. For MERLIN, AIRBUS has done a dedicated study on preliminary design for this calibration path and estimated SNR Eon/Eoff around 43 mainly from detector and fiber entrance face [personal communication]. For CHARM-F, the SNR Eon/Eoff could be estimated from [40] around 59 mainly from fiber entrance face. Table 4 summarizes SNR estimates obtained with Equation (4) and data from Table 3 for the different fluxes in the case of MERLIN and CHARM-F. P on and P off are polarized, P sun is not, E on and E off neither because integrated-spheres depolarize [46]. The maximum number of backscattered photons for MERLIN has been estimated to less than 18000 and for CHARM-F it is about 3550 time more (because the distance between ground and the receiver is 8.5 km instead of 506.3 km). The relative random error (RRE) is the inverse of the SNR and can be expressed in percentage form. Speckle noise is of the same order but smaller than shot noise for MERLIN lidar. Conversely, for CHARM-F lidar, speckle noise is dominant compared to shot noise which becomes negligible with the smaller distance between the receiver and the ground.
Using a mean value of 0.53 for DAOD and a mean value of 1780 ppb for methane content the random uncertainties coming from speckle on XCH 4 estimated with value of Table 4 and Equation (7) are summarized in Table 5. Noise for 7 s average corresponds to 50 km average for MERLIN and to 2 km average for CHARM-F. The random error of 5 ppb from speckle is compatible with the specification of MERLIN random error less than 22 ppb [2] but shows how difficult it will be to reach the target user requirement for random error estimated at 8 ppb.

Discussion and Conclusions
For any satellite or airborne IPDA lidar, speckle induces variability of incident energy on the detector both for atmospheric and energy monitoring branches. Subjective speckle has no impact on the energy absorbed by the detector as long as it does not truncate the energy collected. SNR from speckle can be estimated as above (Table 4) on a shot by shot basis. The associated noises are smaller with respect to the electronic noise, for which averaging is needed. There is no difficulty in calculating such an average for simulated data even if geophysical variability should be taken into account to be realistic [47]. But the noise correlation between the various samples describing one shot and between different shots must be included.
For signal path the speckle pattern changes like the target from one shot to another. If the speckle pattern changes there is no correlation between the noises. But for energy monitoring path, the speckle pattern may be highly correlated from one shot to another. And to limit speckle impact on random noise for XCH 4 spatially averaged, it is necessary to either fully stabilize the speckle pattern, or to deliberately change it over time [48,49]. This second solution is easier to do on mobile lidar and the first one may induce systematic bias on the methane estimation. Using moving parts, specific systems have been designed to change speckle pattern between successive pulses. Tests with CHARM-F showed resulting noise exhibiting pure white noise behavior which can be reduced by averaging [40]. The introduction of such a mechanism prevents any correlation between speckle simulated noises for one shot and another [41].
Simulations of IPDA data with realistic speckle noise may be performed after computations similar to those in Section 3, and with the two hypotheses of full correlation for the samples of one shot (because the noise is fully correlated at this time scale as the speckle pattern is constant during the full registration of one pulse) and of full decorrelation between shots.
In summary, nothing has to be done for sun flux. Only Gaussian random noise has to be added, both on each sampling of energy monitoring fluxes and on each sampling of signal fluxes. Only one random draw is needed for energy monitoring path and signal path per shot. However, a new random draw is necessary for each shot because the speckle pattern vary from one pulse to the other according to the use of the decorrelation mechanism acting between two pulses.
For the MERLIN mission, the SNR on optical energy of backscattered signals due to speckle is around 60, always higher than the shot noise SNR which is less than 50. On the contrary for CHARM-F, speckle SNR is smaller than the shot noise SNR. Speckle occurring on the energy monitoring path is associated to pulse by pulse SNR about 40 for MERLIN (60 for CHARM-F) and so it is the dominant speckle impact even with a decorrelation mechanism. A Gaussian model of these noises has to be included in the Merlin simulator.

Conflicts of Interest:
The authors declare no conflicts of interest.

Appendix A. Speckle Theory
The speckle theory has been developed by J. C. Dainty [30] and J. W. Goodman [31] from previous studies on coherent light [50]. For any beam, its intensity (or brightness) is defined from the amplitude of the electromagnetic wave u at location P i and time t i as: and its degree of first order coherence (g (1) ) as the normalized first order correlation function between point P 1 at time t 1 and point P 2 at time t 2 [51]: which vary between 0 for incoherent light and 1 for first order coherent light. In this equation (A2) g (1) is the spatial complex factor of coherence listed as µ (P 1 , P 2 ) at a given time (t 1 = t 2 ) and the temporal complex degree of coherence listed as γt) at a given location (P 1 = P 2 )). It measures the visibility of interference fringes. Similarly, the beam degree of second order coherence g (2) is defined as the normalized second order correlation function: which vary between 1 and +∞. A monochromatic beam is said chaotic when it can be represented as a Gaussian random process [52] resulting of the well-known random walk in the complex amplitude plane. Thermal light and coherent light after full scattering are chaotic lights. Chaotic light interferes with itself creating speckle pattern. The questions are to determine speckle sizes and the energy collected by a detector put in this speckle field.
For chaotic light, negative exponential relationship is found for the distribution of brightness: where I 0 is the statistical mean value of brightness. The most probable brightness for a speckle is zero and there are more dark speckles in the field than speckles of any other brightness, but there are also rare very bright speckles. The mean intensity of illumination W for an area S during a time ∆T results from the integration of I (P; t) taking into account its spatial and temporal correlations is given by: Similarly the second-order moment is given by: A function S(P) with values between 0 and 1 allows to take into account differences in the contribution of the various points P of the area S to the signal e.g., the sensitivity of the detector or Atmosphere 2019, 10, 540 9 of 13 some movement of the detector during a window time ∆T. In that last case, for a movement in the y direction the function is: S(x, y, t) = S(x, y − vt, 0), The mean value W and the standard deviation σ W are computed using M [31,37] defined by: where R s (∆ξ, ∆η) := S(ξ, η)S(ξ − ∆ξ, η − ∆η)dξdη and Λ(x) = 1 − |x| for |x| < 1, zero otherwise. For chaotic light, the correlation functions verify the following relation [51]: then the SNR is related to M by: where P refers to the polarization index (P is defined for a wave as the ratio of the intensity of the polarized component to the total intensity, |P|= 1 if the light is fully polarized, P = 0 if it is not polarized). M is very well approximated as follow [31] where and Using Wiener-Khintchine and Zernike-Van Cittern theorems, τ c and S c may be computed from the statistics properties of the complex amplitude of the incident light. The characteristic time can be computed from the spread of its frequencies [52]. For a Gaussian beam characterized by its standard deviation σ ν , characteristic time τ c is given by: and, for thermal light through a bandwidth filter of size L filter around mean frequency λ, the characteristic time τ c is given by: On a plane parallel to the emitting surface, the speckle dimensions can be computed from the brightness distribution of the emitted light on the scattering area [31]: where z is the distance between the source and the detector and S eff the effective surface of the source corresponding to the emitting surface S in case of uniform illumination. The coherence area increases when the wave packet propagates as a result from the mix of waves coming from different points of the incoherent source. Assuming speckle and effective area are circular, the coherence area diameter is given by: Furthermore, for a uniform illumination over a circular aperture with diameter d, Fourier transform gives the brightness distributed according to the first-order Bessel function of the first kind which first zero provides the size of coherence area: The small difference between factor 1.27 and factor 1.22 is generally neglected, but in the first case the full area of coherence is counted and in the second case only the disc inside the main lobe. So, the first approach seems more accurate for energy estimation and the second for speckle size measurements.
For completeness, speckle sizes can be estimated not only on a plane parallel to the source but also on plane with angle θ from this direction (see Figure A1).
Atmosphere 2019, 10, x FOR PEER REVIEW 10 of 13 when the wave packet propagates as a result from the mix of waves coming from different points of the incoherent source. Assuming speckle and effective area are circular, the coherence area diameter is given by: Furthermore, for a uniform illumination over a circular aperture with diameter d, Fourier transform gives the brightness distributed according to the first-order Bessel function of the first kind which first zero provides the size of coherence area: The small difference between factor 1.27 and factor 1.22 is generally neglected, but in the first case the full area of coherence is counted and in the second case only the disc inside the main lobe. So, the first approach seems more accurate for energy estimation and the second for speckle size measurements.
For completeness, speckle sizes can be estimated not only on a plane parallel to the source but also on plane with angle θ from this direction (see Figure A1).
It is usual to name as "objective speckle" the pattern observed on a surface due to the interference of scattered waves coming from different points of a rough surface illuminated by laser light. When such a speckle pattern is imaged using an optical system, the entrance pupil may be counted as a source of chaotic light and the interference between waves coming from different points of the imaged plane build in the image plane a new speckle pattern named "subjective speckle". The statistical properties of the scattering object determine speckle distribution in the observation plane (the area lightened). However, it is the size of the entrance pupil that determines speckles dimensions (the smaller the entrance pupil, the bigger the speckles). The subjective speckle size is obtained with the same relation than above, with the numerical aperture NA instead of the one-half angular aperture: Thus, using the characteristics of a system which images the surface giving rise to the speckled the aperture diameter, f the focal distance and G the transverse factor of magnification-the coherence size is given by: .
Then for an uniform source of chaotic light speckle sizes are given in the three directions by the following lengths [53]: It is usual to name as "objective speckle" the pattern observed on a surface due to the interference of scattered waves coming from different points of a rough surface illuminated by laser light. When such a speckle pattern is imaged using an optical system, the entrance pupil may be counted as a source of chaotic light and the interference between waves coming from different points of the imaged plane build in the image plane a new speckle pattern named "subjective speckle". The statistical properties of the scattering object determine speckle distribution in the observation plane (the area lightened). However, it is the size of the entrance pupil that determines speckles dimensions (the smaller the entrance pupil, the bigger the speckles). The subjective speckle size is obtained with the same relation than above, with the numerical aperture NA instead of the one-half angular aperture: Thus, using the characteristics of a system which images the surface giving rise to the speckle-d the aperture diameter, f the focal distance and G the transverse factor of magnification-the coherence size is given by: