Improve navigation bits documentation in ROPP PP user guide

[Ian Culverwell]

ROPP-10 DRR reviewer Riccardo Notarpietro complains that

The Navbit handling is not clear (to me). There is no theoretical description
associated to the interpolation/removal of the navbit.

p30: Not clear what the "lost carrier flag profile variable" is. It should be a
variable containing the navbit. Is this included in some input file? Is it time
stamped? It seems that navbits are interpolated somewhere before making this
bitwise variable available. Which is the module which takes care of this?

p30: What are these "alternative navigation bits"? Who provides these?

p30: it is written "If available, excess phase data navigation bits may be read
from a separate text file (path specified by the navbit file config option),
which contains a sequence of two-digit binary value navigation bits as a
function of GNSS time. These bits are interpolated to the occultation
measurement times and stored as variable LCF". Not clear how they are used. Two
digits? Why? Then for GRAS, Navbits are available in a predefined format (ascii,
with one subframe / GNSS epoch encoded into a string of exadecimals when they
were provided by GSN, NeCDF with again one subframe / GNSS epoch, encoded into
an array of integers).

p30: Eq. 4.3-4.5: .... not clear. If you have external bit available, you should
use it directly correcting the I,Q residual phase. At least this is dobe in the
EUM processing.

p 24: -navfile option. It is stated that the default is to use th e internal
correction. What does it mean? How it works?

At the ROPP-10 DRR we agreed to tidy this up a bit for ROPP-10, but to return to it in earnest at ROPP-11.0.

[Ian Culverwell]

ESA report 2012: "GRAS Radio Occultation Performance Study"

[Ian Culverwell]

Section 7.1.1 ("L1 Preprocessing") of the attached 2012 ESA report describes the DMI/ROPP processing of GRAS navigation bits in more detail than the user guide. The key original references appear to be

1. Sokolovskiy, S., Y.-H. Kuo, C. Rocken, W. S. Schreiner, D. Hunt, and R. A. Anthes, Monitoring the atmospheric boundary layer by GPS radio occultation signals recorded in the open loop mode, Geophys. Res. Lett., 33, L12,813, doi:10.1029/2006GL025,955, 2006.
2. Sokolovskiy, S., C. Rocken, W. Schreiner, D. Hunt, and J. Johnson, Postprocessing of L1 GPS radio occultation signals recorded in open-loop mode, Radio Science, 44, RS2002, doi:10.1029/2008RS003,907, 2009.
3. Beyerle, G., M. E. Gorbunov, and C. O. Ao, Simulation studies of GPS radio occultation measurements, Radio Sci., 38, 1084, doi: 10.1029/2002RS002,800, 2003.

[Ian Culverwell]

While writing this up, I did some tests on the ropp_occ_tool core test file ropp_pp/data/ropp_pp_test.nc (attached). ropp_pp_openloop.f90 returned the following results:

----------------------------------------------------------------------
               ROPP Occultation Pre-processor Tool
----------------------------------------------------------------------

INFO:  Reading configuration file ../config/cosmic_pp.cf.

INFO:  Processing profile    1 of      1
INFO:  Reading input data file ../data/ropp_pp_test.nc.

... (from ropp_io_read_ncdf_get):  'start_time' and yr/mo/dy/hr/mn/sc/ms timestamps differ by 2.000 seconds (probably a leap-second issue) - using yr/../ms timestamp
INFO:  (OC_20090330060731_C001_U999_UNKN) 

...:  Occultation date (dd/mm/yyyy): =  30/03/2009
...:  Occultation time (hh:mm:ss): =  06:07:31
...:  Occultation point (lat,lon): =  (   0.31N,  75.95E)
...:  Undulation (m) =     -0.10091E+03
 
WARNING (from ropp_pp_preprocess):  Data from unknown processing centre ... guessing provenance
INFO (from ropp_pp_preprocess):  COSMIC data preprocessing
INFO (from ropp_pp_preprocess):  No lost carrier flag data in input file
INFO (from ropp_pp_preprocess):  No external NDM data read
INFO (from ropp_pp_preprocess):  COSMIC data: openloop preprocessing
 3.1: LCF =  0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ...
 3.1: imin, imax =  7252 1
 3.2: kmin, kmax =  7252 1
 3.4: DDS =  0. 0.029621517300825313 0.024702292021843184 0.030456155476614354 0.031449725948769124 0.02466501371996672 ...
 3.4: DNB =  -1626074700 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ...
 3.4: CNB =  0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.
 3.4: inb =  -10
... (from ropp_pp_openloop):  Navbit correllation area:          144.4           -0.6
... (from ropp_pp_openloop):  Max. navbit correlation:    0.000 at  -10
 3.5 before shifting navbits: LCF =  0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ...
 3.5 after shifting navbits: LCF =  0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ...
 CNBmax =  0.
 CNBmax =  0.5
 6.1 SIZE(time) =  0
 6.2 ropp_pp_internal_navbit: TGM = 
 6.2 ropp_pp_internal_navbit: TGI = 
 6.3 ropp_pp_internal_navbit: NB = 
 6.4 ropp_pp_internal_navbit: NB = 
 3.6 after CALL ropp_pp_internal_navbit: LCF =  0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ...

Thus, because no open loop data have apparently entered the routine (probably because no LCF in input file), LCF is initially zero, which means imin and imax stay at n and 1 resp., which causes all sorts of chaos - zero navbit correlation functions DNB, arrays of length 0 to be passed to ropp_pp_internal_navbit, etc. I'm sure this isn't intended. * To be addressed at ROPP-11.0 *.

[Ian Culverwell]

ropp_pp_test.nc

[Ian Culverwell]

I don't understand Eqn 4.2 of the ROPP PP user guide either: . This isn't dimensionally consistent, as p is an impact parameter ([p] = L). The answer is probably to be found somewhere in ropp_pp_modelphase.f90, but I can't see it right now.

[Ian Culverwell]

Eqn_4.2_PP_UG.png

[Ian Culverwell]

I'm confused about

  DS(1,:) = k(1)*(phase_L1(:) - phase_LM(:))
  DS(2,:) = k(2)*(phase_L2(:) - phase_LM(:))

in ropp_pp_openloop.f90 too. It appears that the same phase_LM generates different phases (in radians) for L1 and L2. I might have expected something more like

  DS(1,:) = k(1)*phase_L1(:) - phase_LM(:)
  DS(2,:) = k(2)*phase_L2(:) - phase_LM(:)

(since we know that phase_L1/2 have already been divided by k). k ~ 30 m^-1, (i.e. not immensely large or small) so that change wouldn't wreck the equations.

To put it another way, if k(1)*phase_L1 = k(1)*phase_L2, the ROPP formula would give different excess phases (in rads) to L1 and L2. Is that correct?

[Ian Culverwell]

Nor do I really understand

   DO i=2,n
      IF (BTEST(LCF(i),0)) THEN  ! In open loop mode
         IF (BTEST(LCF(i-1),bit_NBQ) .AND. BTEST(LCF(i),bit_NBQ)) THEN
            DS(ic,i) = DS(ic,i-1) +      &
               MODULO(DS(ic,i)-DS(ic,i-1)+pi, 2.0_wp*pi)-pi
         ELSE
            DS(ic,i) = DS(ic,i-1) +      &
               MODULO(DS(ic,i)-DS(ic,i-1)+pi/2.0_wp, pi)-pi/2.0_wp
         END IF
       ELSE                               ! In phase locked loop mode
         IF (Restore_PLL) THEN
            DS(ic,i) = DS(ic,i-1) +      &
               MODULO(DS(ic,i)-DS(ic,i-1)+pi/2.0_wp, pi) - pi/2.0_wp
         ELSE
            DS(ic,i) = DS(ic,i-1) +      &
               MODULO(DS(ic,i)-DS(ic,i-1)+pi, 2.0_wp*pi) - pi
         END IF
      END IF
   END DO

 END DO

in ropp_pp_openloop.f90, especially as it is immediately followed by

  DO ic=1,2
    DO i=2,n
      DS(ic,i) = DS(ic,i-1) + MODULO(DS(ic,i)-DS(ic,i-1)+pi, 2.0_wp*pi) - pi
    END DO
  END DO

I suppose doing the phase accumulation twice shouldn't hurt, but it is confusing. Why is the phase difference limited to (-pi/2, pi/2] in the other cases? I'll have a look in the original references.

[Ian Culverwell]

Having slept on it, I think the second objection is probably wrong. phase_LM is a physical distance in m: the difference between the ray path length and the straight-line distance between the satellites, if the intervening atmosphere were filled with the MSIS refractivity field. So it would lead to different phase lags (i.e. number of cycles) at the L1 and L2 frequencies.

[Ian Culverwell]

I think that Eqn 4.2 should be integral of (dp/dt) dt, so I've made that correction.

[Ian Culverwell]

Outstanding questions.

1. What's the point of the first line of this block, in Sec 6 of ropp_pp_preprocess.f90:

       WHERE(ro_data%Lev1a%phase_L2(:) < ropp_MDTV)
         ro_data%Lev1a%phase_L2(:) = -1.0_wp
         ro_data%Lev1a%phase_L2(:) = phase_LM(:)
       ENDWHERE
   

2. I don't really understand these lines in ropp_pp_openloop.f90:

    DO i=2,n
         IF (BTEST(LCF(i),0)) THEN  ! In open loop mode
            IF (BTEST(LCF(i-1),bit_NBQ) .AND. BTEST(LCF(i),bit_NBQ)) THEN
               DS(ic,i) = DS(ic,i-1) +      &
                  MODULO(DS(ic,i)-DS(ic,i-1)+pi, 2.0_wp*pi)-pi
            ELSE
               DS(ic,i) = DS(ic,i-1) +      &
                  MODULO(DS(ic,i)-DS(ic,i-1)+pi/2.0_wp, pi)-pi/2.0_wp
            END IF
          ELSE                               ! In phase locked loop mode
            IF (Restore_PLL) THEN
               DS(ic,i) = DS(ic,i-1) +      &
                  MODULO(DS(ic,i)-DS(ic,i-1)+pi/2.0_wp, pi) - pi/2.0_wp
            ELSE
               DS(ic,i) = DS(ic,i-1) +      &
                  MODULO(DS(ic,i)-DS(ic,i-1)+pi, 2.0_wp*pi) - pi
            END IF
         END IF
      END DO
   
    END DO
   

but how much do I want or need to? (It's the +/- pi/2 bits I don't get, and Sergey's 2009 paper didn't explain them to me.)

3. Why do eum2ropp (and ropp_pp_grasrs2ropp) use the ingested internal navbits to define the 'alternative' navbit elements (bits 4 and 4), rather than bits 7 and 8, which would seem more appropriate? The problem is that when subroutine ropp_pp_internal_navbit derives the internal navbits, it uses them to define bits 7 and 8, as expected.

[Ian Culverwell]

At the ROPP-11 DRR on 12/11/2021, Riccardo accepted that the ROPP documentation of this stuff is OK now, so mark this as done. (This was done in response to an action from the ROPP-10 DRR.)