Assessing the impact of brightness temperature simulation adjustment conditions in correcting Metop-A SST over the Mediterranean Sea

Remote Sensing of Environment 146 (2014) 214–233

Contents lists available at ScienceDirect

Remote Sensing of Environment journal homepage: www.elsevier.com/locate/rse

Assessing the impact of brightness temperature simulation adjustment conditions in correcting Metop-A SST over the Mediterranean Sea Igor Tomažić a,⁎, Pierre Le Borgne b, Hervé Roquet b a b

AGO-GHER, University of Liège, Allée du Six Aout, 17, Sart Tilman, Liège, 4000, Belgium Centre de Météorologie Spatiale, Avenue de Lorraine, 22307 Lannion, France

a r t i c l e

i n f o

Article history: Received 21 February 2013 Received in revised form 1 October 2013 Accepted 6 October 2013 Available online 12 November 2013 Keywords: Sea surface temperature (SST) Infra-red remote sensing Metop-A AVHRR Retrieval errors

a b s t r a c t Multispectral sea surface temperature (SST) algorithms applied to infrared (IR) radiometer data exhibit regional biases due to the intrinsic inability of the SST algorithm to cope with the vast range of atmospheric types, mainly influenced by water vapor and temperature profiles. Deriving a SST correction from simulated brightness temperatures (BTs), obtained by applying a Radiative Transfer Model (RTM) to Numerical Weather Prediction (NWP) atmospheric profiles and first guess SST, is one of the solutions to reduce regional biases. This solution is envisaged in the particular case of Metop-A Advanced Very High Resolution Radiometer (AVHRR) derived SST. Simulated BTs show errors, linked to RTM, atmospheric profiles or guess field errors. We investigated the conditions of adjusting simulated to observed BTs in the particular case of the Mediterranean Sea over almost one year. Our study led to define optimal spatio/temporal averaging parameters of the simulation observation differences, both during day and night, summer and colder season and for two simulation modes: operational (with reduced vertical resolution – 15 levels – NWP atmospheric profiles and two days old analysis used as first guess SST) and delayed (full vertical resolution – 91 levels – and concurrent analysis used as first guess SST). Each BT adjustment has been evaluated by comparing the corresponding corrected AVHRR SST to the AATSR SST that we adopted as validation reference. We obtained an optimized result across all defined conditions and modes for a spatial smoothing of 15 deg and a temporal averaging between 3 and 5 days. Specifically, analyses based on 10 day averages showed that a standard deviation based criterion favors spatial smoothing above 10 deg for all temporal averaging, while a bias based criterion favors shorter temporal averaging during daytime (b 5 days) and higher spatial smoothing (N 10 deg) for nighttime. This study has shown also the impact of diurnal warming both in deriving BT adjustment and in validation results. © 2013 Elsevier Inc. All rights reserved.

1. Introduction Sea surface temperature (SST), defined as one of essential climate variables, is a key geophysical parameter at the ocean–atmosphere boundary that is routinely retrieved using infrared (IR) radiometers aboard satellite platforms for more than 30 years. Nevertheless, satellite derived SST data still exhibit regional biases mainly due to the intrinsic inability of multispectral SST algorithms (McClain, Pichel, & Walton, 1985; Walton, Pichel, Sapper, & May, 1998) applied to infrared radiometer data. SST from Metop-A/Advanced Very High Resolution Radiometer (AVHRR) has been derived operationally by the European Organisation for the Exploitation of Meteorological Satellites (EUMETSAT) Ocean and Sea Ice Satellite Application Facility (OSI-SAF) at Meteo-France/Centre de Meteorologie Spatiale (CMS) in Lannion (EUMETSAT, 2010) by applying non linear (NL) split window algorithm (Walton et al., 1998) to IR window channels 11 and 12 μm during day ⁎ Corresponding author. Tel.: +32 43663351. E-mail address: [email protected] (I. Tomažić). 0034-4257/$ – see front matter © 2013 Elsevier Inc. All rights reserved. http://dx.doi.org/10.1016/j.rse.2013.10.007

(Eq. (1)) and by applying triple (T37) window algorithm (McClain et al., 1985) to IR channels 3.7, 11 and 12 μm during night (Eq. (2)): NL : STT ¼ aT10:8 þ ðbTCLI þ cSϑ ÞðT10:8 −T12:0 Þ þ d þ eSϑ :

ð1Þ

T37 : SST ¼ ða þ bSϑ ÞT3:7 þ ðc þ dSϑ ÞðT10:8 −T12:0 Þ þ e þ fSϑ

ð2Þ

where T3.7, T10.8, T12.0 are the brightness temperatures at 3.7, 10.8 and 12.0 μm, respectively; Sθ = sec(θ) − 1, θ is the satellite zenith angle, TCLI is the mean climatological value (Operational Sea Surface Temperature and Sea Ice Analysis — OSTIA, Donlon et al., 2012 in our analysis) and the algorithm coefficients have been derived from a simulated brightness temperature data base (François, Brisson, Le Borgne, & Marsouin, 2002). Real time simulated Brightness Temperatures (BTs) are increasingly used in operational Sea Surface Temperature (SST) calculations, either in Optimal Estimation (OE) methods (Merchant, Harris et al., 2009; Merchant, Le Borgne, Marsouin, & Roquet, 2008; Merchant, Le Borgne, Roquet, & Marsouin, 2009) or coefficient based methods (Le Borgne,

I. Tomažić et al. / Remote Sensing of Environment 146 (2014) 214–233

Roquet, & Merchant, 2011; Petrenko, Ignatov, Shabanov, & Kihai, 2011). BTs are simulated in the adequate infrared window channels by applying a Radiative Transfer Model (RTM) to Numerical Weather Prediction (NWP) atmospheric profiles, using a guess SST field (guess SST) as surface temperature. These simulations will be referred to as simBTi and the corresponding BTs, measured by the radiometer onboard the satellite as obsBTi. Coefficient based or OE methods result in expressing the final SST as: SST ¼ guessSST þ

X

ai ðobsBT i −simBT i Þ:

ð3Þ

In the case of coefficient based method envisaged here (Le Borgne et al., 2011), a correction to a multispectral algorithms (Eqs. (1) and (2)) used to operationally derive SST is determined by applying the coefficients of this algorithm to simulations in order to produce a simulation derived SST. This “simulated” SST is compared to the “true” guess SST (a priori estimate of the SST) and the difference is used as a correction term to the operational SST derived with this algorithm. In this case, ai in Eq. (3) are the coefficients of the multispectral algorithms (Eqs. (1) and (2)) to which the correction is applied. Simulations differ from observations, because: - There are differences between the surface guess field and the actual SST field (guess errors) - There are differences between the model atmospheric profiles and the actual profiles (NWP errors) - RTM uses erroneous filter functions (filter errors) - RTM may be inaccurate (model errors) - Profile sampling induces errors (profile sampling errors). For instance the original 91 levels of the European Centre for MediumRange Weather Forecasts (ECMWF) outputs are sampled onto 15 pressure levels for some operational applications at Meteo-France. In the case of a coefficient based method applied here, guess errors are accounted for by Eq. (3), but the other sources of errors should be corrected prior to using simulations in the SST correction scheme. In practice adjustments are made by deriving empirical adjustment values from comparisons between simulations and observations. Two approaches have been used: analytical expressions of the simulation errors as a function of water vapor and satellite zenith angle (Merchant et al., 2008; Petrenko et al., 2011) or a dynamic determination of the error geographical distribution (Le Borgne et al., 2011; Merchant, Harris, et al., 2009). In the framework of the EUMETSAT, the OSI-SAF at Météo-France/ CMS in Lannion has developed a geostationary satellite SST chain using simulated BT which became operational in August 2011 (EUMETSAT, 2011). CMS has been running a prototype chain since November 2011 for testing the same approach on Metop-A/AVHRR data. This prototype aims to test methods that will be applied at full resolution in the Suomi National Polar-orbiting partnership (NPP) Visible Infrared Imaging Radiometer Suite (VIIRS) and Metop-A/ AVHRR) future SST processing chains. The prototype processing is done in near-real time, to be as close as possible to operational conditions and uses OSTIA as guess SST, ECMWF forecast profiles and the Radiative Transfer for TIROS Operational Vertical Sounder (RTTOV) model version 10 (Hocking et al., 2011). Considering the CPU time needed to apply the SST correction, main extra CPU time is used for calculating BT simulations. This has already proved to be feasible at full resolution in operational Spinning Enhanced Visible and Infrared Imager (SEVIRI) processing (EUMETSAT, 2011). It is likely that some subsampling will be necessary in the case of Metop-A/AVHRR operational processing. Differences between simulated and observed BTs are routinely monitored at CMS. These routine comparisons are done for pixels of quality level 3 or above, based on the GHRSST quality referencing system (Donlon et al., 2007; EUMETSAT, 2010) (see Data

215

section), i.e. showing no or little cloud contamination. Fig. 1 shows the daily global Metop-A/AVHRR comparison results for the year 2012. Note that a similar monitoring is done at NOAA (Liang & Ignatov, 2011; Liang, Ignatov, & Kihai, 2009), and their results are displayed through http://www.star.nesdis.noaa.gov/sod/sst/micros/. The mean daily biases and standard deviation values are quite stable over the 9 month period considered. The step observed for Metop-A/ AVHRR 3.7 μm simulations on the 2nd of April 2012 (red circle on Fig. 1) is due to a RTTOV version change (from 10.1 to 10.2). On the contrary daily differences are highly variable in space, with a distribution of positive and negative patches a priori which are difficult to interpret. This is shown by Fig. 2a which presents an example of daily simulation minus observation differences at 10.8μm on the 1st of March 2012 at 00:00 UTC, as produced by the prototype monitoring system. 10.8μm has been selected as representing an intermediate case between 3.7 and 12 μm. The 1st of March has been chosen randomly. Timestamp of 00:00 UTC denotes that all data between 18:00 and 06:00 UTC are aggregated in one file and all data between 06:00 and 18:00 UTC are aggregated to the file with timestamp of 12:00 UTC (for more details see Data section and EUMETSAT, 2010). To better understand the simulation minus observation differences, data filtering was reinforced to avoid as much as possible contamination by clouds or by residual diurnal warming. This was done by using SST quality level 5 (see Data section) only, and by excluding data for which the absolute difference between BT simulations and observations is above 1.5 K and wind is below 2 m/s. These filtered data have been averaged over one month (March 2012) and the spatial distribution of simulation minus observation difference shows clearer patterns than the unfiltered data (see Fig. 2b). In a consequence we decided to derive adjustment values from the geographical distribution of the errors averaged over space and time, as previously done in Merchant et al. (2008). The temporal and spatial scales of this averaging have been determined up to now through brief preliminary studies. For the Metop-A/AVHRR prototype chain we adopted a time averaging of 10 days and a space averaging of 10° but they were not studied in detail. In this article we investigate the practical issues of deriving an adjustment to the simulations by comparing simulations to observations. A previous study (Tomažić, 2012) demonstrated, over the Adriatic region, that using a fine scale regional model in place of a global one (ECMWF) with a coarser resolution does not impact the simulations results in the SST retrieval context. Similarly, the same study (Tomažić, 2012) showed that the use of a 1 km resolution analysis brought no significant improvement compared to using the OSTIA analysis. In consequence in the present study we have used the OSTIA analysis associated to ECMWF profiles. They have been used in the CMS prototype chain (an operational like environment) and in a delayed mode to use analysis instead of forecast and 91 ECMWF profile levels instead of 15. We selected the Mediterranean Sea not to redo simulations over the entire domain and because this area offers a good clear sky coverage with extreme conditions of diurnal warming and a good range of atmospheric conditions. This choice was also motivated by the first regional application of the method to VIIRS data acquired in direct readout mode at CMS and covering the North-East Atlantic and the Mediterranean Sea. In what follows the two simulation sets will be referred to as “operational” and “delayed mode” simulations. Various time and space adjustment scales have been tested for each mode. Each time and space combination has been evaluated by validating the SST bias correction deduced from the corresponding simulations by comparison to the Advanced Along-Track Scanning Radiometer (AATSR) data, which were in practice the only validation data available at large scale over the Mediterranean Sea. The next section presents the data used and the main validation results by comparison with the AATSR data obtained over the Mediterranean Sea for the selected period. The brightness temperature calculations are described in Section 3. This section presents also the simulations minus observation statistics recorded for the operational and delayed mode data. Section 4 describes

216

I. Tomažić et al. / Remote Sensing of Environment 146 (2014) 214–233

a

b

c

Fig. 1. Mean daily differences and standard deviation between simulated and observed Metop-A/AVHRR brightness temperatures for year 2012, a) 3.7; b) 10.8 and c) 12.0 μm. Red circle denotes impact on BT simulations for channel 3.7 μm due to change in RTTOV version (from 10.1 to 10.2).

the method we adopted to test the various adjustment options. In Section 5 we describe the results obtained for averaging time scales from 1 to 30 days and space scales from 0.5 to 50° in latitude–longitude degrees. We present our conclusions in Section 6.

2. Data All data span the time range between June 2011 and March 2012. Three different sets of data have been used: Metop-A 12 hourly L3C files, AATSR L3C file and OSTIA. We have chosen a 0.05° regular grid over the Mediterranean Sea domain (see Fig. 3) as a common grid for all used datasets. All SST products contain information about the quality of each pixel (quality flags) that generally follows the principles defined by the GHRSST project (Donlon et al., 2007).

2.1. Data description 2.1.1. AATSR The AATSR sensor onboard ENVISAT has an advanced calibration system and dual-view capability which leads to an improved atmospheric correction procedure (Zavody, Mutlow, & Llewellyn-Jones, 1995) and to reduced sensitivity to stratospheric aerosol (Merchant & Harris, 1999). Accuracy of the AATSR SST retrieval has been confirmed in numerous studies that used in situ measurements for validation (Corlett et al., 2006; O'Carroll, Eyre, & Saunders, 2008; O'Carroll, Watts, Horrocks, Saunders, & Rayner, 2006) or ship based radiometry in dedicated validation campaigns to confirm absolute accuracy stability (Wimmer, Robinson, & Donlon, 2012). Due to the lack of in situ data in Mediterranean Sea and due to the proven quality of AATSR data we used AATSR data as a validation reference.

Fig. 2. a) Simulated − Observed 10.8 μm BT at 00:00 UTC on the 1st of March 2012.; b) Simulated − Observed 10.8 μm BT for nighttime cases averaged over March 2012.

I. Tomažić et al. / Remote Sensing of Environment 146 (2014) 214–233

217

Fig. 3. 10 day (17 March till 26 March 2011) averaged SST (K) difference between Metop-A and AATSR for nighttime over the Mediterranean domain.

AATSR L2P files obtained from ESA in GHRSST format were remapped (to the nearest neighbor) to the common 0.05 deg grid. From all remapped orbits, L3C files were created centered at 00:00 and 12:00 UTC time, with a ±6 hour time window centered at those times. Due to the AATSR narrow swath (512 km) there were not overlapping pixels, and the only criteria in deriving L3C file were actual the existence of the data in the region of interest. Orbit times over the area vary from East to West from 07:30 and 10:30 UTC by day and 19:30 to 22:30 UTC by night. This corresponds to about 10.5 ± 0.1 h Local Solar Time (LST) by day and 21.5 ± 0.5 h LST by night. Within the L2P files there is information about the pixel quality (quality flags) represented by a scale from 0 to 5 (Corlett, 2008) and to ensure highest possible validation quality we used AATSR data with the quality level 5 (“excellent”). 2.1.2. Metop-A/AVHRR Metop-A/AVHRR SST products are derived operationally at OSI-SAF at full resolution over the whole globe. A non linear split window Eq. (1) is used by day (solar zenith angle smaller than 90°) and a triple window Eq. (2) is used by night (solar zenith angles larger than 110°). In the intermediate solar zenith angles (90 to 110°) a weighted mean of daytime and nighttime algorithm is used to insure a smooth transition. Full resolution granules are delivered as operational products and they also contain information about the pixel quality (quality flag) represented by a scale of six levels (0 = “unprocessed”; 1 = “cloudy”; 2 = “bad”; 3 = “suspect”; 4 = “acceptable”; 5 = “excellent”). The detailed explanation is given in EUMETSAT (2010) report. Granules are also aggregated twice daily from 18:00 to 06:00 UTC (centered on 00:00 UTC) and from 06:00 to 18:00 UTC (centered on 12:00 UTC) and delivered on a global grid at 0.05° resolution (EUMETSAT, 2010). In this processing, Metop-A/AVHRR SSTs full resolution pixels allocated to a grid point are spatially averaged providing they have the same quality level. These gridded products are stored internally at CMS as “workfiles” with observed BTs and ECMWF outputs such as wind field and total water vapor content. The Mediterranean sector of these workfiles constitutes the basic material of this study. Orbit times over the area vary from East to West from 07:00 and 10:00 UTC by day and 19:00 to 22:00 UTC by night. This corresponds to about 10 ± 1 h LST by day and 21 ± 0.3 h LST by night. Validation results are routinely made available on the EUMETSAT OSISAF website (http://www.osi-saf.org). They are also available through the NOAA SST QUAlity Monitor (SQUAM, http://www.star.nesdis.noaa.

gov/sod/sst/squam) or published in research articles (e.g. O'Carroll, August, Le Borgne, & Marsouin, 2012). The CMS prototype started on the 1st November 2011 and the AATSR stopped in March 2012. To include a summer season in the studied dataset the prototype has been run from archives for the June to August 2011 period. 2.1.3. OSTIA OSTIA SST data (Donlon et al., 2012), obtained from NASA/PO.DAAC, were used as guess SST field in simulating BTs. They represent “foundation” SST, free from diurnal warming with a spatial resolution of ~6 km (1/20°). In a previous study (Tomažić, 2012) over Adriatic Sea it was shown that using higher resolution (1 km) guess SST field (CNR UHR, (Buongiorno Nardelli, Tronconi, Pisano, & Santoleri, 2013)) did not have significant impact on the improvement of the simulation qualities. A more elaborate global study (Saha, Ignatov, Liang, & Dash, 2012) showed that OSTIA fields perform very well for the period when ENVISAT/AATSR were available (8th April 2012) and degraded afterwards. The availability of OSTIA data in NRT applications is an important issue, as it is the case for current Metop-A operational processing at Meteo-France/CMS. The closest OSTIA field in time of the full coverage of Metop-A satellite overpass is available only with two days of delay, which creates a larger guess error compared to using OSTIA from the same day. 2.2. Data intercomparison Metop-A/AVHRR and OSTIA SST data sets have been compared to AATSR SST fields over the Mediterranean for the June 2011 to March 2012 period, as follows: - AATSR skin values have been converted to subskin values by adding 0.17 K (Donlon et al., 2002) - Saharan dust events are frequent over this area in spring and summer. To eliminate these cases, we have used the SEVIRI derived SDI (Saharan Dust Index) (Merchant, Embury, Le Borgne, & Bellec, 2006) provided with the OSI-SAF SST products in this region. Based on CMS experience we have discarded cases with SDI larger than 0.1. - Only Metop-A/AVHRR and AATSR quality level 5 have been used. - The maximum absolute time difference allowed between AVHRR and AATSR data is 15 min in daytime to minimize diurnal warming effect (see discussion in Section 2.3 below) and 60 min in nighttime.

218

I. Tomažić et al. / Remote Sensing of Environment 146 (2014) 214–233

- Daytime is defined by a Solar zenith angle smaller than 90° and nighttime by a solar zenith angle larger than 110°. - Due to cloud cover and the AATSR narrow swath, AATSR and AVHRR coincident data may sample different regions of the Mediterranean area each day. To obtain more stable results, statistics are produced per 10 day periods (see Fig. 3). Fig. 4 shows the mean differences between the three data sets per 10 day periods recorded over the domain in the conditions described above, as well as the standard deviation and the number of cases. In nighttime conditions, the agreement between the 3 sources of SST is good, with standard deviations below 0.4 K (see Table 1). In early summer (June, July), Metop-A/AVHRR SST is slightly positively biased and OSTIA is slightly negatively biased compared to AATSR, while in winter time Metop-A/AVHRR SST bias is slightly negative. This bias seasonality was observed before over the Adriatic Sea, where comparison were done between different IR sensors (NOAA17/AVHRR, Metop-A/

Table 1 Total validation statistics between Metop-A SST, OSTIA and AATSR SST for day and night. Day

OSTIA-AATSR Metop-A-OSTIA Metop-A-AATSR

Night

bs ± std [K]

N

bs ± std [K]

N

−0.28 ± 0.57 +0.41 ± 0.60 +0.26 ± 0.39

411,866 10,864,847 343,899

−0.08 ± 0.39 +0.04 ± 0.40 −0.00 ± 0.29

1,346,847 9,661,700 1,171,780

AVHRR; MSG/SEVIRI, Terra/MODIS) and ENVISAT/AATSR (Tomažić & Kuzmić, 2011) but also with comparisons to in situ data for the year 2003 (Tomažić, Kuzmić, Notarstefano, Mauri, & Poulain, 2011). The main reason of higher bias seasonality appears to be related to the observed water vapor seasonality (e.g. Tomažić & Kuzmić, 2009) and SST algorithm not capable to accurately resolve the full range of atmospheric profiles.

Fig. 4. Time series of 10 day SST differences between OSTIA and AATSR (black), Metop-A and OSTIA (blue) and Metop-A and AATSR (red) for a) daytime and b) nighttime. First panel in each subfigure represents biases, second panel represents standard deviations and third panel represents number of differences within 10 days.

I. Tomažić et al. / Remote Sensing of Environment 146 (2014) 214–233

219

Fig. 5. Time difference impact on Metop-A/AVHRR − AATSR SST difference in summer daytime for very low winds (below 2 m/s). First panel shows binned average biases based on Metop-A and AATSR differences with colors depicting the average water vapor content, second panel shows binned standard deviations with color depicting standard deviation of water vapor content for that bin. Third panel shows number of points in each bin.

In daytime conditions, Metop-A/AVHRR is positively biased compared to AATSR by about 0.2 K, whatever the season. OSTIA is negatively biased till the end of November. This bias, which can reach −0.5K in the midst of summer is due to diurnal warming and was expected from a comparison between foundation (OSTIA) and skin (AATSR) SST (converted to subskin values). 2.3. Diurnal warming Considering that we are using data in the Mediterranean Sea in summer time, we paid special attention to diurnal warming. We analyze below the Metop-A/AVHRR minus AATSR error from June to August 2011 through dependence on time difference, wind and model integrated total water vapor column to detect a possible diurnal warming effect. First, a trend with time difference between Metop-A and ENVISAT overpass was identified for very low winds (b2 m/s). As shown by Fig. 5, the total biases range between −0.5 K for Metop-A

passing 50–60 min before ENVISAT and +0.2 K for Metop-A passing just after ENVISAT. We also cross compared this trend with a typical DW cycle over the Mediterranean Sea based on SEVIRI data (Fig. 6) averaged over the last 10 days in August. The mean warming rate is about 0.4 K per hour for low wind speed at 10 LST which is similar to the value obtained with direct measurements between Metop-A/ AVHRR and ENVISAT/AATSR. To minimize the impact of this time difference we used a maximum absolute time difference of 15 min as threshold for all daytime studies. Although in winter time the time difference could have been relaxed, for consistency and simplicity we used the same threshold value through the whole period. After having applied this time filtering, we still observe by day a persisting AVHRR-AATSR difference trend with winds below 4 m/s even for orbits 15 min apart (Fig. 7). Zooming over low wind regimes shows that 3.7 m/s seems to be more precisely the limit and this marks also the onset of small waves on sea surface [Donlon et al., 2002].

Fig. 6. Typical summer SST cycle in the Mediterranean as a function of LST and wind speed, averaged over the last 10 days in August, left full cycle, right zoom around midnight LST.

220

I. Tomažić et al. / Remote Sensing of Environment 146 (2014) 214–233

Fig. 7. Summer 2011 daytime (a) and nighttime (b) SST difference between Metop-A and AATSR in dependence of wind (x) and total column water vapor (z). Winds are binned by 0.1 m/s, and total column water vapor is represented by the color (g/cm2). First panel in each subfigure shows biases, second panel shows standard deviations and third panel shows number of points.

Below this limit (that we rounded up to 4m/s) AATSR and AVHRR in warm summer Mediterranean conditions show discrepancies most probably due to a difference in the compared data effective resolutions. AVHRR data have been indeed averaged over a 0.05 grid mesh, a “SEVIRI like” resolution, whereas AATSR data have been remapped to the nearest neighbor. This may have a smoothing effect for large amplitudes (Gentemann, Minnett, Le Borgne, & Merchant, 2008). By night (Fig. 7b), the agreement between AATSR and AVHRR is not sensitive to wind and the main wind impact is visible in the standard deviation values. Indeed at nighttime, the cooling rate around 21 UTC is at most 0.1 K per hour (Fig. 6), so we used 60 min as the nighttime time comparison threshold for the whole period. It should be noted from Fig. 6 that SST at 21 LST are 0.3 K above the nighttime minimum for low winds. This important information has to be taken in consideration when comparing OSTIA based simulations and 21 LST Metop-A brightness temperature observations for adjustment purposes. Fig. 8 confirms the trend with wind of the guess field used in operational simulations (OSTIA two days before) minus AATSR SST difference.

3. Brightness temperature simulations 3.1. “Operational” simulations CMS has been running a prototype chain since November 2011 for testing the bias correction method applied operationally to SEVIRI data (Le Borgne et al., 2011) on Metop-A/AVHRR data. The prototype processing is done in near-real time, to be as close as possible to operational conditions. This choice has several consequences: — The prototype processing chain uses 3-hourly short lead-time atmospheric forecasts (6, 9, 12 and 15 hour forecasts) from ECMWF, initialized from the analyses at 0000 h and 1200 h UTC every day. — ECWMF fields are gridded to L3C grid (0.05 deg) — At each point of the 0.05° grid where a Metop-A SST value is present, the simulation is performed using the corresponding Metop-A/ AVHRR satellite zenith angle and the closest in time ECMWF forecast, interpolated at grid point location.

I. Tomažić et al. / Remote Sensing of Environment 146 (2014) 214–233

221

Fig. 8. Nighttime: OSTIA − AATSR SST as a function of wind speed.

— The Met Office operational analysis (OSTIA) has been used as surface temperature field. To mimic operational conditions the OSTIA field used in the prototype is the one dated two days earlier than the Metop-A data being processed. — The radiative transfer computations, based on RTTOV version 10.1, use a limited number of vertical levels from the model (15 levels at 1000, 950, 925, 900, 850, 800, 700, 600, 500, 400, 300, 250, 200, 150 and 100 hPa), due to limitations of the operational near-real time dissemination of ECMWF products to Météo-France.

3.2. “Delayed mode” simulations To evaluate the impact of operational constraints, simulations have also been done in a delayed mode environment with the following differences: - ECMWF analysis fields at 00, 06, 12 and 18 UTC were used instead of forecasts - simulations are calculated on ECMWF grid (0.1125 deg) and then interpolated to the Metop-A L3C working grid (0.05 deg) - fields are interpolated in time of the mean Metop-A overpass over the Mediterranean Sea for each available orbit - OSTIA corresponding to the same day as the observations - RTTOV version 10.1 (same version) has been applied to 91 profile values (instead of 15 levels) Simulation minus observation difference time series are displayed only for delayed mode simulations (Fig. 9) due to the similarity with operational simulations whereas all statistics is presented in Table 2. As already observed at the global scale (Fig. 1) simulations at 3.7 and 10.8 μm show practically no biases and simulations at 12 μm are positively biased by about 0.2 K. In the Mediterranean Sea channel 12.0 μm simulations are more positively biased for the delayed mode (~0.22 K) compared to operational simulations (~0.15 K). Standard deviation values are smaller in all channels for delayed mode than for operational simulations. The best simulations (without corrections) across both types of simulations are obtained for the 10.8 μm channel

(bias within 0.05 K and SD ~0.5 K), while delayed mode give the best simulations for the 3.7 μm channel simulations (Table 2). Results shown by Fig. 9 look quite stable with time and rather similar for the operational (not shown) and delayed mode. This apparent stability hides significant short scale time and space variability. Fig. 10 shows 10 day averaged simulations − observations for channels 10.8 μm in October. Channel 10.8 μm simulations have a zero mean bias and October has been chosen as an average month. This figure also shows that operational simulations may differ from observations by ±1 K in coherent patches and that these errors are significantly attenuated in delayed mode. Analyzing the simulation to observation daytime differences in summer time with respect to wind speed (Fig. 11) shows that there is a high dependence with the wind speed with amplitudes between − 0.5 K to 1 K for operational mode, and − 0.4 K to 0.4 K for the delayed mode simulations. At high wind speeds (above 10 m/s), simulation to observation differences are much higher for operational compared to delayed mode simulations. A positive trend for higher winds could be attributed to the different representation of the atmosphere in dry conditions (associated with high winds) due to the reduced number of vertical levels. At low wind speeds (below 4 m/s), some wind dependent negative simulation minus observation differences are expected, due to the fact that simulations are derived with the OSTIA field not affected by diurnal warming. In operational mode, the guess SST field is OSTIA from two days ago, and in warming period (from March to the end of August) simulations should be always lower compared to observations, due to the combined effect of seasonal increasing trend and diurnal warming. The difference between operational and delayed mode is a combined effect of different OSTIA and the different vertical resolutions of the profiles. This is visible in the lower wind speed region, where the slope is higher in the operational simulations compared to delayed mode simulations. Nighttime simulations in summer (Fig. 12) show a positive trend of biases and negative trend of water vapor load with the increasing wind speed, which is more pronounced for operational mode simulations (not shown) than for delayed mode. The OSTIA timing problem and a still evident residual diurnal warming effect in the late evening (see previous section), led to select 4 m/s as the nighttime threshold in

222

I. Tomažić et al. / Remote Sensing of Environment 146 (2014) 214–233

Fig. 9. Metop-A/AVHRR delayed mode differences between BT simulations and observations for channel 3.7 μm (black curve; only nighttime), channel 10.8 μm (red) and channel 12.0 μm (blue) for a) daytime and b) nighttime. First panel in each subfigure represents 10 day averaged biases, second panel represents standard deviation from 10 days files and third panel shows number of points in 10 days period.

deriving BT adjustment fields (see next section) as a compromise between having enough number of available points and excluding residual diurnal warming affected measurements . In winter time (not shown), all these effects are not visible, but for the sake of simplicity we maintain the same filtering throughout the whole period. 4. BT adjustments and SST correction 4.1. BT simulation adjustment method Comparison between BT simulations and observations showed differences due to the combined effect of guess errors, NWP errors, filter errors, RTM errors or profile sampling errors, as introduced in Section 1. The bias correction method we use is able to account for differences between guess field and true SST, through Eq. (3), but not for the other sources of errors. An adjustment of simulations is needed to solve this issue. Our simulated BT adjustments are based on temporal

averaging and spatial smoothing of the simulation minus observation differences. We assume a zero mean difference between the guess SST (used in simulations) and the true SST (corresponding to observations), we assume also that there is no systematic profile or RTM error differences between the adjustment conditions, and the conditions when we apply these adjustments. These assumptions are discussed below. Finding optimal combination of space/time criteria is part of the research, therefore we tested a range of averaging sizes. For temporal averaging, we used averaging up to 30 days (1, 3, 5, 10, 15, 20, 30), where 1 day means no averaging. For spatial smoothing we used different window sizes (Table 3) from 0.1 deg to 50 deg and over the whole domain (Inf). In the prototype a single temporal (10 days) and spatial (10 deg) smoothing was used, while in the study phase, both delayed and operational mode simulations were analyzed over the full range of space/time parameters. The most obvious problem induced by using simulations with no adjustment is the overestimation of channel 12 μm temperatures by

I. Tomažić et al. / Remote Sensing of Environment 146 (2014) 214–233 Table 2 Total statistics for the Metop-A BT simulations derived in delayed and operational mode for 3 channels. Delayed

(simBT − obsBT) CH 3.7 μm (simBT − obsBT) CH 10.8 μm (simBT − obsBT) CH 12.0 μm Operational

(simBT − obsBT) CH 3.7 μm (simBT − obsBT) CH 10.8 μm (simBT − obsBT) CH 12.0 μm

Day

Night

bs ± std [K]

N

bs ± std [K]

N

–

–

−0.03 ± 0.37

9,579,784

+0.03 ± 0.51

11,855,997

+0.04 ± 0.48

9,579,784

+0.23 ± 0.55

11,855,997

+0.23 ± 0.55

9,579,784

bs ± std [K]

N

bs ± std [K]

N

–

–

−0.04 ± 0.43

10,564,572

−0.03 ± 0.59

13,223,490

−0.00 ± 0.54

10,564,572

+0.15 ± 0.62

13,223,490

+0.17 ± 0.60

10,564,572

Day

Night

223

about 0.2 K. Since there is near-zero bias for channel T10.8, this would reduce the T10.8–T12.0 simulated term in Eqs. (1) or (2), leading to underestimation of the atmospheric correction term, an underestimation of the simulated SST and hence to an overestimation of the simulation derived bias correction. We nevertheless considered this option in the test set. One assumption of this adjustment is the average zero bias between guess SST and true SST. For this reason we selected nighttime data. However, Fig. 12 shows an increase of simulated − observed BT difference with wind speed. Part of the problem is residual DW effect at the time of Metop-A nighttime orbits, while using OSTIA that represents true foundation temperature. To try to avoid any systematic biases due to residual cooling we discarded all data showing wind below 4 m/s in the adjustment procedure. Adjusting operational and delayed mode BT simulations we used slightly different procedures in the prototype and the study phase: 1) Sampling and concatenation Prototype: Observed and simulated BT are sampled every 20 points and lines (1° longitude/latitude) and concatenated over 10 days, together with latitude, longitude, satellite and solar angles, model

Fig. 10. Nighttime 10 day median of differences between BT simulations and observation for channel 10.8 μm for a) operational and b) delayed mode simulations in the period between 09/10/2011 and 18/10/2011. Title of each figure gives overall statistics (bias and standard deviation).

224

I. Tomažić et al. / Remote Sensing of Environment 146 (2014) 214–233

Fig. 11. Metop-A/AVHRR summer daytime simulation minus observations for channel 10.8 μm for operational (black line) and delayed (red line) mode.

water vapor, OSTIA SST and operational SST. The following filters were applied in constructing differences: • • • •

|OSTIA − operational SST| b1 K |simBT10.8 − obsT10.8| b 1.5 K Metop-A/AVHRR quality flags 3, 4 and 5 Solar zenith angle above 110° (only nighttime data)

• Solar zenith angle above 90° (to ensure enough data to perform analysis) • |OSTIA − SST operational| b 1 K • | simT10.8 − obsT10.8| b 1.5 K 2) Satellite zenith angle normalization.

Study phase: Simulations and observation fields are aggregated over T days without sampling (full resolution 0.05°) and with applying the following filtering:

For a given channel, the simulation minus observation difference is expressed as a function of the satellite zenith angle secant (operational mode processing) or just absolute value of satellite zenith angle (delayed mode processing) by a quadratic formula (4)

• wind speed above 4 m/s • Metop-A/AVHRR quality flags 4 and 5 • SDI b0.1

simBT i −obsBT i ¼ ai Sz þ bi Sz þ ci

2

Fig. 12. Metop-A/AVHRR summer nighttime simulation minus observations for channel 10.8 μm for operational (black line) and delayed (red line) mode.

ð4Þ

I. Tomažić et al. / Remote Sensing of Environment 146 (2014) 214–233

225

Table 3 Spatial smoothing parameters used in analysis. Inf denotes smoothing over the whole domain (Mediterranean Sea) that produces constant adjustment value for each channel. Space (pixels) Space (deg)

1 0.05

3 0.15

5 0.25

9 0.45

15 0.75

21 ~1

31 ~1.5

51 ~2.5

where Sz is either secant of satellite zenith angle or absolute value of satellite zenith angle. The coefficients in Eq. (4) are derived by regression on the concatenated file where i represents channel number. deltaseci = aiSz + biS2z + ci contains the overall mean BT difference and the influence of the satellite zenith angle, considered over the whole domain and over the averaging time period. Each simulation minus observation BT difference is normalized by adding a satellite zenith angle correction to the raw difference as: deltanormi ¼ simBT i −obsBT i −deltaseci

ð5Þ

3) Building mean difference fields Sampled fields are averaged over time after the satellite zenith angle normalization step. At a given location, deltanormsi are averaged over T days (T = 10 for prototype processing). No limitation is imposed on the number of cases actually contributing to the time averaging. After time averaging, a smoothing (mean moving window filter) is done over the averaged fields in S × S pixel boxes, corresponding to S × 0.05° latitude longitude boxes (Table 3). These fields will be defined as a “normalized BT adjustment fields” (deltanorm_meani). In case of operational processing, smoothing is done over sampled fields in 10 × 10 pixel boxes corresponding to 10° latitude longitude boxes, while in the case of an analysis the smoothing is done on the parameters shown in Table 3. Full resolution fields are then reconstructed by expanding the sampled field to full resolution (re location of the sampled values on the full resolution grid), then smoothing in 20 × 20 pixel boxes (only for prototype processing). This procedure is equivalent to choosing S = 200 in delayed mode.

71 ~3.5

101 ~5

201 ~10

301 ~15

351 ~18

401 ~20

501 ~25

601 ~30

999 ~50

Inf Inf

4) Applying the adjustment fields Simulated brightness temperatures are corrected by using the overall mean BT difference and influence of satellite zenith angle (deltaseci) and normalized BT adjustment fields averaged over the whole domain and time period (deltanorm_meani): corrsimB T i ¼ simBT i −deltaseci −deltanormmea ni

ð6Þ

In Eq. (6), first correction (deltaseci) corresponds to the BT adjustment over the whole domain normalized by satellite zenith angle while second one (deltanorm_meani) corresponds to the local BT adjustment based on spatial smoothing parameters. Comparison between Figs. 9b and 13 as well from the statistics represented in Tables 2 and 4 we see that the adjustment procedure reduced biases (mainly for channel 12.0 μm) and standard deviations both for operational and delayed mode simulations. Similar to uncorrected simulation analysis, delayed mode simulations have lower standard deviation compared to operational simulations while the wind dependence is preserved (not shown) after adjustment. 4.2. SST correction After deriving nighttime BT adjustment fields, simulated BT fields are corrected by the derived BT adjustment for each combination of space-time parameters. We assume that the derived nighttime BT adjustment could be applied also to daytime simulations. After that, applying the Metop-A SST algorithm to adjusted BT gives a “simulated” SST, and SST algorithm error is defined as a difference between simulated adjusted SST and the guess SST used in the simulation step (in our case OSTIA). The corrected SST is then obtained by applying

Fig. 13. Nighttime 10 day averaged differences between adjusted delayed mode BT simulations and observations for all channels, where 20 days of time averaging and 15° spatial smoothing were chosen as an example of correction.

226

I. Tomažić et al. / Remote Sensing of Environment 146 (2014) 214–233

Table 4 Total statistics for the Metop-A adjusted BT simulations derived in delayed and operational mode in the same configuration for 3 channels. Delayed (T20S301 configuration)

(simBT − obsBT) CH 3.7 μm (simBT − obsBT) CH 10.8 μm (simBT − obsBT) CH 12.0 μm Operational

(simBT − obsBT) CH 3.7 μm (simBT − obsBT) CH 10.8 μm (simBT − obsBT) CH 12.0 μm

Night bs ± std

N

+0.03 ± 0.36 +0.05 ± 0.47 +0.05 ± 0.53

5,927,782 5,927,782 5,927,782

Night bs ± std

N

+0.03 ± 0.41 +0.02 ± 0.51 +0.03 ± 0.57

5,264,672 5,264,672 5,264,672

smoothing combinations, where bias stability is defined as a standard deviation of 10 day averaged biases. Operational and delayed mode simulation results are shown on Figs. 16 and 17 for the whole period; and Figs. 18 and 19 for the winter season only. These three parameters give almost independent information so they are analyzed and presented separately. Time series and tables are produced separately for day and night and for operational and delayed mode simulation, and due to the specific conditions in summer time in Mediterranean (intensive DW episodes), we analyzed separately the whole period (with excluded data for wind below 4 m/s) and only the colder season (from 11/2011 to 03/2012) with excluded data for wind below 2 m/s.

5.1. Operational daytime — whole period the Metop-A SST algorithm to BT observations and subtracting the derived SST algorithm error for each space-time combination. The main assumption that the nighttime BT adjustment is applicable in a daytime condition is questionable in conditions of intensive summer DW events because DW may induce specific atmospheric profiles, and thus specific DW induced profile errors, which are not accounted for by nighttime adjustment. In addition, filtering out low wind conditions concerns large areas in the Mediterranean Sea (see Fig. 7 — third subplot) and reduces the representativeness of the nighttime BT adjustment applied both during day and night. Maybe the optimum solution will be to use different BT adjustments for daytime and nighttime, but this would also require SST analysis capable of resolving diurnal cycle. This can be done by combining high temporal resolution SST data from geostationary satellites together with hydrodynamic numerical models to produce a diurnal cycle estimate. Marullo et al. (2014) performed such analysis where they combined Mediterranean Forecasting System (MFS) numerical model analysis with hourly SEVIRI SST data to reconstruct hourly SST field over Mediterranean Sea. 5. Results Adjusted simulations are used in a bias correction method as described in the previous section. Each time and space adjustment combination produces a set of simulations that are used to derive correction terms to the CMS operational Metop-A AVHRR algorithm. The efficiency of the adjustment is measured by the validation results of the corrected SST compared to coincident AATSR data. The comparison conditions are identical to those of the validation experiment described in Section 2. The tested solutions are referred to by the corresponding T (time) and S (space) values (defined in previous chapters) where: - T00S00: corresponds to no BT adjustment; - T01S0.05 corresponds to a pixel to pixel BT adjustment where the corrected SST is equal to the SST field used for deriving simulations (OSTIA) - SInf corresponds to using one unique value, derived by averaging over T days and smoothing over the whole domain The result analysis will be based on two types of figures: 1. time series showing the bias and standard deviation of the operational CMS results (with or without bias correction) and the results obtained for typical combinations of simulation adjustment: T (1, 3 and 30 days) and S (0.05, ~15 and ~50 deg). Figs. 14 and 15 represent time series for the operational and delayed mode simulations, respectively. The best compromise solution (3 days; 15 deg, see Section 5.5 below) is also shown on these figures. 2. Aggregated figures presenting the three key parameters, bias, bias stability and standard deviation as a function of all time and space

The most striking feature in Fig. 14 is the difference between summer and the cold season. In summer, the correction is overestimated in most cases and leads to negative corrected SST biases, whereas during the cold season the correction works as expected. The correction is deduced (with opposite sign) from the difference between the simulated and the guess SST. If simulated SST is larger than it should, then the correction is too (negatively) large. A simple cause of overestimating simulated SST is to overestimate the simulated brightness temperatures. Residual nighttime diurnal warming affected AVHRR observations are larger than the foundation SST (OSTIA). Since the nighttime Metop-A/AVHRR measurements are around 21 LST, Fig. 6 suggests that this temperature should be about 0.2–0.3 K warmer than the minimum temperature observed at 5 LST that is closest to foundation SST. Therefore, the resulting BT adjustment is too large and leads to overestimated simulations, and then to overestimated corrections. Different space/time realizations produce prominently different results in the summer season, while in the colder season all realizations produce similar results. Fig. 14 introduces the main contradiction between bias (1st panel in Fig. 14a and b) and standard deviation results (2nd panel in Fig. 14a and b). Small time averaging (1 day, blue lines) produces low bias and high standard deviation results, whereas the opposite is observed for long time averaging. Similarly small space smoothing leads to large standard deviation. The optimum combination will necessarily be a compromise. If we analyze now the results in Fig. 16, the best result, when using bias as the selection criterion (Fig. 16a — 1st panel), is obtained without time averaging (T01) and with lower spatial smoothing (below 1.5°), where the absolute minimum bias value is −0.049 K. The total range of obtained biases is between −0.15 K and −0.05 K, while the operational not corrected bias is 0.25 K. The increase of bias depends on time averaging, with the highest absolute values obtained for the highest time averaging, but nevertheless all obtained values show improvement over the not corrected results. Using bias stability (Fig. 16a — 2nd panel) as criterion favors lower temporal averaging (3 days) and again with lower spatial smoothing (below 1 deg). In this case the minimum values are the same as the value obtained with operational algorithm without correction. Other values are increasing towards the highest time averaging and spatial smoothing window with maximum values of about 0.23 K. Finally, using overall standard deviation (Fig. 16a — 3rd panel) as a criterion shows a spatial smoothing dependence with highest values corresponding to the lowest smoothing and almost no time averaging dependency. Excluding values of the low smoothing (below ~3.5 deg) obtained values fall within the small range from 0.29 to 0.31 K, an improvement compared to standard deviation obtained without correction (0.35 K). Small time averaging window (1–3 days) with higher spatial smoothing (above 5 deg) window is the best compromise suggested by this study.

I. Tomažić et al. / Remote Sensing of Environment 146 (2014) 214–233

227

Fig. 14. Time series of a) daytime and b) nighttime 10 day averaged SST differences between operational mode corrected Metop-A SST and AATSR SST. Metop-A SST is corrected using combination of several spatial (0.05, 15 and 50 deg) and temporal (1, 3 and 30 day) criteria's with spatial criteria's labeled with “S” and temporal with “T” (e.g. combination of spatial smoothing over 15.05 deg and temporal averaging over 3 days is labeled with T03S15.05, …). For comparison additional lines were plotted showing Metop-A SST difference to AATSR without any correction (“Op” — black line), without applying BT adjustment (“T00S00” — dotted) and correction applied in operational environment (“Op corr” — gray line with star). Colored line with black dots represents the results obtained with optimum parameters (3 days; 15 deg). Each subfigure has three subpanels: time series of 1) biases; 2) standard deviations and 3) number of points.

5.2. Operational nighttime — whole period During nighttime, the operational bias is already close to zero (− 0.024 K) and it could be hard to improve the result if we rely only on bias as selection criterion. Nevertheless a bias minimum is obtained when averaging over 3 to 10 days and for higher spatial smoothing (above 10 deg), where the obtained values are below −0.05 K (minimum is −0.036 K). The time series (Fig. 14b), shows a seasonality in the operational bias signal producing positive biases in summer time, and negative in wintertime, giving an overall bias to near-zero. The bias stability value is an interesting criterion in these conditions (Fig. 16b — 2nd panel), and the largest improvement is obtained for 15 to 20 days averaging

with a high spatial smoothing of 50 deg and a bias stability value of 0.06 K. Considering standard deviation (Fig. 16b — 3rd panel), the optimum is obtained in the upper right corner of the table (largest spatial and temporal averaging). As in daytime however, the variability is small (0.26 ± 0.01 K) for all windows above 17.5 deg. Comparing day and night results, an overlapping between day and night optimal conditions is found between 1 and 5 day averaging, and between 10 and 15° smoothing, based on bias and bias stability results. For day and night, optimizing standard deviation would lead to select the largest spatial scales with contradictory options between day (low time averaging) and night (large time averaging). Altogether a good compromise in nighttime conditions would be found for

228

I. Tomažić et al. / Remote Sensing of Environment 146 (2014) 214–233

Fig. 15. Same as previous figure but with Metop-A SST corrected in delayed mode.

time averaging larger than 3 days and smoothing window above 15 deg.

5.3. Operational to delayed mode — whole period Operational (Fig. 14) and delayed (Fig. 15) mode time series show high similarities, with over correction of the positive error in the summer daytime (both for the daytime and nighttime) and error reduction in winter time. No adjustment results (T00S00) produce the best results in summertime, revealing problematic BT simulation adjustment conditions for this season. In wintertime, the no adjustment options leads to large positive errors, particularly at night. It should be stressed again, that in the case of operational mode, AATSR data coincident with Metop-A data are more independent from the guess SST field used in the simulations (OSTIA two days ago), whereas in the case of delayed mode, the same AATSR data used in validation was also used in deriving the guess field (OSTIA). Since both types of simulations produce similar results, this proves that the method is robust even when using not coincident guess SST fields.

The daytime overall bias comparison (Fig. 17a — 1st panel) shows similar patterns as in the operational mode results (favoring small time averaging), but the absolute values of biases in delayed mode were higher (−0.093 K) compared to the operational mode (−0.049 K), due to the period of very negative bias in September and October which is not included in operational processing for technical reasons. Bias stability (Fig. 17a — 2nd panel) shows the largest differences between the two modes, favoring higher time (20days) averaging in delayed mode against lower averaging (3 days) in operational mode. Both bias stability minima have values similar to the results obtained without correction (~0.18 K). Overall standard deviation (Fig. 17a — 3rd panel) show consistently the minimum values for high spatial smoothing (above 5 deg) and small time averaging with similar absolute values between the two modes. For nighttime, bias result (Fig. 17b — 1st panel) are quite similar to the operational ones, favoring similar time averaging criteria's (3–30 days) and high spatial smoothing (between 10 and 50 deg), with very similar absolute value of biases (−0.049 K for delayed mode and −0.036 K for operational mode). The main difference is found for bias stability, where the delayed mode analyses give more stable values

I. Tomažić et al. / Remote Sensing of Environment 146 (2014) 214–233

229

Fig. 16. Aggregated operational mode 10 day averaged differences for all analyzed spatial and temporal parameters and over the whole period by a) day and b) night. First panel in each subfigure shows biases; second panel shows bias stability and third panel shows standard deviation. Each cell represents averaged value of statistical parameter throughout the analyzed period.

throughout the whole period (Fig. 17a — 2nd panel), with improvement of almost 50% (from 0.060 to 0.034 K) compared to the operational mode. In the delayed mode, optimum values are found for large time averaging (30 days) compared to one day in the operational case. he most favorable standard deviation (Fig. 17b — 3rd panel) condition is obtained with no averaging (1 day) and smoothing over the medium window sizes (6–10 deg). Again, all smoothing above 2.5 deg give overall SD values in the very small range (between 0.248 and 0.266 K). Delayed mode does not lead to conclusions significantly different to those of the operational mode providing the most contradictory time averaging solutions (1 day and 30 days) are avoided.

5.4. Colder period — operational and delayed mode Filtering out low wind conditions however is not totally satisfying, because in the Mediterranean Sea, this favors dry atmospheres and may bias the results. To check the stability of our previous conclusions, we have analyzed overall performances only for the colder period.

For that purpose we created aggregated tables based on the five months in the colder season (Figs. 18 and 19), where we also decreased the wind filtering criteria to 2 m/s, mainly to compensate for the more pronounced cloud mask, since diurnal warming is much less frequent and of smaller intensity. Results suggest that although the overall bias and standard deviation have been improved as expected, the patterns remained which are very similar to the patterns obtained for the total period. For daytime conditions, both the operational and delayed mode analysis give near zero biases with similar patterns as for the whole period, favoring lower time averaging (1 or 3 days) and medium to high smoothing window size (10–50 deg). It should be noted that the overall amplitude in the winter time months is also reduced, ranging between 0 and 0.05 K. Bias stabilities are also reduced in the winter daytime (from 0.20 K to 0.14 K for operational and 0.10 K for delayed mode) with a shift towards higher smoothing windows (15 to 18 deg) and smaller time averaging (1 to 3days). Again, the total range of values is within 0.02K. Winter time standard deviations are optimal when used without time averaging (T01) and for the medium smoothing window size (from 10 deg to 15 deg), whereas for the whole period higher

230

I. Tomažić et al. / Remote Sensing of Environment 146 (2014) 214–233

Fig. 17. Aggregated delayed mode results for different spatial and temporal parameters over the whole period by a) day and b) night. See the description from the Fig. 16.

window sizes were preferred in both cases (25 deg for operational and 30 deg for delayed mode). 5.5. Discussion Analyzing both periods (whole and colder season) and both types of simulations (operational and delayed) leads to the following remarks: 1) The most consistent results throughout the day and night and two simulations types are obtained for the standard deviation results, where in all cases high spatial smoothing, above 10 deg is always favorable with very small variability (within 1–2 hundredths of Kelvin) between different time averaging realizations. The main reason for this is that by smoothing over a larger area we decrease the noise introduced in the simulations through errors in the guess field, atmospheric profiles and RT model. 2) Bias results are consistent across the different simulation modes, but diverge between day and night. For the daytime, results are dependent to the time window averaging, improving from high

to low temporal window and the minimum values are obtained between 1 and 5 days averaging. During nighttime, results are essentially dependent on the spatial smoothing window size, while the optimum results are obtained for smoothing above 10 deg and higher time averaging (3–30 days). This difference between day and night suggests that daytime biases are affected by more transient features (mainly diurnal warming) while nighttime results are affected by long range weather patterns that change the atmospheric structure. Again, the difference between day and night is also driven by the difference in the SST algorithm, where during the daytime a split window algorithm is used and during the nighttime more accurate triple window algorithm is used. 3) Bias stability results show the highest variability between both simulations types and between day and night, but also the smallest variability within each analysis. During daytime, the applied correction basically shifts constantly positive results (0.25 K) to near zero values, and therefore bias stability values are similar before and after correction. Still, the minimum values for both

I. Tomažić et al. / Remote Sensing of Environment 146 (2014) 214–233

231

Fig. 18. Aggregated wintertime operational mode results for different spatial and temporal parameters by (a) day and (b) night. See the description from the Fig. 16.

types of simulations and both periods are obtained for 3 days averaging and lower spatial smoothing (below 15 deg). During nighttime, all corrections give improvement while the overall variability (especially for the delayed mode simulations) is lower compared to the daytime results. Still, the optimum parameters for operational mode simulations favor both higher temporal averaging (3–20 days) and higher spatial smoothing (above 15 deg). 4) Considering all different conditions we suggest using temporal averaging of 3 to 5 days and spatial smoothing about 15 deg. For final test purposes, we marked (colored line with black dot) the obtained optimum combination of parameters (3 days; 15 deg) in time-series (Figs. 14 and 15) and we clearly see a balanced improvement between bias and standard deviation values, not obtained for other chosen parameters while for certain periods we got either improvement in bias or in standard deviation. 6. Conclusions and future work We assessed the impact of BT simulation adjustment conditions in correcting Metop-A SST over the Mediterranean Sea and elaborated a

detailed methodology for optimization of estimates of regional biases between simulated and observed BT. We used two different BT simulations: full vertical (91) ECMWF NWP profiles combined with concurrent OSTIA SST field (delayed mode) and reduced vertical (15) ECWMF NWP profiles combined with two days delayed OSTIA field (operational). For both types of simulations, separately for day and night we performed BT adjustments and SST corrections over different temporal and spatial averaging window sizes. To validate the method we used AATSR data and three statistical parameters derived from 10 day averaged differences between Metop-A and AATSR: bias, bias stability and standard deviation. We took special care in analyzing daytime condition in summer time during intensive DW episodes that tends to mask out atmospheric problems we aimed to correct for. The final choice results from compromising between different conclusions according to each of these parameters are used as selection criterion. Based on all different conditions we found optimum temporal averaging between 3 and 5 days and spatial smoothing of 15 deg. Specifically, based on the standard deviation cost function we consistently obtained that smoothing above 10 deg is the most favorable, almost regardless of temporal averaging. Bias based results diverge between day and night, where daytime analysis are mostly time averaging dependent and favors short temporal averaging (1 to

232

I. Tomažić et al. / Remote Sensing of Environment 146 (2014) 214–233

Fig. 19. Aggregated wintertime delayed mode results for different spatial and temporal parameters by (a) day and (b) night. See the description from the Fig. 16.

5 days) while nighttime results are spatial smoothing dependent and favors higher smoothing (above 10 deg) for any temporal averaging (except no averaging). The main improvement of using full resolution atmospheric profiles compared to reduced ones is evident in the lower standard deviation of the simulation to observation differences and for nighttime bias stability results where almost all space/time realizations produce values lower (between 0.034 and 0.062 K) than the results derived without any correction or with operational simulations (between 0.060 and 0.116 K). Converging conclusions derived from operational and for delayed mode simulations, produced in totally independent ways, is an indication of the consistency of the optimal solution retained. Availability of the guess SST field with even the higher spatial resolution, with possibility to resolve finer scale mesoscale features, could probably decrease the spatial smoothing window. Also, high influence of DW in summer daytime conditions in Mediterranean Sea decorrelates nighttime from daytime, leading ideally to more appropriate separate BT adjustment for day and nighttime. In practice this is still not feasible at present since it would require diurnal warming resolving daytime SST analysis. This methodology together with obtained optimized spatio/temporal parameters is foreseen to be applied at CMS globally for Metop-A/AVHRR

but also for SUOMI-NPP/VIIRS and Metop-B/AVHRR SST algorithms upon successful testing on the global scale. Acknowledgements This work was realized in the context of the BESST (Inter-sensor Bias Estimation in Sea Surface Temperature) - SR/12/158 project funded by the Belgian Science Policy (BELSPO) in the frame of the Research Program For Earth Observation “STEREO II” and EUMETSAT Visiting Scientist program (OSI_AS12_02). The data from the EUMETSAT Satellite Application Facility on Ocean and Sea Ice used in this study are accessible through the SAF's homepage: http://www.osi-saf.org. ENVISAT/AATSR data are provided by the European Space Agency (ESA) via the project ID 5802. ECMWF data are obtained through the MARS archive. We would like to thank Sonia Péré for producing Fig. 6. Finally, the authors would like to thank the three anonymous reviewers for their valuable comments. References Buongiorno Nardelli, B., Tronconi, C., Pisano, A., & Santoleri, R. (2013). High and Ultra-High resolution processing of satellite Sea Surface Temperature data over the

I. Tomažić et al. / Remote Sensing of Environment 146 (2014) 214–233 Southern European Seas in the framework of MyOcean project. Remote Sensing of Environment, 129, 1–16. http://dx.doi.org/10.1016/j.rse.2012.10.012. Corlett, G. K. (2008). Single Sensor Error Statistics: Report for (A)ATSR L2P Project. University of Leicester. Corlett, G. K., et al. (2006). The accuracy of SST retrievals from AATSR: An initial assessment through geophysical validation against in situ radiometers, buoys and other SST data sets. Advances in Space Research, 37(4), 764–769. http://dx.doi.org/10.1016/j.asr. 2005.09.037. Donlon, C. J., Martin, M., Stark, J., Roberts-Jones, J., Fiedler, E., & Wimmer, W. (2012). The Operational Sea Surface Temperature and Sea Ice Analysis (OSTIA) system. Remote Sensing of Environment, 116, 140–158. http://dx.doi.org/10.1016/j.rse.2010.10.017. Donlon, C. J., Minnett, P. J., Gentemann, C., Nightingale, T. J., Barton, I. J., Ward, B., et al. (2002). Toward improved validation of satellite sea surface skin temperature measurements for climate research. Journal of Climate, 15(4), 353–369. Donlon, C., et al. (2007). The Global Ocean Data Assimilation Experiment High-Resolution Sea Surface Temperature Pilot Project. Bulletin of the American Meteorological Society, 88(8), 1197–1213. http://dx.doi.org/10.1175/BAMS-88-8-1197. EUMETSAT (2010). Low Earth Orbiter Sea Surface Temperature Product User Manual, SAF/OSI/CDOP/M-F/TEC/MA/127. [online] Available from: http://www.osi-saf.org EUMETSAT (2011). Geostationary Sea Surface Temperature Product User Manual, SAF/OSI/CDOP/M-F/TEC/MA/181. [online] Available from: http://www.osi-saf.org François, C., Brisson, A., Le Borgne, P., & Marsouin, A. (2002). Definition of a radiosounding database for sea surface brightness temperature simulations: Application to sea surface temperature retrieval algorithm determination. Remote Sensing of Environment, 81(2–3), 309–326. http://dx.doi.org/10.1016/S0034-4257(02)00008-1. Gentemann, C. L., Minnett, P. J., Le Borgne, P., & Merchant, C. J. (2008). Multi-satellite measurements of large diurnal warming events. Geophysical Research Letters, 35(22), L22602. http://dx.doi.org/10.1029/2008GL035730. Hocking, J., Rayer, P., Saunders, R., Matricardi, M., Geer, A., & Brunel, P. (2011). RTTOV v10 Users Guide. Lannion, France: EUMETSAT/NWP SAF [[online] Available from: http:// research.metoffice.gov.uk/research/interproj/nwpsaf/rtm/docs_rttov10/] Le Borgne, P., Roquet, H., & Merchant, C. J. (2011). Estimation of sea surface temperature from the spinning enhanced visible and infrared imager, improved using numerical weather prediction. Remote Sensing of Environment, 115(1), 55–65. http://dx.doi.org/10. 1016/j.rse.2010.08.004. Liang, X., & Ignatov, A. (2011). Monitoring of IR Clear-Sky Radiances over Oceans for SST (MICROS). Journal of Atmospheric and Oceanic Technology, 28(10), 1228–1242. http://dx.doi.org/10.1175/JTECH-D-10-05023.1. Liang, X. -M., Ignatov, A., & Kihai, Y. (2009). Implementation of the community radiative transfer model in advanced clear-sky processor for Oceans and validation against nighttime AVHRR radiances. Journal of Geophysical Research, 114(D6), D06112. http://dx.doi.org/10.1029/2008JD010960. Marullo, S., Santoleri, R., Le Borgne, P., Pere, S., Pinardi, N., Tonani, M., et al. (2014). Combining model and geostationary satellite data to reconstruct the hourly SST field over the Mediterranean Sea. Remote Sensing of Environment, 146, 11–23. McClain, E. P., Pichel, W. G., & Walton, C. C. (1985). Comparative performance of AVHRR-based multichannel sea surface temperatures. Journal of Geophysical Research, 90(C6), 11587–11601. Merchant, C. J., Embury, O., Le Borgne, P., & Bellec, B. (2006). Saharan dust in nighttime thermal imagery: Detection and reduction of related biases in retrieved sea surface temperature. Remote Sensing of Environment, 104(1), 15–30. http://dx.doi.org/10. 1016/j.rse.2006.03.007.

233

Merchant, C. J., & Harris, A.R. (1999). Toward the elimination of bias in satellite retrievals of sea surface temperature 2. Comparison with in situ measurements. Journal of Geophysical Research, 104(C10), 23579–23590. http://dx.doi.org/10.1029/1999JC 900106. Merchant, C. J., Harris, A.R., Maturi, F., Embury, O., MacCullum, S. N., Mittaz, J., et al. (2009). Sea surface temperature estimation from the Geostationary Operational Environmental Satellite-12 (GOES-12). Journal of Atmospheric and Oceanic Technology, 26(3), 570–581. Merchant, C. J., Le Borgne, P., Marsouin, A., & Roquet, H. (2008). Optimal estimation of sea surface temperature from split-window observations. Remote Sensing of Environment, 112(5), 2469–2484. http://dx.doi.org/10.1016/j.rse.2007.11.011. Merchant, C. J., Le Borgne, P., Roquet, H., & Marsouin, A. (2009). Sea surface temperature from a geostationary satellite by optimal estimation. Remote Sensing of Environment, 113(2), 445–457. http://dx.doi.org/10.1016/j.rse.2008.10.012. O'Carroll, A. G., August, T., Le Borgne, P., & Marsouin, A. (2012). The accuracy of SST retrievals from Metop-A IASI and AVHRR using the EUMETSAT OSI-SAF matchup dataset. Remote Sensing of Environment, 126, 184–194. http://dx.doi.org/10.1016/ j.rse.2012.08.006. O'Carroll, A. G., Eyre, J. R., & Saunders, R. W. (2008). Three-way error analysis between AATSR, AMSR-E, and in situ sea surface temperature observations. Journal of Atmospheric and Oceanic Technology, 25(7), 1197–1207. http://dx.doi.org/10.1175/ 2007JTECHO542.1. O'Carroll, A. G., Watts, J. G., Horrocks, L. A., Saunders, R. W., & Rayner, N. A. (2006). Validation of the AATSR meteo product sea surface temperature. Journal of Atmospheric and Oceanic Technology, 23(5), 711–726. Petrenko, B., Ignatov, A., Shabanov, N., & Kihai, Y. (2011). Development and evaluation of SST algorithms for GOES-R ABI using MSG SEVIRI as a proxy. Remote Sensing of Environment, 115(12), 3647–3658. http://dx.doi.org/10.1016/j.rse.2011.09.003. Saha, K., Ignatov, A., Liang, X. M., & Dash, P. (2012). Selecting a first-guess sea surface temperature field as input to forward radiative transfer models. Journal of Geophysical Research, 117(C12), C12001. http://dx.doi.org/10.1029/2012JC008384. Tomažić, I. (2012). Using Fine scale NWP model outputs to improve the retrieval of SST from METOP AVHRR brightness temperatures. Part 1: Initial results based on SEVIRI BT, Eumetsat Visiting Scientist Report, No. OSI_AS11_P03. Tomažić, I., & Kuzmić, M. (2009). ERA-40-aided assessment of the atmospheric influence on satellite retrieval of Adriatic Sea surface temperature. Meteorology and Atmospheric Physics, 104(1–2), 37–51. http://dx.doi.org/10.1007/s00703-008-0015-2. Tomažić, I., & Kuzmić, M. (2011). Improving regional AVHRR SST measurements using AATSR SST data. Proceedings of the 2011 EUMETSAT Meteorological Satellite Conference (pp. 59) [Oslo, Norway]. Tomažić, I., Kuzmić, M., Notarstefano, G., Mauri, E., & Poulain, P.M. (2011). A comparative assessment of satellite-derived Adriatic Sea surface temperature. International Journal of Remote Sensing, 32(17), 4871–4892. http://dx.doi.org/10.1080/01431161.2010.492249. Walton, C. C., Pichel, W. G., Sapper, J. F., & May, D. A. (1998). The development and operational application of nonlinear algorithms for the measurement of sea surface temperatures with the NOAA polar-orbiting environmental satellites. Journal of Geophysical Research, 103(C12), 27999–28012. Wimmer, W., Robinson, I. S., & Donlon, C. J. (2012). Long-term validation of AATSR SST data products using shipborne radiometry in the Bay of Biscay and English Channel. Remote Sensing of Environment, 116, 17–31. http://dx.doi.org/10.1016/j.rse.2011.03.022. Zavody, A.M., Mutlow, C. T., & Llewellyn-Jones, D. T. (1995). A radiative transfer model for sea surface temperature retrieval for the along-track scanning radiometer. Journal of Geophysical Research, 100(C1), 937–952.