Assimilation of GPS Radio Occultation Data

1. Assimilation of GPS Radio Occultation Data Ying-Hwa Kuo UCAR COSMIC Office NCAR MMM Division

2. Outline • Characteristics of GPS radio occultation observation • Factors affecting the results of GPS RO assimilation • Practical considerations for the GPS RO assimilation: • Choices of assimilation variables, • Observation operators, • Choices of Data assimilation systems (e.g., 3D-Var, 4D-Var, EnKF) • Determination of observational errors • Data quality control, … etc • Review of GPS RO assimilation and impact studies • Comparison between WRF 3D-Var and WRF/DART • Comparison of local and nonlocal observation operators

3. The velocity of GPS relative to LEO must be estimated to ~0.2 mm/sec (velocity of GPS is ~3 km/sec and velocity of LEO is ~7 km/sec) to determine precise temperature profiles

4. The velocity of GPS relative to LEO must be estimated to ~0.2 mm/sec (20 ppb) to determine precise temperature profiles

5. Characteristics of GPS RO Data • Limb sounding geometry complementary to ground and space nadir viewing instruments • High accuracy • High vertical resolution • All weather-minimally affected by aerosols, clouds or precipitation • Independent height and pressure • Requires no first guess sounding • Independent of radiosonde calibration • No instrument drift • No satellite-to-satellite bias

6. Problems of using GPS RO data in weather models • GPS RO data (e.g., phase, amplitude, bending angles, refractivity) are non-traditional meteorological measurements (e.g., wind, temperature, moisture, pressure). • The long ray-path limb-sounding measurement characteristics are very different from the traditional meteorological measurements (e.g., radiosonde) or the nadir-viewing passive microwave/IR measurements. GPS RO observation is not a point observation like a radiosonde. • The GPS RO measurements are subject to various sources of errors (e.g., uncalibrated ionospheric effects, tracking errors, super-refraction, optimization of bending angle profiles, …etc). [see Kursinski et al. (1997), JGR, 102 (D19), 23429-23465]

7. Assimilation of GPS RO data:The purpose of data assimilation is to extract the maximum information content of the GPS RO data, and to use this information to improve analysis of model state variables (u, v, T, q, P, …etc).

8. GPS RO measurement and data processing procedures Before we consider the assimilation of GPS RO data, we need to understand what are actually measured and the various data processing steps taken to reduce the data.

9. GPS radio occultation measurements & processing Raw measurements of phase and amplitude of L1 and L2 s1, s2, a1, a2 Multi path, Wave Optics Satellites orbits & Spherical Symmetry Assumption Bending angles of L1 and L2 a1, a2 Bending angle Single path Geom. Optics a s1, s2 Refractivity Ionospheric effect cancellation N High altitude Climatology & Abel inversion Raw measurements of phase of L1 and L2 T, e, P Auxiliary meteorological data See Kuo et al. (2000, TAO, 11 (1), 157-186) for details

10. GPS/MET Variables • Raw measurements of L1 & L2 phase (s1, s2) and amplitude (a1, a2) • Raw measurements of L1 & L2 phase (s1, s2) • Bending angles of L1 & L2 (a1, a2) • Bending angle a(corrected for ionosphere) • Refractivity N (through Abel inversion): • Retrieved T, e, and P

11. Which variables should we use for the assimilation of GPS RO data?

12. Choice of Assimilation Variable Should consider the following factors: • Use the raw form of the data, to the extent possible (e.g., the more processing the less accurate the data due to additional assumptions or auxiliary data used in processing). • Ease to model the observables (and its adjoint) • Minimize the need for auxiliary information (before the assimilation of GPS RO data) • Ease to characterize observational errors • Computational cost

13. Most “raw” form of the data. No assumptions needed. Easy to characterize measurement errors. Observation operator need to model wave propagation inside weather models. Require precise GPS and LEO orbits information. Require ionospheric model to account for ionospheric delays (we don’t have very accurate ionospheric model) Computationally very expensive. Assimilation of L1 and L2 phase and amplitude Pros Cons Not Practical

14. Most “raw” form of the data. No spherical symmetry assumption. Easy to characterize measurement errors. Assume single ray propagation Observation operator needs accurate ray tracing (shooting method required) between GPS and LEO Require precise GPS and LEO orbits information. Require ionospheric model to account for ionospheric delays (we don’t have very accurate ionospheric model) Computationally very expensive. Assimilation of L1 and L2 phase Pros Cons Not Practical

15. Second “raw” form of the data. Does not require precise GPS and LEO orbits information. Shooting method not required. Relative easy to characterize measurement errors. Observation operator needs to perform ray tracing with initial conditions. Require ionospheric model to account for ionospheric delays (we don’t have very accurate ionospheric model) Computationally expensive. Assimilation of L1 and L2 bending angles Pros Cons Major difficulty: Hard to remove ionospheric effects

16. Third “raw” form of the data. Does not require precise GPS and LEO orbits information. Does not require ionospheric model. Shooting method not required. Reasonably easy to characterize measurement errors (still challenging for lower troposphere). Observation operator need to perform ray tracing with initial conditions. Uncalibrated ionospheric effects are a source of error (e.g., residual errors associated with ionospheric correction). Still computationally expensive (hard to implement operationally current generation of computers) Assimilation of neutral atmosphere bending angles Pros Cons A possible choice

17. Observation operator easy to develop (local operator on model variables). Does not require precise GPS and LEO orbits information. Computationally inexpensive (operationally feasible). Easy to characterize measurement errors. Assuming Abel-inverted refractivity as the model local refractivity 4th raw form of the data. Requires initialization by climatology (for upper boundary conditions). [Need to create an “optimized” bending angle profile based on observation and climatology, before the retrieval of refractivity.] Uncalibrated ionospheric effects are a source of error. Bias due to super refraction Assimilation of neutral atmosphere refractivity Pros Cons A possible, most popular choice

18. Requires little or no work in the development of observation operator (as they are model state variables). The retrieved T, q, and P can be assimilated by simple analysis or assimilation methods. Computationally inexpensive. Many data processing steps must be taken before T, q, and P are retrieved. Auxiliary information is needed for retrieval, and it can introduce additional errors. Hard to characterize observational errors (as it is mixed with the errors of the auxiliary information). Bias errors due to superrefraction. Least accurate. Assimilation of retrieved T, q and P Pros Cons Not a Good Choice

19. Recent Development -- linearized non-local observation operators

20. Linearized non-local operators • A new class of linearized non-local (LNL) observation operators have been developed recently that have the following features: • It makes use of simplified ray trajectories (can be straight line or curve line) that do not depend on refractivity. • This linearizes the assimilation problem (wrt refractivity) • Bending angle assimilation requires recalculation of ray trajectory at every iteration (since model refractivity is altered after assimilation) • For the LNL operators, this is not necessary, since ray trajectory does not change for each occultation soundings during the iteration steps. • Abel-inverted refractivity is no longer used as local refractivity. Rather, a new modeled observable is defined as a function of refractivity. • The LNL operators are only slightly more expensive than local refractivity operator, but significantly (about 2 order of magnitudes) cheaper than bending angle assimilation. • The LNL operators account for horizontal refractivity gradients and are much more accurate than local refractivity operator (only slightly less accurate than bending angle obs operator).

21. Three possible LNL observation operators The Abel-retrieved (AR) N is a non-local, non-linear function of the 2-D refractivity in occultation plane. Modeling of RO AR N as the local N may result in significant errors, especially, in the troposphere. Accurate modeling of RO bending angle by ray-tracing is computationally expensive. An alternative: to use simple linearized, non-local observation operators: (i) bending angle (Poli 2004); (ii) refractivity (Syndergaard et al. 2005); (iii) phase (Sokolovskiy et al. 2005) LNL observation operators are NOT meant to represent the true GPS observables (s, a, N). Rather they are “modeled observables”, which are functions of refractivity.

22. References for LNL operators: Ahmad, B., and G. L. Tyler, 1998: The two-dimensional resolution kernel asociated with refrieval of ionospheric and atmospheric refractivity profiles by Abelian inversion of radio occultation phase data. Radio Science, 33, 129-142. Syndergaard, S., E. R. Kursinski, B. M. Herman, E. M. Lane, and D. E. Flittner, 2005: A refractive index mapping operator for variational assimilation of occultation data. Mon. Wea. Rev., 133, 2650-2668. Sokolovskiy, S., Y.-H. Kuo, and W. Wang, 2005: Assessing the accuracy of linearized observation operator for assimilation of Abel-retrieved refractivity: Case simulation with a high-resolution weather model. Mon. Wea. Rev., 133, 2200-2012. Sokolovskiy, S., Y.-H. Kuo, and W. Wang, 2004:Validation of the non-local linear observation operator with CHAMP radio occultation data and high-resolution regional analysis. Mon. Wea. Rev., 133, 3053-3059. Poli, P., 2004: Assimilation of global positioning system radio occultation measurements into numerical weather forecast systems. Ph. D. Thesis, U. of Maryland, 127pp. Poli, P., 2004: Effects of horizontal gradients on GPS radio occultation observation operators. II: A Fast atmospheric refractivity gradient operator (FARGO). Q.J. R. Met. Soc. 130, 2807-2825.

23. Choice of observation operators L1, L2 phase and amplitude L1, L2 phase L1, L2 bending angle Neutral atmosphere bending angle Linearized nonlocal observation operator Local refractivity Retrieved T, q, and P Not practical Complexity Possible choices Not accurate enough

24. Comparison between Ngps vs Nlocal • Ngps: refractivity calculated from ray-tracing and Abel transform based on NCEP global analysis. • Nlocal: refractivity calculated T, e, P of NCEP grid point data. • For most soundings, Ngps is very close to Nlocal, suggesting the validity of spherical symmetry assumption. • For some soundings, where gradients of N are large, Ngps can be significantly different from Nlocal. 62 soundings

25. Case 1: Hurricane Isabel (2003) • Developed in the lower Atlantic ocean, tracked northwest and landed at North Carolina coast on Sept 18, 2003 • The hurricane was category 4 or 5 for a period of 6 days. • The WRF simulation covered a period when the hurricane was category 2. • 24-h forecast from 4-km WRF simulation, valid at 0000 UTC 17 September 2003. A Equivalent potential temperature B Radar reflectivity B A

26. Error of observation operator Original, 4 km WRF horizontal resolution Errors in the troposphere: local refractivity >10%; non-local refractivity <2%; phase <1%

27. Choice of Data Assimilation Systems

28. Data Assimilation Systems History of the main data assimilation algorithms used in meteorology and oceanography, roughly classified according to their complexity (and cost) of implementation, and their applicability to real-time problems. Currently, the most commonly used for operational applications are OI, 3D-Var and 4D-Var. From F. Bouttier and P. Courtier Based on ECMWF Training Materials

29. Choices of Assimilation Systems Factors that need to be considered: • Ability to assimilate non-traditional variables (e.g., bending angles, refractivity, or other modeled observables): • Simpler methods (Cressman objective analysis, nudging, OI) cannot assimilate in-direct variables. • 3DVAR, 4DVAR, EnKF can assimilate any variables that can be expressed as functions of the basic model variables. • Ease for the implementation of observation operators: • 3DVAR and 4DVAR require the development of adjoint of observation operator. • EnKF only needs the forward observation operator. • Computational cost: • 3DVAR much cheaper than 4DVAR & EnKF • 4DVAR & EKF compatible in cost • Ability to assimilate data at the time and location when they are taken (4DVAR & EnKF). • Ability to use model (or dynamics) constraints (4DVAR & EnKF). • Ability to consider flow-dependent background errors (4DVAR & EnKF).

30. Variational assimilation of GPS RO data • Assimilation of N or a requires the use of variational data assimilation (or EKF) systems, as N and a are not model predictive variables. • In Variational Analysis (e.g. 3D- or 4D-VAR, we minimize the cost function: • where xo = x(to) is the analysis vector, xb is the background vector, d is the observation vector, O is the observation error covariance matrix and B is the background error covariance matrix. • H is the forward model (observation operator) which transforms the model variables (e.g. T, u, v, q and P) to the observed variable (e.g. bending angle, refractivity, or other modeled observables).

31. COSMIC (Constellation Observing System for Meteorology, Ionosphere and Climate) • 6 Satellites was launched: • 01:40 UTC 15 April 2006 • Three instruments: • GPS receiver, TIP, Tri-band beacon • Weather + Space Weather data • Global observations of: • Pressure, Temperature, Humidity • Refractivity • Ionospheric Electron Density • Ionospheric Scintillation • Demonstrate quasi-operational GPS limb sounding with global • coverage in near-real time • Climate Monitoring A Joint Taiwan-U.S. Mission FORMOSAT-3 in Taiwan

32. 1.7 Million Profiles in Real Time 4/21/06 – 5/6/2009

33. ECMWF SH T Forecast Improvements from COSMICAssimilation of bending angles above 4 km Sean Healy, ECMWF

34. ECMWF Operational implementation of GPSRO on Dec 12, 2006 ↑ Neutral in the troposphere, but some improvement in the stratospheric temperature scores. Obvious improvement in time series for operational ECMWF model. Dec 12, 2006 Operational implementation represented a quite conservative use of data. No measurements assimilated below 4 km, no rising occultations. Nov 6, 2007 Operational assimilation of rising and setting occultations down to surface Sean Healy, ECMWF

35. 100 hPa Temperature vs. radiosondes Sean Healy, ECMWF NH tropics SH

36. NCEP Impact study with COSMIC • 500 hPa geopotential heights anomaly correlation (the higher the better) as a function of forecast day for two different experiments: • PRYnc (assimilation of operational obs ), • PRYc (PRYnc + COSMIC) • Assimilated ~1,000 COSMIC profiles per day • Assimilated operationally at NCEP 1 May 2007 • Assimilating refractivities from rising and setting occultations at all levels (including low level), provided they pass QC • Results with COSMIC “very encouraging” • Lidia Cucurull, JCSDA

37. UKMO Bias and RMS as function of forecast range mean mean Temp, 250 hPa, SH Wind speed, 100 hPa, SH RMS RMS

38. Prediction of Typhoon Shanshan (2006)

39. Typhoon Shanshan (Sept 10-17, 2006) Operational forecasts using variational assimilation failed to predict the curving of the typhoon. Central SLP pressure

40. COSMIC RO soundings (September 13, 2006) RO soundings are randomly distributed over the domain, provide large-scale information.

41. Assimilation Experiments  WRF/DART ensemble assimilation at 45km resolution for 8-14 September 2006. 32 ensemble members. Control/NoGPS run: Assimilate operational datasets including radiosonde, cloud winds, land and ocean surface observations, SATEM thickness, and QuikScat surface winds. •GPS run: Assimilate the above observations + RO refractivity.

42. Impact of RO Refractivity on Ensemble Forecasts • 16 members • with a finer nested grid of 15km • initialized at 00UTC 13 Sept 2006

43. Probability forecast of accumulated rainfall (24hours, 12Z 14-15 Sep., > 60mm/day, %) OBS NOGPS GPS Rainfall forecast is enhanced in Northern Taiwan with COSMIC data. This is closer to the observed rainfall.

44. Ensemble Forecasts of Central Sea Level Pressure NoGPS GPS Ensemble Mean Observed Ensemble mean Observed Intensity is increased with RO data. Ensembles give significance.

45. Comparison of WRF 3DVAR and WRF/DART forecast of Shanshan (2006) • Assimilation for 24 hours starting 00Z 13 September 2006 using both 3DVAR and WRF/DART ensemble system • Assimilation of CWB conventional data with/without RO data DARTNBNG: NO GPS run using DART DARTNB: With GPS run suing DART CYCLNBNG: NO GPS run using 3dvar CYCLNB: With GPS run using 3dvar • Followed by a 3-day forecast on 14 September 00Z.

46. Zonal wind Analysis along the typhoon centers (125.8E) on 00Z 14 September WRF/DART 3DVAR Zonal wind is stronger in WRF/DART analysis

47. Vorticity Analysis along typhoon centers (125.8E) on 00Z 14 September WRF/DART 3DVAR Vortex is stronger in WRF/DART analysis

48. Analysis at 00Z 14 Sept 2006 WRF-3D-Var WRF/DART GPS - NO GPS EnKF - 3D-Var With GPS Without GPS

49. Typhoon intensity (central pressure) 3DVAR WRF/DART OBS

50. Typhoon track error 3DVAR WRF/DART