Implementation of a GPS-RO data processing system for the KIAPS-LETKF data assimilation system

The Korea Institute of Atmospheric Prediction Systems (KIAPS) has been developing a new global numerical weather prediction model and an advanced data assimilation system. As part of the KIAPS package for observation processing (KPOP) system for data assimilation, preprocessing, and quality control modules for bending-angle measurements of global positioning system radio occultation (GPSRO) data have been implemented and examined. The GPSRO data processing system is composed of several steps for checking observation locations, missing values, physical values for Earth radius of curvature, and geoid undulation. An observation-minus-background check is implemented by use of a one-dimensional observational bending-angle operator, and tangent point drift is also considered in the quality control process. We have tested GPS-RO observations utilized by the Korean Meteorological Administration (KMA) within KPOP, based on both the KMA global model and the National Center for Atmospheric Research Community Atmosphere Model with Spectral Element dynamical core (CAMSE) as a model background. Background fields from the CAM-SE model are incorporated for the preparation of assimilation experiments with the KIAPS local ensemble transform Kalman filter (LETKF) data assimilation system, which has been successfully implemented to a cubed-sphere model with unstructured quadrilateral meshes. As a result of data processing, the bending-angle departure statistics between observation and background show significant improvement. Also, the first experiment in assimilating GPS-RO bending angle from KPOP within KIAPS-LETKF shows encouraging results.


Introduction
Global positioning system radio occultation (GPS-RO; Kursinski et al., 1997) is a limb-geometry remote-sensing technique whereby the time delay of GPS radio signals that have passed through the limb of the Earth's atmosphere are used to determine vertical profiles of measurements related to the refractive index.GPS satellites transmit two microwave signals (1.2 and 1.5 GHz) to receivers on low Earth orbit (LEO) satellites.An occultation occurs when the microwave signals transmitted by one of the GPS satellites as it rises or sets pass through the Earth's atmosphere.During an occultation, the ray connecting the GPS and LEO satellites scans the atmosphere, providing vertical information of the atmosphere from the refraction of the GPS radio signals as measured by the receiver in a low Earth orbit.The raw measurements of radio occultations are the phase and amplitude of radio signals transmitted by the GPS satellites.Based on these measurements and knowledge of the precise positions and velocities of the GPS and LEO satellites, vertical profiles of bending angle and atmospheric refractivity are derived by use of the local spherical symmetry assumption and the Abel inversion (Phinney and Anderson, 1968).The observations have high vertical resolution (0.1 km near surface to 1 km tropopause) and global coverage, even though the horizontal resolution is relatively poor (hundreds of kilometers).Also, they show high accuracy (equivalent to < 1 K; average accuracy < 0.1 K) and precision (0.02-0.05 K) (Anthes, 2011) for a temperature in the vertical range of 10 to 40 km and equal accuracy over either land or ocean (Cucurull et al., 2013).The most powerful benefits of the GPS-RO measurements are no satellite bias and minimal effect on the data by Published by Copernicus Publications on behalf of the European Geosciences Union.1260 H. Kwon et al.: Implementation of a GPS-RO data processing system clouds or precipitation compared to other satellite observations.Because of these benefits, many operational numerical weather prediction centers, such as the Met Office of the United Kingdom, ECMWF, NCEP, Météo-France, Environment Canada, and JMA, have started incorporating GPS-RO soundings into their assimilation systems with clear positive impacts on weather forecasting (e.g., Healy, 2008;Buontempo et al., 2008;Cucurull and Derber, 2008;Aparicio et al., 2009;Rennie, 2010).In particular, GPS-RO data assimilation shows strong sensitivity to upper atmosphere temperature structures, an area that is only weakly constrained by other observations in the analysis and that is prone to large model uncertainties (Anlauf et al., 2011).
The Korea Institute of Atmospheric Prediction Systems (KIAPS) is a government-funded, non-profit research and development institute that was established in 2011 by the Korea Meteorological Administration (KMA).The goal of the KIAPS is to develop the next-generation operational global numerical weather prediction (NWP) system, which can be used for global modeling as well as local areas, particularly optimized to topographic and meteorological features of the Korean Peninsula.The KIAPS has been developing an advanced data assimilation system in addition to a global model (KIAPS integrated model with spectral element dynamic core, or KIM).As one of the data assimilation systems, a local ensemble transform Kalman filter (LETKF) (Hunt et al., 2007) has been successfully implemented for the National Center for Atmospheric Research Community Atmosphere Model with Spectral Element dynamical core (CAM-SE) (Dennis et al., 2012), which has the same grid structure on the cubed sphere as KIM.After a successful evaluation of the KIAPS-LETKF data assimilation system with various observing system simulation experiments (OSSEs) (Kang and Park, 2013), assimilation of real observations of surface and rawinsonde data from NCEP preprocessed data has been performed (Jung et al., 2014).In preparation for GPS-RO data assimilation, preprocessing and quality control modules for bending-angle measurements of GPS-RO data are well implemented in the KIAPS Package for Observation Processing (KPOP) to provide optimal observation for the data assimilation.Finally, we have tested GPS-RO bendingangle data assimilation within the KIAPS-LETKF system to see whether our first version of the GPS-RO data assimilation cycle works well in a coupled system of KIAPS-LETKF and KPOP.
In this paper, we describe the GPS-RO data processing system for bending-angle data assimilation and present preliminary results from bending-angle data assimilation experiments with the KIAPS-LETKF system.In Sect.2, we present the GPS-RO processing system in KPOP.In this section, background ingest and spatial interpolation step, observation operator for bending angle, quality control procedure, and results from the quality control process for bending-angle data assimilation are introduced.In Sect.3, description of analysis cycles within the KIAPS-LETKF system and its preliminary result are presented.Section 4 contains a summary and plans for future work.

Background ingests and spatial interpolation
The observation needs to be compared with the model background for quality control and observation monitoring.
Here, operational global model forecasts of the KMA are used as a model background.The global model is the Unified Model (UM; Davies et al., 2004), which was developed by the UK Met Office.It is a non-hydrostatic, gridpoint model, with the Charney-Phillips grid in a vertical direction.The horizontal resolution is approximately 25 km (N512/ ∼ 0.352 • × ∼ 0.234 • ) and it has 70 vertical levels, with the model top at 80 km.
To produce the simulated observation by use of model background fields, the spatial interpolation of model variables to the observation space is required.We used a bi-linear interpolation method to transform model variables into the observation space in the horizontal direction.In the vertical coordinate of the UM, model information is provided on a staggered height grid, with pressure and density on ρ levels and potential temperature and humidity information on the intermediate θ levels.Therefore, the pressure on model levels where humidity information is stored should be known.The linear interpolation in natural log of pressure is applied for the calculation of the pressure on the intermediate vertical levels.

Observation operator
The purpose of data assimilation is to find an optimal analysis state, depending on the difference between observation and model background (i.e., innovations) and their error statistics.To have such innovations, one should have the observation operator that maps atmospheric variables in the model grid space into the observed variables in observation space (Eyre, 1994).Therefore, the observation operator (forward model) is one of the most important components in the data assimilation system that calculates the simulated observation by use of model variables.
GPS-RO measurements have a two-dimensional limb geometry that, ideally, should be accounted for when they are assimilated into numerical weather prediction systems.Generally, the ray-tracing operator (Gorbunov and Kornblueh, 2003) for the bending angle is reported to lead to best results in terms of standard deviation of the difference between observation and model background (Pingel and Rhodin, 2009).However, the use of a ray-tracing operator requires a high computational cost, so most numerical weather prediction centers (e.g., ECMWF, UK Met Office, DWD, and Météo-France) currently use one-dimensional observation operators for the bending-angle assimilation and inflate the observation errors to partially compensate for this source of forward model error (Anlauf et al., 2011).The assumption of spherical symmetry enables calculation of the bending angle from a one-dimensional integration of the background refractivity lapse-rate profile located at the vertical level of the observation point.Recently, ECMWF has tested the two-dimensional bending-angle operator in their numerical weather prediction system and planned to implement it into the next model cycle.For the same reason as most other operational centers, however, we adhere to a onedimensional bending-angle operator that is included in the Radio Occultation Processing Package (ROPP) developed by the European Organization for the Exploitation of Meteorological Satellites (EUMETSAT) Radio Occultation Meteorology Satellite Application Facility (ROM SAF).The onedimensional bending-angle operator in ROPP solves Eq. ( 1) by use of the algorithm of Healy and Thépaut (2006).It calculates the integral of the bending angle (e.g., Kursinski et al., 1997), given the observed impact parameter a, by use of the interpolated model fields (temperature, pressure, and humidity) at the horizontal location of the tangent point of the observation.
where α is the total ionospheric-corrected bending angle, n is the refractive index derived from the model, and x = nr, where r is distance from the local center of curvature.Assuming that refractivity varies exponentially with x between model levels j and j + 1 and x − a, the bending-angle calculation at a single layer is given by α j = 10 −6 2π ak j N j exp k j x j − a (2) where k j and the error function (erf) are defined as Ray bending above the model top is accounted for by extrapolating the uppermost model parameters and evaluating α top = 10 −6 2π ak j N j exp k j x j − a (5) The bending-angle forward model is composed of the following steps.First, the partial water vapor pressure on each background level is derived from the definition of specific humidity for evaluating refractivity.The refractivity (N ) is then calculated on each background model level by application of the equation given by Bevis et al. (1994) and then interpolated to geopotential heights.For the radio occultation soundings, the scattering term is usually assumed to be negligible, and the ionospheric effects in the equation of refractivity are assumed to have been removed during the preprocessing.Therefore, the refractivity can be expressed (Smith and Weintraub, 1953;Bevis et al., 1994) as where P d is the pressure of the dry air, T is the temperature, e is the partial water vapor pressure, and k 1 (77.607),k 2 (71.6), and k 3 (3.747)are the atmospheric refractivity constants.Using the refractive index n, the impact parameter corresponding to each background level is derived and, finally, the bending angle is computed as a function of the observation impact parameter from the refractivity profile on background impact parameter levels.

Preprocessing and quality control
For the assimilation of GPS-RO bending-angle measurements, preprocessing and quality control modules are implemented into KPOP.Figure 1 is a flowchart that shows the structure and data processing steps of the GPS-RO processing system.As the first step, input/output (I/O) modules for GPS-RO data are implemented to process the binary universal format for data representation (BUFR) data stream of the KMA by use of an ECMWF BUFR decoder.Quality control modules in gross quality control steps are implemented for checking observation locations, missing values for refractivity, bending angles and their errors, physical values for Earth radius of curvature, and geoid undulation.The criteria of physical values correspond with those used in Météo-France.Also, this step involves determining whether the impact height is increasing with geometric height and the identical geolocation for each data point in the profile to screen cases out of rule.The impact parameter is converted to impact height, which is defined as the impact parameter minus the local radius of curvature, for convenience.To take into account the tangent points drift, we make bins in the vertical with 1.5 km intervals and then use the local information of the latitude and longitude for each data point.An observation-minus-background (O-B) check is implemented by the one-dimensional observational bending-angle operator, as mentioned in a previous section, to reject spurious data.This check consists of rejecting data per element of the profile if the O-B value exceeds 5 times the assumed observation error.The observation errors are assigned by applying the error model described in Healy and Thépaut (2006), in which the combination of forward model error and observation error is assumed to vary only with impact height.The standard deviation of the error is assumed to decrease linearly from 10 % of the observed value to 1 % for impact heights from 0 to 10 km, and to stay at 1 % above 10 km, with an absolute lower limit of 3 µrad.

Processing results
We processed 30 days of bending-angle data from 1 to 30 November 2011, which were used for operational data assimilation at the KMA.The GPS-RO data provided by the KMA include measurements from the GRAS instrument on METOP (Loiselet et al., 2000), the constellation of six satellites launched by the Constellation Observing System for Meteorology, Ionosphere and Climate (COSMIC) program (Anthes et al., 2000), and the Communications/Navigation Outage Forecasting System (C/NOFS) Occultation Receiver for Ionospheric Sensing and Specification (CORISS) instrument.Figure 2 shows geographical coverage of GPS-RO profiles over 6 h.The number of total GPS-RO events for 6 h is approximately 500, and the data are shown to be distributed globally over land and ocean.
The KMA global model forecasts are used as model background.The operational global forecast model used at the KMA is the UM (Davies et al., 2004), developed by the UK Met Office.The GPS-RO data were assimilated to produce the analysis with other kinds of observations in operational forecast cycle at the KMA. Figure 3 presents the global mean of departure statistics of bending angle between the observation and background calculations for the month of Novem-ber 2012. Figure 3a and b show the global mean and standard deviation of bending-angle innovations normalized by background and the number of observations as a function of impact height before and after quality control, respectively.Also, Fig. 3c shows the difference of departure statistics between Fig. 3a and b.After the quality control, the mean agreement is best over approximately 5 to 40 km altitude, between −1 % and +1 %, except for the observation from the CORISS instrument.The statistics from CORISS observations show large mean differences between observation and background calculations at the lower atmosphere below approximately 8 km.The increased bias above ∼ 45 km is considered to be associated with a specific model temperature bias of the UM.As we mentioned above, UM forecast is utilized to produce the model background field for data processing.Increased bias of bending-angle innovations (O-B) above ∼ 45 km is also addressed by Burrows et al. (2014).As a reason, they specifically noted the high-altitude temperature bias in the UM.The region located above 40 km altitude presents standard deviations increasing with height and this increase is partly caused by ionospheric contamination.The COSMIC satellites show the smallest global mean O-B value with impact height of other satellites (METOP and C/NOFS).Also, all COSMIC satellites show consistent bending-angle departure statistics from one another.This results from consistency in the data processing.
We compared our results of bending-angle departure statistics to those from the established operational quality control system of the Met Office by Burrows et al. (2014) (not shown).The period (about 1 month) of data processing and the types (satellites) of observations used by the Met Office are similar to our data processing even though the year and month of observation are not exactly matched.The comparison of results of O-B statistics in our study with those of Met Office shows a good agreement in terms of global mean and standard deviation.The results from CORISS instruments show the large value of mean difference of O-B statistics under ∼ 8 km in the Met Office's data processing as well, so they do not use the GPS-RO observations under the height of ∼ 8 km for an operational data assimilation (C.Burrows, personal communication, 2013).
The number of observations has been reduced by 1.07 % of the total number of observations for the entire vertical range after the quality control.The mean value of observation rejection for 10 km vertical intervals ranges from 0.65 % as minimum to 1.61 % as maximum, and such a difference depends on the height.Relatively many observations are rejected under ∼ 20 km. Figure 3d and e show the zonal mean and the standard deviation of the bending-angle departures between observations and background after the quality control.The simulated observations show overall good agreement with one another except over the tropics, where the errors become larger because of large humidity fluctuations.The moist areas in the troposphere tend to have standard deviations in excess of 12 %.For other areas, the standard deviation remains below 4 % except near the tropical tropopause.In general, the standard deviations of the bending-angle departures between observations and background states are slightly larger in the southern high latitudes than in the northern high latitudes.This result is consistent with previous researches (Cucurull et al., 2007;Poli et al., 2009).
We also investigated the bending-angle departure statistics when the tangent point drift is taken into account in GPS-RO data processing.Figure 4 shows the difference of the bending-angle departure statistics between the results from data processing with and without considering the tangent point drift.In data processing, which does not consider the tangent point drift, horizontal location of each data point in a vertical profile is assumed to be fixed, with mean value of latitude and longitude extracted from the BUFR file for all vertical levels.Although no significant difference is seen between the two treatments over most of the area, certain differences have been found at altitudes below 5 km over the tropics in latitude bands.Furthermore, the total number of data points remaining after quality control with tangent point drift consideration is slightly smaller than the other results with 0.08 %.
3 Assimilation of bending angle within the KIAPS-LETKF system

The KIAPS-LETKF system
The forecast model used for the analysis cycles in this study is the CAM-SE, which was developed for climate projection rather than weather prediction.It has a relatively coarse vertical resolution with 30 layers, and its top is near 2.25 Pa (∼ 40 km), which may not be sufficiently high for optimal performance of GPS-RO data assimilation.However, the choice of the forecast model is temporary until the early version of KIM is available for the KIAPS-LETKF system.The reason we chose CAM-SE is because it has the same grid structure as KIM, so the KIAPS-LETKF system implemented to CAM-SE can easily switch the model to KIM without major modifications of the codes.Because KIM has been developed with an option of a much higher top, near 0.01 Pa (∼ 80 km), and its main goal is numerical weather forecast, we believe that our prospective operational settings for GPS-RO data assimilation will be much better with our own model in the near future.In this study, we would like to test whether a coupled system of KIAPS-LETKF and KPOP would work well with GPS-RO bending-angle data and whether our assimilation system produces reasonable increments as the first step of KIAPS-LETKF with GPS-RO data.Horizontal resolution of the forecast model is set at "ne16np4" (∼ 2.  background states in our experiments contain model bias, especially over polar regions.This imperfection of the forecast model is embedded on purpose to see how the data assimilation system works with an obviously forced model bias. For the ensemble data assimilation experiments, we had 30 ensemble members and a 6-hourly analysis cycle.Within the 6 h assimilation window, we considered background ensembles and GPS-RO data every hour, based on a 4D-LETKF algorithm (e.g., Fertig et al., 2007;Harlim and Hunt, 2007) that deals with flow-dependent background uncertainty in space as well as in time.Adaptive multiplicative inflation (Miyoshi, 2011) was implemented, which greatly helped in estimating background error covariance accurately (e.g., Miyoshi and Kunii, 2012;Kang et al., 2012).Within the 4D-LETKF data assimilation, the analysis mean x a and its perturbations X a are determined by Eqs. ( 7) and ( 8) at every grid point.Here, X b (and Y b ) and y o to compute w (l) and P a (l) are considered at every hour in the 4D-LETKF system.State vector of x in the analysis cycle includes the variables of U (zonal wind), V (meridional wind), T (temperature), q (specific humidity), and P s (surface pressure).Then those variables are updated by multivariate data assimilation as described below.
where x b is the ensemble mean of background states, X b is a matrix whose column is each background ensemble perturbation (background state minus ensemble mean of background), , y o is a vector of all the observations, y b and Y b correspond to x b and X b at the observation space, respectively, R is the observation error covariance matrix, which is the assumed diagonal matrix here, and ρ is the adaptive multiplicative inflation factor (Miyoshi, 2011).Subscript (l) indicates the information within a local area centered at the analysis point.In this study, we used a Gaussian function with 500 km of a standard deviation for the horizontal localization function while setting 1800 km as the cutoff distance.For the vertical function, we have another Gaussian function with 0.2σ of standard deviation and 0.73σ of cutoff.The LETKF simultaneously assimilates all observations within the distance just described in space and within the 6 h assimilation window (Hunt et al., 2007).
To assimilate GPS-RO bending-angle data within the KIAPS-LETKF system, we ran ROPP implemented in KPOP, as an observation operator, as many times as the ensemble size.Then, background ensembles at the GPS-RO observation space Y b were provided for the computation of Eqs. ( 7) and (8).

Data processing with CAM-SE backgrounds
For the bending-angle data assimilation within the KIAPS-LETKF system, we first processed the GPS-RO data with hourly background ensembles from the CAM-SE forecast within KPOP system.All processing steps described in Sect.2.3 are applied for the GPS-RO data processing with CAM-SE background.The data processing for ensemble data assimilation was done using each member of background ensembles of CAM-SE forecasts from every 6 h analysis for 2 weeks.The observations passed through the KPOP processing with all 30 ensemble backgrounds are used for www.atmos-meas-tech.net/8/1259/2015/ the bending-angle assimilation within the KIAPS-LETKF system.The background ensembles at the initial time of bending-angle data assimilation cycles are the forecasts from the analysis assimilating sonde and surface pressure station data only.Once the the data assimilation cycle with RO data (EXP_RO cycle) is initiated, the background ensembles used for the data processing are CAM-SE forecasts starting from the analysis assimilating bending angle in addition to sonde and surface pressure data.
Figure 5 presents departure statistics between the observation and background of the bending angle from 15 to 30 November 2012.The background of the bending-angle results from an ensemble mean of 30 background states.Figure 5a and b show the global mean of departure statistics normalized by background and the number of observations as a function of impact height before and after quality control, respectively.Also, Fig. 5c shows the difference of departure statistics between Fig. 5a and b. Figure 5d and e show the zonal mean and the standard deviation of bending-angle departures between observations and background after the processing.The data processing clearly improves the bendingangle departure statistics in terms of global mean and standard deviation for the vertical range from 10 to 35 km.The global mean departure statistics show best agreement between observation and background over approximately 5 to 30 km of the impact height.The zonal mean feature of departure statistics shows large value over the tropic regions.About 13 % of the total number of observation is rejected during the data processing.The number of rejected observations widely ranges from 2.9 to 32 % in the vertical profile.Note that the number of observations around ∼ 33 km (impact height) was considerably reduced (∼ 75 % of the total number of observations for the vertical range of 32-33 km) compared to other levels.This is mainly due to the large difference of bending angle between observation and background for that vertical range before the data processing.This could be affected by the bending-angle calculation using the one-dimensional observation operator described in Sect.2.2 around the model top with a large vertical spacing of the CAM-SE.Another possible reason can be ascribed to a quality of the CAM-SE background for those vertical levels.CAM-SE is a climate model that may not be optimal for the purpose of NWP although we chose it for a test model due to the same grid structure as KIM.Therefore, we expect that replacing the forecast model by KIM would solve this issue.
When compared with the results obtained from use of the KMA forecast as a background, as presented in Sect.2.4, the departure statistics show a larger variability with the CAM-SE background than with the KMA background in terms of global and zonal mean features of statistics.Also, the number of observations is significantly reduced compared to results of KMA forecasts after the data processing.This again implies that the CAM-SE background has poorer forecasts of atmospheric conditions than the KMA background.If we recall the lower horizontal and vertical resolution of the CAM-SE model and some limitations of climate models for good weather prediction as mentioned in a previous section, the finding looks reasonable.Another factor is that the KMA backgrounds already include the operational data assimilation effects of GPS-RO and many other observations on the analysis, whereas our CAM-SE background is the forecast from the analysis assimilating sonde, surface pressure station data, and GPS-RO bending angle only.
We also see the ensemble spread of background states in the observation space of bending angle in Fig. 6.Because the initial background ensembles for the RO data assimilation experiment come from background states of analyses assimilating only sonde and surface pressure data, the ensemble spread over the Northern Hemisphere is smaller than over the Southern Hemisphere.That is, meteorological variables are relatively well constrained by sonde data over the Northern Hemisphere, so the ensemble spread (standard deviation of ensemble perturbations), representing background uncertainty, has small values there.Therefore, we expect that the impact of GPS-RO bending-angle data would be significant over the regions with large ensemble spread (e.g., Southern Hemisphere) when the observation error of bending-angle data is relatively small.

Assimilation experiments and results
We conducted two experiments: a 1-month (November 2011) data assimilation experiment that included sonde and surface pressure station data only (hereafter referred to as CTRL_SONDE) and a 2-week data assimilation experiment that starts from ensemble background after 2 weeks of CTRL_SONDE analysis (15 November 2011) but starts assimilating GPS-RO bending-angle data in addition to sonde and surface pressure station data (hereafter referred to as EXP_RO).Thus, we would like to see if EXP_RO gives reasonable analysis increments compared with CTRL_SONDE, for the last 2 weeks of the analysis period.
Figure 7 shows a vertical cross section of zonally averaged analysis increment differences of U and T between CTRL_SONDE and EXP_RO.We found significant increment as a result of adding GPS-RO data above the upper troposphere to the model top in the Southern Hemisphere.Because we made EXP_RO start from ensemble background after a 2-week analysis of CTRL_SONDE, background states have been well constrained near the surface and in the Northern Hemisphere at the beginning of EXP_RO (Fig. 6), so the remarkable increment over the upper level of the Southern Hemisphere is reasonable.Also, analysis increments appear in the lower troposphere, especially over high latitudes, because the forecast model has serious temperature bias introduced by the inactivated sea-ice model, as described at the beginning of Sect.3.1.That is, departure of background bending angle from the observed bending angle is large for the bending-angle data to correct the background, although the data have a relatively large error over those levels.Indeed, Fig. 6 shows a large background spread at the observation space, which makes the analysis reflect observation effectively.Thus, it shows the additional correction caused by RO data, even in the lower troposphere, where the background states exhibit serious forecast bias.Furthermore, states and errors between two experiments are transferred by the ensemble forecast, so the difference can be found in the lower troposphere through repetition of forecast-analysis cycles.
Horizontal analysis increment at the level where many bending-angle data of GPS-RO are assimilated is presented in Fig. 8.As expected, significant increments occur where sonde data do not sufficiently constrain the background.Since the initial ensembles of EXP_RO were already constrained by conventional data for 2 weeks of CN-TRL_SONDE, analysis increments look dominant where the conventional data do not exist, especially for the 2-week analysis period of EXP_RO.Recall that the spinup period usually shows significant analysis increments where observations are newly assimilated.We expect that the system will show more comparable increments even in the Northern Hemisphere as the forecast-analysis cycles are repeated.
To examine a posteriori uncertainties between CTRL_SONDE and EXP_RO, we looked at the relative uncertainty reduction in Fig. 9 through the following equation: where x a sonde and x a ro indicate ensemble analysis states from CTRL_SONDE and EXP_RO, respectively, and x a sonde and x b ro are the ensemble means of x a sonde and x a ro , respectively, and E[ ] indicates ensemble mean.As a result, positive val-ues mean that analysis uncertainty estimated in EXP_RO is smaller than that in CTRL_SONDE and vice versa.Thus, the values of Fig. 9 indicate how much reduction of the estimated analysis uncertainty has been introduced by adding GPS-RO data into the KIAPS-LETKF system.
We found significant uncertainty reduction estimated over the areas where the analysis increments are shown in Fig. 8.This illustrates that GPS-RO data are assimilated and reduce uncertainties of background and analysis over the areas with many data and relative inaccuracy of background states.Because we used the adaptive multiplicative inflation method (Miyoshi, 2011), which computes inflation parameters in a way that has large inflation where and when innovation is large, to avoid underestimation of background uncertainty, background/analysis tends to have greater inflation factors than unity.In contrast, the multiplicative inflation parameter is not adaptively estimated when there is no observation (no O-B information), and a very small inflation factor (only 1 % inflation) is set.Therefore, background states in EXP_RO tend to be inflated more than those in CTRL_SONDE.Despite greater inflation factors multiplied in EXP_RO than in CTRL_SONDE, resultant analysis of EXP_RO shows significant reduction of analysis spread as a result of assimilating additional data.Before verifying our analysis, we confirmed that GPS-RO bending-angle data are assimilated effectively for the first test with the KIAPS-LETKF system in a way that reduces analysis uncertainty where the data are expected to contribute.
Finally, we compared our analysis with ERA interim reanalysis (Dee et al., 2011) using the equation x a sonde − ERA − x a ro − ERA in Fig. 10, which gives positive values when the analysis of EXP_RO is closer to ERA interim than CTRL_SONDE and vice versa.Even though the forecast model used for this study has difficulties in making good weather prediction, as described in Sect.3.1, we still expect a meaningful correction to be introduced by additional GPS-RO bending-angle data if our analysis method works.Thus, we have computed that equation of analysis improvement for the 2-week EXP_RO data assimilation period.Figure 10 shows encouraging results.We have obtained positive global mean values for both U and T , which represent the positive impact of GPS-RO data on our analysis.We found that the variable of T has much better fit to ERA interim data than that of U after assimilating GPS-RO bending-angle data.This is a reasonable result because bending-angle data have been known to show stronger sensitivity to the temperature than wind.Based on this fact, we may be able to try the "localization of the variables" for even better analysis within the LETKF system, which allows the assimilation of GPS-RO bending-angle data to update only temperature (Kang et al., 2011).
We also took a look at the vertical profiles of analysis improvement in a comparison with ERA interim data in Fig. 11.It shows significant error reduction introduced by adding GPS-RO bending-angle data overall for 2 weeks of  EXP_RO.There are considerable corrections of errors in upper level wind and temperature.In addition, we could apparently find positive impact of RO data even in the lower troposphere, especially over the polar region where there is forecast imperfection due to the inactivated sea-ice model.Figure 12 shows much greater improvement caused by GPS-RO bending-angle data at the level of 20 hPa than the level of 100 hPa, and the global mean of error reduction looks remarkable for both variables of U and T .Those results are from the first version of GPS-RO bending-angle data assimilation within the KIAPS-LETKF system coupled with KPOP.Thus, we want to prove our first achievement in this paper.In the meantime, we continue improving our current system (such as tuning parameters within LETKF, replacing the forecast model by a better one, etc.) so that it could give comparable results to those of other centers (e.g., Cucurull et al., 2007;Anlauf et al., 2011) in the near future.

Summary and future work
A preprocessing and quality control system for bendingangle measurements is successfully implemented in the  KPOP system for GPS-RO bending-angle data assimilation.The preprocessing and gross quality control procedure consists of checks for observation locations, missing values, and physically consistent values for Earth radius of curvature and geoid undulation.An observation-minus-background check is implemented by use of a one-dimensional observational bending-angle operator that is included in the ROPP software.With this GPS-RO data processing system, we have processed bending-angle observations with background states produced by both KMA and CAM-SE forecasts to investigate the bending-angle departure statistics between observation and background.After the data processing, the mean departure statistics show better agreement between the observation and background of the bending angle.In particular, it is more apparent in the use of a CAM-SE background.Overall, the processing results show reasonable quality of departure statistics, which is consistent with previous studies (Cucurull et al., 2007;Poli et al., 2009;Rennie, 2010;Burrows et al., 2014).
We have tested GPS-RO bending-angle data processed by KPOP within the KIAPS-LETKF system.We conducted two experiments, one assimilating only sonde and surface pressure station data and the other assimilating GPS-RO bendingangle data in addition to sonde and surface pressure station data.Comparison between the two experiments shows remarkable difference of analysis increment over high altitudes, oceans, and the Southern Hemisphere, which explains the characteristics of additional data.In the KIAPS-LETKF system, we have obtained reasonable reduction of analysis uncertainty, and the positive impact of bending-angle data has been presented in the comparison with ERA interim reanalysis.What we have shown in this paper resulted from the first test assimilating the GPS-RO bending angle within the KIAPS-LETKF system.Based on this work, we will improve the KPOP and KIAPS-LETKF systems in various ways.We  have started implementing the KIAPS-LETKF system in a much more advanced forecast model, KIM, which has much finer horizontal and vertical resolution with a higher top.This change will bring major improvement to the GPS-RO bending-angle data assimilation system.We can try the "localization of the variables" (Kang et al., 2011) to see if the results of bending-angle data assimilation would be improved.Moreover, in the coupled cycle of the KPOP and KIAPS-LETKF systems, we will try proactive quality control and estimate observation error by the advanced method of Hotta (2014), which estimates short-range forecast sensi-tivity to the observation data and to observation error within LETKF.Also, we will try to use the information from the prior ensemble statistics for the quality control of bending angle within the KPOP system.In addition, we plan to apply a two-dimensional bending-angle operator suggested (Healy et al., 2007) for advanced bending-angle preprocessing and data assimilation.

Figure 1 .
Figure 1.Schematic diagram of the GPS-RO data processing system implemented in KPOP.

Figure 2 .
Figure 2. Six-hour coverage of the GPS radio occultation events on 7 November 2012.The total number of profiles is 591.

Figure 3 .
Figure 3. Departure statistics between the observation and background calculations of bending angles for the month of November 2012.(a and b) Global mean bending-angle innovation and number of observations as a function of impact parameter before and after quality control, respectively.(c) Difference of departure statistics between (a and b) (a minus b).Statistics are calculated based on satellite data and displayed with different colors.(d and e) Zonal mean and standard deviation of bending-angle innovation for all satellites after quality control.

Figure 4 .
Figure 4. Difference of the bending-angle departure statistics for global (a) and zonal mean (b) between the results from data processing considering (TPD) and not considering (no_TPD) the tangent point drift (TPD minus no_TPD).

Figure 5 .Figure 6 .
Figure 5. Same as Fig. 3 for the month of November 2012 except for the use of CAM-SE ensemble backgrounds.

Figure 7 .
Figure 7. Difference of zonally averaged analysis increments of U (left, unit: m s −1 ) and T (right, unit: K) between CTRL_SONDE and EXP_RO for a 2-week analysis from 15 November 2011.

Figure 8 .
Figure 8. Difference of horizontal analysis increment of U (left, unit: m s −1 ) and T (right, unit: K) between CTRL_SONDE and EXP_RO at the level of 100 hPa for a 2-week analysis from 15 November 2011.

Figure 9 .
Figure 9. Relative uncertainty reduction of EXP_RO from CTRL_SONDE, computed by Eq. (9) for the variable U (left) and T (right) at the 100 hPa level for a 2-week analysis period.

Figure 10 .
Figure 10.Improvement of EXP_RO analysis from CTRL_SONDE, toward ERA interim, computed by the equation x a sonde − ERA − x a ro − ERA for U (left, unit: m s −1 ) and T (right, unit: K) at the level of 100 hPa.Positive values indicate analysis of EXP_RO closer to ERA interim data than CTRL_SONDE and vice versa.Global mean values are present at the top left of each figure.

Figure 11 .
Figure 11.Vertical cross section of zonally averaged improvement of EXP_RO analysis from CTRL_SONDE, toward ERA interim, computed by the equation x a sonde − ERA − x a ro − ERA for U (left, unit: m s −1 ) and T (right, unit: K).Positive values indicate analysis of EXP_RO closer to ERA interim data than CTRL_SONDE and vice versa.

Figure 12 .
Figure 12.Same as Fig. 10 except for the vertical level of 20 hPa.