Guofang Chao , Xinrong Wu , , Lianxin Zhang , Hongli Fu , Kexiu Liu , Guijun Han

a College of Oceanic and Atmospheric Sciences, Ocean University of China, Qingdao, China

b Key Laboratory of Ministry of Natural Resources for Marine Environmental Information Technology, National Marine Data and Information Service, Ministry of Natural Resources, Tianjin, China

c School of Marine Science and Technology, Tianjin University, Tianjin, China

Keywords:China Ocean ReAnalysis (CORA)Validation Multigrid 3D-Var assimilation

A B S T R A C T China Ocean ReAnalysis (CORA) version 1.0 products for the period 2009—18 have been developed and validated.The model configuration and assimilation algorithm have both been updated compared to those of the 51-year(1958—2008) products. The assimilated observations include temperature and salinity field data, satellite remote sensing sea surface temperature, and merged sea surface height (SSH) anomaly data. The validation includes the following three aspects: (1) Temperature, salinity, and SSH anomaly root-mean-square errors (RMSEs) are computed as a primary evaluation of the reanalysis quality. The 0—2000 m domain-averaged RMSEs of temperature and salinity are 0.61°C and 0.08 psu, respectively. The SSH anomaly RMSE is less than 0.2 m in most regions.(2) The 35°N temperature section is used to evaluate the ability to reproduce the thermocline, mixing layer, and Yellow Sea cold water mass. In summer, the thermocline is reinforced, with the gradient changing from 3°C in May to 10°C in August. The mixing-layer depth reproduced by CORA is consistent with that computed from the observed climatology. The Yellow Sea cold water mass forms at a depth of 50 m. (3) The reanalysis current is examined against the tracks of some drifting buoys. The results show that the reanalysis current can capture the mesoscale eddies near the Kuroshio, which are similar to those described by the drifting buoys. Overall, the 2009—18 CORA reanalysis products are capable of reproducing major oceanic phenomena and processes in the coastal waters of China and adjacent seas.

1. Introduction

Ocean reanalysis utilizes a combination of historical ocean observations, ocean dynamic models, and data assimilation technology to reconstruct the long-term and multiscale spatiotemporal variations in the historical state of the ocean. Ocean reanalysis products can provide historical ocean-environment information for climate change research,marine disaster prevention and mitigation, marine environmental assurance, and so on.

Ocean reanalysis research has attracted great attention in several international marine organizations, including the World Ocean Circulation Experiment, the Climate Variability and Predictability Programme, and the Global Ocean Data Assimilation Experiment. The United States, European Center for Medium-Range Weather Forecasts(ECMWF), France, and Japan have all developed their own ocean reanalysis products (see Table 1 ). Meanwhile, the development and application of ocean reanalysis products have been promoted in China. In 2009, the National Marine Data and Information Service (NMDIS) took the lead in independently developing and publicly releasing 23-year(1986—2008) China Ocean ReAnalysis (CORA) trial versions ( Han et al.,2011 ). Subsequently, the 51-year (1958—2008) Global and NorthwestPacific CORA v1.0 products ( Han et al., 2013a , b ; Fan et al., 2019 ) were developed and released with improved ocean numerical models and data assimilation methods. The South China Sea Institute of Oceanology (SCSIO), Chinese Academy of Sciences, produced the Reanalysis Dataset of the South China Sea (REDOS). The global ocean reanalysis products of the Coupled Ocean Model developed by the First Institute of Oceanography (FIOCOM) had a time range of two years. The Ocean University of China and the Institute of Atmospheric Physics, Chinese Academy of Sciences, have also made great contributions to ocean reanalysis ( Zhang et al., 2020 ). To enhance the quality of ocean reanalysis, NMDIS has continuously advanced ocean numerical models and assimilation algorithms to produce operational ocean reanalysis products (CORA v1.0, 2009—18).

Table 1 Current ocean reanalysis products.

The structure of this paper is as follows: Section 2 introduces the 2009—18 CORA v1.0 ocean reanalysis system, followed by presentation of the validation and conclusions in Sections 3 and 4 , respectively.

2. CORA v1.0 ocean reanalysis system

2.1. Ocean model

The CORA v1.0 Northwest Pacific Ocean reanalysis system adopts the Princeton Ocean Model with a generalized coordinate system(POMgcs), which applies internal and external mode separation technology. According to Courant—Friedrichs—Lewy conditions ( Courant et al.,1967 ), different time steps can be used for internal and external modes.The mixed coordinates of “upper

z

and lower sigma ”are applied in the model, which not only provide a high vertical resolution (

z

coordinate) near the thermocline but also prevent the step effect at the bottom boundary and shallow water. The horizontal grid is variable, in which the resolution in the Kuroshio area is 1/8°and gradually extends outward to 1/2°. There are 35 vertical layers, spanning from 0 m to 5500 m. To improve the simulations in the upper ocean, local mixing processes, including wave and sea spray effects, are introduced into the Mellor—Yamada 2.5 turbulence closure scheme and the air—sea interfacial turbulence flux model ( Zhang et al., 2017 , 2018 , 2019 ). The open boundary conditions of temperature, salinity, current, and sea level are provided by CORA v1.0 global ocean reanalysis products. Sixteen tidal components are added to introduce the tidal mixing effect. The harmonic constants of tidal elevations and 2D tidal currents are provided by the Oregon State University TOPEX/Poseidon Global Inverse Solution tidal model, version 7. The NCEP (version 1.0) meteorological reanalysis field is used as the meteorological forcing. The variables include 2 m air temperature, 10-m wind, sensible heat flux, latent heat flux, net shortwave radiation, and net longwave radiation.

2.2. Assimilation method

The 2009—18 CORA v1.0 still adopts the multigrid three-dimensional variational (3D-Var) ocean data assimilation method developed by NMDIS ( Li et al., 2008 ; Xie et al., 2005 ), and the temperature-salinity (

TS

) consistency adjustment algorithm based on the instant model

T-S

data table ( Troccoli et al., 2002 ). In the multigrid 3D-Var, the model background is first interpolated to the observation position and subtracted from the observational value to compute the observational innovation.Then, the multiscale information (from longwave to shortwave) of the observational innovation is extracted by refining the analysis grid from coarse to fine. For each grid level, the background term of the 3D-Var is replaced by a smoothing term that penalizes the horizontal secondorder gradient of the analysis field. Thus, there is no background error covariance matrix

B

in the multigrid 3D-Var. The observational error covariance matrix

R

is set to the identity matrix, which has the same magnitude as the smoothing matrix. In the configuration of the CORA model, POMgcs and the observing system, the number of grid levels is set to nine. Each grid level assimilates all kinds of observations. Compared to the traditional 3D-Var, the multigrid 3D-Var mainly has two advantages. One is that it can artificially extract multiscale observational information. The other is that its convergence speed is much faster than the traditional 3D-Var.

In terms of the assimilation of satellite remote sensing observations,the Modular Ocean Data Assimilation System (MODAS) ( Fox et al.,2002 ) is used to retrieve the “pseudo observations ”of 3D temperature and salinity instead of the water column adjustment algorithm used in the 1958—2008 CORA v1.0. The idea of MODAS is briefly introduced as follows: Three statistical regressions between sea surface temperature (SST) and subsurface temperature, between the sea surface height(SSH) anomaly and subsurface temperature, and between both the SST and SSH anomaly and the subsurface temperature, are first obtained using all historical temperature and salinity observation profiles. For the latter two regressions, the SSH anomaly is approximated by the steric height anomaly computed by vertically integrating the profile that has both temperature and salinity observations. Then, the relationship between temperature and salinity for each standard depth is built using the historical profiles. Once the models are established, one can input the satellite-observed SST and SSH anomalies to reconstruct the subsurface 3D temperature and salinity fields, which have the same horizontal resolution as the satellite observation (i.e., 1/4°) and are different from other in-situ observations, such as ARGO profiles. The assimilation time window is compressed from seven days to one day to increase the SSH anomaly signal of the reanalysis. The main parts of the assimilation process are introduced as follows:

(1) Multigrid 3D-Var assimilation: With a set of grids from coarse to fine,the observational increment relative to the background field is used to conduct the 3D-Var analysis successively. For each analysis, the analysis field obtained from the previous coarse grid is input into the next grid level as the background field, and the increment of each analysis relative to the previous grid level analysis is computed.Finally, the analysis increments of all grid levels are collected to obtain the ultimate analysis results.

(2) Adjustment of

T-S

consistency: Before data assimilation, the oneto-one relationship between temperature and salinity is established using the model results. That is, one temperature value corresponds to one salinity value. Once observations of temperature profiles are assimilated, the 3D salinity field is diagnosed by the temperature analysis using the

T-S

relationship, acting as the background field of the salinity profile assimilation. This method can not only effectively assimilate the temperature and salinity observations but also retain the

T-S

relationship of the model as much as possible.(3) Satellite altimetry assimilation: This process contains the following two steps. First, satellite-observed SST and sea level anomaly data are used to retrieve the subsurface 3D temperature and salinity using the MODAS system. Then, the retrieved “pseudo observations ”of temperature and salinity are assimilated to the numerical model with the multigrid 3D-Var method and the above

T-S

adjustment.Compared to the water column adjustment algorithm used in the previous version, MODAS is independent from the numerical model;therefore, it avoids introducing model error into the “pseudo observations ”of temperature and salinity.

2.3. Observation and quality control

The assimilated observations contain in-situ and satellite remote sensing data. The in-situ temperature and salinity profile observations come from the World Ocean Database 2013 (WOD13), Global Temperature and Salinity Profile Programme (GTSPP), and the Array for Realtime Geostrophic Oceanography (ARGO). The satellite remote sensing observations include the SST from the Advanced Very High Resolution Radiometer (AVHRR) and the sea level anomaly from Archiving, Validation, and Interpretation of Satellite Oceanographic (AVISO) data.

The process of in-situ observations includes routine quality control,such as duplication removal, landing inspection, climatological boundary examination, spike inspection, and stability testing, as well as salinity drift calibration of the ARGO data ( Wong et al., 2003 ; Owens and Wong, 2009 ).

2.4. CORA v1.0 production procedure

The production procedure of the 2009—18 CORA v1.0 products includes the following steps: First, NCEP meteorological reanalysis and CORA v1.0 global ocean reanalysis products are regridded to the POMgcs model to form the meteorological forcing and boundary field.Then, the POMgcs model is restarted from the 1958—2008 CORA v1.0 on 31 December 2008. Afterwards, observations are assimilated into the POMgcs every day with the aforementioned assimilation method to generate the reanalysis products.

The domain of CORA v1.0 is (10°S—52°N, 99°E—150°E), with a horizontal resolution of 1/8°—1/2°. There are 35 vertical layers from 0 m to 5500 m. The reanalysis variables include sea level, temperature, salinity, and current. The products have been released to CORA( http://cora.nmdis.org.cn ), Odinwestpac ( http://odinwestpac.org.cn ),the National Marine Data Center ( http://mds.nmdis.org.cn ), and Centre for Marine-Meteorological and Oceanographic Climate Data(CMOC)/China ( http://www.cmoc-china.cn ).

3. CORA v1.0 reanalysis product validatio n

Referring to the validation methods in the Ocean Reanalyses Intercomparison Project ( Balmaseda et al., 2015 ), CORA v1.0 products from 2009—18 are validated from the following three aspects: (1) the error distribution and statistics of temperature, salinity, and SSH anomaly data;(2) the vertical temperature structure in the 35°N section, and (3) surface current validation.

3.1. Temperature, salinity, and SSH anomaly

The 2009—18 CORA v1.0 products are interpolated into the observed spatiotemporal position of temperature and salinity above 2000 m to calculate the error statistics. To evaluate the quality of the reanalyzed SSH, the annual climatology of the POMgcs is subtracted from the reanalyzed SSH to obtain the SSH anomaly. Fig. 1 shows the horizontal distribution of root-mean-square errors (RMSEs), where the maximum values are 3.52°C, 0.46 psu, and 0.38 m for temperature, salinity, and the SSH anomaly, respectively. In most regions, the RMSEs are less than 1°C, 0.2 psu, and 0.2 m for temperature, salinity, and the SSH anomaly, respectively. Large errors in temperature and salinity occur in the nearshore and coastal regions, while large SSH anomaly errors are located near the Kuroshio Extension. Fig. 2 shows the corresponding time series of the RMSEs of temperature and salinity above 2000 m, in which most of the RMSEs are less than 0.5°C and 0.05 psu for temperature and salinity,respectively.

Fig. 1. Horizontal distribution of the RMSEs of temperature (units:°C, left panel), salinity (units: psu, middle panel) above 2000 m, and the SSH anomaly (units: m,right panel) of the reanalysis products.

Fig. 2. Time series of the RMSEs of (a) temperature (units:°C) and (b) salinity(units: psu) above 2000 m.

Fig. 3. Time series of the daily CORA and AVISO SSH anomaly (upper panel),and the correlation coefficients between them (lower panel), in which the black line indicates the linear trend of the correlation.

Fig. 4. Vertical distribution of the monthly RMSEs of (a) temperature (units:°C) and (b) salinity (units: psu) in the upper 2000 m.

Fig. 5. Temperature climatology of CORA v1.0 in the 35°N section for February, May, August, and November.

Fig. 6. CORA v1.0 monthly averaged surface current in the Taiwan Strait in 2018. The black dots indicate the trajectories of surface drifting buoys.

To further illustrate the relationship between the reanalysis products and the AVISO observed sea level anomaly, the time series of the daily CORA and AVISO SSH anomaly and the correlation coefficient between them are shown in Fig. 3 . The basic trend of the daily CORA and AVISO SSH anomaly is consistent, while the big difference of the daily domainaveraged SSH anomaly in the upper panel is caused by the nearshore remote sensing observational errors and reanalysis errors. The average correlation coefficient reaches 0.7, with a maximum of 0.9. The linear fitting results show that the correlation steadily increased from 2009 to 2018.

Fig. 4 presents the vertical distribution of temperature and salinity RMSEs for each month. The 0—2000 m averaged RMSEs are 0.61°C and 0.08 psu for temperature and salinity, respectively. The maximum error of temperature is 1.11°C, which occurs at a depth of 150 m. The maximum error of salinity is 0.17 psu, which occurs at the sea surface in summer.

3.2. The 35°N temperature section

Fig. 5 displays the temperature climatology of CORA v1.0 in the 35°N section for February, May, August, and November. The reanalysis results exhibit good agreement with the observational climatology( Chen et al., 1993 ). In February, because of the wind change, the sea surface is cooled, and the temperature structure becomes more vertically uniform with no thermocline. In May, due to the change in wind direction and the weakness of the wind force, the sea surface is warmed,and the thermocline and the cold center at the bottom begin to form.Entering summer (August), the hydrological situation of the Yellow Sea has changed significantly compared with that in winter and spring. The thermocline is reinforced, with the gradient increasing from 3°C in May to 10°C in August. The depth of the mixing layer reaches 10 m, which is consistent with the observation. The Yellow Sea cold water mass forms under 50 m. In November, due to the gradually increasing wind velocity,the upper mixing layer rapidly deepens to 30 m, while the thermocline gradually disappears. In the bottom layer, the cold water mass obviously decreases, and the cold center of the North Yellow Sea disappears.

3.3. Sea surface current

Fig. 6 plots the monthly averaged surface current in the Taiwan Strait in 2018, where the black dots indicate the trajectories of surface drifting buoys. We can see that the Kuroshio near the Taiwan Strait bends slightly northward first and then extends southeast or east, which is consistent with other studies (e.g., Su and Yuan, 2005 ). The structures of sea surface flow, mesoscale eddies, and strong currents near the Kuroshio captured by the reanalysis products are also similar to those reflected by the drifting buoys. The Kuroshio strengthens from January to May,forming mesoscale eddies adjacent to Luzon Island. The drifting buoys follow the sea surface current into the South China Sea through the Luzon Strait when the Kuroshio weakens from June to December.

4. Conclusions

This paper briefly introduces the ocean reanalysis system used to generate the 2009—18 CORA v1.0 Northwest Pacific Ocean reanalysis products, and evaluates the quality of the products. Compared to the 1958—2008 CORA v1.0 products, local mixing processes, including wave and sea spray effects, have been introduced into the 2009—18 CORA v1.0 reanalysis system. To increase the SSH anomaly signal of the reanalysis, MODAS is applied to substitute the water column adjustment algorithm, and the assimilation time window is compressed from seven days to one day. The validation results show that: (1) The RMSEs of temperature, salinity, and the SSH anomaly are mostly less than 1°C, 0.2 psu,and 0.2 m, respectively. The total RMSEs of these three variables are 0.61°C, 0.08 psu, and 0.04 m. (2) The reanalysis can capture the structure and variation in the 35°N temperature section in the North Yellow Sea. (3) The sea surface flow, mesoscale eddies, and strong current near the Kuroshio can also be properly reproduced by the reanalysis products.Overall, the new reanalysis is shown to be realistically representative of the ocean state in the Northwest Pacific.

Funding

This project was supported by grants from the National Key Research and Development Program of China [grant numbers 2016YFC1401800 ,2017YFC1404103 , 2016YFC1401701 , and 2019YFC1510000 ], the National Natural Science Foundation of China [grant number 41976019 ], and the Tianjin Natural Science Foundation [grant number 18JCQNJC01200 ] .

Acknowledgments

The

T-S

profile observation data were obtained from WOD13,maintained by NODC, USA, the GTSPP project (Coriolis Data Center),and the ARGO global data center ( ftp://ftp.ifremer.fr ). The satellite remote sensing temperature data were from AVHRR ( https://www.avl.class.noaa.gov/release/data_available/avhrr/index.htm ). The SSH anomaly data were from AVISO ( https://www.aviso.altimetry.fr/data.html ). The atmospheric forcing data were from the NCEP reanalysis product (ftp.cdc.noaa.gov/pub/datasets/necp. reanalysis2).The surface drifter buoy data were from the Global Drifter Program( http://osmc.noaa.gov/erddap/tabledap/gdp_interpolated_drifter.html ).