The biogeochemical structuring role of horizontal stirring: Lagrangian perspectives on iron delivery downstream of the Kerguelen Plateau

. Field campaigns are instrumental in providing ground truth for understanding and modeling global ocean biogeochemical budgets. A survey however can only inspect a fraction of the global oceans, typically a region hundreds of kilometers wide for a temporal window of the order of (at most) several weeks. This spatiotemporal domain is also the one in which the mesoscale activity induces through horizontal stirring a strong variability in the biogeochemical trac-ers, with ephemeral, local contrasts which can easily mask the regional and seasonal gradients. Therefore, whenever local in situ measures are used to infer larger-scale budgets, one faces the challenge of identifying the mesoscale structuring effect, if not simply to ﬁlter it out. In the case of the KEOPS2 investigation of biogeochemical responses to natural iron fertilization, this problem was tackled by designing an adaptive sampling strategy based on regionally optimized multisatellite products analyzed in real time by speciﬁcally designed Lagrangian diagnostics. This strategy identiﬁed the different mesoscale and stirring structures present in the region and tracked the dynamical frontiers among them. It also enabled back trajectories for the ship-sampled stations to be estimated, providing important insights into the timing and pathways of iron supply, which were explored further using a model based on ﬁrst-order iron removal. This context was essential for the interpretation of the ﬁeld results. The mesoscale circulation-based strategy was also validated post-cruise by comparing the Lagrangian maps derived from satel-lites with the patterns of more than one hundred drifters, including some adaptively released during KEOPS2 and a subsequent research voyage. The KEOPS2 strategy was adapted to the speciﬁc biogeochemical characteristics of the region, but its principles are general and will be useful for future in situ biogeochemical surveys.


Introduction
The role of iron as a key limiting micronutrient for large phytoplankton in high nutrient-low chlorophyll (HNLC) waters was brought to prominence by Martin (1990) and motivated a series of bottle incubation experiments. Difficulties in unambiguously interpreting the results of these experiments led to the design and the implementation of more ambitious field studies. They were conducted in regions fertilized artificially or naturally with iron. One of the most striking difference between these types of study is the role played by the horizontal transport of iron. Most of the artificial iron fertilization experiments were conducted using specific strategies aimed at minimizing horizontal effects. Some artificial experiments have targeted quasi-isolated eddy cores, so that the water patches which are trapped by the mesoscale circulation are only marginally mixed with the environmental waters on the time scale of the induced bloom. This was the case for the EIFEX, LOHAFEX and SAGE experiments (Smetacek et al., 2012;Harvey et al., 2011;Martin et al., 2013) and the FeCycle natural iron response study (Boyd et al., 2005). Natural fertilization studies have attempted to target homogeneous regions with weak horizontal transport and mesoscale activity, with the hope that the system could be approximated as a one-dimensional water column reactor, with varying degrees of uncertainty for the conclusions (e.g., Blain et al., 2007;Pollard et al., 2009).
However, these quasi-isolated or homogeneous regions are not typical or representative of the vast majority of the ocean. In general, water parcels are stirred by the mesoscale field in a non-local way, experiencing the cumulative effects of several mesoscale structures from their iron enrichment event to the moment of the phytoplanktonic bloom. Horizontal transport does not only modulate the extension of a fertilized region with respect to its sources (e.g., Mongin et al., 2009), but may also indirectly affect nutrient concentrations through physical and chemical processes like mixing and scavenging, during the advection of a water parcel. This is evident for iron-fertilized waters in the open ocean, which present complex patterns and strong contrasts in satellite images of phytoplankton biomass, reflecting the pathways from the iron sources to the wider ocean, coupled with the complexities of biological responses. Thus, for the purposes of quantifying large-scale ecosystem responses to iron inputs, the importance of horizontal transport in redistributing iron-rich waters in complex patterns has become dramatically apparent. The KEOPS2 campaign aimed at addressing the nature of the response to iron fertilization in deep open-ocean waters downstream from the Kerguelen Plateau, hence including both the spatial dimension with respect to fertilization and its modulation by the mesoscale circulation. Together with traditional interdisciplinary stations with long occupation time, with some sites repeated to capture the temporal dynamics of the bloom, the campaign also needed to cover enough sites to document the spatial contrasts occurring in the region.
These requirements obviously put a strain on the allocation of ship time. The problem was exacerbated by the intrinsic spatiotemporal variability of the phytoplanktonic bloom: the time window allowed for characterizing the region could not be expanded by requesting more days for the campaign, but was upper-limited by the duration of the bloom itself. Moreover, a further requirement was that the different sites had to be occupied in the shortest possible time in order to disentangle space from time in the biogeochemical variability.
This paper first presents our efforts to use the prior and real-time satellite information to understand bloom dynamics, and thereby define an optimal sampling strategy. In doing so, we examine both Eulerian eddy fields and Lagrangian maps of water mass origins. The paper then validates these initial perspectives on the circulation against drifters released during the KEOPS2 project and a subsequent voyage. It provides an iron supply and removal budget based on the integration of satellite, ship-based and drifter information. It illustrates how the context provided by these circulation perspectives informed the interpretation of the observed ecological states (with reference to other works in this volume). Finally, it offers perspectives for biogeochemical process study planning.

Satellite products
Altimeter and ocean color products used in this study were specifically produced for the Kerguelen region by Ssalto/Duacs and CLS with support from CNES from mid-May 2011 to July 2012. Altimetry maps were generated from Jason-1, Jason-2 and Envisat along-track observations with procedures analogous to those used for AVISO data (SSALTO/DUACS User Handbook, 2014). Some of the parameters (i.e., mean dynamic topography, correlation scales and measurements errors) were specifically tuned for the region. These products are based on altimetric measurements at 10 day frequency, and interpolated to provide higher temporal frequency information.
In particular, a high resolution (1/8 • ) regional mean dynamic topography (MDT) was computed , applying the same two-step methodology used for the calculation of the global CNES-CLS09 MDT available at that time (Rio, 2012), but with an updated input data set and refined in situ measurements processing. First, a geoid model based on 1 year of GOCE data was used (Pail et al., 2011) together with the CNES-CLS11 Mean Sea Surface  to derive an MDT at a spatial scale larger than 125 km. This first guess was then improved by using synthetic mean geostrophic velocities and synthetic mean heights, calculated by subtracting the ocean temporal variability measured by altimetry from instantaneous in situ measurements of the ocean current velocities and dynamic heights. We used velocities  . Also, a 3 day low pass filter was applied to get rid of the other ageostrophic high-frequency currents (stokes drift, inertial oscillations, tidal currents). The data set of in situ dynamic heights was calculated using the temperature and salinity profiles from CTD (conductivitytemperature-depth) casts (Cora3.2) and Argo floats for the period 1993-2010, together with 75 CTD profiles measured from 15 October 2011 to 19 November 2011 during the KEOPS2 campaign.
Maps of absolute dynamic topography at a resolution of 1/8 • were produced daily in real time (RT-MADT) by CLS/CNES. Near-real-time (NRT-MADT) maps were produced after 6 days by further correcting the real-time products. Only NRT-MADT products are shown in this study.
Geostrophic currents were computed from the maps of absolute dynamic topography by solving the equation for geostrophic equilibrium (Pedlosky, 1987) with centered finite differences. Additionally, wind-induced surface Ekman currents were computed from a 1/8 • wind analysis from the European Center for Medium Range Weather Forecasting (ECMWF). These were added to the altimetry-based geostrophic currents to obtain regional maps of total velocity.
Composite maps of surface chlorophyll concentration with a resolution of 1/25 • were produced every 3 days, also by CLS/CNES. Each map was constructed by a 10 day weighted mean of MODIS and MERIS measurements. We hereafter refer to these as SCHL values and images.

Eulerian diagnostics
Maps of geostrophic and total velocities have been used to compute the Okubo-Weiss parameter, OW (Okubo, 1970;Weiss, 1991). OW is defined as where s sh and s st are the shear and strain deformations, and ξ is the vertical component of vorticity, respectively. Thus, OW quantifies the relative importance of deformation with respect to rotation. Negative values of OW can be used to quantify the strength of mesoscale eddies (e.g., d 'Ovidio et al., 2009;Chelton et al., 2011), since their velocity field is usually dominated by rotation.
In this study, mesoscale eddies have also been identified from satellite velocities using the vector geometry-based eddy detection algorithm developed by Nencioli et al. (2010). The method identifies the center of an eddy as a local minimum of velocity within a region characterized by rotating ve-locity vectors. Eddy boundaries are defined as the outermost closed streamlines across which velocity magnitude is still radially increasing. As in Yu Liu et al. (2012) and Amores et al. (2013), the satellite velocity fields have been linearly interpolated to 1/16 • resolution to improve the algorithm performance. The method has been applied with parameters a = 4 (defining the number of grid points across which the rotation of velocity vectors is inspected) and b = 3 (defining the dimensions in grid points of the area used to define the minimum of velocity).

Lagrangian and hybrid diagnostics
Several approaches were used to identify the effect of the eddy field on iron redistribution. Altimetry-based finite-size Lyapunov exponents (FSLE, Aurell et al., 1997) were computed with the method proposed by (d 'Ovidio et al., 2004) which can be outlined as in the following. In order to construct a map of Lyapunov exponents, at each grid point, trajectories are initialized and constructed by integrating in time the velocity field derived from altimetry. The finite-size Lyapunov exponent provides the maximal exponential rate of separation (averaged in time) among trajectories initialized nearby. Calling δ 0 the initial separation, δ a prescribed final the separation and τ the time at which the separation δ is reached, the Lyapunov exponent λ is defined as Trajectories were derived by applying a Runge-Kutta fourth-order scheme with a time step of 6 h. Velocity fields have been linearly interpolated in both space and time. Initial and final separation distances were set to 0.05 • (defining also the resolution of the Lyapunov map) and to 0.6 • , respectively. For the real-time analysis during the campaign, only backward-integrated FSLE were computed, providing a measure of the intensity of the convergent regions within the surface velocity field. The robustness of the altimetryderived Lyapunov exponent to small-scale errors and fluctuations in the velocity field has been studied in (Cotté et al., 2011), pages 224-225. In practice, the Lyapunov exponent computed backward in time provides the exponential rate at which horizontal stirring has brought water parcels to a close distance (δ 0 ) that were initially located a larger distance (δ) apart. Maxima (ridges) of Lyapunov exponents are used as a way to extract from altimetry data candidates for frontal regions, as a largescale tracer gradient can be amplified by this mechanism.
The same Runge-Kutta scheme used to compute the backward FSLE was applied to retrieve additional diagnostics fundamental to the real-time analysis. Firstly, the ability of mesoscale eddy cores to retain water parcels was quantified by backward advection of particles released within regions with negative OW values, and their duration within these regions was measured. In effect this yields an estimate of how F. d'Ovidio et al.: Lagrangian perspectives on iron delivery downstream of Kerguelen long ago a water parcel entered a given eddy, yielding a retention value expressed in days. Similarly, backward advection was used to identify water parcels originally above the Kerguelen Plateau, defined as the region where ocean depths are shallower than 700 m. This was achieved by finding the trajectories of virtual drifters which had touched the 700 m isobath in the past, in analogy to the method developed for the Crozet experiment (Sanial et al., 2014, where a validation of the method with lithogenic radioisotopes was also performed; see also Sanial et al., 2015 in this issue). Particle trajectories have also been used to estimate the origin of water parcels by horizontal transport, finding the position (longitude and latitude) of a particle 30 days before and mapping this information at its current position.

SVP drifters
More than 200 WOCE-SVP drifters were used to validate the altimetry-based estimation of stirring patterns. 48 drifters have been deployed according to an adaptive strategy during the KEOPS2 campaign and 24 during the MYCTO-3D campaign that took place in the Kerguelen region in January 2014. Furthermore, we collected 120 trajectories measured between 2011 and 2014 in the same region from the Global Drifter Program (GDP -http://www.aoml.noaa.gov/ phod/dac/index.php) historical archive. In order to compare the stirring patterns estimated from altimetry, we computed the Lagrangian diagnostics detailed in Sect. 2.3 using real drifters' trajectories instead of the altimetry-derived (i.e., simulated) ones. For each location along a drifter's trajectory, we computed how much time before that very drifter had crossed the 700 m bathymetry line ("age of the water parcel") and we recorded at which latitude it left the Kerguelen Plateau ("origin of the water parcel"). For drifters that did not touched the 700 m isobath but were close to it (within 100 km) we used their closest position. Fewer than 10 drifters in the study region had left the plateau more than 90 days earlier: this time exceeds the time extent of the backward computation from altimetry so we excluded those trajectories from our analysis. In order to perform a visual and quantitative pixel by pixel comparison with the maps computed from altimetry, we averaged the values along the trajectories on 0.25 • disks, obtaining climatological maps of the Lagrangian diagnostics.

Mesoscale activity and stirring patterns
The circulation of the Kerguelen region (see Fig. 1 for the bathymetry and a sketch of its current system) has been the object of many studies in recent years, and it is now known in great detail (see in particular Park et al. (2014) and the references therein). When attempting to identify the relation between mesoscale physics and patterns of primary produc- tion, standard Eulerian and Lagrangian diagnostics are however quite deceptive. Indeed, there is no obvious association between mesoscale activity and extension of the SCHL plume. Figures 2 and 3 show the situation at the beginning of the cruise while Fig. 4 displays the instantaneous and climatological extension of the bloom onset (respectively, 11 November 2011 and mean of November 2000-2010). The ocean color composite depicts a plume of blooming water extending from the Kerguelen-Heard shelf break (72-75 • E) eastward to 85-90 • E and contained latitudinally between 47 • : 53 • S. This biologically active region is not characterized by anomalous eddy activity: indeed, eddies populate a much larger region, with hotspots of total kinetic energy located equally inside the plume (e.g., 87 • E, 51 • S) or outside (e.g., 85 • E, 54 • S). Some of these eddies appear to be strongly retentive, but their location does not indicate the position of the plume. Even more puzzling, the fronts and local transport effect induced by the mesoscale activity has no obvious relation with the plume either. The Lyapunov exponent technique -which is typically used as a proxy of transport fronts (d' Ovidio et al., 2009) -yields an indication of stronger activity off the shelf, but has no clear large-scale latitudinal gradient. Instead, a region with a much larger latitudinal extent than the plume displays intense stirring without suggesting where the seasonal chlorophyll plume may be constrained.
Some improvement in the outlook for understanding and predicting the SCHL pattern emerges when altimetry is used to reconstruct the 30 day longitudinal and latitudinal origin of surface waters (Fig. 5). The northern flank of the SCHL plume appears to be bounded by the water stream coming from 47 • S (color-coded in orange in Fig. 5b), which interestingly corresponds to the water which peels off from the north-east corner of the shelf break (located at about 72 • E, 47 • S). Similarly, the southern flank seems to correspond to the stream which hugs Heard Plateau from the south, detaching cyclonically around 75 • E, 72 • S. This observation in fact is coherent with previous work (Mongin et al., 2009), which has shown that the large-scale pattern of the satellite-derived SCHL plume is determined by the wintertime advection of water emanating from the plateau: in this regard, the two departure points over Kerguelen and Heard plateaux mark the edges which enclose the Antarctic Circumpolar Current (ACC) branch which proceeds eastward after having passed over low bathymetry, and after being enriched in iron. Does this large-scale relation between horizontal advection and shape of the plume extend down to the mesoscale? This is an intriguing question, because it suggests the possibility of reproducing the fine scale (lobes and filaments) of the SCHL plume by a simple advection scheme based on altimetric velocity. Note that the rationale behind this approach is that horizontal advection preconditions the horizontal spread of iron east of Kerguelen, which then fuels chlorophyll production when restratification increases average light levels. An interesting consequence of this scenario is that this approach permits the shape of the spring bloom before the onset of the bloom to be estimated from altimetry data, as most of the advection occurs during the wintertime. In order to test this possibility, we first developed a new Lagrangian diagnostic, aimed at identifying along backward-in-time trajectories whether a water particle has been in contact with the plateau during its past advective history (see Methods for details). The diagnostic provides both the position (latitude and longitude) and time of the point at which the water parcels departed from the plateau (Fig. 6). These Lagrangian images provide the putative extension of the plume, including the positions of fronts which create contrasts in terms of origin from the shelf break and of time since last contact with the plateau. The temporal information is especially relevant for biogeochemical interpretations, as at a first approximation, the dissolved iron concentration can be modeled as decaying exponentially in time due to scavenging (Mongin et al., 2009).

Biogeosciences
These maps were computed starting from winter 2011 and were then updated daily using near-real-time altimetry. When the ship headed to the region for measuring pre-bloom conditions, the altimetry-derived calculations were used to forecast the extension of the plume and refine the position of the stations.

Validation
The altimetry-based estimation of the plume was validated during the cruise with ocean color images of the bloom, and post-cruise by analyzing the trajectories of the drifters released during the KEOPS2 and a following cruise 2 years later (MYCTO-3D-MAP), as well as from historical trajectories of Lagrangian drifters that had crossed the region.

Ocean color
The first cloud-free chlorophyll image (Fig. 4a) showing the full extension of the plume was received on the ship on 13 November 2011 and corresponded to the 11 November, i.e., after about 3 weeks from the beginning of the field operations. There was a lot of relief on board by noting the overall good correspondence of the entire plume with the altimetrybased calculation, because the calculation had been used to identify the sampling sites where the bloom was indeed expected to occur. Table 1 provides the extension and degree of superposition (congruence) for the altimetry-based forecast and the SCHL plume as seen in ocean color data. In order to define the boundary of the SCHL plume, we used a threshold value of 0.5 µg L −1 , which is a threshold in between typical springtime HNLC values and blooming waters (around 0.1 and above 1 µg L −1 , respectively). Consistently, the extension of the plume obtained with this threshold from ocean color data is approximately the same as the extension of the plume computed by the Lagrangian analysis, setting a threshold of 3 months for the age diagnostic (see Table 1). Although the objective of the altimetry-based diagnostic is not the prediction of the bloom intensity, note that there is at least a visually qualitative correspondence between "old" water parcels (i.e., the parcels that have left the plateau for comparatively longer times) and low SCHL regions, in particular in the region centered at 73 • E, 49 • S. This observation suggests that iron removal after leaving the plateau is a more important determinant of spatial iron distributions than variations in iron sources along the shelf break. This does not preclude variations in Fe supply over the plateau, if flow along the shelf break homogenizes the off-plateau supply (as discussed further in Sect. 3.5). This perspective that age since leaving

Drifters
In order to validate the altimetry-based estimation of stirring patterns, and in particular the origin from, and time since leaving the apparently iron-rich Kerguelen shelf break, 48 SVP Lagrangian drifters were released during the KEOPS2 cruise. An ideal plan would have released the drifters on a  regularly spaced linear array which would have followed the shelf break. As this would have consumed too much ship time, drifters have been released instead on transit from one station to another -when approaching the shelf break -or adaptively, when crossing key dynamical features like fronts. In order to palliate to this sub-optimal release scheme, a second opportunistic release experiment was performed by the MYCTO cruise during January-February 2014. This second experiment targeted particular regions which remained undersampled by the KEOPS2 scheme. Other 2012-2014 SVP historical drifter trajectories were also included in the analysis.
The time since having left the plateau and the latitude of departure from the shelf break also were estimated for each SVP drifter trajectory in analogy with the altimetry-based calculation. Results and comparison are shown in Fig. 7 (first column: altimetry; second row: SVP drifters). As SVP drifters are anchored to a fixed, shallow layer (nominally 15 m) under the effect of Ekman drift, the comparison is performed by including the Ekman component in the altimetrybased advection scheme (see Data and Methods). Not surprisingly, the SVP plot is much noisier than the altimetryderived one due to its asynopticity. However, qualitative analogies can be readily seen. In particular, the presence of a "fast lane" stemming from the northern part of the Kerguelen Plateau and contouring the polar front, a recirculation region centered at 74 • E, 49 • S enclosing older water, and a latitudinal mixing frontier at 85 • E, where old water parcel encounters are dispersed in thin filaments by strong mesoscale activity. These analogies can be quantified by plotting the ages and origin in scatter plots ( Fig. 7e and 7f) and noting their compact relationships. We also averaged the ages in a band contained between 48 and 50 • S (Fig. 8) along the latitude. This plot shows that the "age" directly computed from the drifters has a longitudinal dependence of about 5 days per longitude degree, corresponding to an eastward mean drifting speed of 0.15 m s −1 . On this constant drift, local deformations are visible, reflecting the presence of the retentive (hence aging) mesoscale recirculation structures of Figs. 6 and 7. In particular, the drifter-based age has three local maxima centered at 74, 77 and 79 • E. We repeat the same calculation for different altimetry-based products. The regional and global geostrophic products greatly overestimate the retentive effect of the first recirculation region and place its center closer to Kerguelen (at 73 • E). They both reproduce the mean eastward drift correctly, with the regional product shadowing (e-f) Scatter plots of drifter-derived vs. altimetry-derived data for the age and the origin from the plateau, respectively. Table 1. The table provides the extent and the degree of overlap (congruence) for the altimetry-based forecast -defined through the "water age" diagnostic -and the plume visible from remote-sensed chlorophyll (Chl.). The boundaries of the plume are defined, setting threshold values. The thresholds for the forecasted plumes were chosen to match (within 10 %) the extension of the chlorophyll plume with typical moderate and high values (0.5 and 1 µg L −1 , respectively). The overlap is 36 % for the lower chlorophyll threshold (hence for the larger plume) and 28 % for the higher chlorophyll threshold (smaller plume).

Threshold
Extent of plume Overlap Age Chl. Age Chl. Extent % (days) (µg L −1 ) (km 2 ) (km 2 ) (km 2 ) Plume min 30 1 1.26 × 10 5 1.29 × 10 5 3.59 × 10 4 28 max 120 0.5 2.96 × 10 5 3.17 × 10 5 9.5 × 10 4 36 the drifter-derived calculation closer than the global onewith lower ages, hence stronger and more realistic velocities (black vs. green line). Interestingly, the overestimation of the retentive region does not seem to be originated by errors in the geostrophic velocities, but by the fact that drifters also react to Ekman velocities. When Ekman velocities are included into the regional product (blue line), the effect of the retentive region over the age is predicted remarkably well by the satellite product.

Choice of the KEOPS2 stations
A great advantage stemming from the ability of the multisatellite maps to reveal the presence and locations of regions with contrasting origins and histories is that the smallscale variability in the age since leaving the plateau can be exploited to examine the large-scale contrasts in physicalbiological coupling. Based on the Lagrangian analysis, the chlorophyll plume can be zoned in terms of expected biogeochemical contrasts: a recirculation region; a jet, on the north flank of the polar front; a cold water tongue propagating northward along the age longitude drifter regional product -Ekman regional product -geostrophic global product -geostrophic eastern shelf break and unfertilized, HNLC waters (the "reference" case). For KEOPS2, this perspective was developed in real time. The initial station sampling -consisting of a north-south transect to capture latitudinal variations, with an east-west transect to examine variations with distance from the plateau -was modified opportunistically to examine the development of a bloom within the "fast-lane" jet to the north of the polar front, and to cover features evolving in a quasi-Langrangian time series within the recirculation feature (Blain et al., 2015, this volume). In addition, these perspectives have provided insight with respect to post-voyage analysis and clustering of the results to provide improved insights into the links between iron supply duration and ecosystem responses (e.g., Lasbleiz et al., 2014;Planchon et al., 2015;Trull et al., 2015, this volume).

Estimation of the mesoscale iron field
Assuming that the iron sources at the shelf break are homogenized by local mixing process and that iron dynamics can be modeled by first-order removal, the "age" field can be converted into an iron field (see Appendix). The model has only one free parameter, the scavenging constant λ, which relies both on in situ measurements of DFe concentrations during KEOPS2 and on the determination of τ by satellite altimetry. The stations off the plateau were sorted into two groups, young stations with an estimated age around of 20 days, and old stations with an estimated age of around 60 days (Fig. 9a).
We first estimated the removal constant for abiotic conditions (λ abio ) using the changes in age of the winter surface mixed layer-integrated DFe values. Assuming a firstorder reaction we derive λ abio = 0.041 ± 0.006 d −1 (r 2 = 0.869), for those stations located in the northern part of the plateau. Decreasing the age of the young stations to 10 days does not significantly change the estimate of λ abio (0.038 ± 0.005 d −1 ), but decreases the coefficient of correlation (r 2 = 0.724). The estimation of the mesoscale iron field is displayed in Fig. 9b.
We then determined the removal constant λ bio during bloom conditions, from DFe-integrated values (0-150 m) at the young stations where a large bloom was already well developed at the time of sampling. Combining these values with the integrated values of the source yields λ bio = 0.058 ± 0.005 d −1 .

Estimation of the exported iron flux
Having estimated Q source (Appendix A), λ and τ , we can estimate the exported flux at any station characterized by its age. For example if τ = 20 days (young water), F exp is in the range 2.6 − 4.3 10 −6 mol m −2 d −1 and in the range 0.2−0.8 10 −6 mol m −2 d −1 for τ = 60 days. Finally, the time scale t s , by which this approximation holds, is given by the time at which λt s 1 or t s 1/λ. Thus the estimate of the vertical flux is valid for comparison with short-term deployment of sediment traps (duration of the deployment of 4-5 days).

Discussion
Biogeochemical field studies are becoming more and more interdisciplinary, so that an increasingly large number of parameters have to be collected and only a very limited number of stations can be occupied during the same day. Nowadays the big challenge of a biogeochemical campaign is how to address this trade-off between the biogeochemical analytical resolution and the spatiotemporal coverage with a limited number of multidisciplinary stations. In-cruise knowledge of the biogeochemical provinces present in a region is therefore essential information to avoid wasted efforts, for instance sampling the same conditions multiple times. In this regard, remote sensing is an unavoidable tool, because it is the only observation capable of snapshots of regional variability.
The time scales characteristic of a bloom (days to weeks) are also the ones of the (sub)mesoscale and in particular of horizontal stirring. This coupling can create very complex biogeochemical contrasts which evolve in time during the campaign itself. One obvious source of information for mesoscale transport is satellite altimetry. Our study shows that these data can be mapped into biogeochemical regions by dedicated diagnostics. The model we developed for the KEOPS2 cruise may be considered as an attempt to translate altimetric sea surface height (SSH) patterns into patterns of primary production, maximizing the spatiotemporal information. Indeed, the model provides a pre-condition to the bloom, but does not inform on the intensity (or it does in a qualitative way only), nor on the timing. It is important to notice that this ability to estimate the plume extension at high precision by such a simple model is that the biogeochemical drivers that control the onset of the bloom are relatively simple: high nutrient waters with only one (spring-time) limiting micronutrient; a fixed source of the limiting micronutrient constrained by the topography (the Kerguelen-Heard Plateau), possibly homogenized along the shelf by a boundary current, and hence stable in time and space; deep winter convection all over the open-ocean region, which dissipates possible nutrient vertical inhomogeneity in the upper layer of the water column. Moreover, the circulation in the region is well captured by altimetry, because transport is dominated there by a barotropic current with a strong geostrophic signal (the ACC). These conditions make the Kerguelen region an ideal large-scale laboratory for studying fertilization events occurring in the Southern Ocean HNLC systems.

Future developments of the iron dispersion model
It is tempting to extend the altimetry-based model into the vertical for predicting the phytoplanktonic bloom more quantitatively. Indeed, without explicitly accounting for vertical dynamics, one may be surprised by the possibility of studying nutrient dynamics at all. However, in the Southern Ocean, geostrophic horizontal velocities are typically strongly correlated with the flow within the thick surface layer (of the order of 100 m, see Vivier et al., 2014). Our model therefore is a representation of iron transport in the ocean upper layer, whether it is occurring directly in the mixed layer or inside a deeper reservoir which is upwelled by wintertime convection or other vertical processes. The impossibility of re-solving the vertical mechanisms in space and time therefore does not hinder the capacity of locating the extension of the iron layer which preconditions primary production. On the other hand, the model does not provide any information on the intensity of the bloom and not even on its phenology (i.e., timing). Adding vertical dynamics to this model, or at least a parameterization, preserving at the same time the fine details directly observed by altimetry, is challenging. The most obvious solution, the use of a three-dimensional circulation model, would only maintain the mesoscale dynamics in a statistical sense at best, even when constrained by assimilation. A full three-dimensional biogeochemical model would therefore be a natural choice for estimating the regional biogeochemical budgets, but not useful for a (sub)mesoscale adaptive strategy, which requires the individual fine-scale structures to be tracked at the best possible precision. For future adaptive sampling, a more promising solution may come from methods which attempt to estimate the vertical dynamics by propagating surface information (e.g., SSH or SST, Isern-Fontanet et al., 2008) in the ocean interior. These vertical dynamics would require in situ calibration and validation like an array of density casts. This strategy is more complex, costly and time-consuming than the deployment of Lagrangian drifters, and should probably require at least two ships: one for high-resolution mapping (for instance with a towed vehicle), the second for traditional in-depth biogeochemical and physical stations (e.g., CTD/rosette). In terms of future biogeochemical campaigns, it may be instructive to ask which aspects were the most informative, which ones redundant and which ones should be improved. The most important information for planning the position of the stations was the possibility of sampling waters with contrasted iron concentrations. What would have happened to guiding the adaptive sampling strategy we used without the Lagrangian model? Considering (i) that iron is scavenged during its advection, (ii) that the main source of iron is the Kerguelen Plateau and (iii) that the mean circulation is eastward, it would have probably been tempting to choose the eastward distance from the plateau as a first-order proxy of iron scavenging. However, this choice would have been quite misleading due to the presence of mesoscale activity, which creates meandering "high-speed lanes" and recirculation features so that the iron-to-distance relation may become quite distorted locally. In particular, we showed that in the few hundreds of kilometers east of Kerguelen (which was the region in the range of the ship), the recirculation region discussed above creates a reservoir of "old" -hence iron-depletedwater in the vicinity of the source, and in turn a locally inverted iron-to-distance relation (i.e., iron locally increasing with distance, and not decreasing).
Although the Lagrangian diagnostic discussed here is specific to the Kerguelen situation, the possibility of analyzing multisatellite data with diagnostics based on mechanistic models is a general approach which looks to be an appealing companion of multiplatform campaigns. In this study we have performed part of the validation of the diagnostics after the cruise, because the drifters that we used were released during the campaign. In the case of a multiplatform campaign, a better and safer strategy would be to deploy some of the instruments (like Lagrangian drifters and profilers) before the cruise, in order to have better information on the physical drivers active in the region and to integrate in situ information with remote sensing. The possibility of releasing some instruments before the campaign would also avoid the scenario in which a diagnostic is first used for guiding a campaign, and then discovered to perform poorly.

Conclusions
Biogeochemical campaigns are under heavy pressure to resolve biogeochemical processes at higher precision and to include an increasingly larger number of coupling mechanisms. This need is changing our conceptual view of a biogeochemical process from an approximated zerodimensional or one-dimensional water column to a spatial extended system with full three-dimensional physical drivers on one side, and complex top-down ecological controls on the other one. This paradigmatic change is blurring the boundaries between biogeochemistry, physics and ecology, pushing towards end-to-end studies. In this perspective, the use of Lagrangian tools for dynamically zoning a region, mapping the transport structures present in a region in nearreal time will become even more critical for disentangling the temporal from the spatial biogeochemical variability, and for optimally choosing sampling stations in terms of their representativeness.

A1 Stirring pathways and iron supply
Along a trajectory, the iron content of a water parcel is not affected by transport, because altimetry-derived velocities are almost divergence-free. For a water parcel with coordinate x, y, t, we designate by Q 0 the integrated iron concentration over the plateau and by τ (x, y, t) the time since when the particle has left the plateau. We assume a first-order scavenging process so that the iron content in a given water parcel is given by Q(x, y, t) = Q source e λτ (x,y,t) . (A1) From this equation we see that the spatiotemporal variability of the field Q is completely inherited by the spatiotemporal variability of the field τ . The advantage of rewriting the advection scheme in this form is that the field τ is a purely advective quantity, i.e., it depends only on the velocity field. In physical terms, the field τ can be interpreted as the "age" of the water parcels since leaving the plateau. Calling the time when a particle has left the plateau t 0 (the "birth" of the trajectory), the "age" of that parcel at time t (and hence the length of time for which scavenging has occurred) is τ = t −t 0 . In practice, the time t 0 can be found by integrating the velocity field backward in time from x, y, and calculating the time when the backward-in-time trajectory hits the plateau. Note that τ is not defined for particles which do not come from the Kerguelen Plateau. For these cases, we impose τ = ∞. A map of τ (Fig. 9a) provides a snapshot of the stirring pathways east of the Kerguelen Plateau. Uncertainties in the estimation of the outflow of water from the plateau, initial iron concentration, scavenging rates and plume extension can be expected to be large, and it is difficult to make a complete inventory of the major processes that can affect iron dynamics. Nevertheless, we develop here a simple model for estimating the order of magnitude of the flux of iron associated to the Kerguelen plume.
We start by assuming that the iron concentration in the upper layer (0-150 m), considered as the winter mixed layer in the plume over the season, is entirely controlled by the offplateau iron flux H and by vertical removal processes Fe exp . Thus, the rate of change of iron content Q (in mol m −2 ) can be written as from which at equilibrium (dQ/dt = 0): H = Fe exp . We estimated the horizontal iron supply (denoted H in the Appendix) using altimetry-derived geostrophic surface currents. The vectors of the surface currents were decomposed into along-plateau and cross-plateau components pointing, respectively, outward and inward with respect to plateau bathymetric contours. The computed outflow water flux out averaged over the period 2006-2013 was 60.1 ± 4.510 3 m 2 s −1 . This flow can be partitioned into 45 % coming from the northern part of the plateau (north of 49.5 • N, i.e., north of the polar front), and 55 % coming from the southern part of the plateau. The integrated concentration of iron (Q source ) in the outflow (top 150 m of the water column) was derived from direct measurements of iron concentrations at stations located on the plateau. These values have large variability van der Merwe et al., 2015, this issue). The source from the northern part of the plateau (i.e., Kerguelen Plateau) was estimated to be 153-236 and 138-158 µmol m −2 in winter and spring respectively, whereas the source in the south (Heard Plateau) was estimated to be ∼ 77 and 54 µmol m −2 in winter and spring, respectively. It is likely that this difference results from shallower bathymetry and the proximity of the island in the north (high-end estimate). Because the travel time to reach the downstream plume edge is circa 3 months, the horizontal flux in November-December was calculated using the wintertime Q source estimate (based on the remnant winter water temperature minimum Fe data), and the horizontal flux in January-February was calculated with the springtime Q source estimate (based on the October to November surface water Fe data).
Using the altimetric mean velocities and these Q source estimates, we obtained the total amount of DFe injected into the 0-150 m layer per day. To express this as an areal average flux, we derived the size of the plume from satellite ocean color images. Using a threshold of 0.3-0.4 µg L −1 Chl, the area is 2.5-3.8 × 10 11 m 2 . This is in good agreement with the estimates based on altimetry (2.7-3.3 × 10 11 m 2 ), assuming an age of the water parcel in the plume not older than 90 and 120 days, respectively. Thus the average supply flux for the whole plume is 2.4 ± 0.6 × 10 −6 mol m −2 d −1 in October to November and 1.7 ± 0.4 × 10 −6 mol m −2 d −1 in January to February. These values are obtained propagating the uncertainties on the iron estimation indicated above and have to be considered as orders of magnitudes. Indeed, if the calculation is repeated by including the uncertainties in the altimetry-derived water fluxes and plume extension, the uncertainty over the iron supply flux grows to 60 %. The calculation does consider different iron concentrations for the Kerguelen Plateau and for the Heard Plateau, but it neglects possible smaller-scale iron variability along the shelf break of these two plateaux. Although the along-shelf jet associated with the meander of the polar front, east of Kerguelen Plateau (dashed line in Fig. 1), mitigates this issue providing a homogenization effect, we are not able to exclude that some contrasts in iron concentration still remain. Future campaigns aimed at pinpointing iron sources will undoubtedly help in constraining this uncertainty. The idealized model Eq. (A2) also allows estimation of the vertical flux of iron that is exported below the plume into deeper waters. Substituting in Eq. (A2) a removal term F exp which follows a first-order law: the solution at equilibrium is The amount of iron that is exported (Q exp ) at 150 m during the time t s at each location of the plume is Q exp = Q source e −λτ (1 − e −λt s ). (A5) If λt s 1 (i.e., for removal periods shorter than a month), then the exponential can be expanded at first order: Q exp = Q source e −λτ (1 − 1 + λt s ) = Q source λt s e −λτ . (A6) Thus the vertically exported flux during t s is F exp = Q source λ e −λτ . (A7)