Evaluation of the Tropical Variability from the Beijing Climate Center's Real-Time Operational Global Ocean Data Assimilation System

Ocean initialization plays a critical role in short-term climate prediction because most predictive skill comes from the initial conditions of the upper ocean, particularly that associated with large-scale variability, such as ENSO and the Indian Ocean Dipole (IOD) (Alves et al., 2011). A common strategy to obtain the optimal initialization of the ocean is to assimilate the available ocean observations into ocean models forced by atmospheric fluxes (Xue et al., 2011; Zheng and Zhu, 2015). Various assimilation methods have been applied for the initialization of operational or quasi-operational coupled forecast systems. For example, a three-dimensional variational (3DVAR) analysis with a temporally and spatially varying background error covariance is used in the global ocean data assimilation system (GODAS) of the National Centers for Environmental Prediction (NCEP) (Behringer et al., 1998). In version 3 of the ECMWF (European Centre for Medium-Range Weather Forecasts) ocean analysis system, the incremental 3DVAR method is applied (Balmaseda et al., 2013). A more advanced data assimilation method, the EnKF (Ensemble Kalman Filter), has been employed in TOPAZ4 (the latest version of TOPAZ-a coupled ocean-sea ice data assimilation system for the North Atlantic Ocean and Arctic), which can make background error covariance vary with time (Sakov et al., 2012). Noticeably, the data assimilation in these systems is performed separately for ocean and atmosphere models in order to achieve optimal initialization (Zheng et al., 2009; Alves et al., 2011).

High-quality in-situ and remote sensing data with fine spatial and temporal resolution are crucial for the development of a global operational ocean data assimilation system (Han et al., 2013; Ratheesh et al., 2014). Satellite observations of the ocean, which provide global coverage and real-time measurements of sea level, SST, sea ice, waves, and winds, have become a primary data source for operational oceanography. The launch of TOPEX/Poseidon (T/P) has provided an accurate description of the large-scale sea level and ocean circulation for the first time since 1992. The altimetry products at higher resolution have been supplied following the launch of ERS-1/2, Janson-1 and ENVISAT (Hunt et al., 2007). SST is an important variable for operational oceanography and for assimilation into global ocean models because SST is strongly related to air-sea interactions and can be used to correct errors in forcing fields (e.g., wind and surface heat fluxes). Satellite-derived SST products are available in two types, according to the sensors equipped on satellites. Optical/infrared data, including data from the Advanced Very High Resolution Radiometer (AVHRR), are affected by clouds and volcanic aerosols in the atmosphere (Reynolds and Marsico, 1993). Microwave data can be obtained through non-precipitating clouds and are very beneficial in terms of geographical coverage (Ratheesh et al., 2014). However, satellite observations are mainly confined to surface information. In-situ observations can complement satellite observations by providing measurements below the ocean surface. This is especially true of the Array for Real-time Geostrophic Oceanography (ARGO) project, which has been established since 2000 and now provides important in-situ observations to validate satellite data and improve the initialization of ocean models. ARGO floats measure the temperature and salinity in the upper 2000 m of the ocean and bring unquestionable advances to ocean forecast models (Huang et al., 2008). The joint use of high resolution satellite observations and high accuracy in-situ observations fulfills the operational requirements for temporal and spatial coverage and near real-time access.

A new generation of short-term climate forecast systems is being developed at the Beijing Climate Center (BCC_ CSM1.1) (Wu et al., 2013, 2014), and the access to initial oceanic conditions depends on the development of a global ocean data assimilation system. The first generation of the global ocean data assimilation system of the Beijing Climate Center (BCC_GODAS1.0) has been developed during the period from January 1996 to December 2000, has been operational since 2002 (Liu et al., 2005), and plays an important role in climate forecast services. With the development of satellite technology and the data assimilation method, some limitations of BCC_GODAS1.0 have gradually emerged, especially with regard to not assimilating satellite observations. The physical processes in the BCC_GODAS1.0 ocean model and the ocean component used in BCC_CSM1.1 are different. A new method for assimilating altimeter and SST data under one dynamic constraint based on 3DVAR was developed to solve the problem of multivariate assimilation

(Wang et al., 2012a, 2012b). The new method has been implemented into BCC_GODAS2.0, which now has the ability to assimilate real-time observations and to run automatically. The physical processes and grid of the ocean model used in BCC_GODAS2.0 are the same as those used in BCC_ CSM1.1. Therefore, the ocean initialization field generated by BCC_GODAS2.0 can be used in BCC_CSM1.1 directly. The global ocean reanalysis products of BCC_GODAS2.0 are being applied in BCC_CSM1.1.

This paper is organized as follows: Section 2 provides an overview of the configurations of the ocean model and the 3DVAR analysis scheme, as well as the ocean observations used in the assimilation. Section 3 discusses the model performance and interannual-to-decadal variability in the tropical Pacific. A summary and conclusions are given in section 4.

The ocean model used in BCC_GODAS2.0 is version 4 of the Modular Ocean Model (MOM4), developed at the Geophysical Fluid Dynamics Laboratory (Griffies et al., 2003). Its zonal resolution is 1^°, while the resolution in the meridional direction is 1/3^° within 10^° of the equator, smoothly increasing to 1^° in the poleward direction of 30^°. There are 50 z-levels in the vertical direction, with a 10 m resolution from the surface to a depth of 225 m, gradually increasing to approximately 366 m in the abyssal zone. The maximum depth is approximately 5.5 km. The tripolar grid is adopted to avoid a singularity at the North Pole (Griffies et al., 2005). The two north poles of the curvilinear grid are situated over land areas in North America and Eurasia, respectively. The model uses identical physical parameterization schemes as described in (Griffies et al., 2005), including the isopycnal tracer mixing and diffusion scheme, the Laplace horizontal friction scheme, and the K-Profile Parameterization vertical mixing scheme. The model is forced with 6-h averaged 10-m wind and 2-m air temperature data from the NCEP/NCAR (National Center for Atmospheric Research) Reanalysis I dataset (http://www.esrl.noaa.gov/psd/). The mean climatological river runoff is specified at the coastlines of the model. Surface temperature and salinity are relaxed to monthly Levitus climatology, with restoring time scales of 90 and 120 days, respectively.

A 3DVAR analysis scheme, adapted from (Yan et al., 2004), has been improved by adopting a recursive filter to explicitly avoid the inverse calculation of the background error covariance matrix (Xiao et al., 2008) and save computing resources. The analysis variables are temperature and salinity. There are two steps carried out in the 3DVAR. First, the temperature and salinity profiles along the satellite altimeter track are estimated by minimizing the cost function, defined as follows: \begin{eqnarray*} J&=&({T}-{T}_{ b})^{ T}{B}_{ T}^{-1}({T}-{T}_{ b})+\\ &&[{S}-{g}(T)]^{ T}{B}_{ S}^{-1}[{S}-{g}(T)] +\dfrac{1}{2\sigma^2}[h({T},{S})-h_0]^2 , \end{eqnarray*} where T and S are the column vectors containing the state variables of temperature and salinity, respectively; T _b is the temperature background vector; and B _T and B _S are the background error covariance matrices in the vertical direction for temperature and salinity, respectively. B _T and B _S are both assumed to be diagonal matrices, and the diagonal element is given by the empirical functions used by (Behringer et al., 1998) and (Yan et al., 2004). The variance of the background error in temperature at depth z is given by $$ a_{vT}\dfrac{(d{T}/dz)^{\frac{1}{2}}}{[(d{T}/dz)^{\frac{1}{2}}]_{{ max}}} , $$ where the constant a_vT is determined empirically by tuning the assimilation results and setting to 2.7; dT/dz is the local vertical temperature gradient at depth z; and \([(dT/dz)^\frac12]_ max\) is the maximum value of \((dT/dz)^\frac12\) in the water column. The determination of B _S is analogous to B _T, except the corresponding constant a_vS is set equal to 0.6. g(T) is a nonlinear T-S relationship function, which is determined from datasets of historical observations. h _o is the observed value of the sea surface height. σ is observational error of the sea surface height, which is set to 6 cm here. The function h(T,S) denotes an observation operator that transforms T and S to the surface dynamic height, which is defined as follows: $$ h(T,S)=-\int_0^{z_m}\dfrac{\rho(T,{S},p)-\rho_0(p)}{\rho_0(p)}dz , $$ where z_m is the reference depth, set to 1000 m here. z denotes the vertical coordinate and p denotes the pressure. The function ρ(T,S,p) denotes the sea water state equation for calculating density. ρ₀(p)=ρ(0,35,p) is the reference density. More details can be found in (Fu et al., 2009a) and (Xiao et al., 2008).

Second, the estimated temperature and salinity calculated from the first step are taken as synthetic observations, which are assimilated with the in-situ temperature and salinity observations, such as those from ARGO. The cost function is defined as follows: \begin{eqnarray*} J&\!=\!&({T}-{T}_{ b})^{ T}{E}_{ T}^{-1}({T}-{T}_{ b})+({S}-{S}_{ b})^{ T}{E}_{ S}^{-1}({S}-{S}_{ b})+\\ &&({HT}-{T}_0)^{ T}{O}_{ T}^{-1}({HT}\!-\!{T}_0)\!+\!({HS}\!-\!{S}_0)^{ T}{O}_{ S}^{-1}({HS}\!-\!{S}_0) , \end{eqnarray*} where T and S are the temperature and salinity vectors at the model levels, respectively; and O _T and O _S are the observation error covariance matrices of temperature and salinity, respectively. The observation error covariances here are not simple in-situ observation error but include the synthetic observation error and repetitiveness errors of the variability of unresolved smaller scales. Their orders of magnitude are larger than those of the observation errors in the first step. We set the observation errors of temperature and salinity equal to 0.8^°C and 0.1 psu, respectively. H is the bilinear interpolation operator that interpolates the model grid points to the observation location. T _b and E _T are the background temperature and its corresponding background error covariance matrix in the horizontal direction, respectively. S _b and E _S are the background salinity and its corresponding error variance matrix, respectively. It should be noted that E _T and E _S are not diagonal matrices. The background error covariance matrix is defined using the Gaussian function as follows: $$ E=A\exp\left(-\dfrac{\Delta x^2}{L_x^2}-\dfrac{\Delta y^2}{L_y^2}\right) , $$ where A denotes the background error variance of temperature (set to 2.0^°C²) or salinity (set to 0.15 psu²); ∆ x and ∆ y are the distances of two horizontal grid points in the zonal and meridional directions, respectively; and L_x and L_y are correlation scales in the zonal and meridional directions, respectively. We take L_x=450 km and L_y=650 km for the temperature, and L_x=420 km and L_y=510 km for the salinity, by tuning the analysis. The detail of the method used to set appropriate values for L_x and L_y can be found in (Zhou et al., 2004).

Figure 1.  Vertical distribution of RMSE for the temperature and salinity during 1990-2011. The global average statistics are shown for the free run (blue line) and reanalysis (red line).

In-situ temperature and salinity, satellite-SST and altimeter-derived sea-level anomalies (SLAs) are used for assimilation. The SST observations are from the AVHRR merged product (ftp://eclipse.ncdc.noaa.gov/pub/OI-daily/) created through optimal interpolation (Reynolds et al., 2007). SLAs are from all available satellite altimeters, including T/P, Jason-1/2, ERS-1/2 and Envisat (http://www.aviso.altimetry.fr/en/home.html), which is 0.25^°× 0.25^° and a 7-day average gridded dataset. In-situ temperature and salinity are from the Global Temperature and Salinity Profile Project (GTSPP) and ARGO. Quality control is critically important because erroneous values can cause spurious overturning in ocean models. The temperature and salinity profiles used in the study are mainly subject to unified quality control, including a duplicate check (similar position/time check), a depth inversion check, a temperature and salinity valid range check, an excessive gradient check, and a spike value check. The temperature and salinity valid range check applies a gross filter on observed values, which needs to include all of the expected extremes encountered in the oceans (temperature in the range of -2.0^°C to 39.0^°C; salinity in the range of 0.0 to 41.0 psu). The excessive gradient check requires that the ratio of the difference in adjacent values to the difference in the depths is no larger -0.2^°C m^-1 for temperature and 5.0 psu m^-1 for salinity. The spike value check requires the test value (α) to be no larger than some pre-defined value: $$ \alpha=|V_2-(V_3+V_1)/2|-[(V_3-V_1)/2] , $$ where V₂ is the measurement being tested as a spike, and V₁ and V₃ are the values above and below the tested layer. If α exceeds 6.0 for pressures less than 5000 hPa, or α exceeds 2.0 for pressures less than 5000 hPa, the measurement V₂ is removed for temperature and salinity. Further details can be found in the ARGO quality control manual, version 2.8(http://www.argodatamgt.org/content/download/15699/102401/file/argo-quality-control-manual-version2.8.pdf).

The operational BCC_GODAS2.0 has been running on a Sunway4000A high-performance computer system at the China Meteorological Administration since 1 April 2014. It runs on six nodes (36 processors) and takes two hours (wall clock) to process the observations and three hours to finish the analysis. The data assimilation cycle for the system restarts the model 10 days prior to the most recent time of the NCEP forcing field, while the observations are assimilated within this time window (10 days). The shell scripts controlling the operational run have been implemented to minimize the need for manual operation. A typical daily run starts with the creation of forcing fields for the ocean model. Then, the operational system starts to run after the real-time observations of the SLAs, SST, in-depth temperature and salinity are processed. If the forcing fields are not received by the prescribed time, the model will be forced by the climatological forcing fields.

In this study, we evaluate the monthly results from the BCC_GODAS2.0 reanalysis for 1990-2011.

The temperature and salinity observations from World Ocean Atlas 2001 (WOA01) are used to evaluate the agreement of BCC_GODAS2.0 with observations. The vertical profiles of the root-mean-square error (RMSE) of temperature and salinity are shown in Figs. 1a and 1b, respectively. Large temperature errors are present mainly in the thermocline, and the errors are largest at a depth of approximately 200 m in the free run. The RMSE of the temperature is reduced at most vertical levels by assimilation. The strong relaxation of the SST constrains the surface temperature variation in the model so that the difference between the results of the free run and the reanalysis results is not obvious in the top layer. The overall RMSE of the temperature is reduced by 0.81^°C after the cumulative effect of many assimilation cycles during integration over 22 years. Like the temperature, the RMSE of salinity is reduced by 0.24 psu by assimilation. Note that the improvement can be attributed not only to the assimilation of temperature and salinity but also to the assimilation of the SLA through the dynamic constraints (Yan et al., 2004).

The performance of BCC_GODAS2.0 is further evaluated in the tropical Pacific Ocean (30^°S-30^°N, 120^°E-70^°W), Indian Ocean (30^°S-30^°N, 30^°-120^°E) and Atlantic Ocean (30^°S-30^°N, 80^°W-20^°E), separately. Table 1 compares the RMSEs of the free run and the reanalysis. The reduction of RMSE (RRMSE, in %) of the reanalysis relative to the free run is shown as well. The RRMSE is defined in (Counillon et al., 2014) as $$ { RRMSE}_{ reanalysis}=\dfrac{{ RMSE}_{ free}-{ RMSE}_{ reanalysis}}{{ RMSE}_{ free}} , $$ where RMSE _free and RMSE _reanalysis are the RMSEs of the free run and the reanalysis, respectively. The RMSEs and RRMSEs of temperature (T) and salinity (S) are calculated for surface (0-225 m) and deep (225-1007 m) waters. The RRMSEs of T and S for water depths greater than 1007 m are not examined because the variations of T/S errors with depth are relatively small below 1000 m (as shown in Fig. 1). In general, the RRMSE for deep waters is larger than the surface due to less variability than at the surface. The maximum error reductions are achieved in the Atlantic Ocean due to the large errors in the free run. The RRMSEs are greater than 60% at the surface, and greater than 70% in deep waters for both T and S. In the Indian Ocean, the RRMSE is generally larger than 70% for S, while the RRMSE for T is more than 20%. The RRMSE for the Pacific is more than 50%, while for T in the surface layer it is about 33.2%. The relatively low error reduction of RRMSE for T at the surface is mainly caused by the relatively low RMSE in the free run. Overall, the error reductions are significant.

The mixed layer depth (MLD) is one of the most important quantities of the upper ocean, and its variability strongly influences the physics of the upper ocean and marine biological processes (e.g., de Boyer Montégut et al., 2004; D'Ortenzio et al., 2005; Dong et al., 2007, 2009; Ren et al., 2011; Sallée et al., 2008; Moore et al., 2013). We estimate the monthly mean MLD from the temperature and salinity based on the density criterion of 0.125 (kg m^-3), in which the density is supposed to change from that of the ocean surface of 0.125 kg m^-3 (Levitus, 1982). The monthly minimum and maximum MLD for BCC_GODAS2.0 and WOA01 are shown in Fig. 2. In general, the MLD is deep in winter and shallow in summer for both hemispheres. The minimum MLD in BCC_GODAS2.0 is less than 30 m in the midlatitudes. This is consistent with the WOA01 minimum MLD. In the central tropical Pacific, the MLD is 50-60 m, and in the eastern upwelling region the MLD is relatively shallow (<20 m). A comparison of the minimum MLD with WOA01 illustrates the relatively deep MLD in the Southern Ocean for BCC_GODAS2.0. The possible reasons for the differences include deficiencies in model physics and the temporally constant background error covariance in the 3DVAR. It is well known that the Southern Ocean has abundant eddies. The submesoscale eddy-driven restratification of the mixed layer has not been resolved in the model due to the relatively coarse resolution (1^°× 1^°) at high latitudes. This is the main contribution to the deeper MLD than those of observations. Another possible reason is that the temporally constant background error covariance in the 3DVAR cannot fully capture the variability of eddies in the Southern Ocean, and thus also cannot simulate eddies well.

A comparison of the maximum MLD highlights the deep mixing in the North Atlantic, consistent with the robust North Atlantic deep water formation. The maximum MLDs in the North Pacific are relatively shallow (<100 m), also consistent with the lack of deep-water formation there. Homogenous cold and moderately low salinity water is found east of Japan, indicating the formation of the North Pacific Intermediate Water, consistent with (Talley, 1993). In the Southern Ocean, the observed maximum MLD at 50^°S is approximately 720 m, while the model MLD is over 750 m. Therefore, BCC_GODAS2.0 is generally in good agreement with WOA01.

Figure 2.  Climatology of the monthly minimum and maximum MLD based on the density criterion 0.125 kg m$^-3$ (Levitus, 1982) in BCC_GODAS2.0 and WOA01. Note that the color scales for the minimum MLD and maximum MLD are different. Panels (a) and (c) are for BCC_GODAS2.0 and (b) and (d) are for WOA01. (units: m)

Sea level variability is a good indicator of dynamic ocean processes. Some recent studies have noted that the sea level in the tropical Pacific displays significant interannual-to-decadal variations in response to climate variability (e.g., Fu et al., 2009b; Merrifield et al., 2012; Qiu and Chen, 2012; Zhuang et al., 2013). To evaluate the model's performance over the tropical Pacific, empirical orthogonal function (EOF) analysis is applied to low-frequency (period longer than 13 months) SLAs from both altimeter observations and BCC_GODAS2.0 simulations. The EOF patterns and the associated principal components (PCs) reflect the spatial and temporal variations, respectively. The first EOF mode of the observed SLA accounts for 54.9% of the total variation, and its EOF-1 shows oppositely signed SLAs in the eastern and western basin (Fig. 3a). The EOF-1 of modeled SLAs (Fig. 3b), which explains 45.1% of the total variation, generally captures the observed spatial pattern. The modeled PC-1 tracks the observations quite well (Fig. 3c; correlation coefficient = 0.74, significant at the 95% confidence level), and both show a decadal increasing trend over the past two decades with embedded interannual modulations. A close inspection of the PC-1s reveals that strong sea level oscillations existed in 1994-95, 1997-98, 2006-07 and 2009-10, indicative of the impacts from El Niño-Southern Oscillation (ENSO) events.

Figure 4 further shows the ability of BCC_GODAS2.0 to reproduce sea level trends during 1993-2011. In both BCC_GODAS2.0 and observations, the linear sea level trends decrease zonally from west to east, consistent with the results of the first EOF modes in Fig. 3. In general, BCC_GODAS2.0 satisfactorily reproduces the sea level rises in the western Pacific, and falls in the eastern Pacific. The most significant rises occur at the latitudes of 5^°-15^°N in the northwestern Pacific. The zonally opposite sea level trends are mainly caused by the steady intensification of the trade winds across the tropical Pacific that started in the early 1990s (Merrifield and Maltrud, 2011; Nidheesh et al., 2013). Recent studies have linked the notable sea level trends, together with the easterly trade intensification, to decadal climate modes, especially the Pacific Decadal Oscillation and the decadal variability of ENSO (Merrifield et al., 2012; Zhang and Church, 2012). It is noteworthy that the modeled sea level trends in Fig. 4b have a basin-mean rate of 1.9 mm yr^-1, approximately 25% lower than the mean rate of 2.4 mm yr^-1 in satellite observations, which is probably because BCC_GODAS2.0 does not include the increase in ocean mass resulting from the contributions of ice sheets, melting glaciers and continental waters.

Figure 3.  EOF-1 of low-frequency ($>$13 months) SLAs derived from (a) altimeter observations and (b) BCC_GODAS2.0 simulations for the period 1993-2011. (c) PC-1s of observed and modeled SLA.

Figure 4.  Linear trends of (a) altimetric and (b) BCC_GODAS2.0 simulated sea levels during 1993-2011 \small \parbox[t]13cm(units: mm yr$^-1$).

Figure 5.  The STD of DJF SST during 1990-2011. Panels (a-d) represent SSTs in ERSST, HadISST, the free run and BCC_GODAS2.0, respectively. (units: $^\circ C^2$).

It is well known that the SST in the tropical Pacific shows significant interannual variability, which greatly influences global climate. In this part, the simulated SST of BCC_GODAS2.0 in the tropical Pacific is evaluated. The monthly SSTs from the Hadley Centre Sea Ice and SST dataset (HadISST) on a 1^°× 1^° resolution (Rayner et al., 2003), and the monthly Extended Reconstruction of SST, version 3b (ERSST.v3b), on a 2^°× 2^° resolution (Smith et al., 2008), from 1990 to 2011, are used to validate the analyzed SST.

Spatial distributions of the standard deviation (STD) of the SST in the equatorial Pacific during boreal winter (December-January-February, DJF) are shown in Fig. 5. Improvements in the STD of SST by assimilation are mainly located in the tropical Pacific and Atlantic (Figs. 5c and d). In the observations (Figs. 5a and b), the distributions of the STD are symmetric about the equator, and the maximum centers of the STD are approximately 140^°W. In BCC_GODAS2.0, the spatial distributions of the SST STD closely resemble the observations (Fig. 5d). However, differences between the observed and simulated SST are found near the east coast of South America, while the SST STDs from BCC_GODAS2.0 are weaker than observed, which is probably because the upwelling-favored wind stress used in the model is weaker than in reality, and the model's grid resolution cannot resolve the local bathymetry near the coast of America.

In recent years, two types of El Niño, the canonical El Niño and El Niño Modoki, have been suggested (i.e., Larkin and Harrison, 2005; Yu and Kao, 2007; Ashok et al., 2007; Kao and Yu, 2009; Kug et al., 2009; Wang et al., 2009). The latter is also referred to as the Date Line El Niño (Larkin and Harrison, 2005), central Pacific El Niño (Yu and Kao, 2007), and warm pool El Niño (Kug et al., 2009). The performances of canonical El Niño and Niño Modoki in BCC_GODAS2.0 are presented. Here, canonical El Niño is defined by the Niño3 index, which is the area-averaged SST anomalies (SSTA) in the region (5^°S-5^°N, 150^°-90^°W). El Niño Modoki is identified by the El Niño Modoki Index (EMI), defined by (Ashok et al., 2007) as $$ { EMI}={ [SSTA]}_{ C}-0.5\times { [SSTA]}_{ E}-0.5\times{ [SSTA]}_{ W} , $$ where the square brackets with a subscript represent the area-averaged SSTA over the central Pacific region (C) (10^°S-10^°N, 165^°E-140^°W), the eastern Pacific region (E) (15^°S-5^°N, 110^°-70^°W), and the western Pacific region (W) (10^°S-20^°N, 125^°-145^°E), respectively. Figure 6 shows that the changes in both the Niño3 index and the EMI in BCC_GODAS2.0 agree well with those from HadISST and ERSST. The correlation coefficients of the two El Niño indices between BCC_GODAS2.0 and the two observations exceed 0.9, which are higher than those of the free run (Table 2). The amplitudes of the two El Niño indices in BCC_GODAS2.0 show some differences from the observations. The intensities of the Niño3 index in BCC_GODAS2.0 are somewhat underestimated.

We next explore the possible reasons why the intensities of the Niño3 index in BCC_GODAS2.0 are weaker than those from observations. Figure 7 shows the climatology of the annual SSTs during 1990-2011 in both the simulations and observations. The mean climate state of SSTs in the eastern Pacific in BCC_GODAS2.0 (Fig. 7d) is colder than those from observations (Figs. 7a and b). The warmer climate state of SSTs will strengthen the air-sea interaction and increase the amplitudes of the ENSO warm and cold phases (Zhang and Ding, 2001; Zhang et al., 2008). The difference in the climate state of SSTs between BCC_GODAS2.0 and observations could result in weaker intensities of the Niño3 index in BCC_GODAS2.0. As is well known, the two El Niño indices show significant interannual variability. The Niño3 index in the observations has interannual variability with periods of 5 years and 3 years, while the periods of Niño3 index are 4 years in the simulations. The EMI shows the same interannual variability, with periods of 2.5 years in both the observations and simulations (Fig. 8).

Figure 6.  Time series of DJF-mean (a) Niño3 index and (b) EMI during 1990-2011 with ERSST (green line), HadISST (blue line), free run data (black line), and BCC_GODAS2.0 (red line) (units: ^°C).

Figure 7.  Climatology of annual SSTs during 1990-2011. The black box indicates the region used to calculate Niño3 index. Panels (a-d) represent SSTs in ERSST, HadISST, the free run and BCC_GODAS2.0, respectively (units: ^°C).

Figure 8.  Power spectrum of the monthly (a-d) Niño3 index and (e-h) EMI during 1990-2011. The black curves and red curves represent the calculated spectrum and 90% confidence interval, respectively. The first, second, third and fourth rows represents ERSST, HadISST, the free run and BCC_GODAS2.0, respectively.

Figure 9.  Spatial distributions of composited SSTAs for (a-c) El Niño Modoki I and (d-f) El Niño Modoki II in the boreal autumn (September-November). The first, second, and third row represent ERSST, HadISST and BCC_GODAS2.0, respectively. The composites are calculated from three El Niño Modoki I events (1990/91, 1991/92, and 2002/03), and three El Niño Modoki II events (1992/93, 2004/05 and 2009/10) (units: ^°C).

According to their different impacts on rainfall in southern China and typhoon landfall activity, (Wang and Wang, 2013) separated El Niño Modoki into two subtypes, El Niño Modoki I and II, in which the warm SST anomalies originate from the central tropical Pacific and subtropical northeastern Pacific, respectively. In this study, the simulations of El Niño Modoki I and II in BCC_GODAS2.0 are evaluated. Using the classification of (Wang and Wang, 2013), there are three El Niño Modoki I events (1990/91, 1991/92, and 2002/03) and three El Niño Modoki II events (1992/93, 2004/05 and 2009/10). Figure 9 shows the composite SSTAs for El Niño Modoki I and El Niño Modoki II during the developing autumn (September-October-November) from the three datasets of ERSST, HadISST, and BCC_GODAS2.0. Similar to those in ERSST and HadISST, the warm SSTAs of El Niño Modoki I in BCC_GODAS2.0 are in the central tropical Pacific and are symmetric about the equator. Compared with those in ERSST and HadISST, the warm SSTAs of El Niño Modoki I in BCC_GODAS2.0 are weaker, and their area is smaller. The cool SSTAs are observed in the eastern tropical Pacific (Figs. 9a-c). For El Niño Modoki II, although the warm SST anomalies of BCC_GODAS2.0 show an asymmetric distribution extending and tilting from the subtropical northeastern Pacific to the central equatorial Pacific, as do the observations, their amplitudes are weaker than those in ERSST and HadISST (Figs. 9d-f). In addition, the relationships of the Indian Ocean with the two types of El Niño Modoki are reproduced. (Wang and Wang, 2014) suggested that El Niño Modoki I and II are associated with a positive and negative IOD in boreal autumn, respectively. Here, a positive IOD means warm/cool SSTAs in the western/eastern tropical Indian Ocean (Saji et al., 1999). In BCC_GODAS2.0, the positive and negative IODs are associated with El Niño Modoki I and II, respectively, which is consistent with observations. The amplitude of the positive IOD during El Niño Modoki I in BCC_GODAS2.0 is weaker than observed.

The second-generation Global Ocean Data Assimilation System of the Beijing Climate Center (BCC_GODAS2.0) was implemented in April 2014. The ocean model of BCC_GODAS2.0 is MOM4, which adopts tripolar grids to avoid a singularity at the North Pole. An improved 3DVAR scheme is used in BCC_GODAS2.0 by adopting a recursive filter to speed up computation. Satellite SLA and SST data and in-situ observations, including ARGO and GTSPP, are assimilated in real time. BCC_GODAS2.0 has been run daily in pre-operational mode, driven by six-hourly fluxes from the NCEP. The MLD variation of BCC_GODAS2.0 is generally in agreement with that of WOA01. The simulated SST and SSH successfully depict the interannual variability of the tropical Pacific, which has a great impact on the climate of China. It is very important, for the climate monitoring services of the BCC, to evaluate the performance of the tropical Pacific variability in BCC_GODAS2.0. In the tropical Pacific, the first EOF and PC of low-frequency SLAs from BCC_GODAS2.0 are similar to those from satellite altimeter observations. BCC_GODAS2.0 can also satisfactorily reproduce the sea level rises in the western Pacific, and falls in the eastern Pacific. The STD of the SST in BCC_GODAS2.0 is similar to those of ERSST and HadISST in the tropical Pacific. The Niño3 index and EMI of BCC_GODAS2.0 show high correlation (>0.9) with those of ERSST and HadISST. A power spectrum analysis of the monthly averaged SST shows that the Niño3 index and EMI of BCC_GODAS2.0 have the most significant periods at approximately 5 and 2.5 years, respectively, which is consistent with ERSST and HadISST. BCC_GODAS2.0 can capture the main features of El Niño Modoki I and Modoki II. The performance of the IOD in the BCC_GODAS2.0 SST is similar to that observed in the Indian Ocean and generally agrees with that shown in (Wang and Wang, 2013). The global ocean reanalysis products and real-time assimilation results are being applied in the new generation of the short-term forecast system, BCC_CSM1.1.