Predictions for oil slicks detected from satellite images using MyOcean forecasting data

The pan-European capacity for the Ocean Monitoring and Forecasting (MyOcean) Marine Core Service, implementing the Global Monitoring for Environment and Security (GMES) objectives, targets the provision of ocean state observations from various platforms and analysis and forecasting products to assist, among other downscaling activities, the needs of the operational response to marine safety, particularly concerning oil spills. The MEDSLIK oil spill and trajectory prediction system makes use of the MyOcean regional and Cyprus Coastal Ocean Forecasting and Observing System (CYCOFOS) downscaled forecasting products for operational application in the Mediterranean and pre-operational use in the Black Sea. Advanced Synthetic Aperture Radar (ASAR) satellite remote sensing images from the European Space Agency (ESA) and European Maritime Safety Agency CleanSeaNet (EMSA-CSN) provide the means for routine monitoring of the southern European seas for the detection of illegal oil discharges. MEDSLIK offers various ways, to be described in this paper, of coupling the MyOcean forecasting data with ASAR images to provide both forecasts and hindcasts for such remotely observed oil slicks. The main concern will be the drift of the oil slick and also, in the case of the forecast mode, its diffusive spreading, although some attempt is also made to estimate the changes in the state of the oil. The successful link of the satellite-detected oil slicks with their operational predictions using the MyOcean products contributes to the operational response chain and the strengthening of maritime safety for accidental or illegal spills, in implementation of the Mediterranean Decision Support System for Marine Safety (MEDESS-4MS) regarding oil spills.

The MyOcean project, which is carried out under the GMES fast track services, targets the provision of a Marine Core Service (MCS) based on the provision of ocean state variables, derived from in-situ and satellite remote sensing observations and numerical models, to meet the needs of the response agencies. These needs are driven by the EU Directives, regional and international conventions and treaties related to the environment and civil protection, marine safety, policy making, assessment and implementation. Particularly, MyOcean aims to provide products to assist downscaled and downstream services within key areas, such as Marine Safety.
The operational ocean forecasts for response to oil spill and floating objects have proved to be of assistance to the agencies related to marine safety and for reducing 15 the impact on the marine environment that may arise from pollution incidents. The initial operational response in oil spill incidents is to identify the exact location of the spill, to predict where the slick will drift, where and when it will arrive, which resources will be threatened and the state of the oil after a give time interval. Such downstream predictions need validated high resolution downscaled ocean and meteorological forcing 20 data. MyOcean regional products have been used operationally in CYCOFOS from January 2009 onwards, for initial and lateral forcing for its downscaled applications, for on-line validation purposes, for improving the CYCOFOS waves forecasts (Galanis et al., 2012) and for oil spill predictions (Zodiatis et al., 2010). CYCOFOS consists of downscaled and downstream forecasting systems and has been operational since 25 early 2002 (Zodiatis et al., 2003a(Zodiatis et al., , b, 2010. CYCOFOS hydrodynamical forecasting modules are improved, validated and consolidated periodically (Zodiatis et al., 2003c, following the developments in EU projects relevant to the operational ocean forecasting, such as the Mediterranean Forecasting System towards Environmental Introduction Predictions (MFSTEP), European Coastal-shelf sea Operational Observing and Forecasting System (ECOOP) and the Marine Environment and Security for the European Area (MERSEA). The MyOcean as well as the CYCOFOS downscaled products provide the local, regional and the EU response agencies, such as REMPEC and EMSA, with forecast for the predictions on the movement and the weathering of the oil spills. 5 In recent years, satellite remote sensing images have become one of the most important tools for detecting oil slicks on the sea surface (Klemas, 2010;De Dominicis et al., 2012a, b). In CYCOFOS, ASAR satellite data received from European Space Agency Environment Satellite (ESA ENVISAT) on a regular basis and processed at the Cyprus Oceanography Centre for possible oil slick detection in the Eastern Mediterranean. with routine monitoring of European seas, using ASAR satellite images for illegal oil discharges and to assist them in the implementation of the EU Directive 2005/35. The well established MEDSLIK oil spill and trajectory prediction system (Lardner, 2004;Lardner et al., 2006;Zodiatis et al., 2007) is adapted to use operationally the MyOcean regional products in the Mediterranean and pre-operationally in the Black 5 Sea, in order to support the response agencies for marine safety. The MyOcean CY-COFOS downscaled products in the Levantine Basin are used also by an automatic version of MEDSLIK (internally known as AutoMedTrack) without any operator intervention, along with ASAR satellite data for short forward and backwards predictions. These predictions can be superimposed with the Automatic Identification of Ships (AIS) 10 traffic information to assist the response agencies to identify the ship responsible for the detected oil slicks, contributing in this way in the implementation of the EU Directive 2005/35. This paper aims to describe the use of the MyOcean regional and downscaled forecasts by the MEDSLIK oil spill and trajectory prediction system with particular applica-15 tion of ASAR satellite images.

A brief summary of the prediction systems
As with any oil slick, one of the tasks arising from a satellite observation is to forecast where it is likely to move to with the aim of determining if any significant environmental resources are threatened by the oil, and if so deciding where and how best to protect 20 them. A second task, less urgent but equally important, is to estimate where the oil originated, that is, to trace its path backwards in time. general oil spill and pollutant prediction model that offers several tools for doing these tasks (Lardner, 2011). Reliable forecasts and hindcasts of oil slicks depend on the use of atmospheric and ocean forecasts. A primary such tool for MEDSLIK is the CYCOFOS forecasting system nested in the regional MyOcean products. This comprises marine forecasts (of 5 the ocean currents and sea temperature and salinity) for the Levantine Basin east of 30 • with resolution 1.85 km and high resolution (1 km) forecasts of the North-East Levantine Basin east of 31 • and north of 33 • . These forecasts are generated at six-hourly intervals. In addition, CYCOFOS produces 3 hourly forecasts of sea surface waves that cover the whole of the Mediterranean and Black seas with resolution of 10 km which 10 are used by MEDSLIK to compute the wave-induced drift (Stokes' drift) of a floating object.
Outside the Levantine Basin, use is made of a choice of MyOcean regional and downscaled forecasting systems from the Mediterranean Ocean Forecasting System (MFS) forecasts for the whole Mediterranean which have resolution 7 km (Pinardi et al., 15 2003;Tonani et al., 2008) to the Aegean Levantine Regional Model (ALERMO) system that covers the Mediterranean East of 20 • with a resolution of 3.7 km and several other limited area forecasting systems. For atmospheric data, MEDSLIK uses primarily the operational daily high frequency (hourly) SKIRON forecasts of winds (Kallos et al., 1998) to estimate the wind drift of the 20 slick. This forecasting system also covers the whole of the Mediterranean and Black Seas with resolution of 5.5 km.
The analysis of ASAR satellite images for potential oil slicks results in the generation of an ASCII file of type "*.tab" in which the boundary of each patch of oil is given in the file as a polygon digitised to a sufficiently fine resolution. Figure 2 shows two examples Introduction as visualized in the MEDSLIK visual interface module. These particular slicks will be used below to illustrate the oil spill model results. MEDSLIK offers different ways of making forecasts and hindcasts from such satellite data which will be described briefly in the following section. The method used for tracing the movement of any floating object is based on a two-stage modified Euler algorithm. 5 The first stage computes an approximate horizontal displacement during a time interval from t to dt of a particle whose original position on the surface is (x, y) at time t as where (u, v) are the x-and y-components of the surface velocity of the water, (W x , W y ) the corresponding components of the wind velocity and (S x , S y ) the components of the Stokes drift, all interpolated to the original position (x, y) of the object and the mid-time (t + 1 / 2 dt) of the time step. Also α is the wind-drag coefficient for the object and β the wind drift angle (by default in MEDSLIK taken respectively as 0.032 and 0 for a parcel of oil, although they are user-adjustable parameters). In the second stage, revised displacements, dx and dy are computed using the same expressions as above but with (u, v), (W x , W y ) and (S x , S y ) interpolated to the midpoint (x + 1 / 2 dx, y + 1 / 2 dy) and time (t + 1 / 2 dt). The position of the object at the end of the time-step is then taken as (x + dx , y + dy ).

20
For computation of a hindcast, identical expressions to those above are used for the displacement on any time step from t to t− dt except that dt is replaced by dt everywhere.
It is well-known that for this algorithm dx and dy are accurate to second order, that is, the errors in these variables are proportional to dt 3 . It would of course be feasible 25 to use a method of higher accuracy, such as one of the Runge-Kutta algorithms, but the cost of the extra computational time is hardly worth the benefit given the level of OSD 9,2012 Predictions for oil slicks detected from satellite images the errors in the present generation of forecasts. We have checked the above Euler method by following a forecast particle track with a hindcast starting from the end point of that forecast and the hindcast track lies almost on top of the forecast track, at least over a few days. The Stokes velocity (S x , S y ) is computed from the wave model using the expressions where H is the significant wave height (range from trough to peak), ω = 2π T with T being the mean period; k = ω 2 g is the wave number and d a unit vector giving the wave direction.

10
MEDSLIK models an oil spill by dividing the spill into a large number of Lagrangian parcels of equal initial sizes that are subject to translation by water currents and wind as well as to changes in the oil properties with the ageing of each parcel. The number of parcels may be chosen up to half a million, which may be necessary for a spill that lasts several weeks, although for more moderate spills as few as 10 000 parcels 15 suffice to give a reasonable picture of the spill behaviour. For drift of the slick, the size of the parcels is irrelevant but later we shall include estimates of the oil properties in which case the estimated volume of oil is divided equally among the chosen number of parcels (which was 10 000 in the examples to be shown later).
The parcels are translated on each time step according to the same model described 20 above. Some of the parcels become diffused into the water column by wave action and for these the forecast water velocities are interpolated in three dimensions to the depth of the parcel. There is of course no wind drag for such submerged parcels. The Stokes velocity decreases with the depth z of the parcel by a factor exp(-2kz). Advection of submerged parcels by vertical motion of the water is ignored as being very small.

OSD
These displacements are added to those given in (1). Here r(0,1) is a uniform random 5 number on the interval (0,1) and K h and K v are respectively the horizontal and vertical diffusivities. The third, vertical, displacement is applied only to already submerged parcels of course. It can easily be seen that the root mean square values of these diffusive displacements are A cloud of such particles undergoing random walks with these r.m.s. displacements satisfies the convection-diffusion equation with K h and K v the horizontal and vertical diffusion coefficients (Csanady, 1973;Ahlstrom, 1975;Hunter, 1987). 15 MEDSLIK has the capability of predicting the movement of one or more floating objects. This facility was originally developed primarily for Search and Rescue operations but it can be used to compute the movement of the boundary points of any slick or slicks specified in a tab file. The algorithm used is as described in Sect. 2. This approach has been followed to construct forecasts and hindcasts for 48 h of the This algorithm for generating such fore-and hind-casts from satellite observations has been automated and the resulting pictures are placed on the CYCOFOS website (http://www.oceanography.ucy.ac.cy/cycofos/oil-spill-predictions.php) on a regular 20 basis, immediately upon detection of oil slicks from the satellite images. During the MyOCEAN project more than 500 oil slicks were detected in the Levantine Basin using ESA ASAR images (Fig. 1) and after data processing at OC-UCY, short 24 h forecasts and hindcasts were computed and posted on the above web site, while the corresponding processed ESA ASAR images are placed on the CYCOFOS web-25 site http://www.oceanography.ucy.ac.cy/cycofos/remote-sensing.php on a regular basis too.

Forecast and hindcast trajectories of boundary points
MEDSLIK has the capability to plot the paths of shipping provided by AIS and to compute the minimum distance that any ship path approaches to any hindcast position of the slick, thus to assess whether any ship would have been capable of causing the slick, as is illustrated in Fig. 4. Figure 4a shows a set of actual ship tracks between 27-28 September 2008. The names of the ships have been changed for obvious reasons. On the early morning of 28 September 2008 a slick was observed in the position indicated on the Fig. 4a by a white circle in the centre of the brown square. To simplify 5 the subsequent figure, the tab file was simplified to a rectangle surrounding the slick and the hindcast tracks of the four vertices plotted in Fig. 4a and magnified in Fig. 4b. The hindcast was made for 48 h using the CYCOFOS daily high resolution forecasts for the region. The tracks are marked with a white dot every six hours. The green dots in Fig. 4b

Forecast and hindcast slick movements
An alternative approach to predicting the motion of a slick observed by satellite, more informative than simply following the movement of the polygon outlines in the tab file is to replace the polygons by a large number of parcels. The slick in the satellite file is separated into its patches and the patches are then randomly filled with oil parcels. The 20 subsequent forecast of slick position will then include diffusive spreading (modelled by a random walk of the parcels) as well as advective translation. A hindcast can also be made using this method but in this case diffusion cannot be included since diffusion is an unstable process when run backwards in time and cannot be simulated by a random walk model. 25 Figure 5 shows the positions of the hindcast and forecast slicks at intervals of 12 h for the satellite image shown in Fig. 2a. The same forecasts of winds, waves and water currents have been used as in Fig. 3 as black arrows while the wind vector is shown in white for the centre of each of the 12 hourly positions. The different colours in the slicks represent different oil concentrations, maximum oil thickness being represented by deep red and minimum by pale green. It is clear that the hindcast slicks show little change in shape and oil concentration as we move backwards in time and the slick outlines are similar to those in Fig. 3. On the 5 other hand the forecast outlines show significant diffusive spreading when compared to Fig. 3. Again, it is hard to distinguish the slicks for the final 24 h of the forecast due to the small movement.

Estimates of the ageing processes for an observed slick
If a spill is large enough and near enough to valuable resources, it would be desirable 10 to adopt response measures to stop it or to clean it up. The response agencies decision on what to do would depend on the amount and state of the oil. For example, if the oil is highly weathered and has become too viscous it would be pointless to spray surfactants on it.
MEDSLIK incorporates fate algorithms to estimate the ageing of the oil based for the 15 most part on algorithms developed by Mackay and co-workers (Mackay et al., 1977(Mackay et al., , 1979(Mackay et al., , 1980Mackay and Paterson, 1980; summarised in the appendices of Lardner, 2004 andDe Dominicis et al., 2012a, b). The most significant of these cover the processes of evaporation of the volatile fraction of the oil and emulsification of water with the oil, both of which have dramatic effects on the viscosity of the oil, and dispersion 20 of the oil into the water column by wave action. The volume of oil-water mousse on the surface is primarily influenced by these three processes. In order to compute these processes certain data are needed: the type of oil, the initial volume on the surface and the time of spill, or equivalently, the age of the oil in the observed slick. In the future it may become possible to estimate some of these variables from a satellite image but at 25 present we have to use guesswork or median values. is known. At present analysis of satellite images is not sufficiently refined to be able to determine thickness and the best we can do is assume the slick is all sheen of thickness 10 microns, which will give us a lower bound on the volume. Again, unless we already know the type of the oil, the best we can do is to take a median value and take it as medium grade oil of 33.5 API with average physical properties for this grade. 5 For the age of the slick, it may become possible to estimate this from the amount of lateral spreading of the image. At the moment an upper limit of 24 h can be placed on the age of the slick if it did not appear on the previous day's satellite image (which is often the case). From the volume of the observed slick we can calculate the original volume by mak-10 ing a simulation of the ageing processes alone from the presumed time of origin to the present using a guessed value of the original volume. The resulting final volume is compared with the estimated volume in the slick in order to correct the guess of the original volume of spilt oil. Secondly, using this corrected volume, a forecast run is made from the supposed time of origin to the time of observation (ageing processes 15 only) and then on to the required end of the forecast period (including advection and diffusion as well as further ageing). Thus a complete picture of the ageing is obtained as well as the future movement. Finally, in order to find the past motion of the slick, a hindcast can made as described above with no ageing or diffusion from the time of observation to the time of origin. MEDSLIK can perform all these operations automatically 20 with one mouse click.
As an illustration we have assumed the age of the slick shown in Fig. 2b is 24 h. The estimated volume of the observed slick was 521 cu m which led to an estimated initial volume of spilt oil of 234 cu m. Figure 6 shows the resulting forecast of the slick positions for 48 h and hindcast for the presumed age of 24 h at intervals of 12 h. Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | decreases due to dispersion; the viscosity increases asymptotically to its maximum value.
As a final illustration we consider the ESA ASAR image of the western part of the Black Sea from 23 June 2011, shown in Fig. 8. After data processing at OC-UCY, two probable oil slicks were detected in open water off the coasts of Romania and Bulgaria, 5 as marked on the figure. Forecasts and hindcasts for 24 h of these two slicks have been made by MEDSLIK using the Black Sea/MyOcean regional ocean forecasts together with SKIRON high frequency winds and the CYCOFOS wave forecasts. The results are shown in Fig. 9. The coupling of MEDSLIK with ASAR images was initiated following a relevant request from EMSA. 10 An interesting feature of the southernmost of the two slicks shown in Fig. 8 is the shadow to the west of the main slick, which appears to be a trail of oil left behind by the slick as it moves to the east under the influence of the wind 2 while the forecast/hindcast shown in Fig. 9 shows steady movement to the west. The Skiron forecast wind was indeed from the west during most of the 24 h preceding the satellite image, consistent 15 with this interpretation. However the wind was light, just 2-4 m s −1 , giving a wind drift of 6-12 cm s −1 , while the water currents were westerly with speed 14-17 cm s −1 , as shown by the black arrows in Fig. 9, dominating the effect of the wind. About 3 h before the observation, the wind backed to the south with speed 2 m s −1 then later to the east, increasing to 5 m s −1 by the end of the forecast period. The displacement of the slick 20 to the north is apparent from Fig. 9 as well as its faster drift to the west after the time of observation when the wind was no longer opposed to the currents.

Conclusions
ASAR satellite imagery has already become a useful means of detecting oil spills.
We compute forecasts of the future drift and spreading of such an observed slick as well as hindcasts of its past drift. In doing this MEDSLIK can make use of regional and downscaled forecasts of water currents, sea surface temperature that are available under the MyOcean project and winds and sea waves from SKIRON and CYCOFOS wave module respectively. The simplest technique that MEDSLIK can use is to compute the 5 future and past movements of the boundary points of the slick outline on the satellite image. This approach only permits advection of the slick to be modelled. A more sophisticated technique is to fill randomly the outline with a large number of parcels of oil and to follow the movements of the parcels. This method allows diffusion also to be included at least in the estimation of future movement. It also allows admittedly crude 10 estimation of the future and past states of the oil to be made. Oil slick forecasts are useful primarily in determining if any valuable resources may be threatened by the oil and in suggesting which methods of combating it may be appropriate. Slick hindcasts may help find the source of the oil, assisting in this way the response agencies to implement the EU Directive 2005/35. By plotting data on ship 15 movements obtained from local or regional AIS, MEDSLIK can assess whether any ship was in the same place at the same time as any hindcast position of the oil slick.
At the current time, satellite determination of the state of the oil slick is largely a matter of guessing median values for several parameters so any estimate must be treated with caution. Perhaps future techniques of making and analysing satellite images may 20 permit estimates of such parameters which would make predictions of its state more reliable. For example, an estimate of the age of the slick may be made from the degree of lateral spreading. In addition, the thickness of the oil may be in the future be estimated from by combining satellite SAR's with a multi-frequency (X-, C-, L-band) and polarization radar. Use of optical spectrometers may also be used in an integrated 25 effort.
The computation of the initial volume of oil in the slick as described in Sect. 3.3 does lead to discrepancies with the estimated volume of the observed slick, generally of the order of 10-15 %. This is largely due to the fact that the aging processes are nonlinear while the adjustment of initial volume is done linearly. The successful operational coupling of MEDSLIK with MyOcean regional/downscaled forecasting products and ASAR imageries serves as the pre-cursor service to implement the Mediterranean decision support system for marine safety (MEDESS4MS). During previous EU project, such as the multi-model MERSEA-IP in which MEDSLIK was one of the models used, inter-comparison exercises with in-situ drifters trajectories addressed the prediction accuracies for the trajectories of floating objects using regional, sub-regional and high resolution forecasting products (Brostrom et al., 2008(Brostrom et al., , 2010De Dominicis et al., 2010). In addition, during the ECOOP 10 project the comparison between MEDSLIK predictions for trajectories of in-situ drifters showed that the accuracy of predicted trajectory depends on the accuracy of the flow forecasts and their resolution (Zodiatis et al. 2010(Zodiatis et al. , 2012. The 24 h predictions of these virtual floating objects' trajectories in most cases, when the flow forecasts were in agreement with the in-situ flow, both the simulated and in-situ trajectories 15 were in agreement between them, particularly in terms of direction. Diversion of the virtual drifters trajectories was found only in locations where the flow forecast was not able to capture the in-situ one. Finally, the hindcast analysis of the biggest oil spill incident in the Eastern Mediterranean (Coppini et al., 2011), that of the Lebanese oil crisis in summer 2006, showed that the resolution of the flow forecasts determines the 20 accuracy of the floating oil slick trajectories.
To conclude, it worth mentioning that MEDSLIK fulfills the requirements of EMSA CSN regarding oil spill modeling and their coupling with the CleanSeaNet ASAR data (EMSA, 2009). The experience gained through the above mentioned implementations of MEDSLIK will help to deliver an integrated operational multi model oil spill prediction service in the Mediterranean, connected to existing satellite and in-situ monitoring platforms and AIS, using the regional and downscaled data from the MyOCEAN GMES Marine Core Service.  Tonani, M., Pinardi, N., Dobricic, S., Pujol, I., and Fratianni, C.: A high-resolution free-surface model of the Mediterranean Sea, Ocean Sci., 4, 1-14, doi:10.5194/os-4-1-2008. World Bank: Report No. 39787-LB: Republic of Lebanon, Economic Assessment of Environmental Degradation due to July 2006Hostilities, 2007 tional European Global Ocean Observing System for the Eastern Mediterranean Levantine Basin: The Cyprus Coastal Ocean Forecasting and Observing System, J. Mar. Technol. Soc., 37, 115-123, 2003a. Zodiatis, G., Lardner, R., Lascaratos, A., Georgiou, G., Korres, G., and Syrimis, M.: High resolution nested model for the Cyprus, NE Levantine Basin, eastern Mediterranean Sea: im-10 plementation and climatological runs, Ann. Geophys., 21, 221-236, doi:10.5194/angeo-21-221-2003, 2003b. Zodiatis, G., Lardner, R., Hayes, D., Solovyov, D., and The successful application of the Mediterranean oil spill model in assisting the decision makers during the oil pollution crisis of Lebanon in summer 2006, Rapp. Comm. Int. Mer Medit., 38, p. 214, 2007 OSD 9,2012 Predictions for oil slicks detected from satellite images    OSD 9,2012 Predictions for oil slicks detected from satellite images