Application of the Full Spectrum Inversion Algorithm for 1 Airborne GPS Radio Occultation Measurements

10 With a GPS receiver onboard an airplane, the airborne RO (ARO) technique provides dense 11 lower troposphere soundings over target regions. The large variation of water vapor in the 12 troposphere causes strong signal multipath, which could lead to systematic errors in RO 13 retrievals with the geometric optics (GO) method. The spaceborne GPS RO community has 14 successfully applied the Full Spectrum Inversion (FSI) technique to solve the multipath 15 problem. This paper is the first to adapt the FSI technique to the ARO measurement with its 16 unique perspective of having a receiver traveling on a non-circular trajectory inside the 17 atmosphere. 18 An end-to-end simulation system is implemented to test the newly developed FSI retrieval 19 algorithm for ARO. The forward-simulated GPS L1 signal amplitude and phase is used to test 20 the modified FSI algorithm. The ARO FSI method is capable of reconstructing the fine vertical 21 structure of the moist lower troposphere in the presence of severe multipath, which leads to 22 large retrieval errors in the GO retrieval. The sensitivity of the modified FSI retrieved bending 23 angle and refractivity to the errors in signal amplitude and the measured refractivity at the 24 receiver is presented. Accurate bending angle retrievals can be obtained from surface up to 25 ~250 m below the receiver, where retrieved bending angle near the receiver altitude becomes 26 sensitive to the measurement noise. Abrupt changes in the signal amplitude do not produce a 27 systematic bias in the FSI retrievals when the SNR is high. A 1 % Gaussian noise in refractivity 28 Atmos. Meas. Tech. Discuss., doi:10.5194/amt-2015-378, 2016 Manuscript under review for journal Atmos. Meas. Tech. Published: 18 January 2016 c © Author(s) 2016. CC-BY 3.0 License.

at the receiver causes ~ 0.5 % refractivity error near the receiver that reduces to ~0.05 % near the surface.

Introduction
Global Positioning System (GPS) satellites transmit radio signals that undergo refractive bending and Doppler shift due to the variations in the refractive index of the Earth's atmosphere which become pronounced in the limb direction.With a GPS receiver onboard an aircraft, the airborne RO (ARO) receiver tracks the occulting GPS signals traverse progressively lower (or higher) atmospheric layers when the GPS satellite sets behind (or rises above) the local horizon of the receiver (Healy et al., 2002;Xie et al., 2008).Different from the spaceborne RO measurement, the ARO receiver is located inside the atmosphere with considerable atmospheric refraction near the receiver.Moreover, in addition to the recording of occulting signals below the local horizon (similar practice in spaceborne RO), the RO signals from above the local horizon also need to be recorded to allow the retrieval of atmospheric property below the ARO receiver.Similar to the spaceborne RO, the measurement of the raw ARO signal phase and amplitude can be inverted to retrieve the bending angle (the cumulative atmospheric refraction bending along each ray path) as a function of impact parameter.The impact parameter is a conservative quantity for each signal ray, i.e., the product of the radius and the refractive index at the tangent point (Kursinski et al., 2000).The bending angle can then be converted to refractivity through inverse Abel transformation (Fjeldbo, 1971).The refractivity (N) or the refractive index (n) of the neutral atmosphere depends on the atmospheric temperature (T in Kelvin), total pressure (P in hPa) and water vapor pressure (e in hPa) (Kursinski et al., 1997;2000), such as: The fundamental observables during an ARO event are the time series of phase and amplitude or the signal-to-noise ratio (SNR) of the RO signals.After the precise positions of the GPS and the receiver are known (e.g., Muradyan et al., 2010), the excess phase delay due to atmospheric refraction can be derived by differencing the measured signal total phase with the GPS-receiver line-of-sight (LOS) distance.In this study we simulate GPS L1 signal (1575.42MHz) for airborne receivers at around 14 km, and neglect ionospheric effects, which can be removed through linear combination with dual frequency measurements (e.g., Vorobev andAtmos. Meas. Tech. Discuss., doi:10.5194/amt-2015-378, 2016 Manuscript under review for journal Atmos.Meas.Tech.Published: January 2016 c Author(s) 2016.CC-BY 3.0 License.
Krasil 'nikova, 1994;Hajj et al., 2002).However, the ionospheric effects could be negligible for ARO retrievals, as the ARO retrieval requires the differencing between the RO signals originating from below (negative elevation) and above (positive elevation) the local horizon, which will cancel out the ionospheric effect (Xie et al., 2008) (hereafter referred as X08).
Moreover, the emphasis of the ARO measurement are below the aircraft flying at the lower stratosphere or upper troposphere, where the the atmospheric bending is dominant by the neutral atmosphere.This will eliminate the dual GPS frequencies recording requirement for ARO measurements.The bending angle difference, i.e., the partial bending angle, between the negative elevation and the positive elevation components at the same impact parameter, is then used to derive the refractivity through he inverse Abel transformation.The derivative of this excess phase represents the Doppler shift of the carrier signal.The commonly used geometric optics (GO) method uses the measured signal Doppler and the GPS-Receiver position/velocities to retrieve the bending angle.One major limitation of the GO method is its inability to account for signal interference, known as multipath that frequently occurs in the moist lower troposphere due to large water vapor variations.When multipath occurs, the signal at the receiver consists of the superposition of multiple rays each having its unique impact parameter, and the Doppler shift derived from the signal phase no longer corresponds to a unique ray path or one impact parameter.As a result, the GO method can lead to large retrieval errors.
Various radio holographic methods have been proposed to overcome the limitations of the GO method in the spaceborne RO retrievals (Gorbunov et al., 1996;Gorbunov and Gurvich, 1998;Sokolovskiy 2001;Gorbunov 2002).Full spectrum Inversion (FSI) proposed by Jensen et al. (2003) (hereafter referred as J03) has been applied to invert the spaceborne RO signals, which outperforms the GO methods in the presence of multipath.However, its application for the airborne RO has not been performed yet and need to address the unique characteristic of ARO occultation measurements.This ARO application need to consider the asymmetry of an airborne receiver trajectory inside the atmospheric media and the requirement to measure RO signals originating from both above and below the local horizon, and it is further complicated by the likely irregular flight path of the airborne platform.
In this paper, the development and implementation of the FSI retrieval method for ARO measurements are presented.An end-to-end simulation system is developed to carry out sensitivity analysis of the new FSI algorithm and its comparison with the GO methods.Real occultation geometry of ARO events is used in this simulation study adopted from the PRE-Atmos.Meas. Tech. Discuss., doi:10.5194/amt-2015-378, 2016 Manuscript under review for journal Atmos.Meas.Tech.Published: January 2016 c Author(s) 2016.CC-BY 3.0 License.
Depression Investigation of Cloud-systems in the Tropics (PREDICT) field campaign over tropical Atlantic during August-September of 2010 (Montgomery et al., 2012;Haase et al., 2014;Murphy et al., 2015).During the PREDICT campaign, the Global Navigation Satellite System (GNSS) Instrument System for Multistatic and Occultation Sensing (GISMOS) was deployed on the HIAPER GV aircraft (Garrison et al., 2007).Haase et al. (2014) and Murphy et al. (2015) presented very promising initial results of the airborne RO observations based on the GO retrieval.However, significant refractivity retrieval errors in the middle and lower troposphere were detected in the RO retrievals.These retrieval errors in the lower troposphere are likely caused by the multipath problem that plagues the GO method.The FSI retrieval method presented in this paper is expected to solve the multipath problem and offer high-quality lower tropospheric ARO soundings with high vertical resolution near the tropical storms to improve our understanding the complicated hurricane genesis process.High resolution ERA-Interim reanalysis profiles from European Centre for Medium Range Weather Forecasts (ECMWF) over the campaign area are used to represent the atmospheric condition.
This paper is organized as follows: Section 2 describes the key implementation steps for the FSI method for ARO.An end-to-end simulation system is presented in Section 3. Section 4 presents the FSI application for simulated ARO observation under severe multipath condition The real occultation geometry from the research flights during the PREDICT campaign is also used.Sensitivity of the FSI retrievals to the ARO measurements errors in signal amplitude and the refractivity at the receiver are explored in Section 5.The conclusions are summarized in Section 6.

Theoretical derivation of FSI for airborne RO measurements
The FSI method recognizes the RO signal recording as radio waves of different frequencies determined by the refractive index of the media through which they pass, and accounts for interference of waves with different frequencies.Each wave with its unique frequency corresponds to one single ray path in the GO application.The FSI retrieval of J03 is based on the assumption that the Fourier transform of the RO signal, which is computed using the method of stationary phase, can identify unambiguously the multiple frequencies present in the signal at a given time when certain conditions are met.(2) where, k, a, n rec R rec and R GPS are the wave number, impact parameter, refractive index at the receiver, radius of receiver and GPS from the center of the Earth, respectively.The difference between eq. 13 in J03 and eq. 2 above arise from the receiver position within the atmosphere, where the refractive index at the receiver is greater than 1 and is a function of receiver position.
The last two terms in eq. 2 arise due to the non-spherical trajectories of the GPS and the receiver.
When the trajectories are spherical, the equation can be simplified as The impact parameter (a) of a signal ray is defined as (Kursinski et al., 1997) Using eq. 4 and taking n = 1 at the GPS position, the angle φ GPS (see Fig. 1) can be calculated as In the case of ARO, φ rec = π/2 refers to the local horizon or zero elevation, whereas φ rec > π/2 refers to the positive elevation and φ rec < π/2 refers to the negative elevation (see X08 for detailed description of positive and negative elevation angles).Therefore, φ rec for positive and negative elevation angles are given by eq.5b and 5c, respectively as: The bending angle (α) can then be calculated as

Correction for non-spherical trajectory
The FSI retrieval in eq.3-6 is valid only when the GPS and the receiver are both moving in circular trajectories within one occultation plane.When the GPS and the receiver trajectories deviate from circular, radial velocities and acceleration terms are introduced, represented by the two radial terms in eq. 2.
In real airborne occultation measurements, the perfectly circular trajectory assumption is not valid in part because of the oblateness of the Earth, and the local variation of the aircraft altitude.
To take into account the oblateness of the Earth, Syndergaard (1998) showed that the inversion of the RO data should be performed assuming local spherical symmetry tangential to the Earth's ellipsoid.
In our current approach, we account for the oblateness of the Earth by calculating the local center of curvature for each occultation event.After oblateness correction, a correction has to be applied to account for the non-spherical trajectories of the receiver and the transmitter.In the current algorithm, the correction for non-spherical trajectories has been performed by projecting the position of both the receiver and the transmitter at each epoch to a circular trajectory.Figure 1 shows the schematic diagram of the projection of GPS signal from a noncircular receiver trajectory onto a fixed radius circular trajectory in the occultation.Similarly, the method is applied to the GPS orbit to allow its projection onto a circular orbit.The figure shows the receiver at position P with radius R rec relative to the local center of curvature O.The projection is done along the direction vector of the ray at P, which is determined from the occultation geometry and bending angle at the receiver height obtained from the CIRA+Q refractivity climatological model (Kirchengast et al., 1999).The first estimate of the projected position is determined using the triangle formed by the joining the origin, O, with P and the point where the direction vector at P intersects the reference circular trajectory with radius R rec,0 .
However, the refractivity difference between P and the projected position causes the ray to ) .( 7) where, d is the geometric phase, i.e., GPS-receiver line-of-sight distance, and the excess phase (exphs) is the difference between the total phase (phs) and geometric phase (d).
Since the GPS is located outside the Earth's atmosphere, the refractive index at the GPS altitude is 1.Therefore, the corresponding change in phase at the GPS (dph GPS ) is given by After the projection at each epoch, the new trajectories of both GPS and receiver become circular relative to the local center of curvature.These additional phase terms introduced by the projection are included into the new total phase (for the circular trajectories), i.e., the sum of the original total phase (d+exphs) plus the phase addition terms resulting from the projection for both GPS and receiver, such as When these projections are applied, both R rec and R GPS become constants, and the two radial terms in eq. 2 become zero.After the adjustment, FSI is applied to the modified signal phase and the original signal amplitude with both GPS and receiver on circular trajectories.

Estimation of the bending angle at local horizon
In the case of airborne RO measurements, the impact parameter The impact parameter is a constant for each ray path, and is given by The angles φ GPS and φ rec are given by where ∅ D,E and ∅ I1J are angles between the line-of-sight and the radial vectors of the GPS (R GPS ) and the receiver (R rec ), respectively.Solving equations 13-15 yields  D,E and  I1J , which can be used to find the epoch of the local horizon (∅ I1J = ?)and separate the signal into positive and negative elevation angle components.

End-to-end Simulation System for Airborne RO Soundings
An end-to-end simulation system (Fig. 2) was developed to investigate the performance of the modified FSI algorithm for airborne RO retrieval.The simulation system consists of two major components; (i) a forward simulator, and (ii) an inverse simulator, i.e., FSI retrieval.The The inverse Fourier transform of F(a) can then be expressed as   = ()  ecd.=.6 . ( Equation 16 can be approximated as   =    cd.=.6  =   . cg(=) .( 18) where,   = −.() .
When the bending angle, α(a), for a given atmospheric profile is known, θ(a) can be calculated as follows: The integral in eq.19 can be computed after the open angle θ(a) is derived from a given α(a) using eq.20.The complex signal F(a) can then be calculated by assuming B(a) is a constant.
The total phase and amplitude of the simulated GPS signal can be obtained from the complex phase function u(θ).The excess phase of the signal can be derived by subtracting the GPSreceiver LOS distance from the simulated total phase.The excess Doppler can be further derived by taking the derivative of the excess phase.In this FSF forward model, the input atmospheric condition is represented by a bending angle profile.In this paper, the atmospheric temperature, pressure and water vapor mixing ratio was first used to calculate the atmospheric refractivity profile followed by the forward Abel transformation (e.g., Xie et al., 2008) to derive the bending angle profile as a function of impact parameter.

Inverse ARO simulator (Full-Spectrum Inversion, FSI)
The inverse simulators in the end-to-end simulation comprise both the GO retrieval (e.g., Xie et al., 2008) and the newly developed FSI.Both the GO and FSI retrievals derive the bending angle profiles as function of impact parameter from the input excess Doppler (for GO) or the combination of both the excess phase and amplitude (for FSI).The inverse Abel transform is then applied to retrieve refractivity from the bending angle.
In the following section, the input atmospheric refractivity and/or bending angle profiles to the forward simulator are directly compared to the output from the inverse simulators to assess the performance of the inversion technique, and to quantify the sensitivities of the ARO FSI to the potential errors presented in several key input variables such as the SNR and refractivity at the receiver.

Application of the FSI retrieval for ARO measurements
To assess the performance of ARO FSI retrieval algorithm, we used real occultation geometry and the atmospheric profile of temperature and water vapor from ERA-Interim reanalysis.One specific occultation geometry involves the GPS satellite (Pseudo Random Number, PRN 24) and the airborne receiver were obtained from the ARO measurements made during the PREDICT flight from 1820 -1900 Z on September 14, 2010 (Research Flight No. 19).The grid profiles of ERA-I temperature and water vapor mixing ratio and the calculated refractivity profiles from the ARO sounding region are shown in Fig. 3(a) and 3(e), respectively.Very moist atmosphere with high mixing ratio ~20 g/kg is seen near the surface that decreases rapidly at higher altitude.Above 10 km, the temperature reduces to around 250 K (-23.15 °C), which leads to a very dry atmosphere and the contribution of water vapor to atmospheric refractivity becomes negligible in comparison with that of temperature.Windowing and tapering at the edges using sinusoidal function are then applied to each of the positive and negative elevation components to avoid adding spurious components during the FFT.This tapering creates artificially low amplitude at 600s where the separation point of the two tapered segments locates.
Similarly, the occultation phase and amplitude time series are also divided into positive and negative elevation parts for the inversion retrieval, followed by a similar sinusoidal tapering at the edges, so this amplitude variation near the zero elevation angle does not affect the retrievals.
The time epoch of the local horizon (zero elevation) is estimated by the GO ray-tracing simulation with the CIRA+Q bending angle model and the given occultation geometry.Given the bending angle, the refractivity below the aircraft is obtained through the inverse Abel transform, by integrating the partial bending angle (e.g., difference in bending angle between the negative and positive elevation at each impact height) from the tangent point height up to the receiver height (Fig. 3e).
Near the receiver, the partial bending angle from GO retrieval has alternating positive and negative errors near the receiver, which leads to small overall refractivity errors below 14 km.
But still relatively large errors are observed above 14 km.For FSI retrievals, the partial bending error is positive due to a positive bending angle spike at negative elevation near the zero elevation, which leads to a positive error in refractivity retrieval and propagates downward to lower levels.As the bending angle increases exponentially downwards, the refractivity errors caused by the bending angle error near the receiver height also decrease exponentially downward (Fig. 3f solid lines), which is consistent with the GO simulation study in X08.
Note that during an airborne RO observation, the temperature, humidity and pressure at the aircraft can be precisely recorded with the in-situ sensors.The refractivity at the receiver can then be obtained from these in-situ measurements.At high altitudes with limited availability of water vapor (e.g., above 10 km), the refractivity generally decays exponentially at higher altitudes with a rather constant scale height to be around 7 km.With a further constraint of the in-situ refractivity observation from the aircraft, the simple exponential refractivity model can be used to estimate the bending angle near the receiver at both positive and negative elevations through the forward Abel calculation (e.g., Xie et al., 2008).Figure 3d shows only the top ~250 m of the bending angle retrieval are noisy, due to the high sensitivity of the both GO and FSI retrieval to the measurement noise in Doppler or phase near the zero elevation.By replacing the noisy retrieved bending angle (e.g.top 250 m below the receiver height) with the simple bending angle model constrained by the in-situ refractivity at the receiver, the refractivity errors near the receiver are almost completely removed for both GO and FSI retrievals.This correction of bending angle near the receiver is necessary for ARO retrievals because the errors at the top propagate downwards during calculation of refractivity using the Abel transformation, which adds artificial bias to the retrieved refractivity.
In the lower troposphere, on the other hand, the large refractivity errors in the GO retrieval in the lowest 1 km are due to the bending angle retrieval error in the presence of multipath.The FSI, on the other hand, resolves the fine vertical structure of both bending angle and refractivity in the presence of the multipath in the moist atmosphere near the surface without introducing retrieval biases.

Sensitivity to signal amplitude and refractivity at the receiver
The accuracy of FSI retrieval depends on the accuracy of the measured phase and amplitude of the signal, the occultation geometry, and the refractivity observation at the receiver.Note the sensitivity of the ARO retrieval to the excess phase or Doppler has been explored in the GO retrieval system in Xie et al. (2008), and we do not expect much difference in terms of sensitivity for the FSI retrievals.In this section, we will quantify the sensitivity of FSI retrievals to the errors in signal amplitude (not used in GO retrieval), and the refractivity at the receiver.
The signal amplitude, cannot be measured with the same level of accuracy as the phase of the signal (Kursinski et al., 2000).In the ARO measurements, the amplitude is affected by the aircraft heading and the attitude relative to the line of sight because of the focused antenna gain pattern.Sharp amplitude jumps could be introduced by the changing aircraft direction during the ARO measurement (Wang et al., 2015).Therefore, it is important to know how sensitive the FSI retrieval is to the less accurate signal amplitude.Besides, under low signal-to-noise ratios (SNR) condition, when the signal amplitude is comparable or lower than the noise, phase measurements could have greater uncertainties.To test the sensitivity of the FSI retrieval algorithm to the variations in signal amplitude, we need to account for the possible phase errors that may arise under low SNR condition.To accomplish this task, signal phase and amplitude were first simulated by using the ray tracing method.Then a sinusoidal amplitude function (e.g., in eq.21) was added to simulate the sharp amplitude jumps produced by the changing aircraft direction (e.g., Murphy et al., 2015), and finally Gaussian noise was added to the amplitude to represent the variations in the ARO amplitude measurements.
=  W  +  k   +  ? (). ( where, SNR 0 is the amplitude simulated using the ray tracing method, and k 1 , k 2 , a and b are constants that determine the shape of the resulting amplitude.Figure 4(a) shows the amplitude simulated by using the ray tracing method (blue) and the modified noise-added amplitude (red).
In the simulation, the Gaussian noise power is assumed to be 1% of the peak signal power.The .This leads to an estimation of the noise power to be ~0.56% of the peak signal power.The variance of the noise, i.e., the square root of the ratio of the noise power to the signal power, is then estimated to be 0.75%.For simplicity, in the following simulation study, the noise power variance of 1% , such that the noise power is 1% of the peak signal power, will be used to represent the upper bound magnitude of the noise.
Wang et al, (2015) have shown that at low SNR, increased phase variance results in large errors in the unwrapped phase of the signal.Therefore, to test the impact of signal amplitude errors on the FSI retrievals, it is important to assess its impact on the measured signal phase.To where, the noise in the two components: I n and Q n are assumed to be independent, normal distributed with zero-mean and 1% in variance.The two modified signal components in eq.22-23 were then used to reconstruct the noise-loaded residual phase ( > ) and amplitudes ( > ) as  It is worth noting that in the very low SNR condition, the amplitude error could potentially lead to cycle slip and/or unwrapping errors, which could lead to systematic bias in signal phase or Doppler observation (e.g., Wang et al., 2015).Such biased phase or Doppler will lead to biases in both bending and refractivity retrieval.In a simulation using noise to be 10% of the peak signal power (not shown), the reconstructed residual phase starts deviating from original residual phase at 3.5 -4 km height range due to the large unwrapping errors.The retrieved bending angle errors start to increase at this height, and the refractivity error exceed ± 2 %, causing large uncertainties of the retrieved quantities below 4 km.However, such errors are caused by the degraded observation in the signal tracking state but not introduced by the retrieval process (e.g., the FSI retrieval), and therefore is out of the scope of this paper.
The refractivity at the receiver can be obtained from the in-situ temperature, pressure and water vapor mixing ratio measurement at the aircraft.It is also one of the key parameter used in FSI retrieval.With an ARO receiver flying at ~14 km, a typical flight level during the PREDICT campaign, the water vapor contribution is negligible, so the refractivity at the receiver can be assumed to be only a function of temperature and pressure.To quantify the sensitivity of the FSI retrieval to the refractivity measurement error at the receiver, a Gaussian noise of 1% in the refractivity (~2K error in temperature) at the receiver is added.Fifty realization of random Gaussian noise were added to the base refractivity at the receiver and the statistics of the FSI retrieval errors are then compiled.maximum of ~0.5% refractivity error near the receiver height and decreasing to ~0.05% near surface.However, no systematic bias is introduced by such random in-situ measurement error.

Conclusions and Discussions
In this study, a Full-Spectrum-Inversion (FSI) algorithm is developed and successfully applied on Airborne GNSS RO (ARO) measurements for the first time.Simulation study demonstrate the capability of the FSI method to retrieve the fine atmospheric vertical structure in the lower moist troposphere where frequent multipath occurs.
In the FSI retrieval process, oblateness correction is applied on the original occultation geometry to fulfill the local spherical symmetric assumption.Then the non-spherical trajectories of both the ARO receiver and the GPS satellite are projected onto circular trajectories relative to local center of curvature.The additional phase terms as a result of projection are then added to the measured phase.Afterward, the occultation signal time series is divided into positive and negative elevation components.The separation point, i.e., the local horizon can be estimated with the aid of the CIRA+Q climatological bending angle at the receiver height.The FSI algorithm is then applied to the amplitude and the modified signal phase for the new circular GPS/receiver trajectories.An end-to-end simulation system is used to test the FSI retrieval using the realistic airborne occultation geometry obtained from the PREDICT field campaign.The simulation study revealed that both FSI and GO inversion methods introduce bending angle errors near the receiver height, which lead to >1% refractivity error near the receiver height but decreasing quickly downward at lower altitudes.The bending angle error near the receiver height can be easily corrected with the aid of the in-situ refractivity measurement.
The end-to-end simulation system was also used to quantify the sensitivity of the FSI bending and refractivity retrievals to the noise in two key parameters including the signal amplitude and the refractivity at the receiver.The FSI retrieval showed a weaker sensitivity to signal amplitude errors as compared to the refractivity errors at the receiver.Even the abrupt changes in signal amplitude do not introduce systematic bias to the retrieval, as long as the SNR is high.The 1 % in-situ refractivity errors at the receiver height could introduce a maximum refractivity retrieval error of 0.5 % (1 K) near the receiver but decreasing gradually to ~0.05% (0.1 K) near the Atmos.Meas.Tech.Discuss., doi:10.5194/amt-2015-378,2016   Manuscript under review for journal Atmos.Meas.Tech.Published: January 2016 c Author(s) 2016.CC-BY 3.0 License.Similar to the spaceborne RO in J03, the derivative of the phase of the Fourier transform of the ARO GPS signal (ϕ q ) with respect to open angle (θ, see Fig.1) can be expressed as

)
Atmos.Meas.Tech.Discuss., doi:10.5194/amt-2015-378,2016   Manuscript under review for journal Atmos.Meas.Tech.Published: January 2016 c Author(s) 2016.CC-BY 3.0 License.The bending angle profile (a function of the impact parameter) in eq.6 is then used to calculate the refractivity profile as a function of geometric height using inverse Abel transformation (seeX08 and Healy et al., 2002 for details).
bend, which is approximated by multiplying the straight line projection by the mean refractive index ( I1J ) between the original position and the projected position.As shown in the figure, this projection from position P to P 0 leads to a change in the total phase (phs), zenith angle (φ), Atmos.Meas.Tech.Discuss., doi:10.5194/amt-2015-378,2016 Manuscript under review for journal Atmos.Meas.Tech.Published: January 2016 c Author(s) 2016.CC-BY 3.0 License.and open angle (θ).The changes in dph and θ, dph rec and dθ rec , respectively, can be calculated as ℎ I1J =  I1J − I1J  I1J ±  I1J ? ? I1J + ( I1J,W ?−  I1J ? (a) is not a unique function of open angle (θ).The same a occurs twice, once each at negative and positive elevation angle relative to the local horizon at the receiver where ∅ I1J = ?in Fig. 1 (Zuffada et al., 1999; Healy Atmos.Meas.Tech.Discuss., doi:10.5194/amt-2015-378,2016 Manuscript under review for journal Atmos.Meas.Tech.Published: January 2016 c Author(s) 2016.CC-BY 3.0 License.et al., 2002).To avoid the non-unique relation between a and θ, the GPS signal time series is split into two parts to separate the positive and negative elevation angle measurements.The Fourier Transform is then applied to each part separately.The time epoch, when the occulting GPS is at the local horizon or the ARO receiver (i.e., the separation point between positive and negative elevation angles) is estimated with the aid of the CIRA+Q climatological bending angle at the receiver altitude.Although the exact local horizon depends on the bending due to the real atmosphere, the CIRA+Q climatological model provides a reasonable estimation at the aircraft altitude in the upper troposphere or lower stratosphere.The total bending angle (α), which is estimated from CIRA+Q climatological model, can be expressed as the sum of the bending accumulated from the ray tangent point to both the receiver (α rec ) and the GPS (α GPS ), i.e.,  =  I1J +  D,E .
forward simulator is used to simulate the phase and the amplitude for an ARO signal given an atmospheric refractivity model and occultation geometry.The inverse simulator (or retrieval Atmos.Meas.Tech.Discuss., doi:10.5194/amt-2015-378,2016   Manuscript under review for journal Atmos.Meas.Tech.Published: January 2016 c Author(s) 2016.CC-BY 3.0 License.component) uses the simulated ARO signal to retrieve the atmospheric bending angle and refractivity profiles.3.1 Forward ARO simulator (Full-Spectrum Forward, FSF)Two different types of forward simulators were used in the study.The first option is a ray tracer that simulates the ARO signals as geometric optics rays (e.g.,Xie et al., 2008).Forward simulation using the ray tracing technique becomes problematic for atmospheres with sharp refractivity gradients.To simulate the GPS signal in the presence of sharp refractivity gradients, another forward model called the full-spectrum forward (FSF) simulator was developed based on J03.The complex GPS signal, F(a), can be represented as a function of the phase function u(θ), wave number (k), impact parameter (a) and the open angle (θ) as   = () cd.=.6 .
Figure 3(b) shows the excess phase Atmos.Meas.Tech.Discuss., doi:10.5194/amt-2015-378,2016 Manuscript under review for journal Atmos.Meas.Tech.Published: January 2016 c Author(s) 2016.CC-BY 3.0 License.and excess Doppler obtained from the FSF forward simulator.The excess phase increases monotonically as a function of time, whereas, its derivative, the excess Doppler, becomes nonmonotonic function of the time starting at ~1500 s (see inset figure).Such behavior in Doppler is a strong indication of signal interference due to multipath.The multipath is further illustrated by the times series of the signal amplitude in Fig 3(c), which shows large variations around ~1500 s.This signal amplitude variation is caused by superposition of multiple signals with varying frequencies.When the phase and amplitude are calculated from the Fast Fourier Transform (FFT) of the bending angle profile in eq. 18, the function u(θ) is divided into two separate segments for positive and negative elevation angles based on the local horizon.

Figure 3
Figure3(d)-(f) show the bending and refractivity retrievals from GO and FSI and the refractivity differences from the truth (i.e., the input profiles).The bending angle retrievals fromGO and    FSI are plotted in Fig 3(d)  along with the input ("true") bending angle profile.Two distinct and important features of the ARO retrievals are shown.The first is the large error in the retrieved bending angle near zero elevation when the tangent point is near the receiver.This feature is present for both the GO and FSI methods.In the case of FSI retrievals, these errors arise from the large phase correction resulting from projecting the non-spherical trajectory to the spherical trajectory in eq. 7. The second feature is the large errors in the GO bending angle retrieval associated with multipath in the lower troposphere.The inset in Fig.3(d)shows that the GO retrieved bending angle below impact height of 3 km (corresponding to geometric height ~1 km) deviates significantly from the input bending angle, whereas the FSI retrieval follows closely to the true (input) bending angle.The FSI is capable of resolving the sharp bending angle structure in the presence of multipath as a result of a significant changes in moisture and/or temperature gradients near the surface.
ARO measurements from the PREDICT campaign show the peak SNR (SNR is the amplitude scaled by the noise level) of ~200 v/v and the low SNR of ~15 v/v near the noise floor, when Atmos.Meas.Tech.Discuss., doi:10.5194/amt-2015-378,2016   Manuscript under review for journal Atmos.Meas.Tech.Published: January 2016 c Author(s) 2016.CC-BY 3.0 License. the weak GPS occultation signal is dominated by the noise(Wang et al., 2015).Based on the analysis (e.g., eq.50-51) inWang et al. (2015), the SNR and the amplitude ratio 20 and 1.5, respectively.Assuming that the amplitude at the peak power and the noise floor are dominated by the signal and the noise respectively, the ratio of peak power to the noise power can then be estimated as the square of the amplitude ratio, account for this impact on the measured phase, a realistic simulation following the ARO OL signal processing(Wang et al., 2015)  is carried out.Two different model atmospheric profiles were used.One ERA-Interim profile (12Z, Sep 13, 2010 at 15°N, 77°W) is used to represent the true atmospheric state, and a CIRA+Q climatological model profile is used to provide the initial prediction of the excess phase and Doppler of the expected ARO signals through GO ray-tracing.Given realistic ARO geometry, the excess phases are simulated based on the two model profiles through the ray tracing model, respectively.The phase difference between the two, i.e., the residual phase, is then generated (Fig.4bin blue).In the presence of measurement noise, the residual phase (φ) and amplitude A(t) of the received signal can be expressed as the in-phase (I) and quadrature (Q) components, such as >  = () ?+ () ?.(25)Atmos.Meas.Tech.Discuss., doi:10.5194/amt-2015-378,2016   Manuscript under review for journal Atmos.Meas.Tech.Published: January 2016 c Author(s) 2016.CC-BY 3.0 License.The new residual phase ϕ n was then added to the simulated model phase (i.e., derived from CIRA+Q) to represent the noise-added signal phase as shown in Fig.4(b, in red).

Figure 4
Figure 4(c) shows the difference between the FSI retrieved bending angle and the "true" bending angle, calculated by Forward-Abel integration of the ERA-I refractivity profile.Similarly, Fig. 4(d) shows the percentage error of the FSI retrieved refractivity compared to the input refractivity profile.Both bending and refractivity errors show near zero mean with small variations, which indicate that the large variation in the amplitude measurement does not introduce systematic bias in the FSI bending and refractivity retrievals when the SNR is high.

Figure 5
Figure5(a) and 5(b) show the absolute bending angle error and fractional refractivity error of the FSI retrieval.The bending angle exhibits a small 0.02° errors with near zero mean across all altitudes whereas the refractivity errors also show near zero mean and decrease downward from a maximum near the receiver.A 2K in-situ temperature error at the receiver introduces a Atmos.Meas.Tech.Discuss., doi:10.5194/amt-2015-378,2016   Manuscript under review for journal Atmos.Meas.Tech.Published: January 2016 c Author(s) 2016.CC-BY 3.0 License.surface.In the very low SNR condition, large amplitude variations could introduce large errors in the phase measurement and result in large uncertainties in the ARO retrieved bending angle and refractivity profiles in the lower troposphere.However, such errors are caused by the degraded observation in the signal tracking state that is out of the scope of this paper, but will be worth of further investigation when applying FSI retrieval on the real ARO measurements.

Figure 1 .
Figure 1.Projection of the receiver from its original position P with radius R Rec onto a new position P 0 on a circular trajectory with radius R Rec,0 relative to local center of curvature O.

Figure 2 .
Figure 2. End-to-end simulation system for airborne RO data processing.(Note: derivation of bending from refractivity in Forward-Abel doesn't need the occultation geometry information but only the height and refractivity at the receiver)

Figure 4 .
Figure 4. (a) The amplitude of the received signal simulated by the ray tracing model (blue) and the noise-added amplitude (red), (b) residual phase with (red) and without (blue) the amplitude noise, (c) FSI-retrieved bending angle error, and (d) fractional refractivity error.The radius of the earth has been subtracted from the impact parameter in (c), where the earth's surface is at ~2.5 km.