Broadband spectroscopy of astrophysical ice analogs. I. Direct measurement of complex refractive index of CO ice using terahertz time-domain spectroscopy

Context: Reliable, directly measured optical properties of astrophysical ice analogs in the infrared (IR) and terahertz (THz) range are missing. These parameters are of great importance to model the dust continuum radiative transfer in dense and cold regions, here thick ice mantles are present, and are necessary for the interpretation of future observations planned in the far-IR region. Aims: Coherent THz radiation allows direct measurement of the complex dielectric function (refractive index) of astrophysically relevant ice species in the THz range. Methods: The time-domain waveforms and the frequency-domain spectra of reference samples of CO ice, deposited at a temperature of 28.5 K and annealed to 33 K at different thicknesses, have been recorded. A new algorithm is developed to reconstruct the real and imaginary parts of the refractive index from the time-domain THz data. Results: The complex refractive index in the wavelength range of 1 mm - 150 ${\mu}$m (0.3 - 2.0 THz) has been determined for the studied ice samples, and compared with available data found in the literature. Conclusions: The developed algorithm of reconstructing the real and imaginary parts of the refractive index from the time-domain THz data enables, for the first time, the determination of optical properties of astrophysical ice analogs without using the Kramers-Kronig relations. The obtained data provide a benchmark to interpret the observational data from current ground based facilities as well as future space telescope missions, and have been used to estimate the opacities of the dust grains in presence of CO ice mantles.


Introduction
One of the main problems in unraveling the chemical and physical properties of molecular clouds, in which the star and planet formation process takes place, is to estimate correctly the amount of gas contained.The difficulties in the direct observation of molecular hydrogen constrain the possibility to calculate the total mass of a cloud.An easy alternative could be to use carbon monoxide as a tracer of molecular gas, but in dense and cold regions of the interstellar medium and protoplanetary disks, CO is not a good tracer of gas mass because CO molecules preferentially resides onto dust grains, forming thick icy mantles (e.g Dutrey et al. 1998;Caselli et al. 1999).
In alternative, the dust continuum emission is the best available tool to compute a molecular cloud mass, if dust opacities are known.
The advent of ALMA and NOEMA facilities offer the possibility to observe the dust continuum emission in the millimeter and sub-millimeter part of the electromagnetic spectrum, with very high angular resolution and sensitivity.However, to properly model the dust continuum emission it is necessary to have Send offprint requests to: Barbara M. Giuliano, Alexei V. Ivlev information about its grain size distribution as well as its chemical composition, since the dust opacity depends directly on these paramenters.If we take into account that the dust grains can be covered by ice mantles, at the center of prestellar cores or in protoplanetary disk midplanes, we also need to investigate how the presence of ices is changing the dust opacities.
Unfortunately, no experimental data are available for these cases, and the interpretation of the dust continuum emission measurements relais on calculated opacity values, as the ones tablulated in Ossenkopf & Henning (1994).The goal of our study is to provide laboratory data on the optical properties of CO ice and utilize these data to calculate the opacities of dust grains covered by CO ice mantles.We will compare the opacity values obtained by our study to those available in the literature.
FIR studies on spectral properties of molecular solids, without deriving optical constants, started early, to deepen the understanding of the infrared-active lattice vibrations of simple species.Anderson & Leroi (1966) provided studies on frequencies of CO and N 2 in the range 40-100 cm −1 , and Ron & Schnepp (1967) complemented the available information with CO, N 2 and CO 2 intensity studies in the same frequency range.Moore and Hudson published in 1994 a comprehensive study of FIR spectra of cosmic type ices, including ice mixtures.These data include the analysis of amorphous and crystalline phases of the pure molecular ices, and the authors discuss the implications of the results on the identification based on astronomical observation.An estimation of the band strengths in the FIR region for pure ices and ice mixtures relevant for astrophysical environments can be found in Giuliano et al. (2014;2016).
Recently, the investigation of the THz spectroscopic properties of ice mantles analogs has gained a considerable interest (Allodi et al. 2014;Ioppolo et al. 2014;McGuire et al. 2016).This technique allows to measure directly the intermolecular vibrations in the ice samples, related to the lattice structure, which can be connected to their large-scale structural changes and finally to their thermal history.On the contrary, the spectroscopic features measured in the MIR frequency range are indicative of the intramolecular vibrations of the sample, which can provide a wealth of information on the molecular identification and chemical reactivity.A comparison of our THz experimental data with that observed in the MIR range, could help us to reveal intra and inter molecular vibrations.However, this study is beyond the scope of this paper and will be addressed in our future investigations.
Nowadays, numerous spectroscopic methods are extensively used for dielectric measurements at THz frequencies (Lee 2009); particularly, we would mention the following ones: Fourier Transform Infrared (FTIR) spectroscopy (Griffiths & de Haseth 1986), Backward-Wave Oscillator (BWO) spectroscopy (Komandin et al. 2013), spectroscopy based on photomixing (Preu et al. 2011) or parametric conversion (Kawase et al. 1996;Kiessling et al. 2013), and, finally, THz time-domain spectroscopy (THz-TDS) (Auston 1975;Van Exter et al. 1989).These methods exploit either continuous-wave or broadband sources, operate in different spectral ranges, and are characterized with different sensitivity and performance.Among them, the THz-TDS seems to be the most appropriate for studying laboratory analogs of circumstellar and interstellar ices.In contrast to other approaches, the THz-TDS yields detection of both amplitude and phase of sub-picosecond THz pulses in a wide spectral range as a result of a single measurement; thus, reconstruction of the dielectric response of a sample might be performed without using the Kramers-Kronig relations (Martin 1967) and involving additional physical assumptions.Furthermore, the THz-TDS yields analysis of separate wavelets forming the time-domain response of a sample; thus, it is a powerful method for the characterization of multilayer samples.Thereby, we have selected the THz-TDS as a spectroscopic technique for our experiments.
We aim at the extension of the laboratory data in the far-IR/THz region, and we will show how the employment of the THz-TDS is able to provide the direct measurement of the real and imaginary part of the refractive index of the ice sample.The experimental and theoretical methods employed will be explained in Section 2, the results obtained and how these data are relevant for astrophysical application will be presented in sections Section 3.3 and Section 4, respectively, the conclusions will be illustrated in Section 5.

Experimental and theoretical methods
For this series of experiments a dedicated setup has been designed and developed in the laboratories of the Center for Astrochemical Studies (CAS) located at the Max Planck Institute for Extraterrestrial Physics in Garching (Germany).The setup is composed of a closed-cycle He cryocooler coupled to a THz time-domain spectrometer.The cryocooler vacuum chamber is small enough to be hosted in the sample compartment of the THz spectrometer, and it is mounted on a motor controlled translational stage which ensures the tuning of the cryostat position with respect to the THz beam.The details of the main components of the apparatus and the ice growing procedure are given in the following subsections.

Cryogenic setup
The high power cryocooler has been purchased from Advanced Research Systems.The model chosen is designed to handle high heat loads ensuring a fast cooling.It is equipped with a special interface capable of reducing the vibration transmitted from the cold head to the sample holder at the nm level.This requirement is important in case of spectroscopic measurements in the THz frequency region, where the induced vibration of the sample can cause the increase of the noise level of the recorded spectra.
The cryostat is placed in a 15 cm diameter vacuum chamber, equipped with four ports for optical access as well as for the gas inlet.The optical arrangement is designed to work in the transmission configuration.The optical windows chosen for the measurements at the desired frequency range are made of high-resistivity float-zone silicon (HRFZ-Si), purchased by Tydex.This material features high refractive index of n Si = 3.415, negligible dispersion and impressive transparency in the desired frequency range.The same material has been chosen as a substrate for the ice growing.In order to suppress the Fabry-Perot resonances in the THz spectra, caused by multiple reflection of the THz pulse within the windows, the thickness of the three Si windows must be different from each other.For this purpose we have chosen to use the 1 mm thick Si window as substrate for the ice deposition, and use two windows of 2 mm and 3 mm thickness as optical access to the THz beam.
We performed measurements with a slightly focused THz beam configuration in order to mitigate vignetting and diffraction of the beam at the metal components of the vacuum chamber.A schematic overview of the chamber arrangement is sketched in Fig. 1.The pumping station is composed of a turbomolecular pump (84 ls −1 nitrogen pumping speed) combined with a backing rotary pump (5 m 3 h −1 pumping speed) providing a base pressure of about 10 −7 mbar.The minimum temperature measured at the sample holder in normal operation mode is 5 K.

THz time-domain spectrometer
The THz-TDS setup used for this work has been purchased from the company Batop GmbH.The model chosen is the TDS-1008, with a customized sample compartment in order to allocate the cryostat.It is based on two photoconductive antennas made of low-temperature-grown gallium arsenide (LT-GaAs), which constitute the emitter and the detector of the THz pulse (Lee 2009).The antennas are triggered by a femtosecond laser (TOPTICA, 95 fs, 780 nm) with a pulse repetition rate of 100 MHz and an average input power of 65 mW.Further details on the setup are provided in the Appendix A. Fig. 2 (a) shows the optical path of the laser beam into the optical bench of the spectrometer.Panel (b) in Fig. 2 shows standard broadband Fourier spectra.The THz pulse registered with the beam path empty is then converted in the blue spectrum spanning the frequency range of 0.05 up to 3.5 THz with maximal spectral amplitude centered at about 1.0 THz.In TDS, the timedomain THz waveform is converted to the frequency domain using the direct Fourier transform, which yields the frequencydependent amplitude and phase of the THz wavelet.Since the frequency-domain data is calculated via the direct Fourier transform, the spectral resolution of measurements is determined as ∆ν = 1/∆T , where ∆T is a size of the of the time-domain apodization filter, chosen to avoid the edge effects (i.e. the Gibbs effect) in the frequency-domain.In our experiments we used the 35-ps Tukey apodization filter (see Appendix B), which yields the spectral resolution of about 0.03 THz.

Vacuum
The green spectrum has been recorded with the cryostat placed in the sample compartment.It is converted from a waveform which contains both a first THz pulse (i.e. a ballistic one) and a train of satellite pulses originating from to the multiple THz wave reflections within the windows.The ballistic THz pulse is delayed in the input/output windows and the substrate of the vacuum chamber.The spectrum of this waveform is slightly suppressed due to the Fresnel losses and modulated due to the interference of the ballistic pulse and the satellites.
In Fig. 2 (b), we show the shaded area at lower frequencies, where we expect growing distortions of the experimental data caused by the THz beam diffraction on the 20-mm-diameter aperture of the substrate.Assuming that the THz beam spot formed at the substrate is diffraction-limited, the lateral intensity distribution in the spot is defined by the Bessel function of the first kind (Born & Wolf 1980).The resulting width of the first intensity peak is approximately (3.8/π)( f /D)(c/ν), where D = 25 mm and f = 67 mm stand for the diameter and the back focal distance of the focusing lens, respectively.From this model, we deduce the critical frequency of 0.3 THz, below which less than 95% of the beam energy passes through the substrate aperture.Thereby, considering both the spectral sensitivity of our THz-TDS setup and the diffraction limits, the spectral operation range of our experimental setup is approximately limited within 0.3 to 2.0 THz.The THz-TDS housing is kept under purging with cold nitrogen gas during the entire experiment, to mitigate the absorption features due to the presence of atmospheric water in the THz beam path.The residual humidity measured at the sensor was less then 10 −3 %.

Ice preparation
The ices are prepared using standard technique in which the molecular sample in its gaseous form is allowed to enter the vacuum chamber through a stainless steel 6 mm pipe.The gas flux is controlled by a metering valve.Once the gas is expanding inside the vacuum chamber it condensates on the substrate.
For this set of experiments an ice thickness of the order of mm is required, to fulfill the sensitivity characteristic of the THz-TDS setup.This value is orders of magnitude higher than the usual ice thickness reached using this deposition technique, which is of the order of µm.To deposit such a thick ice in a reasonable amount of time, we have chosen fast deposition conditions, in which a considerable amount of gas is introduced in the vacuum chamber.In these conditions it is very difficult to obtain an ice sample which is homogeneous enough to determine its optical properties.To overcome this problem the gas inlet characteristics must be set accordingly.
We decided to remove any directionality from the gas inlet, keeping the pipe end cut at the vacuum chamber wall, at a distance of approximately 7 cm from the substrate (see Fig. 1).This configuration creates an ice layer on each side of the substrate.During the deposition, the pressure measured inside the chamber is approximately 10 −2 mbar.The ice deposition was divided in steps of 4, 5, and 6 min duration, up to a total of 30 min deposition time, in three different experiments.The final temperature was up to 28.5, 31.2 and 33.1 K for each step, at 4, 5, and 6 min deposition time, respectively.This increase is due to the condensation of the gas onto the cold surfaces of the cryostat, which is producing a heating rate too fast to be dynamically removed from the cooling system during the deposition.After each step, the THz spectrum has been recorded.This procedure was performed in order to rule out possible effects on the ice structure due to different deposition temperatures.As interstellar ices can be commonly found at temperatures as low as 10 K, the temperature of the cold substrate has been kept at 14 K, which is the lowest temperature achievable in the setup in this configuration, due to the fact that the radiation shield of the sample holder must be removed to ensure that no directionality of the molecular beam is present.Before moving to the next deposition step, the system was allowed to thermalize and spectra recorded after each deposition step have been taken at a temperature of 14 K.
As stated in Urso et al. (2016) the analysis of the Raman and IR spectra of experiments performed at increasing tempertaures from 17 to 32 K show no profile variation in the band at 2140 cm −1 which could be ascribed to a structural change in the ice morphology.
The waveform recorded in the time domain is compared, as reference for the measurements, with the waveform recorded for the substrate without ice, kept at a reference temperature of 14 K as well.After the deposition was completed, we have measured the spectra in different regions of the sample, to ensure that the ice morphology is spatially homogeneous.The results obtained from the spectra measured on a grid of 11 points spaced by 2 mm are in agreement within 10%, indicating a uniform ice formation over the substrate.

Derivation of the optical constants
In order to determine the optical constants, the ice thickness must be known.The laser interference technique is a well-established method to estimate the thickness of an ice sample deposited on a substrate as a function of the time.The absolute accuracy of this method is approximately within 5%, but the maximum CO ice thickness that we can measure with this technique is limited to 5 µm before the reflected laser signal becomes too weak to be detected, due to scattering losses occurring both in the bulk and on surface of the film.Thus, this technique is well suited for studying thin layers, when the ice thickness in total is below 10 µm (5 µm on each side of the substrate), but is not appropriate for experiments on thick ice samples.In turn, for the thickness estimation of the mm-size sample of ice featuring rather low THz absorption, we developed a model to perform an initial estimation of the optical properties and of the ice thickness directly from the recorded THz spectra, as described in the following subsections.

Ice parameters modeling
Our model aims to reconstruct the optical properties of ices, which are defined as follows where n and n are the real and imaginary parts of the complex refractive index n, c is the light speed, and α is the amplitude absorption coefficient, which is defined as half of the value of the power absorption coefficient.Equivalently, we can write where ε and ε are the real and imaginary parts of the complex dielectric permittivity ε.
The model is describing the THz wave propagation through the substrate with the ice deposited on both surfaces.The reconstruction of the ice parameters proceeds following three main steps.
The first task is modeling the reference and sample waveforms.Because of the focused arrangement of the THz beams, the electromagnetic wave is assumed to be planar and to interact with the sample interfaces at the normal angle of incidence.This is a common and conventional assumption widely applied in dielectric spectroscopy (Pupeza et al. 2007;Zaytsev et al. 2014).It allows us to describe all the features of the THz pulse interaction with the multilayer sample using the Fresnel formulas, which define the THz wave amplitude reflection at (and transmission through) the interface between the media, and the Bouguer-Lambert-Beer law, which defines the absorption and the phase delay of the THz wave in a bulk medium.Further details on these assumptions are given in Appendix B.
Figure 3 represents the THz wave propagation through the three layers structure -the first ice film, the HRFZ-Si substrate, and the second ice film, where the symbols 0 to 3 and N correspond to different components of the plane wave passing through the multilayered structure.As shown in Fig. 3, for the sample spectrum, we took into account the contribution of the ballistic THz pulse (1) and the satellite pulses (2 and 3), caused by the multiple THz wave reflection in the ice films.The mathematical description of the wave propagation can be found in Appendix B.
The second step consists in estimating the initial thickness l CO,I , l CO,II and the initial complex refractive index n Init of the two ice layers as shown in Fig. 4. The thickness estimation can be derived from the time delay δt 01 , between the ballistic pulse of the reference and sample waveforms (0 and 1 in Fig. 4), the first satellite pulse (2) and the ballistic pulse (1) of the sample waveform δt 12 , the second satellite pulse (3) and the ballistic pulse (1) of the sample waveform δt 13 , in Fig. 4 (a).Since the HRFZ-Si has a very high refractive index, we consider that the refractive indexes of ice n satisfies the inequality n 0 = 1.0 < n < n Si = 3.415, where n 0 is the refractive index in vacuum.
Then, by neglecting imaginary parts in the complex refractive indexes of media, the first assumptions for the real part of the complex refractive index of ice films and both their thicknesses l CO,I and l CO,II are described as follows: Lines 0 to 3 illustrate the ballistic pulse and satellite pulses transmitted in the direction of the antenna-detector, N stands for unaccounted satellites with larger time delays.Solid lines represent the pulses which have been used for the analysis, while dotted lines correspond to neglected pulses.
Eq. ( 3) is obtained from a mathematical model of sample and reference waveform.Futher details on the derivation can be found in Appendix B. From Eq. ( 4) it is possible to obtain information only on the thicknesses of ice films, while the identification of the specific layer (I or II) is not allowed.We observe a linear increase of the ice thickness with the total deposition time.A first assumption for the real part of the complex refractive index of CO ice is in the range of n init = 1.230 to 1.255 ± 0.035 for all the considered deposition intervals; here, the error accounts for an accuracy of the THz pulse peak position estimation.The first assumption for the imaginary part of the complex refractive index of ice has been done considering α init = 0; thus, n init = 0.
Finally, from the first estimation of the ice thickness and complex refractive index, it is possible to reconstruct the THz dielectric response of ice.The reconstruction procedure is reported in Appendix B, while the results are summarized in Fig. 4. Panels (b) and (c) show the growth of the two CO ice layers in time (t), considering different deposition steps of ∆t = 4, 5, and 6 min.
In order to validate the calculation of the ice thickness using this model, we can compare the results obtained from the THz spectral data to the thickness calculation performed with the well established laser interference techniques described in Section 3.2.The good agreement between the two methodologies confirm the validity of the present analysis.

Laser interference technique
In the adopted experimental configuration, a He-Ne laser beam (λ = 632.8nm) is directed toward the sample and reflected at near normal incidence both by the vacuum-sample and samplesubstrate interfaces.The reflected beam is detected by an external diode detector.It is possible to follow the accretion of the ice film by looking at the interference curve (intensity vs. time) of the reflected laser beam.Further details on the laser interference technique can be found in Urso et al. (2016). 1  The results obtained with the two techniques have been compared.The data on the accretion of the ice vs. time obtained from the analysis of the interference curve, measured at the early stage of the deposition process, are in good agreement with the data obtained from the THz spectra.In addition, by using the laser interference technique, we obtain for the CO ice a refractive index n CO = 1.27 that is close to the value obtained by using the THz technique.The good agreement between the two methodologies confirm the validity of the present analysis.

Reconstruction of the THz response
Finally, from the first estimation of the ice thickness and complex refractive index, it is possible to reconstruct the THz dielectric response of ice.The reconstruction procedure is reported in Appendix B, while the results are summarized in Fig. 4. Panels (b) and (c) show the growth of the two CO ice layers versus the deposition time, considering different deposition steps of ∆t dep = 4, 5, and 6 min.
The recorded THz waveforms and their Fourier spectra are presented in Fig. 5 for the optical substrate (used as a reference) and CO ice samples, at increasing thicknesses.In panel (a) the waveform E (t) is shown for the reference (green) and five subsequent deposition steps (black to light red) of approximately 0.45 mm total thickness for each step.The thickness reached after the total deposition time is approximately 2.3 mm, split in two ice layers of ≈ 0.85 mm and ≈ 1.45 mm on top of each side of the substrate.A small source of inhomogeneity can be ascribed to the position of the pipe connected to the pumping system, which is located in a lateral position of the vacuum chamber with respect to the cold substrate.
In the THz spectrum of CO ice, we observe a Lorenz-like resonant peak centered near 1.5 THz (50 cm −1 ) and a second blurred feature close to 2.5 THz (83 cm −1 ), masked by the sharp bands produced by the atmospheric water contamination in the spectrometer sample compartment.The highest frequency accessible in our setup is presently limited by the strong absorption features of residual water and carbon dioxide in the spectrometer case.We plan to change the current setup to an evacuated case in order to get rid of the contamination from the residual atmosphere, and we expect to extend the accessible frequency range up to 4 THz. 1 A publicly-available interface can be found at http://www.oact.inaf.it/spess/ to calculate the refractive index of the ice sample and derive the theoretical interference curve from the amplitude of the experimental curve.The sample thickness is obtained by comparing the two curves and using a procedure described in a document available at this web page.The estimated deposition rate for these deposition conditions is ≈ 0.05 mm/min for layer I and ≈ 0.03 mm/min for layer II.These values are in a reasonable agreement with the results ob-tained employing the laser technique, where the deposition rate is calculated to be 0.02 mm/min.This agreement validates our hypothesis that the ice structure of thick ices does not differ significantly from the structure of thin ices, growing homogeneously over time during the deposition.The calculated optical properties are independent from the total thickness of the ice sample, allowing us to relate the laboratory data to the astrophysical ice conditions.
Figure 6 shows the determination of CO ice parameters, which are the refractive index (a), the amplitude absorption coefficient (b), and both the real (c) and imaginary (d) parts of the complex dielectric permittivity -see Eqs. ( 1) and (2).

Discussion
A benefit of the direct reconstruction of the optical properties of the ices, provided by THz-TDS, is the detection of the frequency dependent amplitude and phase of the waveform in a broad frequency range as a result of a single measurement.These data eliminate the need of using the Kramers-Kronig relations for the reconstruction of the optical properties, excluding additional distortion of the experimental data by edge effects, which frequently appears as a result of the Hilbert integral transformation.This is of particular importance when dealing with broadband spectral kernels which are usually present even when operating at low temperature.
We could also compare the refractive index of CO ice at THz frequencies with that previously calculated in the mid-IR range, from Hudgins et al. (1993), Ehrenfreund et al. (1997), and Baratta & Palumbo (1998) for a CO ice deposited at 10 K.However, a direct comparison between refrative index values found by different authors is not straightforward.As discussed by Loeffler et al. (2005) and Baratta & Palumbo (2017), the density and in turn the refractive index of an ice sample could strongly depend on the experimental conditions such as temperature, growth angle, and deposition rate.We plan to test the effect of the change in the deposition conditions on the ice structure and check if this change will affect the optical constants.The results will be published in a forthcoming paper.
We compare the spectroscopic signature of the CO ice in our experiments with previous studies available in the literature.Data on far-IR spectra of solid CO were reported by Anderson & Leroi (1966) and Ron & Schnepp (1967).These studies investigated the absorption spectra of amorphous CO, deposited at 10 K on a crystalline quartz substrate, between 250 and 30 cm −1 .Two bands are visible at 50 and 83 cm −1 (1.5 and 2.5 THz).The spectral features observed in our experiments are in excellent agreement with these data, even though the 2.5 THz feature is masked by atmospheric water bands in our setup.The THz-TDS technique has a low sensitivity, requiring very thick ice layers to be detectable.Using our fast deposition rate the MIR vibrational bands of CO were strongly saturated within the first 30 seconds of deposition.We did some preliminary test on the ice growing using the vibrational bands in the NIR range and we were able to follow the ice growing up to ca. 160 µm thickness.Also in this case the ice growing is constant and linear in time during deposition, but the error associated with the thickness estimate using the NIR band is large (20-30%).
It is not surprising, then, that in the data reported by Ioppolo et al. (2014), on THz and mid-IR spectroscopy of interstellar ice analogs, the CO ice absorption band in the THz region has not been observed.The ice thickness in these experiments, in fact, has been estimated to be in all cases less than 10 µm.In our case, the minimum thickness required to be able to observe the 1.5 THz feature is estimated to be of the order of hundreds µm.
The data obtained are then employed to calculate the dust opacity for a given grain size distribution, as reported for exam-  and 2).For all deposition intervals, the dielectric curves demonstrate the existence of a Lorenz-like absorption peak, centered near 1.5 THz and featuring similar bandwidth.Distortions of the results seen at frequencies below 0.3 THz (such as an oscillatory character of n and ε as well as an increase of ε with decreasing frequency) are due to diffraction effects (see Sec. 2.2). ple in Ossenkopf & Henning (1994) (see Appendix C for additional details on the method we developed to reproduce their results).Following their approach, we report in Fig. 7 the calculated dust opacities assuming different ice coatings and different experimental data for the optical constants of ices.Data from Ossenkopf & Henning (1994) are labelled OH94 and reported for bare grains (green dotted line), thin (blue dotted), and thick ice mantles (orange dotted).The labels V = 0, 0.5, and 4.5 indicate the volume ratio between the core of refractory material and the ice mantle (see Appendix C).Conversely to the present work, Ossenkopf & Henning (1994) assume a H 2 O:CH 3 OH:CO:NH 3 = 100:10:1:1 mixture ice mantle composition, i.e. water-based, with a minor amount of methanol, carbon monoxide, and ammonia.As described in detail in Appendix C, following the same procedure as in Ossenkopf & Henning (1994), we extend the real and the imaginary parts of the refractive index from Hudgins et al. (1993) to longer wavelengths and we include spherical carbonaceous impurities.
The dielectric functions found by our experiments refer to pure CO ice, that explains the differences in the opacities shown in Fig. 7 calculated in the 100 -1000 µm range (solid lines).To have a more relevant comparison, we calculated the opacity using the same grain distribution and refractory materials as in Ossenkopf & Henning (1994) but with the optical values of CO ice coating extrapolated from Baratta & Palumbo (1998) (labelled BP98), where the refractive index of CO ice deposited at 12 K is calculated from the spectrum recorded in the 4400 -400 cm −1 (2.27 -25 µm) infrared spectral range.Their opacities (Fig. 7, dashed lines) show an agreement with our results, except for the contribution due to the presence of a CO ice absorbing feature at approximately 200 µm, which is absent in the extrapolated data, as expected.
When we compare the refractive index, we note that the real part by BP98 (n = 1.28) is reasonable close to our data (n = 1.24), while for the imaginary part the discrepancy is very large, since n is completely determined by the absorption feature at 200 microns that is out of the range investigated by Baratta & Palumbo (1998) (see Figure 6).It is worth to mention that Anderson & Leroi (1966) and Ron & Schnepp (1967) report a CO absorption feature at 2.5 THz, which is not visible in our spectra because it is masked by the atmospheric water contamination.This feature should further decrease the actual value of n at higher frequencies (above 2.5 THz).If we could measure the 2.5 THz feature and extend the calculation of the optical properties to this value, we would probably obtain a lower value of n , slightly increasing the discrepancy with the data by Baratta & Palumbo (1998).The reconstruction of the THz dielectric response of ices without the use of Kramers-Kronig relations, which is provided by THz pulsed spectroscopy, can provide an independent methodology to determine the optical properties of ice samples and validate the previous studies.
Since to compute the opacity both the real and the imaginary part of the dielectric constant are employed, we can infer from the calculated opacity curve that the imaginary part, which shows the biggest difference from the data presented by Baratta & Palumbo (1998), is not playing a major role in the determination of the opacity.This conclusion might be different for other absolute values or in a different spectral range.
Values of opacities for some selected wavelengths are also given in Table 1.The motivation for the choice of CO ice arises from the need of investigating the ice mantle properties for sources in which drastic CO depletion is expected, such as prestellar cores or protoplanetary disks midplanes (e.g.Caselli et al. 1999;Pontoppidan et al. 2003).In these cases, an ice mantle rich in CO can be formed, and it will influence the optical properties of the dust grains.Thus, it is interesting to compare how the opacities change when ice mantles with diverse chemical composition are present.The common moleculat species in astrophysical ices, such as H 2 O, CO 2 , NH 3 and possibly N 2 , present absorption features in the FIR range.Therefore, the study of the influence of the spectroscopic features on the opacity of the ice mantles is important and we plan to extend this study to other κ, cm g 2 10 Fig. 7. Calculated and reference opacities of astrophysical dust with CO ice and ice mixtures as a function of the wavelength.Dotted lines labelled with OH94 refer to bare grains and ice mixtures by Ossenkopf & Henning (1994), dashed lines with BP98 to CO ice by Baratta & Palumbo (1998), solid lines to the CO data by the present work.V indicates the volume ratio between refractory core and ice mantle, for which we follow Ossenkopf & Henning (1994), where V = 0 (black) is the bare grain, V = 0.5 (blue) thin ice, and V = 4.5 (orange) thick ice.See text for additional details.

Conclusions
In this work we have presented a study on the optical properties of solid CO at temperature and pressure conditions significant for astrophysical applications.While previous data in the mid-IR frequency range are available in the literature, to the best of our knowledge, this study is the first to provide the complex refractive index and complex dielectric permittivity of CO ice in the THz range.
We have shown that the ability of THz-TDS to measure both the amplitude and the phase information about the transmitted pulse provides direct reconstruction of the complex dielectric function of ices without the use of the Kramers-Kronig relations.The THz spectral features of ices can have a large bandwidth, such as the CO absorption line at 1.5 THz.In this case, an implementation of the Kramers-Kronig relations (i.e., the Gilbert B. M. Giuliano et al.: THz complex refractive index of CO ice transform), relying only on the power transmission/reflection spectrum, could lead to edge effects and resulting distortion of the dielectric response.Such distortions are of particular importance when a spectral feature of the sample is located near the border of the spectral operation range.Our results justify that THz-TDS setup is an appropriate instrument for accurate measurements of dielectric properties of ices at THz frequencies.
These results have been used to calculate the opacities of the dust grains covered by a CO ice layer.Discrepancies with currently available opacities suggest that measurements like those presented here are needed to provide a better interpretation of dust continuum emission, including dust and gas mass estimates.In addition, they will provide further insight on the radiative transfer processes based on ice analogs optical and physical properties.
Acknowledgements.The authors acknowledge the anonymous referee for providing very useful comments which significantly improved the quality of the manuscript.Also, the authors acknowledge Mr. Christian Deysenroth for the very valuable contribution in the designing and development of the experimental setup.We would like to thank Volker Ossenkopf for fruitful discussion to offering insight into the dust modelling methodology.The work of Arsenii A. Gavdush on alignment of the THz-TDS setup and digital processing of the THz waveforms was supported by the Russian Foundation for Basic Research (RFBR), Project #18-32-00816.The work of Gennady A. Komandin on solving the inverse problem of THz time-domain spectroscopy was supported by the Russian Science Foundation (RSF), Project #18-12-00328.Tommaso Grassi acknowledges the support by the DFG clusterof excellence "Origin and Structure of the Universe" (http://www.universe-cluster.de/).This work was funded by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) -Ref no.FOR 2634/1 ER685/11-1.This work has been partially supported by the project PRIN-INAF 2016, The Cradle of Life -GENESIS-SKA (General Conditions in Early Planetary Systems for the rise of life with SKA).

Appendix A: THz-TDS optics
In this section we provide a detailed description of the configuration of the laser beam inside the optical compartment of the setup, already shown in Fig. 2. The power of the laser beam is attenuated and successively divided into equivalent channels.Thus, the antenna-emitter is pumped and the antenna-detector is probed with an equal average power of about 20 mW.The optical delay between the pump and probe beams is varied using a double-pass linear mechanical delay stage from Zaber with the travel range of 101.6 mm and the positioning accuracy of < 3 µm.The THz radiation undergoes 10 kHz electrical modulation in order to synchronously detect the THz amplitude using a lock-in detection principle.
The THz beam emitted by the photoconductive antenna is collimated by an integrated HRFZ-Si hemispherical lens and then focused on a substrate window using a polymethylpentene (TPX) lens with the focal length of 67 mm and the diameter of 25 mm.After passing through the cryostat vacuum chamber, it is collimated by an equal TPX lens in the direction of the antennadetector.Finally, the THz beam is focused onto the photoconductive gap of the antenna-detector by an equal integrated HRFZ-Si hemispherical lens.In our measurements, during the waveform detection, we used a time-domain stride of 50 fs (it allows for satisfying the Whittaker-Nyquist-Kotelnikov-Shannon sampling theorem (Nyquist 1928)), at time-domain window size of 100 ps and an averaging time of 0.1 sec at each time-domain step with no waveform averaging.
The signal measured at the antenna-detector is recorded in the time domain E (t) and converted into its Fourier-spectrum E (ν) via where t and ν stand for time and frequency.

Appendix B: Reconstruction of the terahertz dielectric permittivity
Eqs.
(3) and ( 4) are obtained as follows.We defined the complex amplitude of electromagnetic wave E 0 = |E 0 | exp (iϕ 0 ), which interacts either with a bare substrate (which we use as reference), or with two ice layers (when detecting the sample waveform).
The pulses complex amplitudes depend on the initial complex amplitude of electromagnetic wave E 0 , as well as on ices and reference layers, which define the THz-wave reflection and transmission at the interfaces (R-operators and T -operators, respectively) and its absorption and phase delays in a bulk material (P-operators).Then, we define the amplitudes of the ballistic reference pulse as follows: where l is the thickness of the medium, and the symbols 0 and Si indicates respectively the vacuum and the HRFZ-Si, the ballistic sample pulse as follows: E S,b = E 0 P 0 l 0 − l Si − l CO,I − l CO,II T 0,CO P CO l CO,I × T CO,Si P Si (l Si ) T Si,CO P CO l CO,II T CO,0 where l CO,I and l CO,II are defined as in Section 3.1, and the two satellite sample pulses as follows: E S,1s = E 0 P 0 l 0 − l Si − l CO,I − l CO,II T 0,CO P CO l CO,I × T CO,Si P Si (l Si ) T Si,CO P CO l CO,II T CO,0 × R CO,Si P 2 CO l CO,I R CO,0 = E S,1s exp iϕ S,1s , E S,2s = E 0 P 0 l 0 − l Si − l CO,I − l CO,II T 0,CO P CO l CO,I × T CO,Si P Si (l Si ) T Si,CO P CO l CO,II T CO,0 × R CO,0 P 2 CO l CO,II R CO,Si = E S,2s exp iϕ S,2s , (B.3) which are clearly observed in reference and sample TDS waveforms in Fig. 4 (a).Then, by neglecting the phase changes during the reflection at the interfaces of absorbing media (these changes are weak since we study rather low-absorbing dielectric materials) as well as the distortion of optical pulses due to dispersion of material parameters (the dispersion of HRFZ-Si is very low, while the dispersion in ice is negligible due to its small thickness), we calculated the phases of these pulses: where n 0 and n Si are defined as in Section 3.1 and n CO is the refractive index of the CO ice.These phases are used to calculate the time delays between the pulses, which are indicated in Fig. 4 (a): Solving this system of equations yields Eqs. ( 3) and (4).
Figure 5 (b) shows the Fourier spectra E (ν) of the reference and sample TDS waveforms.In order to filter out the contribution of the satellite THz pulses (caused by the interference in the input and output windows and the substrate), and to improve the analysis of the frequency-domain data, we apply equal apodization procedure (window filtering) to all waveforms, where E (t) and E filt (t) stand for the initial and filtered waveforms, H (t) defines the apodization function, t 0 defines a position of the apodization filter towards the THz waveform, is a Tukey apodization filter (Tukey et al. 1986) with the width τ, and ω stands for a parameter of the filter smoothness (for ω = 0 the window has a rectangular form, for ω = 1 it is the Hann window Harris (1978)).As shown in Fig. 5(a), we use the Tukey window centered at the maximum of the reference THz waveform, with the smoothness parameter of ω = 0.1 and the width of 40 ps (which yields the frequency-domain resolution of 0.025 THz).
Let us consider the Fresnel formulas (Born & Wolf 1980), defining the THz wave amplitude reflection at (and transmission through) the interface between media m and k: where R m,k (ν) and T m,k (ν) stand for coefficients of the complex amplitude reflection and transmission, respectively, while n m (ν) + in k (ν) is the complex refractive index of the media.The relation between complex amplitudes of the THz wave right after the emitter (z = 0), E 0 (ν), and at the position z along the beam axis, E (ν, z), is given by If the thicknesses and the refractive indexes of all layers are known, Eqs.(B.8) and (B.9) yield description of all peculiarities of the THz pulse interacting with multilayer structures (Zaytsev et al. 2014(Zaytsev et al. , 2015;;Gavdush et al. 2019).
We derive the equations describing the complex amplitudes of the reference E R (ν) and sample E S (ν) spectra, assuming only wavelets inside the Tukey apodization.
For the reference spectrum, we obtain where the indexes 0 and Si correspond to the free space and the HRFZ-Si medium; l 0 and l Si are the total length of the THz beam path and the thickness of the HRFZ-Si substrate, respectively; T 0,Si (ν) and T Si,0 (ν) are the transmission coefficients for the respective interfaces (between the free space and the HRFZ-Si), defined by Eq. (B.8); P 0 (ν, z) and P Si (ν, z) are operators describing the THz wave propagation in the free space and the HRFZ-Si, respectively, as given by the exponential factor in Eq. (B.9).
As shown in Fig. 3, for the sample spectrum we take into account the contribution of the ballistic THz pulse (1) and the satellite pulses (2 and 3), caused by the multiple THz wave reflection in the ice films.This yields the following equation E S /E 0 = P 0 l 0 − l Si − l CO,I − l CO,II T 0,CO P CO l CO,I × T CO,Si P Si (l Si ) T Si,CO P CO l CO,II T CO,0 × 1 + R CO,Si P 2 CO l CO,I R CO,0 +R CO,0 P 2 CO l CO,II R CO,Si , (B.11) where the summation terms in the brackets correspond to the wavelets 1, 2 and 3; l CO,I and l CO,II stand for thicknesses of the first and second ice films; R CO,Si (ν) and R CO,0 (ν) are the transmission coefficients for the respective interfaces, see Eq. Considering all reflection and transmission operators in Eq. (B.9), we note that T Th depends only on the refractive index of the HRFZ-Si substrate n Si , which is known a priori, as well as on the parameters of the CO ice to be determined; it excludes contribution of several factors, such as the unknown complex amplitude of the TDS source E 0 (ν), the unknown total length of the THz beam path l 0 , and, finally, the thickness of the HRFZ-Si substrate l Si , which is known, too, but can slightly vary owing to angular deviations of the substrate during the vacuum chamber assembling.
The experimental transfer function T Exp is calculated in a similar manner, relying on the Fourier spectra of the experimental sample E S and reference E R waveforms (after applying the Tukey apodization), .15)bhmie.pyfor the bare grains 4 and with bhcoat.pyfor the coated grains 5 , both from (Bohren & Huffman 1983).In Ossenkopf & Henning (1994) there are three cases, bare grains without ice, thin ice, that has a volume ratio of V = 0.5, and thick ice, where V = 4.5.The radius of the refractory core a and the fraction V determine the radius of the mantle a coat = a(V + 1) 1/3 .

Appendix C.3: Opacities
With Q abs it is possible to retrieve the opacities averaged on the grain size distribution ϕ(a) as where a min to a max is the range where the size distribution is valid, but considering only the refractory material, i.e. silicates or carbonaceous dust.Ossenkopf & Henning (1994) assume a min = 5 × 10 −7 cm, a max = 2.5 × 10 −5 cm, ρ 0 = 2.9 g cm −3 (silicates) and ρ 0 = 2 g cm −3 (amorphous carbon), and ϕ(a) = a −3.5 .
The total opacity is κ tot (λ) = 0.678κ S i (λ) + 0.322κ AC (λ), where the two terms are respectively the silicates and carbonaceous opacities including ice coating, and the two coefficients are calculated from the volume ratio discussed in Sect.3.1 of Ossenkopf & Henning (1994).In particular, Ossenkopf & Henning assume a volume ratio of the refractory components V AC /V Si = 0.69 (their Sect.3.1), that can be converted into the corresponding mass ratio M AC /M Si = 0.4758 by using the relation M i = ρ i V i , where ρ i are the bulk densities of the two refractory components (being in Ossenkopf & Henning the opacity defined per unit mass of the refractory material), so that κ tot = (M AC κ AC + M S i κ S i )/(M AC + M S i ).

Fig. 2 .
Fig.2.THz-TDS setup for the spectroscopy of ices.(a) A schematic of the setup, where FS laser stands for the femtosecond laser, M stands for the optical mirrors, ATT stands for the attenuator of the laser beam intensity, BS stands for the optical beams splitter, DS stands for the mechanical double-pass delay stage, PCAE and PCAD stand for the photoconductive antenna-emitter and antenna-detector, respectively, L stand for the TPX lenses, S and W stand for the HRFZ-Si substrate and windows, respectively.(b) Spectra of THz waveforms E (ν) transmitted through the empty THz beam path or the THz beam path with the cryostat; the shaded region below ≈ 0.3 THz indicates the spectral range where distortions due to the THz beam diffraction on the aperture of the substrate are expected.

l
Fig.3.A time-distance diagram illustrating the THz wave propagation through the HRFZ-Si substrate with ice films deposited on its surfaces.Lines 0 to 3 illustrate the ballistic pulse and satellite pulses transmitted in the direction of the antenna-detector, N stands for unaccounted satellites with larger time delays.Solid lines represent the pulses which have been used for the analysis, while dotted lines correspond to neglected pulses.

Fig. 4 .
Fig. 4. Calculation of the thicknesses l CO,I , l CO,II and the complex refractive index n of the CO ice films (see Fig.3).(a) Time delays between the ballistic THz pulses of the reference (0) and sample (1) waveforms, δt 01 ; the first satellite pulse (2) and the ballistic pulse (1) of the sample waveform, δt 12 ; and the second satellite pulse (3) and the ballistic pulse (1) of the sample waveform, δt 13 (the pulses are marked as in Fig.3).(b,c) Estimates for the thicknesses of the two ice films as a function of the total deposition time t dep for the different deposition intervals of ∆t dep = 4, 5 and 6 min; the first assumption for the real part of the refractive index of ice is between n init = 1.230 and 1.255 ± 0.035, for the imaginary part we first set n init = 0.

Fig. 5 .
Fig. 5. Evolution of the THz pulse and its spectra during the CO ice deposition.(a) Reference waveform E (t) transmitted through the cryostat with the empty substrate (green), and sample waveforms (black to light red) transmitted through the substrate with the CO ice deposited on its surfaces.(b) Fourier spectra |E (ν) | of the reference and sample THz waveforms calculated with the use of Tukey apodization.The waveforms in (a) and spectra in (b) correspond to different values of the total deposition time t dep (indicated), with the deposition intervals of ∆t dep = 6 min.

Fig. 6 .
Fig. 6.Optical properties of CO ice.(a) Real part of the refractive index, (b) amplitude absorption coefficient, (c) real and (d) imaginary parts of the dielectric permittivity (see Eqs.1 and 2).For all deposition intervals, the dielectric curves demonstrate the existence of a Lorenz-like absorption peak, centered near 1.5 THz and featuring similar bandwidth.Distortions of the results seen at frequencies below 0.3 THz (such as an oscillatory character of n and ε as well as an increase of ε with decreasing frequency) are due to diffraction effects (see Sec. 2.2).

|T
(B.8); the definition of the remaining factors (P and T) is similar to that in Eqs.(B.10).We point out that the complex amplitudes (E) and all factors (T , R, and P) in Eq. (B.10) and (B.11) are frequencydependent.Equations (B.10) and (B.11) form a basis for the reconstruction of the THz dielectric response of ices.The reconstruction is performed via the following minimization procedure:n (ν) = argmin n(ν) [Φ (n (ν) , ν)] , (B.12)where argmin is an operator which determines the minimum argument of the vector error functional Φ.The latter is formed from the complex theoretical T Th and experimentalT Exp transfer functions, Th (n, ν)| − T Exp (ν) φ [T Th (n, ν)] − φ T Exp (ν) ..| and φ [...] are, respectively, the absolute values and the phases of the complex functions.We define the theoretical transfer function T Th as the sample spectrum (Eq.B.11) normalized by the reference spectrum (Eq.B.10), T Th = T 0,CO T CO,Si T Si,CO T CO,0 T 0,Si T Si,0 × P CO l CO,I + l CO,II P 0 −l CO,I − l CO,II × 1 + R CO,Si P 2 CO l CO,I R CO,0 +R CO,0 P 2 CO l CO,II R CO,Si .(B.14)