# Triple oxygen isotope systematics of evaporation and mixing processes in a dynamic desert lake system

Claudia Voigt, Daniel Herwartz, Cristina Dorador, Michael Staubwasser

<figure><img src="image_1.png" /></figure>

## To cite this version:

Claudia Voigt, Daniel Herwartz, Cristina Dorador, Michael Staubwasser. Triple oxygen isotope systematics of evaporation and mixing processes in a dynamic desert lake system. Hydrology and Earth System Sciences, 2021, 10.5194/hess-25-1211-2021. [https://doi.org/10.5194/hess-25-1211-2021]

HAL Id: hal-03000434

https://hal.science/hal-03000434

Submitted on 7 Jun 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_3.png" /></figure>

# Triple oxygen isotope systematics of evaporation and mixing processes in a dynamic desert lake system

Claudia Voigt<fnref n="1" />,<fnref n="2" />, Daniel Herwartz<fnref n="1" />, Cristina Dorador<fnref n="3" />, and Michael Staubwasser<fnref n="1" />

**Correspondence:** Claudia Voigt (voigt@cerege.fr)

Received: 26 May 2020 – Discussion started: 15 June 2020

Revised: 25 January 2021 – Accepted: 3 February 2021 – Published: 10 March 2021

**Abstract.** This study investigates the combined hydrogen deuterium and triple oxygen isotope hydrology of the Salar del Huasco, an endorheic salt flat with shallow lakes at its centre that is located on the Altiplano Plateau, N Chile. This lacustrine system is hydrologically dynamic and complex because it receives inflow from multiple surface and groundwater sources. It undergoes seasonal flooding, followed by rapid shrinking of the water body at the prevailing arid climate with very high evaporation rates. At any given point in time, ponds, lakes, and recharge sources capture a large range of evaporation degrees. Samples taken between 2017 and 2019 show a range of δ¹⁸O between -13.3‰ and 14.5‰, d-excess between 7‰ and -100‰, and ¹⁷O-excess between 19 and -108 per meg. A pan evaporation experiment conducted on-site was used to derive the turbulence coefficient of the Craig-Gordon isotope evaporation model for the local wind regime. This, along with sampling of atmospheric vapour at the salar (–21.0 ± 3.3‰ for δ¹⁸O, 34 ± 6‰ for d-excess and 23 ± 9 per meg for ¹⁷O-excess), enabled the accurate reproduction of measured ponds and lake isotope data by the Craig-Gordon model. In contrast to classic δ²H–δ¹⁸O studies, the ¹⁷O-excess data not only allow one to distinguish two different types of evaporation – evaporation with and without recharge – but also to identify mixing processes between evaporated lake water and fresh flood water. Multiple generations of infiltration events can also be inferred from the triple oxygen isotope composition of inflow water, indicating mixing of sources with different evaporation histories. These processes cannot be resolved using classic δ²H–δ¹⁸O data alone. Adding triple oxygen isotope measurements to isotope hydrology studies may therefore significantly im-

prove the accuracy of a lake’s hydrological balance – i.e. the evaporation-to-inflow ratio (E/I) – estimated by water isotope data and application of the Craig-Gordon isotope evaporation model.

## 1 Introduction

The majority of water in the hydrologic cycle on Earth represents a dynamic isotopic equilibrium between evaporation from the ocean, precipitation, and continental runoff, defining a linear relationship between δ¹⁷O and δ¹⁸O or δ²H and δ¹⁸O – i.e. the global meteoric water line (GMWL; Clark and Fritz, 1997; Luz and Barkan, 2010). However, a number of continental water reservoirs, e.g. lakes in semi-arid or arid environments, may deviate from that state as a result of evaporation – which imparts an enhanced kinetic isotope effect on these waters – or because of mixing with pre-evaporated water. Such deviations are described by the secondary parameters ¹⁷O-excess [= δ'¹⁷O – 0.528 · δ'¹⁸O with δ' = 1000 · (δ/1000 + 1)] and d-excess [= δ²H – 8 · δ¹⁸O], which assume values of 33 per meg and 10‰, respectively, for waters plotting on the GMWL. The progress of evaporation causes a systematic decrease of ¹⁷O-excess and d-excess that is principally predictable by the Craig-Gordon (C-G) isotope evaporation model. The C-G model is the foundation of assessing the hydrological balance, i.e. the evaporation-to-inflow ratio (E/I), of lakes using water isotope measurements. The model variables of the C-G equation – relative humidity, temperature, the isotopic composition of atmospheric vapour, continuous groundwater recharge, wind

<footnote n="1">¹Institute of Geology and Mineralogy, University of Cologne, Zülpicher Str. 49b, 50674 Cologne, Germany</footnote>

<footnote n="2">²Aix Marseille University, Centre National de la Recherche Scientifique (CNRS), CEREGE, Aix-en-Provence, France</footnote>

<footnote n="3">³Centro de Biotecnología, Universidad de Antofagasta, Angamos 601, 1270300 Antofagasta, Chile</footnote>
turbulence above the water surface, and salinity – depend on local climatic conditions and require quantification in order to apply the model (Gonfiantini, 1986; Horita et al., 2008; Gonfiantini et al., 2018). Of these variables, wind turbulence remains beyond determination by direct measurement in the field. Importantly, the mixing of different lake water sources, particularly if these represent precipitation from different seasons or different evaporation histories, presents a significant challenge in deriving the hydrologic balance of lakes from water isotopes in the C–G model that cannot be resolved using classic ²H / ¹H and ¹⁸O / ¹⁶O data (Gonfiantini, 1986). Without the means to address wind turbulence and mixing of sources directly, isotope hydrology cannot provide an unambiguous estimate of a lake’s hydrological balance.

The analysis of the ¹⁷O/¹⁶O ratio in H₂O in addition to the commonly investigated ²H/¹H and ¹⁸O/¹⁶O ratios has become an increasingly applied tool in hydrogeological research, because of some favourable properties of the triple oxygen isotope system. For example, while the secondarily derived d-excess parameter is largely dependent on relative humidity and temperature, the ¹⁷O-excess parameter has shown to be temperature insensitive in the natural range between 11–41 °C (Barkan and Luz, 2005; Cao and Liu, 2011). Also, the effect of high salinity on the evaporation trajectory of the C–G model is less significant (Surma et al., 2018). Hence, the combined analysis of δ¹⁷O, δ¹⁸O, and δ²H provides information on fundamental processes of the hydrological cycle such as humidity, moisture sources, evaporation conditions, and mixing, which would be masked by temperature variability if δ¹⁸O and δ²H were analysed alone (Angert et al., 2004; Barkan and Luz, 2005, 2007; Landais et al., 2006, 2008, 2010; Surma et al., 2015, 2018; Herwartz et al., 2017). As such, the general potential of the ¹⁷O-excess parameter as a tool to quantitatively reconstruct paleo-humidity from plant silica particles (Alexandre et al., 2018, 2019) and lake sediments (Evans et al., 2018; Gázquez et al., 2018) could be demonstrated.

The seasonal dynamics in climate variables, flooding, and their respective effect of mixing on the triple oxygen isotope composition of lakes have not yet been investigated. The present study's objectives are to (1) test the potential of triple oxygen isotope analyses to resolve the fundamental hydrologic processes of evaporation, recharge and mixing of sources that cannot be resolved by the classical δ²H-δ¹⁸O analyses; (2) test the robustness of the C–G model as a reliable monitoring tool for the lake hydrologic balance in a hydrologically complex and seasonally dynamic desert lake system; and (3) demonstrate the potential of triple oxygen isotope analyses to derive the hydrological balance of lakes from water isotope and climate monitoring. For this purpose, we measured the isotopic data from springs, ponds, and lake water samples from the Salar del Huasco, N Chile, covering a 300‰ salinity range, and compared them with evaporation trajectories modelled for local seasonal conditions. All variables of the C–G model were constrained by meteorological

<figure><img src="image_4.png" /><figcaption>**Figure 1.** Study area. (a) Catchment of the Salar del Huasco (Salar del Huasco basin) with drainage. The blue triangle marks the sampling location of the Collacagua river. (b) Overview map (DEM derived from SRTM data, created using ArcGIS 10.5.1).</figcaption></figure>

data from an on-site weather station, a pan evaporation experiment, isotopic measurements of atmospheric water vapour samples, and inflowing water sources. Due to the extreme dynamics in climatic and hydrological parameters, the Salar del Huasco provides an ideal environment to test the potential of triple oxygen isotope analyses in resolving different hydrological processes like permanent recharge and highly episodic mixing with flood water and to identify changes in the hydrological balance of lakes.

## 2 Study area

The Salar del Huasco is an endorheic salt flat located in the south of a longitudinal volcano-tectonic depression at the western margin of the Chilean Altiplano at about 3770 m (Fig. 1). It covers an area of about 50 km², but only 5% of the surface is permanently covered by water (Risacher et al., 2003). The hydrological balance of the salar is controlled by a shallow groundwater table, perennial streams, rare and highly seasonal precipitation, and episodic injection of runoff water after heavy rainfalls or snowmelt. The environment is characterized by exceptionally high evaporation rates and high variability in relative humidity and temperature throughout the year with considerable diurnal amplitude.

An on-site weather station (20°15.42' S, 68°52.38' W) has been operated by the Centro de Estudios Avanzados en Zonas Áridas (CEAZA) at the north-western margin of the salar since September 2015 (CEAZA, 2020). Over this period, mean annual relative humidity and temperature are 41% and 5°C, respectively. Both relative humidity and temperature are highly variable throughout the year with
mean seasonal values of 54% and 8°C in austral summer (December–March) and 33% and 2°C in austral winter (June–September). Day–night fluctuations range between 0% and 99% and from –19°C to 23°C. Mean annual wind speed is approximately 2.5 m s⁻¹, most times coming from S to SW direction. Typically, it is calm in the morning (<1 m s⁻¹) but very windy in the afternoon with 5 m s⁻¹ average wind speed and gusts up to 24 m s⁻¹.

Long-term mean annual precipitation is 150 mm a⁻¹, but interannual variability is high due to their dependence on wind patterns (Garreaud and Aceituno, 2001; Risacher et al., 2003). Precipitation occurs mainly in austral summer (December–March) and is related to easterly airflow from the Atlantic Ocean and the Amazon Basin (e.g. Aravena et al., 1999; Garreaud et al., 2003; Houston, 2006). In contrast, less frequent winter rains and snow are associated with the interaction of cold air masses from the Pacific with tropical air masses from the Amazon Basin (Aravena et al., 1999). The intensity of convective storms and variability of the con tinental effect in the Amazon Basin lead to significant iso topic variability in summer rains (Aravena et al., 1999). In contrast, winter rains comprise generally higher δ¹⁸O values and are isotopically less variable (Aravena et al., 1999). This may be attributed to lower variability in the contribution of pre-evaporated water from the Amazon Basin in the dry sea son and the non-convective nature of winter storms.

The hydrogeological system of the Salar del Huasco con sists of three aquifers, where the upper and the intermedi ate aquifer are located in the basin's sediment infill, sepa rated by a thin aquitard, and the lower aquifer is formed in volcanic bedrock (Acosta and Custodio, 2008). The upper aquifer is recharged by the Collacagua river that drains the northern part of the basin and completely infiltrates 10 km before reaching the salar (Acosta and Custodio, 2008). In pe riods of heavy rainfall, the river directly flows into the salar leading to widespread flooding of the northern area (Fig. 2e). Several springs and creeks around the salar and the shal low groundwater table contribute to the formation of peren nial lakes. Furthermore, ephemeral ponds emerge episod ically due to precipitation, surface and subsurface runoff, and a rising groundwater table in the rainy season (austral summer). Episodic precipitation, surface runoff, ephemeral flooding, and variable contributions from different ground water and stream sources result in a highly dynamic system with strongly fluctuating lake level.

3 Sampling

Natural water samples from the Salar del Huasco were taken
during field campaigns in September 2017, September 2018,
and March 2019 (Fig. 2 and Table S2). The sample set in-
cludes springs, perennial lakes (Laguna Grande in the W,
and Laguna Jasure in the SE) and ephemeral ponds of which
several were apparently recharged and some apparently stag-

nant. Lakes and ponds are all very shallow (5–30 cm). Sam-
ples were taken from the water surface, and temperature, pH,
and conductivity were measured on-site using a digital pre-
cision meter Multi 3620 IDS.

Based on satellite images and field observations, five ma-
jor hydrological subsystems could be identified (Fig. 2). The
Laguna Grande in the western part of the salt flat is interme-
diate saline with values varying seasonally between approxi-
mately 20 and 35 g L⁻¹.

A channel originating from the southern end of this la-
goon extends along the southern margin of the salt flat. This
channel contributes to flooding of the southern area during
the rainy season, while lowering of the water table over the
course of the year leads to isolation of ponds (Fig. 2e and
g). The large range in salinity (1 to >100 g L⁻¹) observed in
shallow lakes and ponds at the south to south-western mar-
gin of the salar may indicate that their hydrological balance
alternates between recharge and isolation.

In the south-eastern area, a perennial spring fed a chain of
lakes with generally low salinity of <2 g L⁻¹. Some ponds
were visibly connected by streams and creeks. Others located
closer to the eastern margin of the salt flat were topograph-
ically elevated and therefore isolated from the main south-
eastern inflow. Satellite images indicate that ponds closer to
the salar's centre may represent a mixture of source water
from the south-eastern inflow and a further subsystem in the
north (Fig. 2e).

The northern subsystem is probably fed by subsurface in-
flow of the Collacagua river. Its area may become widely
flooded in austral summer (Fig. 2e) but dries up rapidly af-
ter the rainy season (Fig. 2g). Ponds from the northern area
sampled in September 2017 were characterized by an ex-
treme salinity gradient from 2 to 343 g L⁻¹, broadly decreas-
ing from east to west (Table S2). These ponds were often
salt-encrusted and showed black, sulfur-reducing microbial
mats at the bottom. In March 2019, this entire area was cov-
ered by a large low-salinity lake.

In the north-western area, a spring originating from a
small, vegetated hill fed a series of ponds with salin-
ities <1 g L⁻¹. An adjacent chain of ponds sampled in
September 2017 about 750 m to the south-east with salini-
ties between 1 and 6 g L⁻¹ was not visibly connected to this
spring and must have been sustained by the shallow ground-
water table. Three topographically elevated ponds close to
the margin had higher salinity of 8 to 40 g L⁻¹. The area
around these ponds was more vegetated, indicating that they
were older, i.e. represented an earlier flooding but became
isolated from recharge. Despite extensive rainfall-induced
flooding in February 2019, the north-western part of the salt
flat was almost dry in March 2019.

Besides water samples from surface waters, we sampled
atmospheric vapour using a Stirling cooler cycle system (Le-
Tehnika, Kranj, Slovenia) built after Peters and Yakir (2010).
Sampling was carried out on three days during field cam-
paigns in September 2017 and March 2019 (Table 1). At-
<figure><img src="image_5.png" /></figure>

**Figure 2.** Illustration of the hydrological situation at the Salar del Huasco over the sampling period with sample locations for field campaigns in September 2017 (a, b), September 2018 (c, d), and March 2019 (e, f). Panel (g) reflects the hydrological situation at the Salar del Huasco 6 months after the last sampling campaign. Several hydrological subsystems (coloured areas) were identified based on satellite images (Copernicus Sentinel data, 2017, 2018, 2019) and field observations. Different symbols of sample locations refer to the corresponding hydrological subsystem (see text for details; DEM derived from SRTM data, created using ArcGIS 10.5.1).
<figure><img src="image_6.png" /><figcaption>**Figure 3.** Conceptual comparison of isotope effects associated with the three principal hydrological processes – simple evaporation without recharge, evaporation with recharge, and isotopic mixing – in the diagram of (a) $^{17}$O-excess over $\delta^{18}$O and (b) d-excess over $\delta^{18}$O. In the case of simple evaporation, the isotopic composition of an evaporating water body evolves along the black line. Continuous recharge drives the isotopic composition of the evaporating water body below the simple evaporation trajectory onto the green line. The isotopic composition of a recharged lake is determined by the evaporation-to-inflow (E/I) ratio, where increasing E/I leads to higher $\delta^{18}$O and lower $^{17}$O-excess values. In the case of a terminal lake (E/I = 1) all inflow is balanced by evaporation. The dashed grey line exemplifies admixture of fresh water similar in isotopic composition to the inflowing water, e.g. flood water, to the evaporated brine of a terminal lake.</figcaption></figure>

mospheric vapour was sampled at ground level. The water vapour was captured by streaming air with a flow rate of 600–800 mL min$^{-1}$ through a 4 mL vial, which was attached to the cold finger of the Stirling cooler and insulated. Extremely low absolute humidity values at the Salar del Huasco required sampling over several hours to yield 0.1–0.2 µL of atmospheric vapour.

Additionally, pan evaporation experiments with 600, 800, and 1000 mL of fresh water (0.805 mS cm$^{-1}$) filled in stainless-steel evaporation pans (Ø 20 cm) were carried out on-site over a period of three days (20°15.4' S, 68°52.4' E). The rationale behind this experimental design was to account for the extreme diurnal and significant average day-to-day variability of all variables in the C–G model by letting evaporation progress over several days. By varying the initial volume, we were able to achieve a wide spread of fractional evaporative water loss. Samples were taken every day around 18:00 local time (CLT) and additionally after 13:00 CLT on the third day. The 600 mL pan dried up before the end of the experiment so that no sample could be taken in the evening of the third day. Air temperature, relative humidity, and wind speed were monitored locally at the experiment at about 1.5 m above ground using a Kestrel 5500 weather meter (Fig. S2). Mean temperature and relative humidity were 4 °C and 35 % over the whole period of the experiment. As temperatures dropped below 0 °C in the night, a considerable fraction, if not all, of the water in the pans froze during the night. This led to the formation of either a solid ice block or a thick ice layer above the remaining liquid. Consequently, evaporation from pans was mostly restricted to daytime when samples were thawed. The effective evaporation time interval was assumed to correspond to the period when $T>0$ °C, resulting on average air temperature and rela-

tive humidity values of 10 °C and 22 % (Fig. S2). Water temperatures measured during sampling were found to be up to 5 °C lower than ambient air temperatures. However, during midday, temperatures of water may exceed air temperatures by several degrees due to solar heating of the pans. Winds were very strong between 12:00 and 19:00 with an average wind speed of 5 m s$^{-1}$ and gusts up to 14 m s$^{-1}$ coming from S to W direction.

# 4 Methods

## 4.1 Isotope and chemical analyses

The hydrogen and triple oxygen isotope composition of water samples were analysed by isotope ratio mass spectrometry (IRMS). Complementary concentration data of $Na^+$, $K^+$, $Ca^{2+}$, $Mg^{2+}$, $Cl^-$, and $SO_4^{2-}$ in natural samples were determined by ICP-OES (Table S2). A short discussion of the chemical data is provided in the Supplement (Sect. S2, Fig. S7).

Hydrogen isotope ratios are measured by continuous-flow IRMS of H$_2$. Water samples are injected in a silicon carbide reactor (HEKAtech, Wegberg, Germany) that is filled with glassy carbon and heated to 1550 °C, where they are reduced to H$_2$ and CO. The produced gases are carried in a helium gas stream (100–130 mL min$^{-1}$), separated by gas chromatography (GC) and finally introduced in a Thermo Scientific MAT 253 mass spectrometer for hydrogen isotope analysis. The long-term external reproducibility (SD) is 0.9 ‰ and 1 ‰ for $\delta^2$H and d-excess, respectively.

For triple oxygen isotope analysis, water samples are fluorinated, followed by dual-inlet IRMS of O$_2$. The method is described in detail in Surma et al. (2015) and Herwartz et
<table><caption>**Table 1.** Sampling intervals and isotopic data of atmospheric vapour.</caption><thead><tr><th>Day</th><th>Time</th><th>δ<sup>17</sup>O<br>(‰)</th><th>δ<sup>18</sup>O<br>(‰)</th><th><sup>17</sup>O-excess<br>(per meg)</th><th><i>n</i></th><th>δ<sup>2</sup>H<br>(‰)</th><th>d-excess<br>(‰)</th><th><i>n</i></th></tr></thead><tbody><tr><td>20 September 2017</td><td>18:30–21:30</td><td>-11.1</td><td>-21.0</td><td>19</td><td>2</td><td>-136</td><td>32</td><td>1</td></tr><tr><td>21 September 2017</td><td>18:30–21:30</td><td>-9.4</td><td>-17.8</td><td>16</td><td>2</td><td>-114</td><td>28</td><td>1</td></tr><tr><td>mean</td><td></td><td>-10.3</td><td>-19.4</td><td>18</td><td></td><td>-125</td><td>30</td><td></td></tr><tr><td>SD</td><td></td><td>1.2</td><td>2.3</td><td>2.1</td><td></td><td>15</td><td>2.8</td><td></td></tr><tr><td>20 March 2019</td><td>10:00–13:00</td><td>-12.9</td><td>-24.3</td><td>33</td><td>2</td><td>-154</td><td>41</td><td>1</td></tr><tr><td>overall average</td><td></td><td>-11.1</td><td>-21.0</td><td>23</td><td></td><td>-135</td><td>34</td><td></td></tr><tr><td>SD</td><td></td><td>1.7</td><td>3.2</td><td>9</td><td></td><td>20</td><td>6</td><td></td></tr></tbody></table>

al. (2017). In brief, 2.8 µL of water is injected in a heated CoF₃ reactor (370 °C) that is continuously flushed with helium (30 mL min⁻¹). The produced oxygen gas is cryogenically purified and trapped in 1 of 12 sample tubes of a manifold. The manifold is connected to a Thermo Scientific MAT 253 for dual-inlet IRMS analysis. The long-term external reproducibility (SD) is 0.12‰, 0.25‰, and 8 per meg for δ¹⁷O, δ¹⁸O, and ¹⁷O-excess, respectively. All isotope data herein are reported on SMOW-SLAP scale (Schoenemann et al., 2013). The scale is usually contracted using the setup described herein. This can partly be attributed to blank contribution (Herwartz et al., 2017). We observed an increase in scale contraction over the usage period of a CoF₃ reactor filling and a reduction of precision and accuracy of isotopic data, indicating that the blank contribution increases with time. To account for this effect, SMOW-SLAP scaling was performed daily using internal standards. Isotope measurements with anomalous high scaling factors or standard deviations were discarded.

## 4.2 The Craig-Gordon isotope evaporation model at the Salar del Huasco

The Craig-Gordon (C-G) isotope evaporation model forms the basis for numerous models describing the isotopic composition of natural lakes (e.g. Gonfiantini, 1986; Gat and Bowser, 1991; Horita et al., 2008). There are two principal evaporation scenarios – one without recharge (hereinafter termed “simple evaporation”) and one with continuous recharge (hereinafter termed “recharge evaporation”) (Craig and Gordon, 1965; Criss, 1999; Horita et al., 2008). In the case of simple evaporation, the isotopic composition of the lake is only controlled by the degree of evaporation (Fig. 3; Gonfiantini et al., 2018):

$$R_W = f^B \cdot \left(R_{WI} - \frac{A}{B} \cdot R_V\right) + \frac{A}{B} \cdot R_V,$$

where $R_{WI}$ denotes the initial isotopic ratio (²H/¹H, ¹⁷O/¹⁶O, ¹⁸O/¹⁶O) in the water body, $R_V$ is the isotopic ratio in atmospheric vapour, and $f$ is the fraction of residual

water. Note that $R$ can be inferred from the $δ$ notation by $R = δ/1000 + 1$. The parameters $A$ and $B$ describe the isotopic fractionation associated with evaporation dependent on the relative humidity, $h$, normalized to the water surface temperature:

$$A = - \frac{h}{\alpha_{\text{diff, l-v}} \cdot (1 - h)} \\ B = \frac{1}{\alpha_{\text{eq, l-v}} \cdot \alpha_{\text{diff, l-v}} \cdot (1 - h)} - 1.$$

Diffusive ($α_{diff, l-v}$) and equilibrium ($α_{eq, l-v}$) isotope fractionation cause a systematic increase in δ¹⁸O and decrease in d-excess and ¹⁷O-excess with progressive evaporation. A detailed description of all variables used herein is given in the Supplement (Sect. S1, Table S1).

In the case of recharge evaporation, the isotopic effect of continuous inflow is accounted for by the evaporation-to-inflow ratio (E/I) (Fig. 3; Criss, 1999):

$$R_{WS} = \frac{\alpha_{\text{eq, l-v}} \cdot \alpha_{\text{diff, l-v}} \cdot (1 - h) \cdot R_{\text{WI}} + \alpha_{\text{eq, l-v}} \cdot h \cdot E/I \cdot R_{\text{V}}}{E/I + \alpha_{\text{eq, l-v}} \cdot \alpha_{\text{diff, l-v}} \cdot (1 - h) \cdot (1 - E/I)}.$$

Here, $R_{WI}$ refers to the isotopic ratio in the inflowing water. Under steady-state conditions $E/I \le 1$, whereas the lake begins to desiccate when evaporation exceeds the inflow ($E/I > 1$).

The effect of wind turbulence – which leads to variable proportions of diffusive and turbulent isotope fractionation in dependence on wind speed – is accounted for by inserting a correcting exponent to the diffusive fractionation factor, $α_{diff, l-v}^n$ (Dongmann et al., 1974). The turbulence coefficient $n$ can vary between 0 (fully turbulent atmosphere) and 1 (calm atmosphere) but typically assumes values of $n \ge 0.5$ under natural conditions (Gonfiantini, 1986; Mathieu and Bariac, 1996; Surma et al., 2018). Recent laboratory experiments suggest that part of the evaporating water is removed by spray and microdroplet vaporization without any isotope fractionation when wind speed exceeds ~0.5 m s⁻¹, which
may require a modification of the above-described C–G approach (Gonfiantini et al., 2020). We evaluated the effect of a possible partial evaporative water loss without fractionation on the simple evaporation trajectory by introducing a “virtual outlet” $f_{out}$ in the isotopic mass balance equation, which finally results in a modification of the parameter B:

$$B_{out} = \frac{1 - f_{out}}{\alpha_{eq, l-v} \cdot \alpha_{diff, l-v}^{n} \cdot (1 - h)} + \frac{f_{out}}{(1 - h)} - 1.$$

Salinity affects isotope activities and increases fluid viscosity thereby decreasing the vapour pressure above the water body. In the C–G model, this may be accounted for by correcting equilibrium fractionation factors for the classic salt effect (Horita, 1989, 2005; Horita et al., 1993) and using effective rather than actual relative humidity. The salinity effect requires consideration within the $\delta^2\text{H}-\delta^{18}\text{O}$ system at salinities >100 g L$^{-1}$ and is much less significant for the $\delta^{17}\text{O}-\delta^{18}\text{O}$ system (Sofer and Gat, 1972, 1975; Horita, 1989, 2005; Surma et al., 2018). To avoid unnecessary complication of the more principal approach of this study, salinity effects were neglected in our model calculations.

The C–G model does not account for mixing processes, e.g. as a result of flooding or snowmelt, but can be used to calculate such effect by applying mass balance: $\delta X_{mix} = f \cdot \delta X_1 + (1 - f) \cdot \delta X_2$, where $\delta X_1$ and $\delta X_2$ represent the isotopic composition of two different water sources and $\delta X_{mix}$ is the isotopic composition of the resulting mixed water body with X denoting $^{17}\text{O}$, $^{18}\text{O}$, or $^{2}\text{H}$, respectively. Mixing curves are inversely shaped to the evaporation trajectories in triple oxygen isotope space (Fig. 3). At the Salar del Huasco mixing may occur episodically due to flooding and fluctuations in the groundwater level as detectable in satellite images (Fig. 2). Such mixing processes are likely transient due to the rarity of flooding events. In our model approach, we assume that the isotopic composition of the admixed water is similar to that of spring water, which approximates the isotopic composition of groundwater reasonably well (Uribe et al., 2015).

Simple evaporation, recharge evaporation, and mixing are principally well resolvable in triple oxygen isotope space, whereas in $\delta^2\text{H}-\delta^{18}\text{O}$ space, they tend to merge within data uncertainty (Fig. 3). All three trajectories can be expected in a dynamic arid hydrological setting such as the Salar del Huasco. Water affected exclusively by evaporation must progress along either of the two principal evaporation trajectories defined by the C–G model, depending on the recharge conditions. On the other hand, episodic flooding events initiate mixing leading to transient deviations from these general evaporation trends. The relative isotopic difference between inflowing water ($\delta_{WI}$) and atmospheric vapour ($\delta_{V}$) mainly determines the resolution of evaporation trajectories for different relative humidity in the triple oxygen isotope plot (Surma et al., 2018). For the given boundary conditions at the Salar del Huasco, evaporation trajectories show low sensitivity to changes in relative humidity, temperature, and

the turbulence coefficient in the diagram of $^{17}\text{O}$-excess over $\delta^{18}\text{O}$ (Fig. S3). Thus, seasonal and diurnal variability in these climate variables should have a low impact on the triple oxygen isotope composition of lakes and ponds at the Salar del Huasco, which is favourable for the purpose of this study to resolve different hydrological processes of evaporation, recharge, and mixing.

### 4.3 The isotopic composition of atmospheric vapour

A few atmospheric vapour samples were extracted on-site using a battery-powered Stirling cooler (Peters and Yakir, 2010). Due to the low relative humidity at the study site, extractions lasted 2–3 h, and only three samples could be taken.

In the absence of sufficiently abundant direct measurements, the isotopic composition of atmospheric vapour may be inferred from precipitation data, e.g. derived from the Online Isotopes in Precipitation Calculator (OIPC) (Bowen, 2020), assuming isotopic equilibrium ($\delta X_{V} = (\delta X_{P} + 1000)/\alpha_{eq, l-v} - 1000$, where $\delta X_{P}$ refers to the isotopic composition of precipitation and X denotes $^{17}\text{O}$, $^{18}\text{O}$, or $^{2}\text{H}$, respectively). However, the applicability of the equilibrium assumption is questionable for the Salar del Huasco, where precipitation mostly occurs in rare thunderstorm events originating from easterly sources. More importantly, rain samples cannot be representative of the local vapour because throughout the majority of the year air masses originate from westerly Pacific sources, which generally do not lead to precipitation (Garreaud et al., 2003). To confirm the predominant western origin of atmospheric vapour for the time of this study, the HYSPLIT Lagrangian model (Stein et al., 2015) was used to calculate 7-day air mass back-trajectories in daily resolution (12:00) for the month prior to each of the sampling campaigns for ground level (3800 m above sea level, m a.s.l.) and 1500 m above ground (5300 m a.s.l.) (cf. Aravena et al., 1999). Additionally, air mass back-trajectories were modelled in hourly resolution for the time of vapour sampling to confirm the Pacific origin of sampled water vapour. Because direct vapour extraction in the arid environment requires several hours of attendance, we were only able to extract a limited number of point-in-time samples. To augment the scarce database, we verified our $\delta^{18}\text{O}_{V}$ indirectly from the on-site evaporation experiments.

### 4.4 Isotope turnover time of ponds

The isotopic composition of a pond represents an integrated signal over the turnover time of the isotopes, which needs to be estimated to accurately apply the C–G model. This isotope turnover time is specific for each pond and depends on its surface-to-volume ratio, recharge, and the evaporation rate. Only very small water volumes may capture diurnal variations as shown in pan evaporation experiments (Surma et al., 2018). In contrast, deeper water bodies integrate over longer time intervals. Diurnal cycles and changing conditions over
<figure><img src="image_7.png" /></figure>

**Figure 4.** Oxygen and hydrogen isotope data of the Collacagua river, springs, lakes, and ponds sampled in the Salar del Huasco basin during field campaigns in September 2017 (**a**, **b**), September 2018 (**c**, **d**), and March 2019 (**e**, **f**). Colour coding refers to different hydrological subsystems (see legend and cf. Fig. 2). Note that the symbol size can be larger than the error bars. Trajectories for simple evaporation (black) and evaporation with recharge (green) were modelled using input parameters as summarized in Table 3. To calculate temperature and relative humidity, daytime (06:00–18:00 CLT) conditions of the 10 d (March 2019) or 20 d (September 2017, September 2018) prior to sampling were considered (cf. Sect. 6.3). Additionally, admixture of fresh groundwater (dashed grey line) as it might occur during flooding or snowmelt is exemplified for the case of a terminal lake ($E/I=1$). The global meteoric water line (GMWL) serves as reference. Note that the apparent lack of samples on a recharge evaporation trajectory in September 2018 and March 2019 is simply the result of a sampling focus on sources and on ponds cut off from recharge during these campaigns.
several days or weeks will be smoothed out, particularly if local conditions, e.g. a rough wind regime, lead to a well-mixed lake. The isotope turnover time may be approximated considering the turnover rate of ponds. For a terminal lake, the turnover rate (TR) equals the vaporization rate ($\varphi_{vap}$), which is a function of the potential evaporation ($E_{pot}$) and relative humidity (h): $\text{TR}_{E/I=1} = \varphi_{vap} = E_{pot} / (1-h)$ (Gonfiantini et al., 2018). Mean annual potential evaporation at the Salar del Huasco of 2290 mm (DGA, 1987) together with the mean annual relative humidity of 40 % results in a mean turnover rate of 318 mm per month (11 mm per day). As isotopic steady-state conditions are reached asymptotically, the isotope turnover time is better approximated by the isotope "half-turnover time": $t_{1/2} = \ln(2)/(d \cdot \text{TR})$, with d denoting the depth of the water body. Ponds in the Salar del Huasco are very shallow with depths up to 30 cm, which results in an isotope half-turnover time of <20 d. Considering seasonal variability of the potential evaporation rates and relative humidity (150 mm per month and 33 % in austral winter, 200 mm per month and 54 % in austral summer; DGA, 1987; CEAZA, 2020), the isotope half-turnover time varies between 14 and 26 d for a 30 cm deep pond. To account for the isotope half-turnover time described above, evaporation trajectories were modelled with temperature and relative humidity values averaged over 10 or 20 d prior to sampling in March and September, respectively (Fig. 4). Diurnal variations in the evaporation rate were accounted for using daytime (06:00–18:00 CLT) values.

# 5 Results

## 5.1 Natural waters in the Salar del Huasco basin

Springs sampled in September 2017 comprise average values of $-12.5 \pm 0.6\%_o$ for $\delta^{18}$O, $2 \pm 3\%_o$ for d-excess, and $11 \pm 7$ per meg for $^{17}$O-excess, which are in good agreement with published isotopic data of springs and wells in the Salar del Huasco basin ($-12.56 \pm 1.36\%_o$ in $\delta^{18}$O and $3.2 \pm 5.7\%_o$ in d-excess; data from Fritz et al., 1981; Uribe et al., 2015; Jayne et al., 2016). The springs' isotopic composition shows slight intra- and interannual variability comprising average $\delta^{18}$O, d-excess and $^{17}$O-excess values of $-12.3 \pm 0.5\%_o$, $1 \pm 2\%_o$, and $6 \pm 7$ per meg in September 2018 and $-12.6 \pm 0.7\%_o$, $6 \pm 2\%_o$, and $2 \pm 6$ per meg in March 2019, respectively. The Collacagua river, which was sampled at about 15 km distance from the salar, revealed similar isotopic composition with $\delta^{18}$O, d-excess, and $^{17}$O-excess values of $-12.4\%_o$, $2\%_o$, and $5$ per meg in September 2018 and $-12.5\%_o$, $6\%_o$, and $5$ per meg in March 2019, respectively.

The ponds in the Salar del Huasco show typical evaporation trends of increasing $\delta^{18}$O values with decreasing d-excess and $^{17}$O-excess values ranging from $-11.2\%_o$ to $14.5\%_o$ in $\delta^{18}$O, from $-1\%_o$ to $-100\%_o$ in d-excess, and from 19 to $-108$ per meg in $^{17}$O-excess (Fig. 4). Ponds from

the north-western and south-eastern subsystem are generally less enriched in $\delta^{18}$O compared to ponds from the northern subsystem. Ponds in the south-western area span a wide range in isotopic composition from $-5.5\%_o$ to approximately $14.5\%_o$, including the highest $\delta^{18}$O values observed in our sample set. The permanent Laguna Grande has a composition in the intermediate range of all investigated ponds showing considerable seasonal and interannual variability between $-4.2\%_o$ and $3.2\%_o$ in $\delta^{18}$O, $-2\%_o$ and $-42\%_o$ in d-excess, and $-15$ and $-39$ per meg in $^{17}$O-excess.

## 5.2 Atmospheric vapour

Two atmospheric vapour samples taken during the field campaign in September 2017 yield on average $-19.4 \pm 2.3\%_o$ for $\delta^{18}$O$_V$, $30 \pm 3\%_o$ for d-excess$_V$, and $18 \pm 2$ per meg for $^{17}$O-excess$_V$. An additional vapour sample taken in March 2019 comprises a slightly more depleted $\delta^{18}$O$_V$ value of $-24.3\%_o$ and d-excess$_V$ and $^{17}$O-excess$_V$ values of $41\%_o$ and 33 per meg, respectively. The overall average of all three measurements is $-21.0 \pm 3.3\%_o$, $34 \pm 6\%_o$, and $23 \pm 9$ per meg for $\delta^{18}$O$_V$, d-excess$_V$, and $^{17}$O-excess$_V$, respectively.

Seven-day air mass back-trajectories calculated using the HYSPLIT model suggest that air masses during vapour sampling were mainly derived from Pacific sources (Fig. 5). A predominantly western origin of air masses for the month prior to sampling is confirmed by 7d air mass back-trajectories (Fig. 5). These model results support the assumption that direct vapour measurements reflect the mean annual isotopic composition of atmospheric vapour at the Salar del Huasco.

The measured mean $\delta^{18}$O$_V$ value of $-21.0 \pm 3.3\%_o$ is comparable to the mean annual $\delta^{18}$O$_V$ value of $-21.8\%_o$ calculated based on monthly precipitation data derived from the OIPC database (Bowen et al., 2005). The apparent similarity between the OIPC model and our measurements may be coincident resulting from the seasonality of precipitation sources with a major contribution of the depleted Amazon moisture source and a minor contribution of the relatively enriched winter snow moisture source (cf. Aravena et al., 1999).

## 5.3 Pan evaporation experiments and the determination of the local turbulence coefficient

Water samples from the evaporation experiments comprise increasing $\delta^{18}$O and decreasing $^{17}$O-excess and d-excess values with increasing degree of evaporation (Fig. 6). The observed trends are consistent between experiments carried out with different initial volume. The isotopic composition of the initial water was $-11.0\%_o$, $4\%_o$, and 16 per meg, for $\delta^{18}$O, d-excess, and $^{17}$O-excess, respectively. Over the 3d experimental period, the largest enrichment in $\delta^{18}$O ($\Delta = 31.4\%_o$) and depletion in d-excess ($\Delta = -120\%_o$) and $^{17}$O-excess ($\Delta = -155$ per meg) is observed for the evaporation pan
<figure><img src="image_8.png" /><figcaption>**Figure 5.** HYSPLIT 7 d air mass back trajectories modelled for ground level at the Salar del Huasco (~3800 m above sea level, m a.s.l.) (black) and 1500 m above ground level (5300 m a.s.l.) (red) in daily resolution for the period from 23 August to 22 September 2017. The thick red and black trajectories represent modelled back trajectories for the time of vapour sampling on 20 September 2017, 21 September 2017, and on 20 March 2020.</figcaption></figure>

<table><caption>**Table 2.** Model input parameters of the simple evaporation model for pan evaporation experiments. The isotopic composition of initial water and atmospheric vapour was determined by measurement. The effective evaporation time interval was assumed to correspond to the period when T > 0 °C to account for the freezing of water in the pans during the night.</caption><thead><tr><td>Parameter</td><td>Value</td></tr></thead><tbody><tr><td>δ<sup>18</sup>O<sub>WI</sub> (‰)</td><td>-11.0</td></tr><tr><td>d-excess<sub>WI</sub> (‰)</td><td>4</td></tr><tr><td><sup>17</sup>O-excess<sub>WI</sub> (per meg)</td><td>16</td></tr><tr><td>δ<sup>18</sup>O<sub>V</sub> (‰)</td><td>-19.4</td></tr><tr><td>d-excess<sub>V</sub> (‰)</td><td>30</td></tr><tr><td><sup>17</sup>O-excess<sub>V</sub> (per meg)</td><td>18</td></tr><tr><td>Temperature (°C)</td><td>10</td></tr><tr><td>Relative humidity (%)</td><td>22</td></tr></tbody></table>

with the lowest initial volume (600 mL), which reflects the maximum fractional loss by evaporation relative to the initial volume.

The pan evaporation experiment data were used to empirically determine the turbulence coefficient *n*. The turbulence coefficient may be derived by the best fit of the C–G evaporation trajectory through all experimental data in plots of δ¹⁸O and δ²H over the fraction of remaining water. Alternatively, the turbulence coefficient may be determined by fitting the evaporation trajectory through the experimental data in the plot of d-excess over the residual fraction. The latter method is advantageous, because here the evaporation trajectory is predominantly controlled by the magnitude of the turbulence coefficient (Fig. 7a) and only barely sensitive to the other

<table><caption>**Table 3.** Model input parameters of the simple and the recharge evaporation model for natural ponds and lakes from the Salar del Huasco. The isotopic composition of inflowing water (δ<sub>WI</sub>) was inferred from the spring with the most depleted δ¹⁸O value. The isotopic composition of atmospheric vapour (δ<sub>V</sub>) was estimated by direct measurement. Evaporation trajectories were modelled with temperature and relative humidity values averaged over 10 d (March 2019) or 20 d (September 2017, September 2018) prior to sampling (Fig. 4), corresponding to the mean seasonal isotopic “turnover time” (cf. Sect. 6.1). Diurnal variations in the evaporation rate were accounted for using daytime (06:00–18:00 CLT) values. The value of the turbulence coefficient was derived from pan evaporation experiment data (cf. Sect. 5.3).</caption><thead><tr><td>Parameter</td><td>2017</td><td>2018</td><td>2019</td></tr></thead><tbody><tr><td>δ<sup>18</sup>O<sub>WI</sub> (‰)</td><td>-13.3</td><td>-13.0</td><td>-13.2</td></tr><tr><td>d-excess<sub>WI</sub> (‰)</td><td>6</td><td>5</td><td>6</td></tr><tr><td><sup>17</sup>O-excess<sub>WI</sub> (per meg)</td><td>15</td><td>7</td><td>9</td></tr><tr><td>δ<sup>18</sup>O<sub>V</sub> (‰)</td><td>-21.0</td><td>-21.0</td><td>-21.0</td></tr><tr><td>d-excess<sub>V</sub> (‰)</td><td>34</td><td>34</td><td>34</td></tr><tr><td><sup>17</sup>O-excess<sub>V</sub> (per meg)</td><td>23</td><td>23</td><td>23</td></tr><tr><td>Temperature (°C)</td><td>7</td><td>6</td><td>12</td></tr><tr><td>Relative humidity (%)</td><td>27</td><td>23</td><td>36</td></tr><tr><td>Turbulence coefficient</td><td>0.55</td><td>0.55</td><td>0.55</td></tr></tbody></table>

variables of the C–G equation (Fig. S3). Importantly, this approach is insensitive to δ¹⁸O<sub>V</sub>, the variable that is least constrained in our study (Fig. 7b).

Using model input parameters summarized in Table 2, fitting the evaporation model through all experimental data in the plot of d-excess over the fraction of remaining water implies *n* = 0.44±0.04. This value falls within error in the low-
<figure><img src="image_9.png" /><figcaption>**Figure 6.** Isotopic data of pan evaporation experiments in the diagrams of (a) $^{17}$O-excess over $\delta^{18}$O and (b) d-excess over $\delta^{18}$O. Pan evaporation experiments were carried out with different initial volume of 600 mL (diamonds), 800 mL (circles), and 1000 mL (squares). Note that the symbol size can be larger than the error bars. The solid black line in panel (b) shows the best fitting evaporation trajectory ($n = 0.55$) for the simple evaporation model. The red lines represent the evaporation trajectories modelled using the virtual outlet model with $n = 0.59$ and an outlet of 20% (see text). Dashed lines illustrate the evaporation trajectories modelled with the simple evaporation model using the same turbulence coefficient as for the virtual outlet model. Model input parameters are summarized in Table 2.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>**Figure 7.** Sensitivity of evaporation trajectories in a plot of d-excess over residual fraction to (a) the turbulence coefficient ($n$) and (b) the isotopic composition of atmospheric vapour ($\delta^{18}O_V$). Isotopic data of pan evaporation experiments with initial volume of 600 mL (diamonds), 800 mL (circles), and 1000 mL (squares) are shown for comparison. Note that the symbol size is in most cases larger than the error bars. The solid line represents the modelled evaporation trajectory for simple evaporation using input parameters as summarized in Table 2 and a turbulence coefficient of 0.44. Dashed lines show model results for (a) different turbulence coefficients in intervals of 0.1 and (b) varying $\delta^{18}O_V$ in intervals of 5‰, keeping other parameters constant.</figcaption></figure>

ermost range of turbulence coefficients ≥ 0.5 that are typi-
cally observed under natural conditions worldwide (Merlivat
and Jouzel, 1979; Gonfiantini, 1986; Mathieu and Bariac,
1996; Surma et al., 2018; Gázquez et al., 2018) and appar-
ently reflects excessive evaporation during midday at higher
than average turbulence caused by prevailing strong winds
with gusts up to 24 m s⁻¹. However, fitting the evaporation
trajectory independently to δ¹⁸O and δ²H data results in a
discrepancy between the individually derived turbulence co-
efficients. The δ¹⁸O and δ²H data of our evaporation exper-
iment may only fit the C–G model at an unrealistically low
value of n = 0.34 for δ¹⁸O, while no fit can be achieved for
any value of 0 ≤ n ≤ 1 in the diagram of δ²H vs. the residual
fraction (Fig. 8). Similar discrepancies were observed in lab-

oratory experiments at wind speeds above 0.5 m s⁻¹ carried
out by Gonfiantini et al. (2020).

In this particular experiment, the freezing at night-time
might have biased the relationship between d-excess and the
fraction of remaining water. When the ice begins to melt in
the morning, a considerable fraction, if not the whole result-
ing water film on top of the pan's frozen surface layer, may
evaporate in isolation from the bulk of water underneath the
ice. The fraction of total water in the pan would thus be re-
duced without affecting the isotopic composition of the bulk
ice. Sublimation of ice at night would add to this effect. A
similar effect would be introduced by microdroplet vaporiza-
tion without isotope fractionation during daytime at strong
winds as suggested by Gonfiantini et al. (2020). Applying
<figure><img src="image_11.png" /></figure>

**Figure 8.** Evolution of (a) $\delta^2$H, (b) $\delta^{18}$O, and (c) d-excess of remaining water in the evaporation experiments. The solid lines show evaporation trajectories obtained by fitting the simple evaporation model independently to $\delta^{18}$O ($n = 0.34$), $\delta^2$H ($n = 0$), and d-excess ($n = 0.44$) vs. residual fraction. Note that the $\delta^2$H data cannot be fitted for any reasonable turbulence coefficient. The dashed lines in (a) and (b) are evaporation trajectories modelled using the same turbulence coefficient as for the d-excess data.

a “virtual outlet” introduced by these authors – i.e. a fraction of water that is lost without fractionation – in the isotope mass balance equation yields an excellent fit in all three plots of d-excess, $\delta^{18}$O, and $\delta^2$H vs. residual fraction for $n = 0.59 \pm 0.06$ with 20±4% loss of water without fractionation (Fig. 9). The hypothesized loss of water without isotope fractionation – either during thawing, due to sublimation or microdroplet vaporization – would have a negligible effect on isotopic data in a diagram of d-excess over $\delta^{18}$O (Fig. 6b). At given boundary conditions, the best fit for this purely isotopic correlation is obtained for a value of $n = 0.55 \pm 0.09$,

<figure><img src="image_12.png" /></figure>

**Figure 9.** Isotopic evolution of remaining water in evaporation experiments, similar to Fig. 8. Here, the red lines show evaporation trajectories obtained by fitting the virtual outlet model ($n = 0.59$, outlet = 20%). The dashed lines represent evaporation trajectories for the simple evaporation model obtained using the same turbulence coefficient as for the virtual outlet model.

which is within error identical to the turbulence coefficient $n = 0.59 \pm 0.06$ derived from the virtual outlet model.

The disadvantage of the latter model is the requirement of a $\delta^{18}$O$_V$ value, which in this case is only poorly constrained by a few measurements and in turn translates into a higher uncertainty. For the above calculation, $\delta^{18}$O$_V = -19.4 \pm 3‰$ was used, as derived from our own vapour measurements carried out during the period of the experiment. The obtained value of $n = 0.55 \pm 0.09$ is in good agreement with the global range of reported turbulence coefficients (Merlivat and Jouzel, 1979; Gonfiantini, 1986; Mathieu and Bariac, 1996; Surma et al., 2018) and thus indirectly supports the accuracy of the measured $\delta^{18}$O$_V$ value.
In a diagram of ¹⁷O-excess vs. δ¹⁸O, the isotopic data fall below the predicted evaporation trend, regardless of which model or turbulence coefficient is used (Fig. 6a). This mismatch may be caused by two effects: (1) partial melting during the thawing of our experiment in the morning may result in uncontrollable mixing effects that are only observable in a ¹⁷O-excess over δ¹⁸O diagram (Fig. S5), and (2) diurnal variations of temperature and relative humidity generate corresponding changes of the theoretical isotopic end point of evaporation leading to a diurnal evolution of the simple evaporation trajectory (Fig. S6; Surma et al., 2018). Both effects become increasingly important, as the experiment progresses to smaller residual water volumes, consistent with our data. High-resolution sampling of similar evaporation experiments may aid in resolving these two effects in the future.

In conclusion, regarding this experiment, we suggest that the complementary analysis of all water isotopes in an evaporation experiment principally allows the determination of the turbulence coefficient. The effect of freezing and thawing on the total water loss was not fully resolvable from our experimental data. Hence, we used the turbulence coefficient of *n* = 0.55 for the following modelling calculations, which was derived from the isotopic correlation between d-excess and δ¹⁸O, independently from the fraction of remaining water.

# 6 Discussion: isotopic hydrology in the Salar del Huasco basin

## 6.1 Mixing in tributaries and groundwater aquifers

The Collacagua river and its tributaries originate from springs at different altitudes in the Salar del Huasco basin (cf. Fig. 1) and could be expected to comprise a somewhat lower δ²H and δ¹⁸O composition than evaporated salar water due to the altitude effect (Uribe et al., 2015). Samples from the Collacagua river fall below the local meteoric water line (LMWL)/GMWL (Fig. 10), reflecting a significant evaporation history. Likewise, the δ¹⁸O values of springs in the Salar del Huasco basin are slightly higher than local precipitation (δ¹⁸O = -17‰ to -13‰; Fritz et al., 1981; Aravena et al., 1999; Uribe et al., 2015; Scheihing et al., 2018); the d-excess and ¹⁷O-excess generally fall below the LMWL/GMWL (Fig. 10). This is a common observation for groundwater in the Atacama Desert and other desert environments (Aravena, 1995; Surma et al., 2015, 2018). Based on diagrams of d-excess over δ¹⁸O, this offset has been attributed to evaporation of precipitation during infiltration into the soil in previous studies (Aravena, 1995; Fig. 10b). However, all spring and river samples fall below any reasonable evaporation trend in triple oxygen isotope space (Fig. 10a). As such, evaporation along the river path or within the aquifer cannot be solely responsible for the observed isotopic composition. Additional modification of river and groundwa-

ter by mixing of different sources must be invoked to explain the triple oxygen isotope composition. Mixing most likely occurs during infiltration into the soil between precipitation and older, evaporated connate water in the vadose zone. Mixing may also take place with water adsorbed onto salts, e.g. halite (NaCl), or structurally bonded water of minerals, e.g. mirabilite (Na₂SO₄·10H₂O), which we observed frequently during sampling in the Salar del Huasco environment.

## 6.2 Hydrological processes in the salar

The isotopic composition of lakes and ponds from the Salar del Huasco sampled in the period from September 2017 to March 2019 principally reflects the evaporation trend predicted by the C–G model for given boundary conditions (Table 3; Fig. 4). Most of the ponds fall on the recharge evaporation trajectory indicating that evaporation occurs while recharge takes place by surface inflow from springs, streams, and creeks as well as subsurface inflow by groundwater. Ponds with high δ¹⁸O values sampled in the southern/southern area of the salar in September 2017 and September 2018 fall within the envelope spanned by the trajectories for recharge evaporation and simple evaporation at the same boundary conditions (Fig. 4a, c). These ponds may have been isolated from recharge a short time before sampling due to the general lowering of the water table during the dry season and evolved from the recharge evaporation trajectory towards the simple evaporation trajectory without a general change in climatic boundary conditions.

Other seemingly isolated ponds from the overflow region connecting the Laguna Jasure with the northern subsystem, sampled in March 2019, clearly fall on the simple evaporation trajectory (Fig. 4e). These ponds were probably remnants of a major flooding in the previous austral summer and cut off from recharge for long enough prior to sampling to evolve all the way from the recharge evaporation trajectory or a mixing curve (see below) to the simple evaporation trajectory. In general, both evaporation trends – recharge evaporation and simple evaporation – are well resolved in the plot of ¹⁷O-excess over δ¹⁸O (Fig. 4a, c, e) but are indistinguishable in a diagram of d-excess over δ¹⁸O (Fig. 4b, d, f).

Ponds sampled in September 2017 in the northern and eastern part of the salar fall generally below the predicted recharge evaporation trajectory but instead within the envelope for mixing defined by the inflow and terminal lake endmembers (Fig. 4a). Thus, in addition to evaporation and recharge, mixing of younger, less evaporated floodwater and older, more evaporated lake water is at least of temporary importance. At the Salar del Huasco, the seasonality of precipitation and the impact of evaporation cause fluctuations in the groundwater table and occasional flooding. An increase in the groundwater table after heavy rainfalls or snowmelt may lead to admixture of fresh water to the pre-evaporated shallow subsurface flow or pond water. The episodic nature of precipitation should result in rather transient mixing events
<figure><img src="image_13.png" /></figure>

**Figure 10.** Oxygen (a) and hydrogen (b) isotopic data of the Collacagua river and springs from the Salar del Huasco basin. Colour coding refers to different hydrological subsystems (cf. Fig. 2). The isotopic range of precipitation (−17‰ to −13‰) is derived from previously published data (Scheihing et al., 2018, and references therein). The local meteoric water line (LMWL: δ²H = 7.93 · δ¹⁸O + 12.3; Boschetti et al., 2019) was derived from precipitation data of northern Chile. Due to the lack of a comprehensive dataset of ¹⁷O in precipitation, the global meteoric water line (GMWL) was used as reference in the triple oxygen isotope plot (a). The shaded grey area indicates the isotopic evolution of precipitation undergoing evaporation. The Collacagua river and all local springs fall below the common evaporation trend in a plot of ¹⁷O-excess over δ¹⁸O (a). This offset can be explained considering mixing of precipitation with evaporated water in the vadose zone during infiltration into the soil. The effect of mixing is schematically illustrated by the dashed grey lines.

at the Salar del Huasco. This episodic mixing model differs from the continuous mixing hypothesis proposed by Herwartz et al. (2017) at the Salar de Llamará in the Central Depression of the Atacama Desert, where continuous admixture of isotopically light groundwater to a pre-evaporated subsurface flow may lead to systematic variations in the isotopic composition of inflowing water along the flow path.

## 6.3 Impact of climatic dynamics

Considerable seasonal and diurnal variability in environmental conditions at the Salar del Huasco impose variations in climatic boundary conditions during the isotope turnover time of a few weeks to a few months. Moreover, the observed mixing of different generations of precipitation during infiltration and possibly longer residence time of water within the aquifers feeding the spring sources of the lakes and ponds in the Salar del Huasco should result in some variability of the starting point of an evaporation trajectory. Consequently, triple oxygen isotope evaporation trajectories of the C–G model are subject to shifts over time. To evaluate the impact of such variability on the C–G model’s resolving power of hydrological processes and detection of mixing from model mismatch in particular, different scenarios were simulated within the range of observed source water and climate variability. The following simulations were focussed on our main field campaign in September 2017, where mixing was identified to have affected the isotopic composition of ponds.

The effect of source variability is most visible in a broader range of evaporation trajectories in the diagram of ¹⁷O-

excess over δ¹⁸O for throughflow ponds and lakes with E/I < 0.5 (Fig. 11a). At the Salar del Huasco, source variability may account for the offset of many of the ponds, but a number of those with seemingly high E/I ratios fall well below the simulated envelope.

The impact of variability in the isotopic composition of ambient vapour (δ¹⁸O~V~) that could potentially affect evaporation trajectories in the C–G model was estimated using the standard deviation (1SD = 3.3‰) of the three direct measurements distributed over the time interval of this study (Fig. 11c, d). Again, many of the ponds are plotted inside the predicted evaporation envelope defined by variability of δ¹⁸O~V~ except for those ponds with high E/I that again fall outside the range. To explain these ponds with an evaporation trajectory shifted by variability in the ambient vapour composition, δ¹⁸O~V~ values of about −30‰ or even lower would be necessary. This seems to be unlikely considering our measurements (−21.0±3.3‰). Most importantly, however, the combined ²H/¹H, ¹⁷O/¹⁶O, and ¹⁸O/¹⁶O data cannot be reproduced by invoking a very low δ¹⁸O~V~ value. With an extreme δ¹⁸O~V~ of −30‰, the evaporation trajectory may be forced through those ponds below the envelope in ¹⁷O-excess, but that would result in an evaporation trajectory falling well above the majority of ponds in the diagram of d-excess over δ¹⁸O (Fig. 11c, d). The model outcome is very well constrained in terms of variability in δ¹⁸O~V~ because trends move in opposite directions between the two plots.

Sensitivity tests suggest a low impact of seasonal and diurnal variability in relative humidity, temperature, and wind
<figure><img src="image_14.png" /><figcaption>**Figure 11.** Illustration of the impact of variability and uncertainty in the model input parameters on the evaporation trajectory exemplified for lakes and ponds sampled in September 2017. The green line shows the modelled recharge evaporation trajectory for mean daytime conditions (06:00–18:00 CLT) at the Salar del Huasco. The shaded green areas in (a) and (b) illustrate the uncertainty introduced to the modelled trajectory by spatial variability in the isotopic composition of inflowing water ($\delta_{WI}$). Panels (c) and (d) show the impact of variable $\delta^{18}O_V$ on the evaporation trajectory illustrated by multiples intervals of the SD (1 SD = 3.3‰; shaded areas). The dark solid lines represent corresponding evaporation trajectories in (c) and (d). Panels (e) and (f) illustrate the effect of uncertainty in the assumed effective evaporation time interval. Recharge evaporation trajectories have been modelled for (1) daily mean (00:00–24:00 CLT; $h = 33\%$, $T = 4$ °C) (dark green curve), (2) daily mean $h = 20\%$ and $T = 10$ °C weighted for the average diurnal distribution of temperature and wind speed (brown curve), and (3) daytime mean (06:00–18:00 CLT; $h = 27\%$, $T = 7$ °C) conditions (light green curve).</figcaption></figure>
conditions on the isotopic composition of lakes and ponds in the Salar del Huasco (Fig. S4). To verify this, we compared recharge evaporation trajectories for daily mean (00:00–24:00 CLT; *h* = 33%, *T* = 4 °C) and daytime mean conditions (06:00–18:00 CLT; *h* = 27%, *T* = 7%). Additionally, we modelled the recharge evaporation trajectory for daily mean relative humidity and temperature (*h* = 20%, *T* = 10 °C) weighted for the mean diurnal distribution of temperature and wind speed. High temperature and high wind speed amplify evaporation and were thus more strongly weighted. Evaporation trajectories of all three scenarios fall close to each other in triple oxygen isotope space (Fig. 11e), demonstrating that variability in relative humidity and temperature cannot account for the larger deviations of ponds from the recharge evaporation trajectory in the triple oxygen isotope plot that we attributed to either mixing or cut-off from recharge and evolution along a simple evaporation trajectory.

Collectively, these sensitivity analyses support our conclusion of mixing with flood water. This conclusion is further corroborated by the relationship between δ¹⁸O and salinity of the ponds (Fig. S8). Most of the ponds show a significantly higher salinity than expected from their degree of evaporation inferred from their isotopic composition. This may be partly caused by redissolution of previously precipitated salts. However, the fact that ponds particularly from the northern area fall on distinct lines indicates mixing of highly evaporated water with fresh flood water. Flooding in this case does not necessarily invoke inundation originating from the springs and inflow channels along the salar's margin but may as well be the result of a rising groundwater table, which will mostly affect the more central, low-lying ponds that are also more evaporated. In fact, all samples that indicate mixing show advanced evaporation, i.e. plot rather to the right in the diagram. These results demonstrate the capability of the ¹⁷O-excess parameter to resolve mixing processes. In contrast, the d-excess parameter cannot distinguish mixing from evaporation at variable climatic boundary conditions at the Salar del Huasco (Fig. 11f).

## 7 Conclusions

When applied to the triple oxygen isotope system, the classic Craig-Gordon isotope evaporation model reliably predicts the recharge evaporation trend of throughflow ponds within the hydrologically complex and seasonally dynamic groundwater-recharged Salar del Huasco, N Chile, at ambient climatic boundary conditions and thus, by inference, from probably most lacustrine systems in general. The different hydrologic processes of evaporation and mixing of waters with different evaporation histories are demonstrably resolved in the triple oxygen isotope data. Cut-off from recharge during drought and subsequent isotopic evolution on a simple pan evaporation trend is also observable in individual ponds. Likewise, flooding, i.e. the rapid mixing of

evaporated water from the lake with fresh flood water of a different isotopic composition, may be identified using triple oxygen isotope data. In the classic δ²H-δ¹⁸O system, the two types of evaporation trajectories as well as mixing curves overlap. Therefore, these processes are generally poorly resolvable with conventional δ²H-δ¹⁸O measurements and the principle variable defining the hydrologic balance of a lake – i.e. the evaporation-to-inflow ratio (E/I) – cannot be determined unambiguously by the C-G model using the classic δ²H-δ¹⁸O system alone.

A fundamental requirement for the application of the C-G model is an on-site estimate of the wind turbulence coefficient. Although the turbulence coefficient may affect the C-G model in a non-linear fashion at high wind speed (Gonfiantini et al., 2020), our on-site experiment suggests that a turbulence coefficient representing average ambient wind conditions may be estimated reasonably well by using the empirical relationship between d-excess and δ¹⁸O.

The C-G model – and thus, the critical E/I parameter – may thus be well constrained based on the analysis of all three isotope ratios – ²H/¹H, ¹⁷O/¹⁶O, and ¹⁸O/¹⁶O. This is due to the fact that evaporation trajectories in plots of d-excess and ¹⁷O-excess over δ¹⁸O are differently affected when input variables change. Triple oxygen isotope data are also capable of identifying mixing of precipitation with older, pre-evaporated connate water during infiltration in a groundwater recharge area. Combining the advantages of triple oxygen isotope data in terms of resolving fundamental hydrological processes – specifically mixing and the different types of evaporation with and without recharge – with the advantage of the classic δ²H-δ¹⁸O system regarding the estimation of the critical wind turbulence exponent, a high level of accuracy in the hydrologic balance estimate of a lake – or an aquifer – may be achieved.

*Data availability.* All data reported herein are provided in the supplement and are freely available at the Collaborative Research Centre 1211 database at https://doi.org/10.5880/CRC1211DB.36 (Voigt et al., 2021). It comprises the results of the ion concentration (ICP-OES) and isotopic analyses of water samples from ponds and lakes at the Salar del Huasco investigated in this study, as well as the isotopic data of the pan evaporation experiments carried out on-site.

*Supplement.* The supplement related to this article is available online at: https://doi.org/10.5194/hess-25-1211-2021-supplement.

*Author contributions.* MS, CV and DH designed the study. CV, MS and CD were responsible for field work. CV conducted data analysis. CV, MS and DH evaluated the data. CV, MS and DH wrote the manuscript with contributions from CD.
*Competing interests.* The authors declare that they have no conflict of interest.

*Acknowledgements.* We thank Franc Megušar and Le-Tehnika (Kranj, Slovenia) for producing the Stirling cooler cycle system. Further, we thank the two anonymous reviewers, whose constructive comments significantly improved the manuscript.

*Financial support.* This research has been supported by the German Research Foundation (DFG) (grant no. 268236062).

*Review statement.* This paper was edited by Genevieve Ali and reviewed by two anonymous referees.

## References

Acosta, O. and Custodio, E.: Impactos ambientales de las extracciones de agua subterránea en el Salar del Huasco (norte de Chile), Bol. Geol. Min., 119, 33–50, 2008.

Alexandre, A., Landais, A., Vallet-Coulomb, C., Piel, C., Devidal, S., Pauchet, S., Sonzogni, C., Couapel, M., Pasturel, M., Cornuault, P., Xin, J., Mazur, J.-C., Prié, F., Bentaleb, I., Webb, E., Chalié, F., and Roy, J.: The triple oxygen isotope composition of phytoliths as a proxy of continental atmospheric humidity: insights from climate chamber and climate transect calibrations, Biogeosciences, 15, 3223–3241, https://doi.org/10.5194/bg-15-3223-2018, 2018.

Alexandre, A., Webb, E., Landais, A., Piel, C., Devidal, S., Sonzogni, C., Couapel, M., Mazur, J.-C., Pierre, M., Prié, F., Vallet-Coulomb, C., Outrequin, C., and Roy, J.: Effects of leaf length and development stage on the triple oxygen isotope signature of grass leaf water and phytoliths: insights for a proxy of continental atmospheric humidity, Biogeosciences, 16, 4613–4625, https://doi.org/10.5194/bg-16-4613-2019, 2019.

Angert, A., Cappa, C. D., and DePaolo, D. J.: Kinetic $^{17}$O effects in the hydrologic cycle: Indirect evidence and implications, Geochim. Cosmochim. Ac., 68, 3487–3495, https://doi.org/10.1016/j.gca.2004.02.010, 2004.

Aravena, R.: Isotope hydrology and geochemistry of northern Chile groundwaters, Bull. Ins. Fr. Etudes Andin., 24, 495–503, 1995.

Aravena, R., Suzuki, O., Peña, H., Pollastri, A., Fuenzalida, H., and Grilli, A.: Isotopic composition and origin of the precipitation in Northern Chile, Appl. Geochem., 14, 411–422, https://doi.org/10.1016/S0883-2927(98)00067-5, 1999.

Barkan, E. and Luz, B.: High precision measurements of $^{17}$O/$^{16}$O and $^{18}$O/$^{16}$O ratios in H₂O, Rapid Commun. Mass Sp., 19, 3737–3742, https://doi.org/10.1002/rcm.2250, 2005.

Barkan, E. and Luz, B.: Diffusivity fractionations of H₂¹⁶O/H₂¹⁷O and H₂¹⁶O/H₂¹⁸O in air and their implications for isotope hydrology, Rapid Commun. Mass Sp., 21, 2999–3005, https://doi.org/10.1002/rcm.3180, 2007.

Boschetti, T., Cifuentes, J., Iacumin, P., and Selmo, E.: Local meteoric water line of northern Chile (18°S–30°S): An application of error-in-variables regression to the oxygen and hydrogen stable isotope ratio of precipitation, Water, 11, 1–16, https://doi.org/10.3390/w11040791, 2019.

Bowen, G. J.: The Online Isotopes in Precipitation Calculator, available at: http://www.waterisotopes.org, last access: 23 September 2020.

Bowen, G. J., Wassenaar, L. I., and Hobson, K. A.: Global application of stable hydrogen and oxygen isotopes to wildlife forensics, Oecologia, 143, 337–348, https://doi.org/10.1007/s00442-004-1813-y, 2005.

Cao, X. and Liu, Y.: Equilibrium mass-dependent fractionation relationships for triple oxygen isotopes, Geochim. Cosmochim. Ac., 75, 7435–7445, https://doi.org/10.1016/j.gca.2011.09.048, 2011.

CEAZA, (Centro de Estudios Avanzados en Zonas Áridas): Estación Salar de Huasco, available at: http://www.ceazamet.cl, last access: 23 September 2020.

Clark, I. and Fritz, P.: Environmental Isotopes in Hydrogeology, Lewis Publishers, New York, USA, 1997.

Craig, H. and Gordon, L. I.: Deuterium and oxygen 18 variations in the ocean and the marine atmosphere, in: Stable Isotopes in Oceanographic Studies and Paleotemperatures, edited by: Tongiorgi, E., Laboratorio di Geologia Nucleare, Pisa, Italy, 9–130, 1965.

Criss, R. E.: Principles of Stable Isotope Distribution, Oxford University Press, Oxford, UK, 1999.

DGA: Balance Hídrico de Chile, MINISTERIO DE OBRAS PÚBLICAS, Santiago, Chile, 1987.

Dongmann, G., Nürnberg, H. W., Förstel, H., and WGENER, K.: On the enrichment of H₂¹⁸O in the leaves of transpiring plants, Radiat. Environ. Bioph., 11, 41–52, https://doi.org/10.1007/BF01323099, 1974.

Evans, N. P., Bauska, T. K., Gázquez, F., Brenner, M., Curtis, J. H., and Hodell, D. A.: Quantification of drought during the collapse of the classic Maya civilization, Science, 361, 498–501, https://doi.org/10.1126/science.aas9871, 2018.

Fritz, P., Suzuki, O., Silva, C., and Salati, E.: Isotope hydrology of groundwaters in the Pampa del Tamarugal, Chile, J. Hydrol., 53, 161–184, https://doi.org/10.1016/0022-1694(81)90043-3, 1981.

Garreaud, R. D. and Aceituno, P.: Interannual Rainfall Variability over the South American Altiplano, J. Climate, 14, 2779–2789, https://doi.org/10.1175/1520-0442(2001)014<2779:IRVOTS>2.0.CO;2, 2001.

Garreaud, R. D., Vuille, M., and Clement, A. C.: The climate of the Altiplano: Observed current conditions and mechanisms of past changes, Palaeogeogr. Palaeocl., 194, 5–22, https://doi.org/10.1016/S0031-0182(03)00269-4, 2003.

Gat, J. R. and Bowser, C.: The heavy isotope enrichment of water in coupled evaporative systems, in: Stable Isotope Geochemistry: A Tribute to Samuel Epstein, edited by: Taylor, H. P., O’Neil, J. R., and Kaplan, I. R., The Geochemical Society, St Louis, MO, USA, 159–168, 1991.

Gázquez, F., Morellón, M., Bauska, T., Herwartz, D., Surma, J., Moreno, A., Staubwasser, M., Valero-Garcés, B., Delgado-Huertas, A., and Hodell, D. A.: Triple oxygen and hydrogen isotopes of gypsum hydration water for quantitative paleohumidity reconstruction, Earth Planet. Sci. Lett., 481, 177–188, https://doi.org/10.1016/j.epsl.2017.10.020, 2018.

Gonfiantini, R.: Environmental isotopes in lake studies, in: Handbook of Environmental Isotope Geochemistry, edited by: Fontes, A.
D. and Fritz, P., Elsevier Science, New York, NY, USA, 119–168, 1986.

Gonfiantini, R., Wassenaar, L. I., Araguas-Araguas, L. J., and Aggarwal, P. K.: A unified Craig-Gordon isotope model of stable hydrogen and oxygen isotope fractionation during fresh or salt-water evaporation, Geochim. Cosmochim. Ac., 235, 224–236, https://doi.org/10.1016/j.gca.2018.05.020, 2018.

Gonfiantini, R., Wassenaar, L. I., and Araguas-Araguas, L. J.: Stable isotope fractionations in the evaporation of water: The wind effect, Hydrol. Process., 34, 3596–3607, https://doi.org/10.1002/hyp.13804, 2020.

Herwartz, D., Surma, J., Voigt, C., Assonov, S., and Staubwasser, M.: Triple oxygen isotope systematics of structurally bonded water in gypsum, Geochim. Cosmochim. Ac., 209, 254–266, https://doi.org/10.1016/j.gca.2017.04.026, 2017.

Horita, J.: Stable isotope fractionation factors of water in hydrated saline mineral-brine systems, Earth Planet. Sci. Lett., 95, 173–179, https://doi.org/10.1016/0012-821X(89)90175-1, 1989.

Horita, J.: Saline waters, in: Isotopes in the Water Cycle: Present and Past Future of a Developing Science, edited by: Aggarwal, P. K., Gat, J. R., and Froehlich, K. F. O., Springer, Dordrecht, the Netherlands, 271–287, 2005.

Horita, J., Cole, D. R., and Wesolowski, D. J.: The activity-composition relationship of oxygen and hydrogen isotopes in aqueous salt solutions: II. Vapor-liquid water equilibration of mixed salt solutions from 50 to 100°C and geochemical implications, Geochim. Cosmochim. Ac., 57, 4703–4711, 1993.

Horita, J., Rozanski, K., Cohen, S.: Isotope effects in the evaporation of water: a status report of the Craig-Gordon model, Isot. Environ. Healt. S., 44, 23–49, https://doi.org/10.1080/10256010801887174, 2008.

Houston, J.: Variability of precipitation in the Atacama Desert: Its causes and hydrological impact, Int. J. Climatol., 26, 2181–2198, https://doi.org/10.1002/joc.1359, 2006.

Jayne, R. S., Pollyea, R. M., Dodd, J. P., Olson, E. J., and Swanson, S. K.: Spatial and temporal constraints on regional-scale groundwater flow in the Pampa del Tamarugal Basin, Atacama Desert, Chile, Hydrogeol. J., 24, 1921–1937, https://doi.org/10.1007/s10040-016-1454-3, 2016.

Landais, A., Barkan, E., Yakir, D., and Luz, B.: The triple isotopic composition of oxygen in leaf water, Geochim. Cosmochim. Ac., 70, 4105–4115, https://doi.org/10.1016/j.gca.2006.06.1545, 2006.

Landais, A., Barkan, E., and Luz, B.: Record of $\delta^{18}\text{O}$ and $^{17}\text{O}$-excess in ice from Vostok Antarctica during the last 150 000 years, Geophys. Res. Lett., 35, L02709, 2008.

Landais, A., Risi, C., Bony, S., Vimeux, F., Descroix, L., Falourd, S., and Bouygues, A.: Combined measurements of $^{17}\text{O}$excess and d-excess in African monsoon precipitation: Implications for evaluating convective parameterizations, Earth Planet. Sci. Lett., 298, 104–112, https://doi.org/10.1016/j.epsl.2010.07.033, 2010.

Luz, B. and Barkan, E.: Variations of $^{17}\text{O}/^{16}\text{O}$ and $^{18}\text{O}/^{16}\text{O}$ in meteoric waters, Geochim. Cosmochim. Ac., 74, 6276–6286, https://doi.org/10.1016/j.gca.2010.08.016, 2010.

Mathieu, R. and Bariac, T.: A numerical model for the simulation of stable isotope profiles in drying soils, J. Geophys. Res.-Atmos., 101, 12685–12696, https://doi.org/10.1029/96JD00223, 1996.

Merlivat, L. and Jouzel, J.: Global climatic interpretation of the deuterium-oxygen 18 relationship for precipitation, J. Geophys. Res., 84, 5029–5033, https://doi.org/10.1029/JC084iC08p05029, 1979.

Peters, L. I. and Yakir, D.: A rapid method for the sampling of atmospheric water vapour for isotopic analysis, Rapid Commun. Mass Sp., 24, 103–108, https://doi.org/10.1002/rcm.4359, 2010.

Risacher, F., Alonso, H., and Salazar, C.: The origin of brines and salts in Chilean salars: A hydrochemical review, Earth-Sci. Rev., 63, 249–293, https://doi.org/10.1016/S0012-8252(03)00037-0, 2003.

Scheihing, K., Moya, C., Struck, U., Lictevout, E., and Tröger, U.: Reassessing Hydrological Processes That Control Stable Isotope Tracers in Groundwater of the Atacama Desert (Northern Chile), Hydrology, 5, 3, https://doi.org/10.3390/hydrology5010003, 2018.

Schoenemann, S. W., Schauer, A. J., and Steig, E. J.: Measurement of SLAP2 and GISP $\delta^{17}\text{O}$ and proposed VSMOW-SLAP normalization for $\delta^{17}\text{O}$ and $^{17}\text{O}$excess, Rapid Commun. Mass Sp., 27, 582–590, https://doi.org/10.1002/rcm.6486, 2013.

Sofer, Z. and Gat, J. R.: Activities and concentrations of oxygen-18 in concentrated aqueous salt solutions: Analytical and geo-physical implications, Earth Planet. Sci. Lett., 15, 232–238, https://doi.org/10.1016/0012-821X(72)90168-9, 1972.

Sofer, Z. and Gat, J. R.: The isotope composition of evaporating brines: Effect of the isotopic activity ratio in saline solutions, Earth Planet. Sci. Lett., 26, 179–186, https://doi.org/10.1016/0012-821X(75)90085-0, 1975.

Stein, A. F., Draxler, R. R., Rolph, G. D., Stunder, B. J. B., Cohen, M. D., and Ngan, F.: NOAA's HYSPLIT atmospheric transport and dispersion modeling system, B. Am. Meteorol. Soc., 96, 2059–2077, https://doi.org/10.1175/BAMS-D-14-00110.1, 2015.

Surma, J., Assonov, S., Bolourchi, M. J., and Staubwasser, M.: Triple oxygen isotope signatures in evaporated water bodies from the Sistan Oasis, Iran, Geophys. Res. Lett., 42, 8456–8462, https://doi.org/10.1002/2015GL066475, 2015.

Surma, J., Assonov, S., Herwartz, D., Voigt, C., and Staubwasser, M.: The evolution of $^{17}\text{O}$-excess in surface water of the arid environment during recharge and evaporation, Sci. Rep.-UK, 8, 4972, https://doi.org/10.1038/s41598-018-23151-6, 2018.

Voigt, C., Herwartz, D., Dorador, C., and Staubwasser, M.: Isotopic and ion concentration data of ponds and pan evaporation experiments conducted at the Salar del Huasco, CRC1211 Database (CRC1211DB), https://doi.org/10.5880/CRC1211DB.36, 2021.

Uribe, J., Muñoz, J. F., Gironás, J., Oyarzún, R., Aguirre, E., and Aravena, R.: Assessing groundwater recharge in an Andean closed basin using isotopic characterization and a rainfall-runoff model: Salar del Huasco basin, Chile, Hydrogeol. J., 23, 1535–1551, https://doi.org/10.1007/s10040-015-1300-z, 2015.