The thermal emission of Saturn’s icy moons-Effects of topography and regolith properties

Context. The effects of space weathering and other alteration processes on the upper surface of Saturn’s icy moons are yet to be explored. Aims. We present a thermophysical model parametrised by way of regolith properties such as porosity, grain size, and composition, as well as the local topography. The modelled surface temperature and apparent emissivity are intended to be compared to measurements taken by Cassini’s Composite Infrared Spectrometer (CIRS), using its focal plane FP1. We study how they are impacted by the topographic model and the regolith properties. Methods. As an example, we coupled the topography of the Dione moon with our model. Simulations provide the thermal history of the surface elements of the shape model included in the FP1 footprints at the viewing geometries along one CIRS observation. The heat transfer in the regolith may occur through conduction or radiation. Its bolometric albedo, A, and hemispherical emissivity, εh, are expressed as a function of grain properties. Results. The model roughly reproduces the observed variations of surface temperature, TF, and apparent emissivity, εF, in the chosen example, while assuming uniform regolith properties. The dispersion of temperatures within the footprints due to the difference in local time of the surface elements explains most of the directionality of the apparent emissivity, εF (Em), at emission angles of Em ≥ 30◦. Adding topography at the 8-km scale amplifies this effect by a few percent. Refining the scale to 1 km increases it again by a single percent but at a high computational cost. This particular anisotropy of εF (Em) cannot be explained by the directional emissivity, εd, of the regolith. The temperature TF is less affected by this dispersion or by the topographic resolution. Adding regional variations of grain size significantly improves the agreement between the model and observations. Conclusions. This model demonstrated its good performance and, thus, it is ready for testing current hypotheses on regolith processing by space weathering on Saturn’s icy moons, such as regional changes in grain size.


Introduction
Beyond large-scale meteorite impacts or geological activity, there are many processes at work that effectively alter the surfaces of Saturn's atmosphere-less icy moons. Irradiation by solar ultraviolet photons or bombardment by magneto-spheric plasma particles are among the processes known to alter the regolith properties through chemical alteration, amorphisation, sputtering, or molecular desorption (Johnson et al. 2013). The regolith may also be polluted by exogenous elements originating from Saturn's E ring or Phoebe ring . The structural alteration of the surface over time and the intricacy of these processes are difficult to constrain (Baragiola 2003). However, they may be observed through surface features or radiometric asymmetries.
Thanks to numerous multi-wavelength observations, the Cassini mission (2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014)(2015)(2016)(2017) has provided a much more detailed picture of the surface of these moons. Asymmetric patterns between their leading and trailing hemispheres have been confirmed or detected, either in terms of colour and Bond albedo at visible and near-infrared (VNIR) wavelengths or in terms of thermal inertia at thermal infrared wavelengths (Hendrix et al. 2012;Howett et al. 2010). The brightness increase observed in the VNIR suggests that the E ring icy grains could impact these surfaces by sandblasting or by deposition (Verbiscer et al. 2007;Kempf et al. 2018;Howett et al. 2018 for a review). The Cassini VIMS (Visible and Infrared Mapping spectrometer) data provide evidence that the leading hemisphere of Dione, Tethys, and Rhea exhibit deeper absorption bands of water ice, suggesting the presence of fresh E ring particles (Scipioni et al. 2013(Scipioni et al. , 2014Stephan et al. 2010Stephan et al. , 2012Stephan et al. , 2016. Hendrix et al. (2012) suggest that the regolith grains are smaller on the trailing hemisphere of Mimas than on its leading hemisphere, which is consistent with a bombardment with E ring particles.
The alteration of icy surfaces by high energy electrons produces lens-shaped patterns at the equator of Mimas, Tethys, and Dione (Paranicas et al. 2001(Paranicas et al. , 2014. In Saturn's system, keV electrons are expected to bombard their trailing hemisphere while MeV electrons should collide with their leading hemisphere. The extent of the impacted areas varies according to the moon (Paranicas et al. 2014;Nordheim et al. 2017;Howett et al. 2018). Lens-shaped patterns have been discovered in the leading hemispheres of icy moons that look darker in the visible and brighter in the IR/UV colour ratio maps . They are correlated with thermal anomalies observed by CIRS and interpreted as higher thermal inertia in the leading hemispheres (Howett et al. 2011(Howett et al. , 2012. No such anomalies were discovered on Rhea (Howett et al. 2014). It is indeed expected that MeV A8, page 1 of 11 electrons alter the surface down to a depth of a few centimetres, where the diurnal thermal waves penetrate and the thermal inertia anomalies are detected (Paranicas et al. 2014;Schaible et al. 2017;Nordheim et al. 2017). However, the predicted lens-shaped features in the trailing hemisphere, induced by below-MeV electrons at shallower depths, have not been observed and may be hidden by some competing process such as E ring deposition or cold plasma bombardment (Nordheim et al. 2017;. The darkening and reddening observed on the trailing hemispheres of Tethys, Dione, and Rhea  may result from radiolysis induced by colliding protons or ions within this cold plasma (Johnson et al. 2013;Schenk et al. 2018;Howett et al. 2018). The VIMS near-infrared spectra did not reveal any change in size or composition correlated with the lenses (Hendrix et al. 2012).
If the observed lens-shaped patterns can be explained by an alteration of the surface by MeV electrons, the mechanism of grain sintering and the reason for thermal inertia increases in the icy regolith remain poorly understood (Ferrari & Lucas 2016;Schaible et al. 2017). To explore this direction further, we propose a thermophysical model that relates thermal inertia to regolith properties, namely: grain size and composition, layer structure, porosity, and the quality of grain contacts.
The thermal emission of Saturn's icy moons has been thoroughly observed by the CIRS instrument during Cassini's flybys. Analysis of this voluminous dataset by Howett et al. (2010Howett et al. ( , 2011Howett et al. ( , 2012Howett et al. ( , 2014Howett et al. ( , 2016 yielded major advances in the detection of regional thermal patterns. It brought new constraints on the Bond albedo and the thermal inertia, and helped build scenarios on their origin. Thermal inertia were found to be very weak, around 10-20 J m −2 K −1 s −1/2 (hereafter noted as MKS). These constraints were derived assuming smooth surfaces, unit emissivity (ε obs = 1), and vertical homogeneity of the icy regolith. The hypothesis made on emissivity arises from the fact that most of the analyses were performed on daytime data of the CIRS mid-infrared focal plane FP3 to maximise the spatial resolution, but at the cost of a degeneracy between temperature and emissivity, FP1 measurements being used to determine night-time temperatures. High porosity of the regolith is generally cited to explain such a low thermal inertia. Ferrari & Lucas (2016) related such low values to the presence either of amorphous ice grains, whose thermal conductivity is much lower than that of the crystalline form, or by the presence of rough grains with irregular faces that are likely to reduce the solid conduction through grain contact. They showed that such low values can then be obtained with normal porosity of particle beds. Thermal emission largely depends on regolith albedo and emissivity that balance the absorbed and emitted radiation on the surface and that depend on both grain size and composition. Therefore, the unit emissivity assumption may lead to an underestimation of the measured surface temperature and impact the derivation of thermal inertia.
The topography, mesoscale roughness, and porosity of the upper regolith layers also impact the energy deposition on the surface and the local thermal history, mainly due to the various slopes and shadowing effects (Spencer 1990;Lagerros 1996;Ferrari 2018). Davidsson et al. (2015) investigated the effect of topography on thermal emission with artificial topographies. They showed that it strongly depends on the roughness type, even when the degrees of roughness are identical. They outlined the bias introduced by these topographic models on the determination of thermal inertia and surface roughness. So it is also expected that the temperature observed in the CIRS field-of-view (FOV) depends on the local topography.
In this paper, we put together the thermal model developed by Ferrari & Lucas (2016), which is controlled by regolith grains and layer properties, and the known topography of Saturn's icy moons at various spatial resolutions ranging from 1 km to 8 km. In Sect. 2, we couple heat transfer in the near-surface with topography to predict the thermal emission as a function of regolith properties, topography, the area covered by the footprints of the FOV on the surface, and the epoch. The effect of these parameters on the directional and temporal variations of Dione thermal emission is investigated in Sect. 3. Finally, the ability of the model to mimic the measurements along a latitudinal scan on the surface of this moon around noon local time is evaluated.

Thermal emission of icy regoliths
The thermal emission of Saturn's icy moons was captured by the CIRS focal plane FP1 in the wavelength range from 17 to 1000 µm (Flasar et al. 2004). The footprints corresponding to the 3.9 mrad-wide FOV of the instrument vary from a few tens to a few thousand kilometres (Flasar et al. 2004). They comprise different local slopes and orientations relative to the sun and the instrument, so that various thermal histories and temperatures are captured at the time of observation. Thermal emission is a mixture of black/gray body thermal emissions. The temperature, T F , and the apparent emissivity, ε F , measured in the FOV may vary as well with the viewing geometry. We describe the articulation between surface models and a thermal model to account for the effect of topography and regolith properties on the thermal emission. We take advantage of an observation of Dione to study the impact of the complexity of the shape model on the measured T F and ε F .

Topography of icy moons
Thanks to the lifespan of the Cassini mission and the multiangular viewing geometries during its flybys, shape models have been derived for some icy moons (Gaskell et al. 2008). Specific formats (containers) known as DSK (Digital Shape Kernel) have been generated to use these models within the NASA's Navigation and Ancillary Information Facility (NAIF) SPICE system environment of observation ephemeris and geometry (Acton et al. 2018). These shape models are made of triangular plates and vertices positioned in the inertial body reference frame. The DSK environment provides routines to compute the position, visibility and illumination angles of any plate at a given epoch. The Python wrapper (SpicePy, Annex et al. 2020) of the NAIF/SPICE software package N066 was used to link the shape file with the thermal model and the CIRS pointing vectors. The surfaces have be tessellated at various ground sample distances (GSD) and vectors are stored in an implicitly connected quadrilateral format (ICQ, Gaskell et al. 2008). The GSD can be approximated by GSD ∼ πR M /2q √ 3 with R M the radius of the moon and q a parameter that reflects the grid sampling that is usually taken as a power of 2. It corresponds to a latitudinal resolution of 180/2q √ 3 deg per plate. According to Gaskell's nomenclature, the shape model at a resolution of q = 64 has to be considered as a global topography model (GTM) of good spatial resolution. Shape models can be found up to q = 512. This topographic model is referred as the q-topo case (Fig. 1). The surface of Dione has been also approximated by a tessellated ellipsoid of variable resolution, parametrised by a latitudinal sampling of 180/n deg per plate, that is n ∼ 2q √ 3. The ellipsoid is generated and tessellated using MATLAB, with a A8, page 2 of 11 longitudinal sampling twice as large, since the number of samples per axis should be constant. It is then encapsulated in a DSK (Digital Shape Kernel) format to be exploited in the same way as the topographic model. The tessellated ellipsoid model is referred as the n-ell case. This surface model allows for studying the effect of the dispersion of local times of plates in the footprint on the thermal emission, independently of the effect of the topography, which itself introduces an additional dispersion due to the different slopes, orientations, shadows, and visibility conditions (taken into account by the q-topo model).
The viewing and illumination angles are calculated at the barycentre of each plate within the footprints of the observation. The relative distance from the plate to the center of the footprint is also stored to estimate its relative contribution to the total heat flux in the footprint, taking into account the FOV spatial sensitivity (Sect. 2.2.1).
Finally, we have the 'flat' case, where the surface is approximated by a single plate whose normal is that of the underlying ellipsoid at the FOV intercept point. This is the simplest surface model, where neither the effect of dispersion in local times in the footprint nor that of the topography are taken into account. The illumination and viewing geometries are calculated at this point at the time of acquisition and the model calculates the thermal history of this virtual surface element.

Thermal emission of the surface elements
The model for the thermal emission of planetary regolith inherits the formalism developed by Ferrari & Lucas (2016). It links the effective conductivity and the thermal inertia of regolith with its compositional and structural properties. Thermal inertia can be highly time and space dependent as the effective thermal conductivity may have a significant temperature and grain size dependent radiative component, even in icy cold regoliths (Ferrari & Lucas 2016). Heat transfer by conduction through contacts may be indeed very limited due to grain roughness and supplanted by raditive heat transfer through the pores. Here, this formalism is coupled with the heat diffusion equation so as to account for the thermal history of any surface element of the GTM.

Heat diffusion and boundary conditions
Temperature distribution in the regolith, T (z, t), is a function of depth, z, and time, t. It is determined by solving the heat diffusion equation: with ρ(z) as the effective volume density, k E as the effective thermal conductivity, and C P (T ) as the specific heat capacity. The value of ρ(z) depends on the regolith porosity, p, which may vary with depth but which we keep constant here ( Table 1). The effective thermal conductivity can be either constant (k E = k E0 ) or may depend on temperature, porosity, and grain size (k E = k E (T, z)). Radiative heat transfer may indeed occur in parallel with conduction heat transfer, so that k E (T, z) = k C + k R . Ferrari & Lucas (2016) studied the variations of the effective thermal conductivity with grain size, ice phase, the type of grain contact (tight or loose), and the model used for the radiative conductivity k R . The latter is proportional to the emissivity, the grain size, R, the temperature, T 3 , and the porosity, p (Table 1, Breitbach & Barthels 1980;Gundlach & Blum 2012). Ferrari & Lucas (2016) also showed that: (1) Thermal inertia values as low as 10-20 MKS at Saturn temperatures (80 ± 30 K) could be obtained with regolith of normal porosity (p = 0.3-0.8) but with loose grain contact or amorphous ice phase; (2) The thermal inertia of the regolith is minimum for grain sizes in the range of 100 µm-1 mm. For lower sizes, it increases as contacts between particles increase and radiative conductivity decreases. For larger grain size, the radiative conductivity dominates as it scales with the grain size; (3) The thermal inertia of amorphous ice at 80 K is minimum for grain size in the range 10-100 µm, while it can reach up to 100 MKS for grains of several centimetres in size. In such a model, thermal inertia is therefore expected to vary spatially depending on local time, temperature, latitude, and regolith properties (porosity, grain size or ice phase). Any plate radiates into space from its skin-layer (z = 0) with the boundary condition: where µ 0 stands for the cosine of the illumination angle relative to the local normal of the plate (q-topo or n-ell shape models) or to the normal of the ellipsoid at the center of the FOV intercept (flat surface model). In addition, T S is the surface temperature, A its bolometric Bond albedo,ε h its average hemispherical emissivity, S the solar constant, and D AU the heliocentric distance measured in AU at the epoch of acquisition. The bottom layer at depth may receive an internal up-welling thermal flux φ 0 , which is null here. The diffusion equation is solved with a Crank-Nicolson scheme. Vertical sampling is fixed to ten samples per skin depth δ p , which is estimated for the average temperature at the epoch of acquisition for a diurnal cycle. Time sampling is adapted to ensure convergence on temperatures within 0.1-0.2 K in a reasonable calculation time, that is, up to 0.3 sec per A8, page 3 of 11 A&A 655, A8 (2021) Average hemispherical emissivity (Eq. (11)) µ 0 , µ, µ 0, j , µ j Incidence and emission angle cosines at center of FOV, of plate j Baragiola (2003).
plate over a full diurnal cycle of Dione and for very low thermal inertia. This is below the error bar of ∼0.5 K on observed temperatures. In cases where the regolith does not conduct heat, either because the thermal conductivity is very low (k E ∼ 0) or the porosity is very large (p ∼ 1), the thermal inertia is near zero and the surface temperature is calculated by: At this point, we can check that the surface temperature, T S , which is a solution of Eq.

Bond albedo and emissivity
The surface of Saturn's icy moons is covered with water ice, most likely crystalline ice (e.g. Cruikshank et al. 2005). Both crystalline and amorphous phases may coexist at various depths depending on their thermal and space weathering history (Ferrari 2018). Their VNIR reflectance spectra typically exhibit a reddening below 1 µm and a quenching of strong water ice absorption above 3.6 µm. This indicates the presence of tholins as no other darkening components can explain such effects . Pitman et al. (2010) constrained the bolometric Bond albedo A of the leading and trailing hemispheres of icy moons using VIMS data: A ranges from 0.37 ± 0.08 to 0.93 ± 0.11 for Dione and Enceladus trailing hemispheres, respectively. A decrease in A may be due either to an increased grain size or to an increased fraction of water ice contaminants (Emery et al. 2005). The fraction of water ice contaminants varies from a few to up to tens of percent and the superficial grain size ranges from 1 µm to a few hundred µm, depending on the photometric model Scipioni et al. 2013).
Modelling the VNIR data is beyond the scope of this paper, however, the bolometric Bond albedo, A, and the average hemispherical emissivity,ε h , which depend on grain size, composition, and temperature, are primordial factors controlling the energy balance (Eq. (2)). Thus, we developed a coherent approach ensuring consistency between the solar reflectance and the infrared emittance. For that purpose, we introduce Hapke's formalism (Hapke 2012).
Let consider an intimate mixture of water ice and tholins, the optical constants of which are available over a large wavelength domain from the UV to the millimetre waves (Fig. 2). Hapke (2012) expressed the plane albedo, A h , and the hemispherical emissivity h of a regolith made of grains with size, R, as a function of their single scattering albedo, ω 0 , their cosine asymmetry factor, ξ, and the regolith porosity. Emissivity may be also directional, d (µ), depending on the emission angle (Table 1). For a quasi-continuous medium made of arbitrary single-scatterers, we can use Hapke's Isotropic Multiple Scattering Approximation (IMSA) model to write the hemispherical reflectance r h and the hemispherical emissivity h : and with κ P the porosity factor set to 1, γ λ = 1 − ω 0,λ and r 0 = (1 − γ λ )/(1 + γ λ ). The directional hemispherical reflectance r h (µ 0 , λ), also called plane albedo, is the total fraction of incident collimated radiation per unit area scattered into the upward hemisphere. In case of anisotropic scattering within the slabs, the single scattering albedo ω 0 can be replaced by ω * 0 : and This approximation works quite well with the hemispherical reflectance, however, the performance is reduced when considering the directional reflectance (Hapke 2012). The spherical albedo A S (λ) is the total fraction of incident irradiance scattered by a body in all directions. In the case of isotropic scatterers, it is written as: The bolometric Bond albedo A is the average of the spherical albedo weighted by the solar spectral irradiance B λ (T ): The value of B λ (T ) is estimated by the Standard Zero Air Mass Solar Spectral Irradiance ASTM E490 over the 0.2-1000 µm range. The corresponding solar constant S = 1366 W m −2 . The radiance emitted in the direction µ from an isothermal layer at temperature T is I(λ, µ) = d (µ, λ)B λ (T ), where the directional emissivity d (µ, λ) = γ λ H(µ) can be approximated by: The average hemispherical emissivity over the thermal infrared domain at Saturn distance for a surface element covered with regolith is: with B λ (T S ) as the blackbody radiation of the plate at surface temperature, T S . The single scattering albedo ω * 0 of a grain with size R is calculated by the Mie theory using the complex refractive index of crystalline and amorphous water ice. When appropriate, the geometrical optics approximation is used to speed up the calculations. Bertie & Whalley (1967) data were complemented by Hudgins optical constants at longer wavelengths for crystalline water ice around 100 K (Hudgins et al. 1993). In case of pure amorphous water ice, we assumed that k ∼ λ −0.5 in the far infrared region (Fig. 2). The optical constants of tholins are taken from Khare et al. (Khare et al. 1984). We verified that an increase in the tholins fraction reddens and darkens the VNIR spectrum as the imaginary part of the complex refractive index (the extinction coefficient) of tholins decreases and is much larger than that of water ice between 0.2 and 0.7 µm (Fig. 2).
As expected, the bolometric Bond albedo A decreases with the grain size, R, and the tholins fraction up to 20% (Fig. 3a). It is close to unity for pure ice micrometre grains as suspected on the Enceladus trailing hemisphere (Pitman et al. 2010). For a tholins fraction of a few percent, A ∼ 0.6, which is a common value on icy moons, for grain sizes between 10 µm and a few hundred µm. The minimum observed albedo A ∼ 0.4 is obtained for a fraction of 20%.
The model shows that the average hemispherical emissivity does not change whatever the tholin fraction (Fig. 3b). It highly depends on grain size:ε h = 0.5-0.6 for sub-micro-metric grains; it is lowest at 0.3 for micro-metric grains and grains ≤ 0.1 µm; it steeply increases for grain sizes between 10 µm and 40 µm; finally,ε h = 0.90-0.97 for larger grains. Given a grain size, it exhibits some dispersion when the surface temperature varies from 50 to 110 K (Fig. 3b). Similar curves are obtained with amorphous water ice. Whileε h depends on T S , we approximatē ε h (T S ) ∼ε h (80 K) in the boundary condition (Eq. (2)) to simplify the resolution of the transient thermal cycle. It is therefore assumed thatε h depends on grain size but not on temperature.
As a consequence, the regolith temperature is also expected to vary with grain size, since it is scaling with ((1 − A)/ε h ) 1/4 : grains ≤ 0.2 µm are poor emitters and absorbers; grains ∼1 µm are poor emitters and good absorbers; and larger ones are good absorbers and emitters. Figure 3c displays the characteristic temperature T slow = T 0 /2 1/4 of a slowly rotating moon with low thermal inertia and covered with icy grains polluted by a variable fraction of tholins (0%, 5%, 20%). For pure ice, T slow increases with R ≥ 0.3 µm. As the albedo decreases with increased tholins fraction, the temperature of µm-to-10-µm-sized grains becomes prominent. The temperature of grains ≥60 µm is almost constant for tholins fractions greater than a few percent.
This indicates that the bolometric Bond albedo, A, is mainly controlled by the grain size and pollutant fraction in ice, while the average hemispherical emissivity,ε h , in the infrared domain is mainly controlled by the grain size, the icy composition and marginally the temperature. We may keep the tholins fraction as a model parameter that fixes A on the basis of the ice phase and the grain size. However, in the following, A will be set up a priori, taken from Pitman et al. (2010). Observations in the visible domain are indeed the most appropriate to constrain the bolometric Bond albedo. It is also useful to limit the number of parameters to perform the inversion of the model within acceptable time frames.  (11)). The error bars figure the dispersion ofε h when T S varies from 50 K to 110 K for a tholins fraction of 5% in crystalline water ice grains. The curves are similar in cases of amorphous water ice (Eq. (11)). (c) Typical surface temperature T slow of a plate covered with regolith of grain size, R, for various tholins fractions (0%, 5%, 20%) in the case of a slowly rotating moon.

Diurnal thermal cycles
of Dione was about 10 MKS. A constant thermal conductivity k E0 = 10 −3 W m −1 K −1 yields similar thermal inertia with the values of C P (T ), ρ 0 , and p provided in Table 1. For a constant effective thermal conductivity k E0 (Fig. 4a), the regolith average temperature over the diurnal cycle increases conversely to the average hemispherical emissivity of grains: micro-metric grains (0.2-2 µm) with lower emissivity are hotter (Fig. 3b). Regolith covered with grains larger than a few hundred µm exhibit lower average temperature because of their higher emissivity. The thermal inertia decreases with higher porosity, so that the amplitude of the diurnal thermal gradient expectedly increases.
When the thermal conductivity is both radiative and conductive (Fig. 4b), the thermal cycle depends more on grain size, temperature, and bulk conductivity of the ice phase. It also varies with the type of contacts between grains. Ferrari & Lucas (2016) showed that thermal inertia decreases conversely to grain size in the range of µm-to-mm for a regolith porosity p ≥ 0.3. It increases again for cm-sized grains. In case of amorphous ice with loose contacts, the thermal inertia for grains size ≤ 1 µm reaches 60 MKS as the contact conductivity scales with 1/R and decreases below 22 MKS for larger ones. The thermal gradient significantly increases. For cm-sized grains the thermal inertia increases again and the gradient is reduced. The thermal gradient and average temperature are all the more dependent on the grain size and bulk conductivity of the ice phase. Because their bulk conductivity indeed differs by almost a factor of ten, the thermal inertia of amorphous ice is lower and the diurnal thermal gradient is larger (Fig. 4b).

Temperature and apparent emissivity in CIRS FP1 FOV
The temperature estimated from the measured spectrum in the FP1 footprint depends on the thermal history of all the surface elements the footprint contains. Considering N visible elements, the modelled spectrum I F (λ) can be calculated as the weighted sum of the elementary spectra as with µ j the cosine of the local emission angle of plate j, S j its surface, and I j (λ, µ j ) its modelled emission spectrum. This latter can be either for an isotropic hemispherical emissivityε h, j or for a directional emissivity ε d, j . The emission spectrum I j (λ, µ j ) of each surface element is weighted by the spatial sensitivity w j of the FP1 detector that is a function of its distance to the center of the FOV (Sect. 2.1). The sensitivity to the distance is assumed to be a 2D Gaussian with a FWHM of 2.42 mrad and an aperture diameter of 3.9 mrad (Flasar et al. 2004). The apparent temperature, T F , and emissivity, ε F , in the FOV are derived from the best fit of I F (λ) = ε F B λ (T F ) over the wavenumber (w n ) range from 20 to 350 cm −1 (wavelength range from 29 to 500 µm) around the observed emission peak. This is the same range wherein the CIRS FP1 emission spectra I obs are fitted in order to optimise the signal-to-noise ratio and to minimise the errors on T obs and ε obs . The noise equivalent spectral radiance (NESR) may be prominent above this wavenumber range depending on the spectral resolution used and the measured temperature and emissivity (Flasar et al. 2004). Figure 5 shows simulated and observed emission spectra within two footprints selected among the 350 obtained during the observation (Sect. 3, Fig. 6). The temperature and emissivity are lower in southern footprint #5 than in footprint #3 (Sect. 3, Fig. 7) and the data are noisy above 400 cm −1 . For the chosen model parameters -which do not result from any inversion procedure -it can be seen that the inclusion of spatial sensitivity in the FP1 FOV has variable impact on the modelled emission spectrum I F and is most detectable at southern latitudes (footprint #5). Increasing the thermal conductivity K E0 from 0 to 0.001 W m −1 K −1 reduces the apparent temperature (the peak shifts to lower wavenumbers) in these footprints which have been acquired close to 12:30 local time.

Results
The impact of topography and regolith properties on the directional variations of T F and ε F is discussed for the observation CIRS_222DI_DIONE001_PIE. Its design follows a south-north latitudinal scan, from about −50 • S to +16 • N, at an average west longitude of 123 • W (Fig. 6). The spatial extent of FP1 footprints is about 270 km at the equator. N = 3413 surface elements of the GTM are flown over during that sequence, with an average of 600 elements per footprint. Local time is almost constant during the scan, between 12:30 and 12:50. For the GTM (q = 64), the GSD is about 8 km (given R M = 561 km), which fits the latitudinal sampling of the n-ell ellipsoidal shape model with n = 222. The incidence angle at the FOV intercept with the ellipsoidal surface of the moon varies from 80 • to 25 • along the scan while the emission angle varies from 55 • to a few degrees. The observed temperature T obs decreases from 95 K to 78 K with the incidence angle (Figs. 7a and 8b). The apparent emissivity ε obs highly varies along this scan, from 0.5 at southern latitudes (high emission angles) to 0.92 close to the normal viewing angle (Fig. 7b). At ε obs ≤ 0.75, the measured temperature T obs is constant around 80 K. Releasing the constraint of unit emissivity (ε obs = 1) highlights the directional anisotropy of the apparent emissivity as a function of the emission angle Em.

Effect of topography
Simulations of apparent temperature and emissivity within the footprint are performed assuming a regolith with a bolometric Bond albedo A = 0.63 and a grain size, R = 2 mm (Fig. 7). Their average hemispherical emissivityε h ∼ 0.92 is indeed comparable to the observed value ε obs near nadir (Figs. 3b and 7b). The thermal conductivity is either null (by default) or constant with k E0 = 0.001 W m −1 K −1 (when specified), so that the thermal inertia is either Γ = 0 or Γ = 15 MKS for a porosity p = 0.5.
Considering a 'flat' surface model, the apparent emissivity is found to be constant along the scan: ε F ∼ε h . The apparent temperature is consistent with T obs , although slightly too cold for southern latitudes (footprints #4 and #5). The value of ε F varies from 0.93 to 0.86 (∼7.5%) with Em when the hemispherical directional emissivity e d (µ) -Eq. (14) -is included in the model, which is much lower compared to the observed variation.
In the case of a tessellated ellipsoid (n-ell), ε F varies from 0.92 to 0.71 (∼23%) with Em. At the spatial resolution of ∼8 km (n = 222), this model indeed takes into account the dispersion in local time due to the dispersion in longitude of the different surface elements in the field of view. If it is not yet sufficient to fit ε obs (Em), the directional variation introduced here (with zero conductivity and no spatial sensitivity of the FP1 detector) is already significant. Increasing the spatial resolution n does not change the results.
Taking into account the local topography at the equivalent 8 km scale (q-topo with q = 64) slightly increases the slope down to ε F = 0.70 at footprint #5, with null thermal inertia and no spatial sensitivity taken into account. Improving the spatial resolution of the GTM down to the kilometre scale (q = 512) does not significantly change T F (ε F ) and ε F (Em) in this observation. Even at southern latitudes, the variation of emissivity as a function of q ≥ 64 does not exceed 1%.
On the other hand, the introduction of the spatial sensitivity of FP1 reduces the directional emissivity for Em ≥ 35 deg (footprints #4 and #5) as it appears to reduce the radiometric weight (w j ) of hotter plates, reducing the total specific intensity, I F , at larger wavenumbers (Fig. 5). The apparent temperature decreases only slightly, while the apparent emissivity increases more significantly. Setting the thermal inertia to 15 MKS naturally decreases the temperature by about ∼15 K, as expected for a scan performed around noon local time, while ε F (Em) reaches Observed I obs , fitted ε obs I(T obs ) and modelled spectra I F for footprints sampled during the FP1 observation CIRS_222DI_DIONE001_PIE (see Figs. 6,7). The modelled spectra are calculated for amorphous ice with grain size, R = 2 mm, and regolith properties A = 0.63, p = 0.5, and K E0 = 0 (Γ=0) or K E0 = 0.001 W m −1 K −1 (Γ ∼ 15). The effect of the spatial sensitivity of FP1 is included (+SS) or not. The Dione GTM at spatial resolution q = 64 is used. Top: footprint #3 located at latitude and west longitude (−30.89 • , 123.84 • ). Bottom: footprint #5 located at latitude and west longitude (-50.40 • , 130.02 • ).
its highest slope with a 27% variation between normal viewing and Em ∼ 55 deg. Increasing the thermal inertia may help better fit ε obs (Em) but it would even more decrease the temperatures at noon and widen the shift between T F and T obs (Fig. 7a).
Accounting for the dispersion in local times in the footprints with a tessellated ellipsoid (n-ell shape model) definitely impacts the directionality of the apparent emissivity. Including the 8 kmscale topography has a second-order influence on this. Improving this q-resolution by a factor of 8 does not have much impact on the results of the model at Em ∼ 55 deg. Taking into account the spatial sensitivity of the FP1 detector is important. Finally, we notice that with the chosen parameters corresponding to reasonable a priori in terms of Bond albedo and thermal inertia, we can only approximate the data obtained while assuming uniform regolith properties. However, we can globally explain the behaviour of T F (ε F ) and ε F (Em), while introducing a tessellated surface and releasing the condition ε obs = 1.

Effect of thermal conductivity and grain size
The grain size controls the energy balance at the surface and the thermal cycle by affecting the average hemispherical emissivity. The less emissive micro-metric grains are the hottest, given the Bond albedo A = 0.63 (Fig. 3b, c). For a constant thermal conductivity k E0 and a given shape model, when varying R from 0.2 µm to 2 cm (Fig. 8-Top a, b), the modelled temperature T F is too cold, except for micro-metric grains. A regolith made of grains of a few hundred microns to a few millimetres has an apparent emissivity close to that observed. For micro-metric or sub-micro-metric grains, T F agrees with T obs but ε F remains much lower than ε obs (Fig. 8-Top c).
If the thermal conductivity is both radiative and conductive ( Fig. 8-Bottom), the grains must be composed of amorphous water ice for the model to yield thermal inertia ranging between 2 and 23 MKS, in cases of tight contacts and p = 0.7 (Ferrari & Lucas 2016). Introducing size and temperature dependent thermal conductivity directly impacts T F (Fig. 8-Bottom a, b) but only slightly changes the apparent emissivity ε F . This latter is mainly driven by the dispersion in local hour angle (Fig. 8-Bottom c versus Fig. 8-Top c). Such conductivity introduces a large discrepancy in the thermal history among plates covered with grains larger than a few hundred of µm: the thermal conductivity is indeed dominated by radiative conductivity which scales with R: the larger R the larger k E ∼ k R . Therefore, cm-sized grains have larger thermal inertia and their temperature decreases near noon local time (Fig. 8-Bottom a, b). The thermal conductivity is the lowest for grain size in the range 1-100 µm, so that their temperature is the highest at local noon. It  Fig. 6). Inside the oval, the grain size is R = 2 mm and outside of it R = 20 µm. The Bond albedo is uniform, either A = 0.49 or A = 0.63. The thermal conductivity of amorphous ice is constant, either k E0 = 0 or k E0 = 0.001 W m −1 K −1 with porosity p = 0.8 (Γ = 10 MKS) for both regions. The topography of GTM (q = 64) and the spatial sensitivity of FP1 are included. The emissivity of surface elements is assumed to be isotropic (Eq. (13)).
increases again for very small grains as conduction through loose contacts scales with 1/R and dominates the radiative conductivity. In this case, the modelled temperature is too cold compared to the measurements. Finally, the question of regional variations in regolith properties can be also addressed with this model, for example a regional variation in grain size, which may originate from sintering after bombardment by MeV electrons in the lens-shaped region of icy moons leading hemispheres (Fig. 6). We considered here a regolith composed of 20 µm-sized icy grains covering the leading hemisphere outside of the central 70 • -wide lens-shaped region, which is itself covered by 2 mm-sized grains made of amorphous ice. The choice of the size is dictated by the temperature and the emissivity observed at the two extremes of the scan, near the equator and at southern latitudes. At southern latitudes outside the oval shape, T obs ∼ 80 K and ε obs ∼ 0.55, whereas at the equator T obs ∼ 95 K and ε obs ∼ 0.92 ( Fig. 8-Top, Fig. 3). We performed four simulations assuming either A = 0.63 (Pitman et al. 2010;Howett et al. 2010) or A = 0.49 (Howett et al. 2014) with either null or constant thermal inertia (k E0 = 0.001 W m −1 K −1 ) and porosity p = 0.8 (Fig. 9). When taking into account the regional variability and the local topography at 8 km-scale (GTM with q = 64) together with the FP1 spatial sensitivity, the modelled ε F (Em) is close to ε obs (Em) in all cases (Fig. 9b). This is only slightly sensitive to the Bond albedo and the low thermal inertia. The apparent temperature T F is best reproduced by a regolith with either albedo A = 0.49 and Γ = 10 MKS or A = 0.63 and Γ = 0 MKS. Considering regional variations in regolith properties around the lens-shaped region of the leading hemisphere, on the Dione moon particularly, seems to help in reproducing the observed apparent temperatures and emissivity, whereas a single population can hardly accomplish this (Fig. 8). Of course, an estimation of the physical properties of these regions with the model should be carried out for each of Saturn's moons by analysing as many relevant observations as possible.

Conclusion
In this paper, we present a thermophysical model with the aim of analysing the temperatures and the apparent emissivities of Saturn's icy moons observed with CIRS FP1. We include their topography and a parametrisation via the regolith physical properties such as grain size, water ice phase, porosity, and the type of grain contact. We conclude that: 1. The dispersion in local times of plates within the footprints has a definite impact on the directional anisotropy of the apparent emissivity. Including the 8 km-scale topography has a second-order influence on this. Improving this resolution down to the 1 km-scale has small impact on both apparent temperature and emissivity at Em ∼ 55 deg. Taking into account the spatial sensitivity of the FP1 detector is important. The directional anisotropy of the apparent emissivity cannot be explained by the directional emissivity ε d (R, µ) at the scale of the plate as proposed by Hapke (2012). 2. The model can account for regional variations in regolith properties and test current hypotheses on its evolution by space weathering, such as a change in grain size between a lens-shaped region and southern and northern latitudes, which may be covered with smaller grains. 3. Considering the hemispherical emissivity varying with grain size and the topography to analyse CIRS FP1 observations of icy moons with this model will yield constraints on the physical properties of their surface and may lead to some changes in the current estimates of their thermal inertia and their regional variations.
This model is adapted to address the question of the actual properties of atmosphere-less planetary surfaces evolving under diverse weathering processes by analysing their infrared emission. The analysis of CIRS FP1 observations of Mimas, Dione, and Rhea is underway in order to test current hypotheses on regolith processing by space weathering, such as regolith changes in grain size.