The Mediterranean ocean forecasting system: first phase of implementation (1998–2001

Abstract. The Mediterranean Forecasting system Pilot Project has concluded its activities in 2001, achieving the following goals: 1. Realization of the first high-frequency (twice a month) Voluntary Observing Ship (VOS) system for the Mediterranean Sea with XBT profiles for the upper thermocline (0–700 m) and 12 n.m. along track nominal resolution; 2. Realization of the first Mediterranean Multidisciplinary Moored Array (M3A) system for the Near-Real-Time (NRT) acquisition of physical and biochemical observations. The actual observations consists of: air-sea interaction parameters, upper thermocline (0–500 m) temperature, salinity, oxygen and currents, euphotic zone (0–100 m) chlorophyll, nutrients, Photosinthetically Available Radiation (PAR) and turbidity; 3. Analysis and NRT dissemination of high quality along track Sea Level Anomaly (SLA), Sea Surface Temperature (SST) data from satellite sensors to be assimilated into the forecasting model; 4. Assembly and implementation of a multivariate Reduced Order Optimal Interpolation scheme (ROOI) for assimilation in NRT of all available data, in particular, SLA and VOS-XBT profiles; 5. Demonstration of the practical feasibility of NRT ten day forecasts at the Mediterranean basin scale with resolution of 0.125° in latitude and longitude. The analysis or nowcast is done once a week; 6. Development and implementation of nested regional (5 km) and shelf (2–3 km) models to simulate the seasonal variability. Four regional and nine shelf models were implemented successfully, nested within the forecasting model. The implementation exercise was carried out in different region/shelf dynamical regimes and it was demonstrated that one-way nesting is practical and accurate; 7. Validation and calibration of a complex ecosystem model in data reach shelf areas, to prepare for forecasting in a future phase. The same ecosystem model is capable of reproducing the major features of the primary producers’ carbon cycle in different regions and shelf areas. The model simulations were compared with the multidisciplinary M3A buoy observations and assimilation techniques were developed for the biochemical data. This paper overviews the methodological aspects of the research done, from the NRT observing system to the forecasting/modelling components and to the extensive validation/calibration experiments carried out with regional/shelf and ecosystem models. Key words. Oceanography: general (ocean prediction; instruments and techniques) Oceanography: physical (currents)


Introduction
The Mediterranean Forecasting System (Pinardi and Flemming, 1998) has established two major goals: Scientific: to explore, model and quantify the potential predictability of the ecosystem fluctuations at the level of primary producers from the overall basin scale to the coastal/shelf areas and for the time scales of weeks to months through the development and implementation of an automatic monitoring and a nowcasting/forecasting modelling system, the latter called the Mediterranean ocean Forecasting System (MFS) as a whole. Pre-operational: to demonstrate the feasibility of a Mediterranean basin operational system for predictions of currents and biochemical parameters in the overall basin and coastal/shelf areas and to develop interfaces to user communities for dissemination of forecast results.
These goals should be achieved in three phases: 1. First phase (1998)(1999)(2000)(2001): a pilot project for the implementation of the observing system backbone and demonstration of forecasting capabilities at basin scale; 2. Second phase (2002)(2003)(2004)(2005): consolidation and upgrade of the observing system for the physical components, extension of observations to biochemical variables, demonstration of sub-regional forecasting capabilities at the five-day range, three-dimensional ecosystem model implementation; 3. Third phase (2006)(2007)(2008): observing system verification and further extension toward operationality, shelf areas' primary producer forecasts, consolidation of products from forecasts.
This paper describes the methodological approach and the results of the pilot project phase, also kmown as the Mediterranean Forecasting System Pilot Project (MFSPP). The detailed results are described in several technical papers that follow this overview and they are referenced along the text.
Here we try to give the idea of the overall development of the "system", from observations to modelling, and the technological developments that are required to achieve forecasting capabilities.
The shelf areas of the Mediterranean Sea are outlined in Fig. 1. They are composed of extended shelves, such as the Adriatic Sea and the Tunisian shelf areas, and narrow shelf areas limited by steep continental slopes. Wind and thermohaline driven currents characterize the circulation in the deeper part of the basin: these currents are forcing the shelf regions laterally, making the problem of predicting the flow field in the shelf areas fully coupled with the open ocean.
In addition to this lateral coupling, the ocean-atmosphere interaction mechanisms are very intense both at the synoptic and the slowly varying atmospheric variability scales (Korres et al., 2000). The biochemical fluxes in the lower compartment of the trophic food web (bacteria and primary producers) are strongly affected by this varying physical environment (Legendre and Rassoulzadegan, 1995) at least as much as from the land derived inputs.
Thus the coastal environmental prediction problem at short and medium time scales requires the understanding and modelling of the large spatial and long time scales, as well as the local and short scales. A possible methodological approach is to "downscale" the large/long scale processes to the local/short scale, hypothesizing a conceptual model that parameterizes the effects of the large scale at the local level through nesting of different resolution observational networks and models. In Fig. 2 we illustrate the components of the two major system modules, e.g. observing and modelling, with data assimilation being part of modelling. The large/long and local/short system components are different but they need to be designed and developed together to maximize the benefits from each of them.
In general, we may say that the coastal environmental prediction problem can be connected to the design of an interdisciplinary observing system coupled with numerical predictive models of atmospheric, oceanic and marine biochemical state variables. The most important component of the predictive system is the assimilation engine that merges the observations with different model state variables, trying to reduce the uncertainties associated with the knowledge of the initial condition. However, the coastal environmental prediction problem has a multiplicity of system components that could limit the predictability time. They are: (1) the limited predictability of the atmospheric forcing directly influencing the coastal dynamics; (2) the lateral boundary condition uncertainties, considering both the open boundary conditions and the river runoff uncertainties, affecting the long-term memory and short-term variability of the system; (3) the adjustment processes to the downscaled large-scale initial field; (4) the initial condition specification for all the dynamical variables of interest; (5) the flow field nonlinearities.
We have subdivided the MFSPP results into four components: 1. the observing system; 2. the basin scale forecasting system; 3. the regional/shelf scale modelling/forecasting system; 4. the ecosystem modelling/forecasting system. While all the aspects of the first two components were implemented and demonstrated in Real Time or Near Real Time (NRT) during the Pilot Project, the third and fourth components consisted of validation/calibration experiments to serve as a scientific base for the demonstration of forecasting capabilities in the second and third phase of the MFS program.  Pinardi et al., 2002). The red ink words indicate the part of the system that has been implemented during MFSPP. In the following, Sect. 2 overviews the methodological approach and the results of the observing system component. Section 3 describes the basin scale forecasting and regional/shelf scale modelling is analyzed in Sect. 4. Section 5 overviews the findings of the ecosystem modelling and a discussion of future phases is offered in Sect. 6.

The VOS-XBT system
The MFSPP-VOS design was defined in 1999 (Manzella et al., 2001). It consisted of seven tracks to be covered twice a month for a period of nine months with a 12-nautical mile resolution and 460-760 m XBT probes. Although it would be preferable to collect data everywhere down to 760 m or more, the choice of different XBT types was dependent on the ship of opportunity speed. The actual ship tracks network is shown in Fig. 3.
During the first period of measurements (September -December 1999) the long trans-Mediterranean route was broken into two different pieces, one from Haifa to Messina and the other from Palermo to Gibraltar. However, in January 2000 the Ship Company changed the route, resulting in a direct Haifa -Gibraltar track. The track from Spain to North Africa was done using two different alternative transects (from Barcelona to Arzew or from Barcelona to Skikda) to match the biweekly frequency. The track from Sete (France) to Tunis (Tunisia) had a high variability, due to the bad weather conditions associated to the strong Mistral winds.
XBT data were transmitted in NRT via the ARGOS satellite transmission system (Du Penhoat, 1999), with a constant sub-sampling of each profile at 15 inflection points. This is different from the world ocean system that may change the number of sub-sampling points, depending on the number of inflection points of the profile. In addition, the full resolution profiles were also archived and sent to the ENEA collecting center (located in La Spezia, Italy) within one month. We call this second set of measurements delayed mode full resolution VOS-XBT profiles.
Since the forecasting system is intermittent and has a weekly frequency, it is useful to see how much data is delivered every week. In Fig. 4 we present the amount of data collected each week both as NRT decimated profiles and as delayed mode full resolution profiles. As we can notice, the NRT data were only 50-60% of the collected profiles, a low percentage in terms of costs of the VOS system. The AR-GOS transmission failures are large also due to problems in reconstructing the profile with a fixed number of decimation points (Manzella et al., 2003). The positive points about this VOS system are: 1. the data quality and amount is substantial and important for the assimilation/initialization of forecasts (see Raicich and Rampazzo, 2003;Demirov et al., 2003); 2. The NRT quality control management system can be done from a centralized data collection center that releases the data via internet technology.; 3. The XBT data themselves offer the opportunity to monitor the basin scale seasonal changes in an accurate way and contribute in a substantial way to the long-term monitoring of trends in the upper thermocline (Fusco et al., 2003).
The major drawback of the present system is the decimation and transmission system. The decimation does not allow for the proper reconstruction of the subsurface temperature structure due to the limited amount of decimation points (Manzella et al., 2003). The number of decimation points is limited by the low bit rate satellite transmission system. In addition, the software is not flexible enough to adapt to the changing stratification conditions in the Mediterranean Sea. Consequently, many profiles were not transmitted when the stratification was absent. Thus, at the end of the Project, a cellular phone data transmission protocol has been drawn up that offers the possibility of transmitting in NRT the full resolution profiles. More details about a new data management system for VOS XBT is described in Manzella et al. (2003).

The M3A buoy system
The new concept of the Mediterranean Multisensor Multidisciplinary Array is contained in two design features: 1) different moorings are allocated depending on the maintenance needs; 2) the different moorings communicate through a subsurface acoustic transmission system, in order to send data in real time to land. The precise system design is described in detail in Nittis et al. (2003). The salient features are shown in Fig. 5. They are: 1) a surface buoy, equipped with ARGOS transmitter, hosting sensors for meteo-marine (surface atmospheric temperature, wind speed and direction, atmospheric pressure, humidity, wave height and direction) and 1.5 m depth marine variables (temperature, conductivity, turbidity, dissolved oxygen and chlorophyll-a); 2) a central mooring line, called line 1, connected to the surface buoy hosting the hydroacoustic modem, four conductivity and temperature sensors located at fixed depth within the first 500 m of the water column; 3) a biochemical mooring line, called line 2, consisting of another underwater modem that transmits to line 1, and four packages of CTD, turbidity, chlorophyll-a, dissolved oxygen and PAR sensors distributed in the first 100 m, corresponding to the upper part of the euphotic layer, and one nutrient analyzer for nitrates; 4) a third mooring line, called line 3, consisting of an upward looking ADCP device, positioned at 500 m, sampling at several tens of meters interval the water column. The M3A is designed to match the basic needs for multidisciplinary long time series to validate/calibrate biogeochemical and hydrodynamic models. Last, but not least, the chlorophyll profile may serve to find feature models for the subsurface extrapolation of surface chlorophyll data from satellite colour sensors.
The M3A system was deployed in January 2000 northwest of the city of Heraklion, Crete, at a depth of 1030 m. The mooring lines 1 and 3 could be maintained every six months, only while two month period visits were necessary for mooring line 2. The M3A design allows for easy recovery of line 2. Bio-fouling was the primary problem in obtaining reliable data from line 2, especially from the turbidity and PAR sensors.
The data transmitted and/or collected by the buoy system are presented in Fig. 6 for some of the observed parameters. The line 1 sensors captured the seasonal thermocline formation and the intrusion events of salty water, as well as fresher waters in the upper thermocline. Thus, the sensor distribution can give a coarse resolution measure of the water mass transformation processes occurring at seasonal time scales. It is evident that the spring-summer subsurface chlorophyll maximum is well captured by the line 2 sensors. However, malfunctioning of the sensors due to bio-fouling can be detected in January at 115 m and in the period September 2000-January 2001 at 40 m. As shown later, the M3A data have shown to be a formidable observational data set for ecosystem model calibration and testing.

The NRT satellite data
The NRT satellite observations processed were: 1) Sea Level Anomalies (SLA) from ERS-2 and Topex/Poseidon (T/P) and 2) Sea Surface Temperature (SST). Both data sets have been used for the initialization of weekly forecasts. The ocean color data were also analyzed to produce maps of surface productivity and chlorophyll but they were done in delayed mode and they will not be described here.
Every Wednesday (T0) different SLA data sets are released, i.e. along-track T/P and ERS-2 SLA for the nominal period T0-22 days to T0-2 days and a map of SLA (signal and mapping error) centered at T0-7 days. Along-track SLA data are filtered, sub-sampled and corrected for alongtrack biases.
For verification purposes, the NRT maps from T/P and ERS-2 are intercompared with the Delayed Mode (DM) data that are more accurate, allowing for a reduction of orbit errors. The mean error associated with the NRT system could be estimated by calculating the average for all the grid points of the root mean square (rms) of the differences between NRT and DM maps. Table 1 gives the rms error for the combined T/P and ERS-2 maps (TPERS) and for the maps done only with the T/P. For combined TPERS maps the rms is about 3.76 cm and 3.21 cm for T/P maps. This error represents a non-negligible part of the mean variance of SLA (comprised between 6 to 9 cm). The main source of difference between the DM and the NRT data is the orbit error and the different amount of data used in the two analyses. In fact, the NRT system computes a combined map, every week, at time T0-7, using only three weeks of data (two before and one after T0). For the DM maps, four weeks are used instead, equally spaced around the central day of the analysis. The data assimilation system explained in the following section uses, in particular, the along-track SLA data to produce analyses every week utilizing the past one-two weeks of SLA data. In Fig. 7 we present a typical two-week coverage, which was used to produce an analysis. The coverage is high and almost evenly distributed over the Mediterranean Sea, except for the ERS-2 tracks that have a 35-day repeat cycle. Each two weeks, T/P covers 32 and ERS-2 43 tracks.
The NRT SST product is composed of weekly SST maps on the model grid. The analysis is made by the melding of two images collected at different centers from the NOAA-14 and 15 AVHRR sensors. One center is located in Lannion, France and the other in Rome, Italy. The data used come only from night-time passes, and an objective analysis technique interpolates on the cloudy pixels, using a monthly mean climatology as first guess. The data are produced every week on Thursday (T0+2), giving rise to the largest delay in the final release of the forecasting products. In the future, daily SST maps should require less pre-processing, so that the release of the data set will be closer to the start day of the forecast.
The characterization of the variability during the years 1999 and 2000 and the bivariate analysis of the SST and SLA signal is given in Buongiorno Nardelli et al. (2003). Here it will be sufficient to say that the Western Mediterranean was interested by the large variability of eddies in the Algerian basin and around the Balearic Islands, while the Ionian Sea was interested by the large event of cyclonic anomalies which probably started in 2000, and by large positive anomalies in the Pelops gyre area (southwestern corner of the Peloponnesian peninsula) and the Levantine basin.
3 Basin scale forecasting component Every deterministic forecasting system is composed of a numerical model and a data assimilation scheme that produces an optimal estimate of the initial condition from observations and past model output. During the Pilot Project, the basin scale forecasting model, called Ocean General Circulation Model-OGCM, was based upon a modified version of the Modular Ocean Model (Pacanowski et al., 1990). The resolution is 1/8×1/8 • of latitude and longitude and 31 levels in vertical. Such a numerical grid has been calibrated extensively in the past years to simulate the mesoscale, seasonal and interannual variability of the Mediterranean Sea (Pinardi and Masetti, 2000;Demirov and Pinardi, 2002). The model equations and parameterizations are summarized in Appendix A.
In order to make the forecast, we need to couple the OGCM with atmospheric forcing. For the Mediterranean, the asynchronous coupling of the atmospheric surface operational analyses with the OGCM was already established (Castellari et al., 1998(Castellari et al., , 2000. The atmospheric surface variables are input into bulk formulas that describe the radiative and turbulent heat fluxes at the air-sea interface. Such bulk formulas are chosen in order to maintain a negative net surface heat budget on a long-term average, thus allowing the production of deep waters at different rates, depending on the winter heat losses. The surface boundary condition for heat is then: The terms on the right-hand side of Eq. (1) are the net longwave emitted by the surface, Q B , the latent, LE, and the sensible heat flux, H . They depend upon the air temperature at 2 m, T a , the sea surface temperature computed by the model, T o , the total cloudiness, C, the relative humidity computed from dew point temperature at 2 m, rh, the 10 m wind velocity modulus, |ν w |. The different heat bulk expressions for the terms in Eq.
(1) were determined by Castellari et al. (2000) and they are: the Bignami et al. (1995) for Q B , the Gill (1982) for LE and Kondo (1975) for H fluxes. The important concept is that T o comes from the model integration itself, while all the other meteorological parameters are assigned from independent data. This surface flux formulation is called interactive, since the heat fluxes depend upon the state of the ocean directly. In the present phase, the surface meteorological fields are given on a grid of 0.5×0.5 • , every six hours but hopefully, the resolution will increase in the near future.
For the salinity flux, we consider: where S * is a reference salinity, taken to be the monthly mean climatology, z 1 is the first model layer thickness, i.e. 10 m, and T is a relaxation time taken to be 5 days. The climatological salinity field has been computed from the latest hydrological Medatlas data set (Fichaut et al., 1998) and it is presented in Brankart and Pinardi (2001).
The momentum flux is written as: where τ w is the wind stress calculated from the surface winds with the Hellerman and Rosenstein (1983) formula. The model has rigid lid fields, but it computes the sea level following Pinardi et al. (1995). This is purely a diagnostic quantity that does not contain external gravity waves and it is similar to the implicit free surface solutions at large scales, as shown by Dukowicz et al. (1994).
The Reduced-Order Optimal Interpolation scheme developed by De Mey and Benkiran (2002), named SOFA, was implemented, together with the OGCM. The scheme works intermittently, producing an analysis every week and it is described in detail in Demirov et al. (2003) and Sparnocchia et al. (2003). In summary, an optimal estimate of the initial condition or nowcast is computed every week using: where x a is the analysis or nowcast, x f is the first guess, that in our case is a model simulation, K is an approximate Kalman gain decomposed into vertical empirical orthogonal functions (EOFs) and horizontal correlations, y o is the observation and H is the observational operator that interpolates the model first guess into the observational position in (λ, θ) and transforms the model variables into the observed variable. Our implementation of the SOFA scheme is multivariate and x a and x f contain temperature, salinity and barotropic stream function. The Pilot Project produced an analysis and a 10-day forecast every week, starting from January 2000. The time line of events during each operational week is shown in Fig. 8. The analysis is produced at noontime on Tuesday of each week by a smoother-filter assimilation system for SLA and VOS-XBT, described in detail by Demirov et al. (2003). The smoother-filter procedure is as follows: alternatively, every week, the XBT and SLA are assimilated with a smoother mode type of OI, producing an analysis for the previous Tuesday with respect to the present week of forecast. After that, a filter mode OI scheme is carried out to assimilate the remaining data, in order to arrive at the present week Tuesday, the starting day of the forecast. For SST, a simpler assimilation scheme is used, since the data are weekly averages. The SST is used in a flux correction term added to Eq. (1) as follows: where Q T T * = 69 W m 2 • C and T * is the weekly mean SST. This corresponds to a 7-day relaxation constant for the 10 m deep surface layer of the model. This value is quite large if compared to global values computed by Oberhüber (1988), but it can be justified since SST is changing weekly.
The data sets arrive with a maximum delay of two days with respect to Tuesday of each week. Then, on Friday, two successive cycles are run: an analysis cycle that uses ECMWF surface analyses, and the forecast cycle that uses 10-day atmospheric surface fields from ECMWF. Considering that during the first day of the forecast, from Tuesday to Wednesday noontime, the model is adjusting dynamically the assimilated initial condition (the analysis), and the forecast is actually available two days after the usable start time of the forecast, i.e. noon on Friday.
The forecast skill is shown by means of the root mean square basin average error between forecast and analysis in Fig. 9 for the entire 2001 year. In Demirov et al. (2003), the same analysis is carried out for a six-month period in year 2000, during the Targeted Operational Period of the project. The error is also calculated as the difference between the analysis and the nowcast, the so-called persistence forecast error. This implies that the initial conditions are persisted and compared with the actual dynamical forecast. Following the meteorological experience, numerical predictions are thought to be valid only if the forecast error beats the persistence error.
The forecast error is always below 0.8 • C for all the different weekly 10-day forecasts. During most of this year, only SLA was assimilated, and due to the intermittent assimilation scheme and the two weeks assimilation cycle (Demirov et al., 2003), the rms error growth has a biweekly frequency, as shown in Fig. 9 for the 30 m level. The forecast error is lower than the persistence error at all depths, except for a few cases during late summer, where possibly the temperature corrections due to SLA data are not completely correct. It is interesting to notice that during this period the 5-m forecast error is really an estimate of the quality of the numerical forecasting system, since the analysis uses relaxation toward SST satellite data during the whole week of intercomparison. Then, we conclude that the system is clearly accurate and better than persistence at the surface.

The regional/shelf scale modelling/forecasting component
The MFS Pilot Project developed and implemented four regional models at 5-km resolution and nine shelf models at 1.5-3 km resolution nested within the forecasting OGCM. The areas of implementation are shown in Fig. 10. This exercise was preparatory toward the actual forecasting with regional and shelf models that will occur in the next phase. All models, comprehensive of the OGCM, were run with perpetual year forcing. The OGCM was run for seven years while the regional/shelf simulations were three years long, reaching a repeating seasonal cycle and being continuously fed by the coarser model fields. The perpetual year forcing experiments are a well-known practice in physical oceanography of the Mediterranean Sea and they allow one to produce dynamically consistent fields at seasonal time scales that usually compare well with climatological observations (Roussenov et al., 1995;Wu and Haines, 1996;Pinardi and Masetti, 2000). Since the Mediterranean has a large seasonal cycle, this is the first test of the model accuracy and potential for predictions.
This mode of operations consists of forcing the model with climatological monthly mean heat, salt and momentum fluxes for several years. The model reaches then a statistical steady state, composed of a repeating seasonal cycle in temperature, salinity, flow across basin straits, etc. The perpetual year forcing used here is described in Korres and Lascaratos (2003). The heat fluxes were computed using the bulk formulas described in the previous section, but also using the prescribed and observed SST from the global analyses of Reynolds and Smith (1994). In addition, the meteorological surface fields come from the ECMWF re-analysis data set (Gibson, 1997), which covers the period January 1979 to December 1993. A monthly mean net heat flux, defined as: is then computed and used as a forcing of all models. Such heat flux formulation is not interactive and thus it is necessary to correct for the mismatch between the Reynolds SST and the actual SST of the model. Both the OGCM and the regional/shelf models needed then to implement a heat flux The error is computed as the root mean square difference between forecast and analysis ( red curve, forecast error) and the analysis minus the initial conditions (blue curve, persistence error). The error is computed for the 52-weekly forecasts done in the year and it is computed over a one-week interval. correction as follows: where [C 1 ] = W m 2 • C , T o is the model SST and T * is a reference climatological monthly mean temperature field.
While the OGCM used only Eq. (2), the regional/shelf models used a salt flux boundary condition such as: where [C 2 ] = m s , S o is the model surface instantaneous field and S * is the monthly mean reference surface salinity climatological field from Brankart and Pinardi (2001). Here , where E is the same as in Eq. (4), the precipitation, P , is taken from the climatological values of Fig. 10. The regional/shelf models implemented in MFSPP: 4 regional models at 5 km resolution and 9 shelf models at 1.5-2 km resolution. Jaeger (1976) and the river runoff, R, is specified by monthly mean values wherever possible.
Values of C 1 ranged from 10 to 40 W m 2 • C in different model implementations while for C 2 the value was 0.7-1 m/s everywhere.
One way to decrease the value of C 1 in the regional models is to use the monthly mean heat fluxes directly from the OGCM. The latter are changed by the heat flux correction so that the overall heat budget of the basin becomes negative. This is shown in Fig. 11, where we show the monthly mean heat flux and the heat flux correction for the seventh year of perpetual forcing of the OGCM. The annual mean heat flux without flux correction is +18 W/m 2 , the annual mean heat flux correction term is -23 W/m 2 , producing a negative annual mean budget, i.e. -5 W/m 2 , which is consistent with the known negative heat budget of the Mediterranean Sea (Castellari et al., 1998). This procedure allowed the regional models to almost halve the value of C 1 . This is due to the fact that, using the corrected OGCM heat fluxes, the heat flux correction in the regional models will correct only for a local mismatch between lateral heat transport at the open boundaries and the surface heat flux over the specific region. By the same reasoning, the regional models gave the shelf models their corrected surface heat budget and this also decreased the value of C 1 in the shelf models. This kind of heat flux nesting is novel and was experimented successfully in all the MFS models.
A complete set of open boundary specifications for nesting the models at different resolution and with different physics was developed. The regional and shelf models were all derived from the Princeton Ocean Model (Blumberg and Mellor, 1983), except the one that uses the OPA code (Madec et al., 1998) in the northwestern Mediterranean. Thus, the nesting had to consider in general rigid lid fields as boundary conditions for free surface models.
The coupling between the OGCM and the regional/shelf models is realized by specifying the total velocity field interpolated over the finer model grid, imposing the interpolation constraint (see Appendix B) developed during the project. This constraint allows for the total transport to be preserved after the interpolation is carried out from the coarse to the fine grid. In addition, between the OGCM and the regional models we impose: 1. a zero-gradient boundary condition for the free surface (i.e. the free surface is not nested); 2. the upstream advection of T , S is used, such as: used three alternative different conditions. They are: where is equal to ±1 depending on the position of the open boundary: = 1 for an eastern and northern boundary, = −1 for a western and southern boundary. Here H high = h high + η high and H coarse = h coarse + η coarse .
Condition (a) was used both for rigid lid to free surface and free to free surface nesting. Condition (c) was only used for free to free surface nesting. Condition (b) was only used for rigid lid to free surface nesting. These conditions should be enforced only at outflow points. The determination of outflow/inflow points is done on the basis of the U norm coarse direction with respect to the open boundary only; 4. For the tangential barotropic component the coarse and high-resolution fields were simply equalized.
The lateral boundary fields were 10-day averages from the OGCM and regional model simulations.
The regional and shelf models were capable of reproducing the major dynamical features of the known climatological circulation. Several experiments were carried out in each region by varying the initial condition from the coarser model, the nesting lateral boundary conditions time frequency and kind, the surface flux corrections and the heat penetration laws (Appendix A). In Fig. 12 we show a typical example of the intercomparison between monthly mean fields of the OGCM and the nested regional model . The high resolution regional model develops more realistic features such as the Mersa-Matruh and the Shikmona gyre systems (southern and southeastern Levantine area).
Each region/shelf decided the "best" model implementation toward forecasting in the future phases. Results for the different model implementations are discussed in Korres and Lascaratos (2003) for the Levantine and Aegean Sea, Zodiatis et al. (2003) for the Cyprus area, Kourafalou and Barbopoulos (2003) for the Northern Aegean, Triantafyllou et al. (2003a) for the Cretan Sea, Echevin et al. (2003) for the northwestern Mediterranean, Drago et al. (2003) for the Maltese shelf, Sorgente et al. (2003) for the Sicilian Strait, Zavatarelli et al. (2003) for the Adriatic and Northern Adriatic, Brenner (2003) for the Israelis shelf area.

The ecosystem modelling/forecasting system
This part of the MFSPP was dedicated to the calibration/validation of a one-dimensional ecosystem model, formulated on the basis of the ERSEM code (European Regional Seas Ecosystem Model, Baretta et al., 1995) for the biochemical components coupled with a one-dimensional version of POM. The physics considers only vertical diffusion and Coriolis acceleration as basic hydrodynamic processes to simulate the vertical structures of the water column. The vertical diffusion coefficients for temperature and salinity are described by the Mellor and Yamada (1983) turbulence closure sub-model. The physical fields of temperature and salinity are simulated under the vertical boundary conditions of heat and salt fluxes explained in Eqs. (2) and (5). The biological components of the model are grouped into "functional groups", according to their trophic level, and subdivided according to feeding method or size, as indicated in Fig. 13. ERSEM describes the cycling of carbon, nitrogen, phosphorus and silicate through the pelagic food web.
The MFSPP version of this model was modified to consider a variable carbon to chlorophyll ratio . This is due to the fact that this ratio varies between a value of 10 to 200, and it is then difficult to interpret both in situ and surface chlorophyll data from satellites in terms of carbon biomass. Additionally, the light limitation and adaptation of phytoplankton can be now parameterized in terms of chlorophyll variability, which is a natural state variable, instead of optimal light irradiance, as it was used before.
The one-dimensional ecosystem model was validated against historical data at different sites in deep waters around the Mediterranean. The aim was to try to identify site specific parameters in the ecosystem model and to see if the model could shift by itself between different ecosystem regimes. One of the results obtained is shown in Fig. 14. Turley et al. (2000) found a highly significant positive correlation between primary production and bacterial production in both the western and eastern basins of the Mediterranean, indicating that the primary production is a significant source in open ocean waters of DOC for bacterial production. The onedimensional ERSEM model, implemented in seven different sites across the Mediterranean , is capable of capturing this feature, as shown in Fig. 14, where model and data are intercompared. Although simulated production rates are very sensitive to the initial condition upon which those simulations are based, a generic model parameterization is sufficient to represent these differences and produces results in line with the observations. Four different regional model implementations were carried out with approximately the same one-dimensional ERSEM and POM coupled model. The model was implemented with a perpetual year forcing, previously explained for the hydrodynamics. Such forcing allows the capturing of seasonal biomass variability of primary producers in different shelf areas and the comparison with regional/shelf historical data sets. The results of such validation at seasonal time scales are offered in Vichi et al. (2003) for the northern Adriatic Sea, Triantafyllou et al. (2003b) for the Cretan Sea.
In preparation for forecasting, the one-dimensional model was also run with high-frequency atmospheric forcing and data assimilation techniques were developed. The results of high-frequency atmospheric forcing and intercomparison with M3A data are described in Triantafyllou et al. (2003b) and Siddorn and Allen (2003). In order to show the potential, Fig. 15 represents the intercomparison of one such simulation for the M3A station for the period February to July 2000. The model exhibits a good fit with data, except for the period between days 100 and 140, which is coincident with the onset of stratification and the formation of a deep chlorophyll maxima.
The model simulation in the Cretan Sea, at the M3A site, and the newly collected data at high frequency, provided the first calibration/validation exercise for two different data assimilation schemes. The first, based upon the Extended Kalman Filter  and the second based upon the Singular Evolutive Extended Kalman filter (Hoteit et al., 2003). For example, the assimilation of this data, using the EnKF method, results in a marked improvement in the ability of ERSEM to hindcast chlorophyll (Fig. 14). The predictability window of the EnKF appears to be at least 2 days, which  Fig. 13. The trophic structure of the ERSEM model (Baretta et al., 1995). indicates that the methodology might be suitable for future operational data assimilation systems using more complex three-dimensional models.

Discussion and future outlook
In this paper we overviewed the first phase of the Mediterranean Forecasting System implementation results, the development of technologies to demonstrate the practical feasibility of forecasting at basin scales and the preparation of the regional/shelf scale and ecosystem models for the future phases. First of all, the practical feasibility of 10-day basin scale forecasts has been demonstrated. In addition, it has been proved that the forecast can be released with a two-day delay with respect to the actual start day of the forecast, with near- real-time exchange of atmospheric, satellite and in situ data. The exchange is based upon Internet technology, satellite transmission from the ship to land and near-real-time quality control procedures. Two major drawbacks were found: the decimated profiles from the VOS-XBT have a cost/benefit ratio that is too high, due to a missing temperature signal in the deeper layers and the failures of the satellite transmission system. The second is the lengthy computation of weekly SST that should be substituted by daily, almost unprocessed data. The remaining part of the observing system is proven to be reliable and robust, capable of supporting actual forecasting at the weekly time scales. Moreover, new nesting techniques have been developed and implemented in several regional and shelf dynamical regimes, between different resolution and different physics models.
Finally, a generic ecosystem model, both for open ocean and shelf area conditions, has been tested and calibrated. It is shown that such a functional group, and a biomass-based model can reproduce the biomass fluctuations in the euphotic zone, both at seasonal and synoptic time scales.
The second phase of the MFS will start rapidly and is based upon the conclusions from the first phase; it will aim at: 1. Improve and expand the existing near-real-time largescale monitoring system. In particular, it will develop: 1) new XBT and fluorimeter sensors; 2) a multiple launcher for XBT; 3) a multiparametric non-expandable sensor, 4) Jason-1 and Envisat SLA NRT data analysis; 5) daily SST products; 6) a network of M3A stations with new sensors and more reliable data transmission; 2. Add new observing system components in terms of automated technology. This means: 1) a subsurface vertical profiling floats program, i.e. the Mediterranean component of global ARGO (http://www.argo.ucsd.edu); 2) a glider experiment on a VOS-like track across the basin; 3. Improve the 10-day basin scale forecast system and demonstrate the feasibility of near-real-time three to five-day forecasts in different nested regional areas; 4. Develop the asynchronous ocean-atmosphere coupling with high resolution atmospheric forcing over regional areas; 5. Implement the three-dimensional ecosystem models coupled to the forecasting system for future predictions of biochemical elements variability; 6. Consolidate the dissemination of forecasts to a wide user community and develop applications with endusers.
This second phase will start in 2002 and last until 2005. The major scientific issues to study are the combination of all different measurements to improve the global Mediterranean forecasting system, the establishment of useful limits for the forecast skill scores and the understanding of the oceanatmosphere coupling for short and high resolution ocean spatial scales. The second major objective is to demonstrate the practical feasibility and accuracy of real-time nested regional forecasts. Last, but not least, the second phase will prepare for the three-dimensional biochemical flux modelling and assimilation in coastal areas. outflow/inflow system at Gibraltar. They are different from zero only from the middle of the Gibraltar Strait outwards into the Atlantic Ocean, in a box approximately 3 × 3 • of latitude and longitude. The formula used is: λ(x, y, z) = α + αθ x − x 0 e −αx 2 − 1 , where α −1 = 1 day, θ is the step function that becomes unity nine grid points (x 0 ) west of Gibraltar and α − 1 2 = 30 km.

Appendix B Interpolation constraint and correction
We define the Interpolation Constraint (IC) as: where U total coarse is the coarse grid total velocity field, U total int is the interpolated total velocity field, (h coarse , η coarse ) and (h high , η high ) are the bathymetry and the free surface of the coarse and high resolution grids, respectively. The aim is to find a U total high such that (B1) is preserved. In order to do so, 3 steps are necessary: Step 1: Calculate on the coarse grid: which is a transport in m 3 /s. Obviously, if the coarse model is rigid lid, η coarse = 0.
Step 2: Interpolate the into the finer grid and calculate: Step 3: Compute the final velocity field using : U total high (x, y, z, t) = U total int (x, y, z, t) − T corr S F (x, y, z) , where T corr = T int − T coarse and During the Pilot Project, most of the regional/shelf models have applied F = 1.