Evaluation of the search and rescue LEEWAY model in the Tyrrhenian Sea: a new point of view

The trajectories’ prediction of floating objects above the sea surface represents an important task in search and rescue (SAR) operations. In this paper we show how it is possible to estimate the most probable search area by means of a stochastic model, schematizing the shape of the object appropriately and evaluating the forces acting on it. The LEEWAY model, a Monte Carlo-based ensemble trajectory model, has been used; here, both statistical law to calculate the leeway and an almost deterministic law inspired by the boundary layer theory have been considered. The model is nested within the subregional hydrodynamic model TSCRM (Tyrrhenian Sicily Channel Regional Model) developed in the framework of PON-TESSA (Programma Operativo Nazionale; National Operative Program – TEchnology for the Situational Sea Awareness) project. The main objective of the work is to validate a new approach of leeway calculation that relies on a real person in water (PIW) event, which occurred in the Tyrrhenian Sea in July 2013. The results show that by assimilating a human body to a cylinder and estimating both the transition from laminar to turbulent boundary layer and the drag coefficients, it can be possible to solve a force balance equation, which allows the search area to be estimated with good approximation. This new point of view leads to the possibility of also testing the same approach for other different categories of targets, so as to overcome the limitations associated with the calculation of the leeway in the future by means of standard statistical law.


Introduction
Meteocean and environmental forecasting is increasingly being used in operational decision-making in the sea for demographic, geographic and strategic applications.Safety of lives and assets at sea are a shared objective of many countries.Having an efficient ocean forecast system is essential to improve the prediction of sea state and to provide useful environmental ocean information, so as to increase the effectiveness of search operations (Breivik et al., 2013).In the event of an accident, timely search and rescue (SAR) intervention is helpful in significantly lowering the loss of lives and also to contain the damage.A considerable amount of resources is currently invested in maintaining SAR capabilities by major maritime nations.However, on many occasions, the available SAR capabilities prove to be inadequate to provide timely assistance in distress scenarios.Nevertheless, the effectiveness of the available SAR capabilities can be increased if we take advantage of high-quality environmental forecast data for SAR planning.An example of operational ocean forecasting and services coupled with search and rescue (SAR) activities is represented by the Global Ocean Data Assimilation Experiment (GODAE BLUElink) operational ocean prediction system (Brassington et al., 2007), used by the Australian Maritime Safety Authority.A list of global and regional operational ocean forecasting systems, supported by GODAE (Bell et al., 2009), can be found in Davidson et al. (2009), while a recent review of the evolution of the global and regional forecasting system from GODAE into GODAE OceanView (Bell et al., 2015) is described in Tonani et al. (2015).These systems are based on ocean general circulation models and data assimilation techniques that are able to correct the model using information inferred from different types of observations, acquired by various sensors and platforms.
In coastal areas, high-frequency radars are crucial to acquire surface two-dimensional current data (Barrick et al., 1977;Barrick et al., 2012;Cianelli et al., 2015;Paduan and Rosenfeld, 1996).Maps acquired by a coastal highfrequency (HF) radar have been used both for obtaining backtracked trajectories of floating objects (Abascal et al., 2012;Abascal et al., 2009b;Berta et al., 2014) and for SAR applications (Ullman et al., 2006;Ullman et al., 2003).Some European countries have invested important efforts toward the implementation of national HF radar networks (Quentin et al., 2013); (Carrara et al., 2014), but at present just one unified HF coastal radar network has been started in the Mediterranean Sea.It is the Tracking Oil Spills and Coastal Awareness (TOSCA) network, covering the Aegean, Adriatic, Tyrrhenian, Ligurian and Balearic seas; it started as support for the decision-making process related to marine accidents concerning oil spill pollution and search and rescue (SAR) operations (Bellomo et al., 2015).The TOSCA network currently covers sensitive and environmentally relevant areas, affected by intense ship traffic and/or the presence of oil pipelines; nevertheless it is only the first step towards the building of an integrated HF radar system in the Mediterranean regional alliance, and its coverage is still relatively small, so it is neither always involved in the ocean forecasting systems nor in related connecting applications such as SAR planning.
The first step in marine search planning is to determine the most probable area containing the searched object, and that is even more important if the target is a person in water (PIW).The definition of the probable search area is essentially linked to the quantification of some unknowns, such as the last known position (LKP) and typology, shape and dimensions of the objects.Moreover, a PIW or object without propulsion is also subject to drift from ocean currents, wave action and direct wind action (Davidson et al., 2009).Hence, effective SAR planning requires accurate environmental forecasting, especially regarding surface currents, temperature and surface winds on a short-range timescale.Today, such information can be provided in near-real time (NRT) from global and regional operational oceanography systems following on from the successful implementation of the Global Ocean Observing System (GOOS).
In the Mediterranean Sea, operational forecasts starting from the year 2000 are provided by the Mediterranean Forecasting System (MFS) (Tonani et al., 2008).The MFS uses a horizontal grid of resolution of 1/16 • and provides detailed forecasts at a regional scale for the whole Mediterranean Basin.The fundamental part of the system is the assimilation scheme for the blending of the observations into the model in order to provide the most accurate description of the past and the best initial condition for the forecast fields on a large scale (Oddo et al., 2009).
In many cases, SAR activities may also require fields of interest on a high spatial resolution.Subregional forecast systems that can resolve small-scale processes as well as fronts characteristics of the study area are necessary.A numerical technique widely used in operational oceanography to increase the horizontal resolution is the downscaling procedure (Pinardi et al., 2003).It permits the transmission of information at the interconnecting boundaries (temperature, salinity and velocity) from the coarse-resolution grid to the fine-resolution grid (Sorgente et al., 2003).This has been achieved through implementation of a nested high-resolution ocean model for the central Mediterranean Basin developed in the framework of the project Development of TEchnology for Situational Sea Awareness (TESSA).This hydrodynamic model, named Tyrrhenian Sicily Channel Regional Model (TSCRM), is based on the Princeton Ocean Model (POM) (Blumberg and Mellor, 1987) and has a horizontal resolution of 1/48   Here the leeway is calculated by statistical law as in Breivik and Allen (2008) (green line in Fig. 5).The probability of containment (POC) is also computed, and it varies from blue/minimum to red/maximum value in each grid box of 0.6 × 0.6 • .At rescue time (b), the target (magenta point) is outside the particles' cloud; it will be included in the cloud only 7.5 h later (c).
Education, University and Research of Italy.The general aim is to improve products and services of operational oceanography in southern Italy and to integrate them with technological platforms.The latter are set up to disseminate information for situational sea awareness (SSA), and also in support of SAR activities at sea.In Italy these activities are performed by the coastguard, within the competence of the Ministry of Infrastructure and Transport of Italy.
In this paper we present some numerical results demonstrating the prediction of PIW trajectories in the central Tyrrhenian Sea.We use two different approaches for the leeway calculation: the first approach is standard and it consists in leeway calculation by means of the statistical parameters shown in the Allen and Plourde (1999) table.The second one is a variant of the first approach: it consists in leeway calculation by means of the force balance equation based on modifications of the target boundary layer.The aim of this work is to validate the latter approach.We hope that in the near future, the new method will permit a practical formula of leeway calculation to be defined that is helpful not only in similar PIW cases but also in other categories of targets.
The new approach is embedded in the LEEWAY model, a drift model based on a stochastic approach (Breivik and Allen, 2008), combined with environmental forecast data from TSCRM and the European Centre for Medium Weather Forecasts (ECMWF).The following sections of the paper are organized as follows: Sect. 2 describes the models used and the numerical experiments; an analysis of the model results is presented in Sect.3, together with possible future extensions.A summary and concluding remarks are given in Sect. 4.

Material and methods
In this section we describe our drift forecast operational numerical model developed to predict the trajectories of floating objects (Sect.2.1), the characteristics of the LEEWAY model, also including our variation of the calculation of the leeway (Sect.2.2), and the different experiments (Sect.2.3) based on PIW cases.

Ocean model
The current forecasting model used in this work is the Tyrrhenian Sicily Channel Regional Model (TSCRM), an operational, subregional, nested ocean model implemented for the central Mediterranean Sea during the framework of the project TESSA.The TSCRM subregional ocean model covers the area from 8.98 to 16.5 • E in longitude and from 31 to 43 • N in latitude (Fig. 1).It is a free surface threedimensional primitive equation finite difference hydrodynamic model, based on POM (Blumberg and Mellor, 1987).It solves the equations of continuity, motion, conservation of temperature and salinity.The model assumes that the fluid is hydrostatic and the Boussinesq approximation is valid.The density is calculated by adaptation of the UNESCO equation of state revised by Mellor (1991); the horizontal viscosity and vertical mixing coefficients are computed respectively using the Smagorinsky approach (Smagorinsky, 1993) and the Mellor and Yamada (1982) turbulence closure scheme.It is an extension of the Sicily Strait sub-Regional Model im-plemented by Sorgente et al. (2011), and was successively made operational by Fazioli et al. (2016) and embedded into the regional model of the MFS (Tonani et al., 2008) through the downscaling procedure (Sorgente et al., 2003).This modeling technique, widely used to produce numerical weather predictions (Koch and McQueen, 1987), allows the model results to be downscaled from the regional scale to the subregional scale by providing values of temperature, salinity and velocity (three-dimensional fields' daily mean) at the interconnecting boundaries from the coarse-resolution grid to the fine resolution.This permits a more detailed description of the circulation in those areas where the regional model cannot resolve the small mesoscale features.By means of the downscaling technique, TSCRM receives information at lateral open boundaries using an offline one-way asynchronous nesting technique (Zavatarelli et al., 2002), coupled to the regional model with which it is embedded.Surface momentum and buoyancy fluxes are interactively calculated through standard bulk formulae.The latter use sea surface temperature as predicted by the model, and 6 h (00:00, 06:00, 12:00, 18:00 UTC) and 0.25 • atmospheric variables provided by ECMWF operational analyses (Sorgente et al., 2011).TSCRM produces 1 h and 1/48 • fields that are used with the 6 h and 0.25 • wind velocity fields to force the Lagrangian model LEEWAY.

LEEWAY model
The drift of a floating object above the sea surface is the result of the balance between hydrodynamic and aerodynamic forces.This balance induces the object to move with a certain angle relative to downwind direction (Richardson, 1997;Breivik and Allen, 2008); therefore it is fundamental to take the leeway into account in the estimation of the trajectories of drifting objects with good approximation.
The term leeway refers to an object's motion induced by the atmospheric wind (10 m reference height) and waves relative to ambient current (between 0.3 and 1.0 m depth).This definition standardizes the reference levels for the measurements of leeway for SAR objects and provides a practical way to utilize current and wind vectors from numerical models (Allen, 2005).The empirical relation between leeway and wind speed is given by empirical coefficients (linear regression coefficients and their standard deviations) for a total of 63 different targets classes, listed as a result of field campaigns performed by the US coastguard (Allen and Plourde, 1999).These coefficients give an estimate of the relation between the wind speed and the leeway velocity vector, converted into the following components: a downwind component and positive (right of the wind) and negative (left of the wind) crosswind leeway components (Allen, 2005).They are made explicit in the following equations: The wind velocity vector W 10 m is the wind speed measured at 10 m reference height, or simulated by a weather forecast model.The term L d of Eq. ( 1) represents the downwind leeway component, while the terms L c of Eqs. ( 2a) and ( 2b) are the crosswind leeway components, which show divergence of the SAR object from the downwind direction.The parameter a denotes the leeway rate (leeway to wind ratio), b is the regression intercept at zero wind and is the regression residual (Breivik et al., 2012).The determination of the leeway, or rather the relation between the wind speed and the leeway speed and divergence angle, gives a measure of how much the wind directly pushes on an object floating at sea.It is necessary to assume that the empirical coefficients of linear regression of Eqs. ( 1), ( 2a) and (2b) also include the contribution of the Stokes drift (Breivik and Allen, 2008).Following this, we assume that leeway can only be expressed as a function of the wind and we use a practical definition of leeway which does not tell the wind and wave influence apart.
This approximation can also be extended to other small objects and to most sea states, according to Breivik et al. (2012).
By estimating the linear regression coefficients of L d and L c leeway components, it is possible to recreate the expected drift of an object by means of modeled or measured data concerning the wind and sea current.As the LEEWAY model uses empirical formulae and imperfect approximation to the hydrodynamic laws, a Monte Carlo technique works to generate an ensemble through random perturbations; the ensemble allows the uncertainties of the forcings (wind and current), leeway drift properties (draught, length and beam) and the last known position (LKP) of the search object to be taken into account.The trajectory is then obtained by estimating the corresponding probability density function of containment of the object at each time step, and its envelope is the estimation of the search area (Davidson et al., 2009).
An exhaustive description of the stochastic approach used in the LEEWAY suite is given by Hackett et al. (2006).
In this work we present a variation of the calculation of the leeway, replacing the linear regression equation described by Eqs. ( 1), ( 2a) and ( 2b) with an almost (due to uncertainties  in forcing fields) deterministic law which is the balance of aerodynamic and hydrodynamic forces acting on the target, evaluating the modification of its boundary layer.The target is supposed to be approximately symmetrical; therefore the lift can be neglected.We suppose that such a target is tracked by the surface current under the hypothesis of steady motion: where F P is the weight force, F A is the Archimedes force, F ADw is the active drag force in water and F ADa is the active drag force in air.The first two forces give the emerged / submerged ratio, whereas the last two forces are responsible for the transport on the horizontal plane.Equation ( 3) is solved on vertical and horizontal planes.On the vertical plane, the equation is the balance between the weight force and the Archimedes reaction, and its solution means just using the reduced mass of the target.
On the horizontal plane, we solve the equation that distinguishes between the laminar and the turbulent boundary layer.The resistance to motion for a laminar boundary layer is only given by the viscosity friction and so we solve Stokes' law as follows: where µ w and µ A are the water and air viscosity, respectively; R W and R A are the submerged and emerged radii; u c and w are the current and wind velocity, while u B is the object velocity.
If the boundary layer is turbulent, the resistance to motion is given by the balance of the kinetic energy exchange between the moving body and the fluids:  5).Now the particles' cloud includes the target at the rescue time, but the dispersion is much larger than that shown in Fig. 3.
where S S and S E are the submerged and emerged area and C D W and C D A are the corresponding drag coefficients.
Finally we estimate the probability of containment (POC) by overlapping a spatial grid on the geographical one; we arbitrarily set each grid box to 0.06 × 0.06 • and we calculate the particles' percentage for each grid box, assigning a color bar from blue to red to distinguish the corresponding POC values from a minimum to a maximum value.The probability cumulative value is calculated by summing the values of highest probability of the individual cells; this meant that it could correspond to cells that are not geographically consecutive.A final record is kept to store the vertices of the cells that have a high probability.

Numerical experiments
In this work, the LEEWAY model is used to reproduce a real event that occurred in the western Tyrrhenian Sea, along the coast of Sardinia.On 11 July 2013 a man was seen for the last time at approximately 20:30 UTC (start point) on board a ferry that was on duty from Cagliari (Sardinia) to Civitavecchia (Italy), and only shortly before midnight (end point) was the alert given.His body was retrieved at 10:30 UTC on 12 July 2013 at a point of known coordinates (Table 1).
The LKP is a critical step, as the accuracy of this information is critical for the outcome of the search.In this work, we represent it as a line between start and end points (Table 1) recorded from the automatic information system (AIS) data provided by the general headquarters of the Italian coastguard (Fig. 2).
Numerical simulations of PIW trajectories are run using two different methodologies.The first one consists in the leeway calculation by means of the equations described in Eqs. ( 1), ( 2a) and (2b), while in the second one the almost deterministic laws described by Eqs. ( 3) and ( 5) are used.
In both the methodologies we rely on the conclusions of Breivik and Allen (2008) about the forcings' perturbation; therefore they are considered of secondary importance, espe- cially in our stable and relatively homogeneous conditions.We use the mean standard deviation of the current and wind components, calculated on the spatial grid and for the period of the simulation (43 h): the estimate is 0.07 m s −1 for the easterly current component and 0.1 m s −1 for the northerly one, while it is 2.68 m s −1 for the easterly wind component and 2.93 m s −1 for the northerly one.
We perform three different sets of experiments planting a total of 1000 particles in eight time steps from start to end points of last known position, and subsequently estimating the smallest distance between the mass center of each subset of particles and the retrieval point.
The first set of experiments relates to a check of the target configurations (sitting, vertical, horizontal/survival, horizontal/deceased and unknown status) shown in the Allen and Plourde (1999) table, as well as a study of sensitivity related to the change of the drift from a persistent direction to an opposite one (jibing) relative to downwind direction.This change will result in a sign change on the crosswind component.The jibing is set in the model as a frequency fixed throughout the simulation.Such information is very impor- tant: if the object does not jibe, the initial probability distribution could rapidly split into two equal probability distant areas, according to the initial uncertainty of the tack relative to downwind.In contrast, the more frequent the jibing is, the more central the distribution will be (Allen et al., 2010).We set a jibing frequency from 1 to 8 % for each listed position.
The second set of experiments is performed after choosing the best solution from the previous study; on the re-Figure 8.Here the standard error coefficients are set separately on the offset or on the slope of the regression line according to the wind velocity (blue lines in Fig. 5): for wind velocity lower than 10 m s −1 , the leeway perturbation coefficient is 4 times that in Breivik and Allen (2008) only on the offset of the regression line, while for wind velocity greater than 10 m s −1 , only that of the slope was changed.The particles' cloud includes the target at the rescue time (panel (b) in the figure), and the dispersion is reduced with respect to that of the preceding experiment.gression line used to calculate the downwind and crosswind leeway components, we check different coefficients of the time-invariant Gaussian perturbation ( n ), which Breivik and Allen (2008) introduced to include the variance that increases with wind speed.
Finally we carry out the third set of experiments: here the solutions are calculated using the almost deterministic approach.The human body geometry is assimilated to a cylinder with a height / width ratio between 4 and 7, according to Hoerner (1965); then the projected frontal area, which is predominantly responsible for the drag, is calculated.Due to air in the lungs, we suppose that the body is submerged in a standing position up to thorax; therefore we calculate both the submerged and emerged part.To include the error of the calculation of the real surface exposed to the flow, we introduce two different perturbations: a vertical random oscillation between the thorax and the neck and a random rotation around the vertical axis; at the same time we assume a nonsignificant Stokes drift.Many experiments are performed to find the critical Reynolds number, which marks the transition from the laminar to the turbulent boundary layer.Each critical value between 120 000 and 170 000 is coupled with all drag coefficients between 0.98 and 1.12 (the latter estimated according to Zdravkovich (1997) curves), resulting in different runs.Finally, we find the critical Reynolds number and the drag coefficient corresponding to the solution, where the distance between the planted particles' subsets and the retrieval point/time is the shortest.
To check the reliability of the almost deterministic approach, we run the model in a different geographical area, using a target similar to a person and referring to the critical Reynolds number vs. drag coefficient found in the last experiment.The data, provided by the Italian coastguard, were collected during SAR training; on 13 November 2013  8. Now the mass center of the third subset of particles (planted from 21:30 to 22:00 UTC) has the smallest distance at the rescue time: it is about 8.9 km and it is very similar to that of the first experiment.The smallest absolute distance is about 2.8 km, estimated 31.5 h later, and it belongs to the mass center of the fourth group.at 08:30 UTC a dummy was planted in the Sicilian Channel and it was retrieved 13 h later.The geometrically scaled target and the statistics of the new forcings are included in the model; we estimate a standard deviation of 0.09 m s −1 on the forcings for the easterly current component and 0.11 m s −1 for the northerly one, while it is 4.07 m s −1 for the easterly wind component and 4.25 m s −1 for the northerly one.
The data from the latter experiment are shown in Table 3  and 4.

Results
In this section, the main results of each set of experiments are shown.The first set of experiments is about the configuration of the target and the effects of the corresponding jibing; the coefficients of the regression line are set as in Breivik and Allen (2008) (green lines in Fig. 5).Here the unknown status, which has the highest error variance in the statistical calculation of the leeway, gives the best results; the change of the jibing values does not change the results significatively.This is why we opt for the LEEWAY standard value of 4 %.
The final map is shown in Fig. 3 where the target is outside the particles' cloud at the retrieval time.The distance between the mass center of each planting group and the retrieval point/time is shown in Fig. 4: we verify that the third group, planted between the 21:30 and 22:00 UTC, is the closest to the target, with a distance which oscillates around the value of 9 km; the smallest absolute distance (about 2 km) is reached by the fourth group, planted between 22:00 and 22:30 UTC, 31.5 h later than the retrieval time.
The second set of experiments studies the influence of the leeway coefficients' perturbations on the results.The leeway estimated by means of experimental data is heteroscedastic (the variance in fact increases with wind velocity; Allen and Plourde, 1999); therefore a spread around the regression lines in Eqs. ( 1), ( 2a) and (2b) is necessary.It is recreated, perturbing the standard error from a normal distribution N(0,σ ).In Fig. 5 all the regression lines checked are shown.As a first step we have checked the PIW σ values shown in the Allen and Plourde (1999) table (black lines in Fig. 5) and we have verified that they induce a dispersion that surpasses the advection.In theory we could correct this effect if we were able to better calculate the standard error on the regression Figure 10.Here the slope standard error coefficient is 1/5, i.e., 4 times that proposed in Breivik and Allen (2008), while the offset error is 0.5 (unchanged compared to that proposed in Breivik and Allen (2008)).The yellow lines in Fig. 5 are the corresponding regression lines' beams.Now the configuration seems good and the real event seems well reproduced.coefficients' estimators, i.e., if we knew the corresponding covariance matrix.We do not have this information, so we have checked different correction coefficients of the standard errors proposed in Breivik and Allen (2008), varying arbitrarily the initial values 1/20 in the slope standard error and 1/2 in the offset one, until a good final configuration is obtained, i.e., when the real event is reproduced with the smallest possible dispersion.We have verified that the spread is such that the particles' final cloud includes the target at retrieval time (Fig. 6) when the slope is a n = a + n / 5 and the offset is b n = b + 2 • n (red lines in Fig. 5), i.e., when the perturbation coefficients are 4 times larger than the ones proposed in Breivik and Allen (2008).The third group is again the nearest to the retrieval point, and the relative distance oscillates around the value of 5.5 km (Fig. 7); such a group has the absolute smallest distance (about 4.3 km) 5 h later relative to the retrieval time.As a result, the uncertainty of the time of the accident is reduced, and the distance between the mass center of the favored planted group and the retrieval point is almost halved; but the dispersion still seems high.We have improved this result, preserving the initial heteroscedastic but selectively increasing the originary correction coefficients according to wind velocity; if the wind velocity is < (>) 10 m s −1 , only the offset (slope) standard error is increased 4 times (blue lines in Fig. 5).Now the particles' cloud includes the target at the rescue time, and the dispersion of the particles is reduced (Fig. 8), but the distance between the mass center of the third group and the retrieval point/time is very similar to that of the first experiment, resulting in about 8.9 km (Fig. 9).As in the first experiment, the fourth planted group has the smallest absolute distance again (about 2.8 km), 31.5 h later than the retrieval time.
Finally we checked the error standard variations only on the offset or only on the slope, and we verified that a good final configuration is obtained when the offset error perturbation coefficient is the same as that proposed in Breivik and Allen (2008) but the slope error is different, resulting in the value 1/5 (yellow lines in Fig. 5).Now the smallest parti- cles' dispersion is obtained (Fig. 10); the third group has the smallest absolute distance from the retrieval point (about 3 km) and it is recorded at the retrieval time (Fig. 11).In other words, we can say that now the final configuration is good and the real event seems to be well reproduced.
The previous experiments pointed out the different performances of the model when the real statistics on forcings and the statistical regression line to calculate the leeway velocity are used.Not only the heteroscedasticity of the experimental dataset compiled by Allen and Plourde (1999) needs to be correctly accounted for to optimize the model performances, but we also need to choose the best perturbation coefficients to use in the regression line.We do not know them in advance and we only can say that they should be such that the dispersion does not swamp the advection role and at the same time, it should cause the POC to be maximized.
The third set of experiments consists in simulating the PIW drift by means of the almost deterministic law.After checking the projected frontal area by means of the allometric parameters for a "standard human body" according to Herman (2006), we choose a height / width ratio equal to 4.44.The final critical Reynolds number range is estimated between 165 000 and 168 000.The Reynolds numbers' heights are always between O(10 4 ) and O(10 5 ); due to these values and to the role of roughness on the bound-ary layer modifications (Yeo and Jones, 2011), we conclude that the PIW boundary layer in our experiment is often laminar, with transition to turbulence prevalently precritical.That means that the drag crisis never occurs, the boundary layer never separates and the drag coefficient in the turbulent boundary layer is constant for a large Reynolds numbers range.In our experiments the drag coefficient corresponding to the best results is 1.12.All the parameters to solve the force balance equation are shown in Table 2.
The results now show that at retrieval time the cloud includes the target (Fig. 12) and at same time the dispersion of the particles looks low and qualitatively comparable with that shown in Fig. 10, but now the absolute distance of third group of particles from the retrieval point is less than half, resulting in about 1.3 km (Fig. 13).
In Fig. 14, where the snapshots of Fig. 13 are shown without the probability map, it is visible that the particles move according to surface circulation patterns; some mesoscale and sub-mesoscale cyclonic and anticyclonic structures are visible, according to complex circulation of the area as described in Pascual et al. (2013), Rinaldi et al. (2010), Iacono et al. (2013) and Astraldi and Gasparini (1994).In particular, the area where the particles are planted is characterized by three principal surface structures: two cyclonic gyres (named "A" and "C"), separated by an anticyclonic gyre (named Table 3. Data for the second experiment: a dummy planted at a known point and rescued after 12.5 h.• N, are persistent during the period of the simulation and show circulation features according to summertime ones as described in Rinaldi et al. (2010) and Marullo et al. (1994).
Most of particles are planted along the western branch of gyre A while the smallest amount is planted between the cyclonic gyre C and the anticyclonic gyre B. After 18 h, the particles' set is separated into two subsets: the gyre A drives the larger subset westwards near the Sardinia coast while the variability of the current between gyres C and B drives the remaining particles to the open sea.This analysis allows us to conclude that the persistent summer gyre A is responsible for the motion of PIW in our experiment, and we think that in all probability, the person fell into water north of the retrieved point.We checked our approach in a different area (Sicilian Channel), using experimental data (Table 3) and scaling our target geometrically (Table 4).The results show that the target is inside the particles' cloud, of which the minimum distance, estimated to be about 5.8 km, is reached at the rescue time (Fig. 15).The target is located on the external border of the cloud (Fig. 16), then the POC is from medium to low value.Figure 17 shows that the particles' cloud is planted on the external part of an anticyclonic vortex, in an area where the current flowing along the southern Calabria coastline interacts with the current coming from the north, where the Messina Strait, separating Sicily from the Italian Peninsula, is located.It is a narrow passage where the Tyrrhenian and Ionian Sea sub-basins are interconnected; the interaction between its bathymetry and topography and the strong currents induce the generation of inertial eddies and strong horizontal current shears, generally located along both Sicily and Calabria main capes (Cucco et al., 2016).Although the Mediterranean Sea is characterized by very small tidal displacements which do not influence the circulation significatively (Sannino et al., 2015), in the Messina Strait, large gradients of tidal displacement are recorded because the semidiurnal tides in the Tyrrhenian and Ionian Sea are approximately in phase opposition.The tides in the Messina Strait are then the principal forcing of the circulation, which develops mainly along the main axis of the channel; in particular, during the flood the flow is northward, whereas during the ebb the flow is southward (Cucco et al., 2016).
In our experiment the dummy was planted in an area where the northern border of the anticyclonic vortex, located in the area 15.3-16.3• E, 37.2-37.65 • N, mixes with both the surface flow carrying the Ionian waters from the southern Calabria coast to the eastern Sicilian coast and the part of the current coming from the Messina Strait.The TSCRM model does not include the tidal forcing but the effects of the strong current from the strait, visible in Fig. 17 and probably induced by the topography, can be assimilated by magnitude order to the real transport during the experiment.This assumption is supported by results of Cucco et al. (2016), which show that the water fluxes calculated including the tidal contribution change between 800 % (considering the instantaneous flow) and 60 % (considering the residual flow), and the thermohaline properties of the coastal areas are distributed over an external area of about 1500 km 2 extending up to 60 km from the strait.The greatest impact of the tidal forcing is on the surface waters; therefore it plays an important role in the drifting of floating targets in this area.In general, both Fig. 14 and Fig. 17 underline the importance of having an efficient operational ocean prediction system for the SAR activities, not only as forcing of the dispersion Lagrangian model, but also because it gives much information useful to analyze the particles' cloud final evolution.In particular, in Fig. 17, it seems that approximation with the POC calculation is greater than that estimated in the Sardinia experiment, but we think that it can be improved if the set of particles is clustered according not only to planting time, but also to splitting the cloud into as many subsets as there are hydrodynamic structures; in this way each subset will have its own POC, which could be high even if it includes a small number of particles.We are currently working to test this hypothesis.

Conclusions
We have presented the results of numerical experiments performed to estimate the probability of containment (POC) of a person in water (PIW).We have referred to a real event that occurred in the western Tyrrhenian Sea.The LEEWAY model, nested within the subregional hydrodynamic model TSCRM, has been used.The real statistic of the forcings on the observation period was accounted for and the last known position has been based on the trajectory of the ferry during 3.5 h, recorded by the automatic information system; to reproduce the event, a set of 1000 particles planted in eight time steps was used.
The objective of the work has been to validate an almost deterministic law to calculate the leeway, based on the boundary layer theory; these results have been compared with those obtained estimating the leeway by means of the standard statistical approach.In this approach, the covariance www.nat-hazards-earth-syst-sci.net/16/1979/2016/Nat.Hazards Earth Syst.Sci., 16,[1979][1980][1981][1982][1983][1984][1985][1986][1987][1988][1989][1990][1991][1992][1993][1994][1995][1996][1997]2016    matrix of the line regression coefficients' estimators would be required to estimate the leeway correctly; it is not known, so we have checked different arbitrary coefficients of the corresponding standard errors.The best results have been obtained when a new value on the slope was set equal to 1.5, while the offset one was set equal to 0.5.
To calculate the leeway in the almost deterministic approach, the PIW was cylinder-shaped with a height / width ratio equal to 4.44; the critical Reynolds number was found in the range between 165 000 and 168 000 and the drag coefficients were estimated equal to 1.12.Now that the results have shown that at the retrieval time, the particles' cloud includes the target, and the third subset of particles is representative of the event, then we can also estimate the time of the accident with good approximation.These results are comparable with the results estimated by means of the standard statistical approach, only when the coefficients of the regression line standard error are correctly set, but we cannot know them in advance.The important result is that in the new approach, the real event is correctly reproduced and the distance of the mass center of the favored subset of particles is more than halved compared to the best solution obtained through the statistical approach, resulting in 1.3 km.The tests in the Sicilian Channel using a similar target have confirmed the reliability of the method.
A second important result is that the particles' distribution in both experiments is coherent with the hydrodynamic structures, highlighting the importance of having an efficient operational ocean prediction system for SAR activity; the hydrodynamic field is required not only because it forces the Lagrangian model, but also because it allows the final results to be read in a more full way.We think, in fact, that the results of the POC can be improved by splitting the set of particles into as many subsets as there are hydrodynamic structures, so that each subset can significatively contribute with its own POC; it corresponds to the consideration of as many timeevolving probability density functions of the location of the search object as there are hydrodynamic structures.
Finally, the Stokes drift will have also to be included.
The results obtained here encourage other different target categories to be explored.We believe that in the long term it will allow the limitations connected with the standard approach based on the statistical parameters to be overcome, providing leeway-drift formulae for practical use based on appropriate Reynolds critical numbers and drag coefficients for specific targets.This idea is supported by the realistic option to execute tests in naval tanks, simulating different meteomarine conditions, scaling opportunely specific targets and then overcoming the costs and difficulties of acquiring field data.
The data used to carry out this research are free and available on request by writing to the author.

Figure 1 .
Figure 1.Geographical domain of the Tyrrhenian Sicily Channel Regional Model (TSCRM).The color scale represents the sea bottom.

Figure 2 .
Figure2.The last known position, represented by the blue particles' distribution at the alert time; it is based on the automatic information system (AIS) that gives the ferry position every 6 s (green line on the left of the picture).The coordinates of the start and end point are defined in Table1.

Figure 3 .
Figure 3. Particles' cloud estimated by the LEEWAY model overlapped on the hydrodynamic field.Here the leeway is calculated by statistical law as inBreivik and Allen (2008) (green line in Fig.5).The probability of containment (POC) is also computed, and it varies from blue/minimum to red/maximum value in each grid box of 0.6 × 0.6 • .At rescue time (b), the target (magenta point) is outside the particles' cloud; it will be included in the cloud only 7.5 h later (c).

Figure 4 .
Figure 4. Distance of the mass center of each group planted every 30 min from 20:30 UTC (group 1) on 11 July 2013 to midnight the following day (group 8), corresponding to the solution in Fig. 3.It may be possible to verify that the lower distances at the rescue time (values at right and low part of the panel) are reached by the mass center of the third subset of particles (planted from 21:30 to 22:00 UTC).The smallest distance is reached 31.5 h later by the fourth subset (planted from 22:00 to 22:30 UTC) and is about 2 km.On the upper right side of the image, the trajectories of each group's mass center throughout the simulation (43 h) are visible.

Figure 5 .
Figure 5. Downwind (up) and crosswind (down) components of leeway vs. 10 m wind velocity.Each line's beam, distinguishable by color, is generated by perturbing the error standard from a normal distribution N (0,sigma) and arbitrarily changing the corresponding correction coefficient.The black lines are generated by means of random perturbation of the originary sigma value shown in the Allen and Plourde (1999) table.The green lines are obtained by means of the coefficients of the standard errors proposed in Breivik and Allen (2008) (hereafter called reference coefficients).The red lines are generated by multiplying the reference coefficients by 4; therefore here the regression line coefficients are a n = a + n / 5 (slope) and b n = b + 2 • n (offset).The blue lines are obtained by selectively increasing the reference coefficients (according to wind velocity) 4 times.Finally the yellow lines are calculated by only multiplying the slope reference coefficient by 4; therefore now the regression coefficients are a n = a + n / 5 (slope) and b n = b + −0.5 • n (offset).

Figure 6 .
Figure 6.Here the standard error coefficients are 4 times greater than ones in Breivik and Allen (2008) (red line in Fig.5).Now the particles' cloud includes the target at the rescue time, but the dispersion is much larger than that shown in Fig.3.

Figure 7 .
Figure 7. Distance of the mass center of each group planted every 30 min from 20:30 UTC (group 1) on 11 July 2013 to midnight the following day (group 8), corresponding to the solution in Fig. 6.Now the mass center of the third subset (planted from 21:30 to 22:00 UTC) has the smallest distance at the rescue time and it is about 5.6 km.The smallest absolute distance is about 4.3 km and it is reached 5 h later by the same subset of particles.The trajectories of each group's mass center throughout the simulation seem to have changed slowly with respect to the preceding experiment.

Figure 9 .
Figure9.Distance of the mass center of each group planted every 30 min from 20:30 UTC (group 1) on 11 July 2013 to midnight the following day (group 8), corresponding to the solution in Fig.8.Now the mass center of the third subset of particles (planted from 21:30 to 22:00 UTC) has the smallest distance at the rescue time: it is about 8.9 km and it is very similar to that of the first experiment.The smallest absolute distance is about 2.8 km, estimated 31.5 h later, and it belongs to the mass center of the fourth group.

Figure 11 .
Figure 11.Distance of the mass center of each group planted every 30 min from 20:30 UTC (group 1) on 11 July 2013 to midnight the following day (group 8), corresponding to the solution in Fig. 10.Now the mass center of the third subset of particles (planted from 21:30 to 22:00 UTC) has the smallest absolute distance at the rescue time (about 3 km).

Figure 12 .
Figure 12.Particles' cloud simulated by the model overlapped to the hydrodynamic structures.Here the leeway is calculated by means of the force balance equation.The cloud includes the target at the rescue time, and the dispersion does not seem to swamp the advection; in particular the particles' cloud seems to separate during a longer period of time according to hydrodynamic structures (e and f).

Figure 13 .
Figure 13.Distance of the mass center of each group planted every 30 min from 20:30 UTC (group 1) on 11 July 2013 to midnight the following day (group 8), corresponding to the solution in Fig. 12.Here the leeway is calculated by means of the force balance equation.It may be possible verify that now the solution is better than the previous ones: at the rescue time the third group (planted from 21:30 to 22:00 UTC) has the smallest absolute distance and is about 1.3 km.

Figure 14 .
Figure 14.Particles' dispersion in the Bonaria experiment, overlapped on the hydrodynamic field.Three principal surface structures are visible in the experiment: two cyclonic gyres, A and C, separated by the anticyclonic gyre, B. Cyclones A and C are persistent during the period of the simulation.Most particles are planted along the western external part of the gyre A while the smallest amount is planted between the cyclonic gyre C and the anticyclonic gyre B. After 18 h, the particles' set is separated into two subsets: gyre A drives the larger subset westwards near the Sardinia coast, whereas the variability of the current between gyres C and B drives the remaining particles to the open sea.

Figure 15 .
Figure 15.Distance of the mass center of the particles set planted at a known time/point to simulate the trajectory of a dummy in the Sicilian Channel.Again the leeway is calculated by means of the force balance equation.On the upper right side of the picture the trajectory of the mass center throughout the simulation (21 h) is visible.

Figure 16 .
Figure 16.Particles' cloud overlapped on the hydrodynamic field in the experiment executed in the Sicilian Channel.In (d) the configuration at the rescue time is visible.Far from the coast the current is prevalently eastward but the particles' subset influenced by the current near the coast simulate the trajectory from the planting point to the rescue point/time well.

Figure 17 .
Figure 17.Particles' dispersion in the Sicilian Channel experiment, overlapped on the hydrodynamic field.The particles' cloud is planted on the external part of an anticyclonic vortex, where the current flowing along the southern Calabria coast interacts with the strong current coming from the north, where the Messina Strait is located.The particles are driven according to the surface circulation pattern.

Table 1 .
Data describing the incident.