Documentation for RefNLPEP2D forward operator and variable transforms

docs/inside/jedi-components/ufo/obsops.rst

@@ -1352,6 +1352,136 @@ Examples of yaml:
        - name: [atmosphericRefractivity]
        threshold: 3
 
+GNSS RO non-local pseudo excess phase (NLPEP 2D)
+------------------------------------------------
+
+Description:
+^^^^^^^^^^^^
+
+The RefNLPEP2D forward operator is a JEDI UFO implementation of the non-local 
+pseudo excess phase approach of Sokolovskiy, et al, 2005. This is logically similar
+to the GPSEPH forward operator in WRFDA, but shares no code with that earlier 
+implementation.
+
+RefNLPEP2D is a 2D forward operator based integrating excess refractive index, 
+computed from refractivity, over a predefined GNSS RO ray path that does not
+change during the data assimilation process. The same ray path is used to compute a 
+derived observation of non-local pseudo excess phase from refractivity observations 
+and the assumption of spherical symmetry. This operator allows a model's horizontal 
+variations in refractivity to be used when simulating GNSS RO observations.
+
+The RefNLPEP2D operator assumes that observations of nonLocalPseudoExcessPhase
+have been generated prior to entering the forward operator code. These 
+observations must have used the same ray path definition that is used in the
+forward operator. This is most easily accomplished within JEDI using the 
+:code:`NonLocalPseudoExcessPhase` variable transform to convert observed
+atmosphericRefractivity to nonLocalPseudoExcessPhase.
+
+Configuration options (ObsOperator):
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+All obs operator options are under the variable-grouping :code:`obs options`.
+
+* :code:`n_horiz`: the horizontal points the operator integrates along the 2d plane. Default is 11. Has to be an odd number > 0.
+
+* :code:`res`: The horizontal resolution of the 2d plane in km. Default is 11 km. The length of the ray path is :code:`n_horiz` times :code:`res`.
+
+* :code:`ray_length`: An alternate way to specify the ray path length in km. If set, it will override :code:`n_horiz` and :code:`res`. However, the ray length is only approximate, as it may need to be adjusted to ensure :code:`n_horiz` is an odd number. Therefore, use of this parameter is not recommended.
+
+* :code:`top_2d`: the highest geometric height to apply the 2d operator. Default is 90 km. Heights above :code:`top_2d` will use only a single vertical column of model data located at the tangent point and assume a spherically-symmetric atmosphere. Otherwise each ray will use n_horiz columns of model data.
+
+* :code:`ray_path_gen_type`: Subclass of GnssroRayPathGenerator to use to create the pre-defined ray path. Default is :code:`StraightLine`. This is currently the only option, but other path geometries are planned in the future. 
+
+* :code:`refr_algo`: Subclass of RefractivityCalculator to use to compute refractivity and its derivatives from temperature, pressure, and specific humidity. Default is :code:`ReugerBevis`, which uses the Reuger, 2002 model if :code:`use_compress` is true, and the Bevis, 1994 model is :code:`use_compress` is false. No other values are currently supported, but this configuration could be used in the future to test other refractivity models.
+
+* :code:`use_compress`: true if atmospheric compressibility should be used; false otherwise. Defaults is true.
+
+Configuration options (ObsSpace):
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+* :code:`obsgrouping`: Applying sequenceNumber as group_variable can get RO observations to be grouped into vertical profiles. Otherwise RO data will be treated as single observations. Because 2D rays can use information from multiple levels in the same RO profile, grouping by RO profiles is necessary for proper data handling in this operator.
+
+Configuration options (ObsFilter):
+
+* :code:`Variable Transform`: filter = :code:`NonLocalPseudoExcessPhase` is required to compute derived observations of non-local pseudo excess phase from observations of atmosphericRefractivity read from input files. This filter is required unless input observations files contain nonLocalPseudoExcessPhase as a pre-computed ObsValue. Configuration of :code:`res`, :code:`top_2d`, and :code:`n_horiz` within the variable transform must match the configuration of the GnssroRefNLPEP2D forward operator. :code:`SkipWhenNoObs` must be false if any MPI rank has no obs in it or filter FinalCheck will fail.
+
+* :code:`Variable Transform`: filter = :code:`GnssroRefractivityGradient`: This transform is optional and is only needed for diagnostic purposes. It computes the vertical derivative of refractivity with respect to geometric height (dN/dz) and saves this as a derived observation. It optionally computes and saves a ducting category for each vertical level if :code:`calculate ducting flag` is true, and the ducting category of the most extreme level, constant for all levels in a single RO profile, if :code:`calculate profile ducting flag` is true. The ducting category as defined as fractions of the critical ducting value -1/Rearth. The smallest fraction isgiven by the configuration parameter :code:`min super refraction` and must be >= 0 and < 0.6. The other ducting category thresholds are are hard-coded as 0.6, 0.7, 0.8, 0.9, and 1.0. :code:`SkipWhenNoObs` must be false if any MPI rank has no obs in it or filter FinalCheck will fail. 
+
+* :code:`Domain Check`: a generic filter used to control the maximum height one wants to assimilate RO observation.
+
+* :code:`ROobserror`: a RO specific filter. Use generic filter class to apply observation error method.  More information on this filter is found in the :doc:`observation uncertainty documentation <obserrors>`. 
+
+  * options: :code:`filter variables`. The :code:`name` = :code:`nonLocalPseudoExcessPhase` must be specified when used with the RefNLPEP2D forward operator.
+  * options: :code:`errmodel` = :code:`NCEP`, with more to come. The default errmodel is :code:`NBAM`, but this value is not supported for the nonLocalPseudoExcessPhase variable. One must specify errmodel explicitly for the RefNLPEP2D forward operator.
+
+Examples of yaml:
+^^^^^^^^^^^^^^^^^
+
+.. code-block:: yaml
+
+ observations:
+   observers:
+   - obs space:
+       name: GnssroRefNLPEP2D
+       obsdatain:
+         engine:
+           type: H5File
+           obsfile: Data/ufo/testinput_tier_1/gnssro_obs_2025070400_refnlpep2d.nc4
+         obsgrouping:
+           group_variables: 
+           - sequenceNumber
+           sort_variable: "height"
+           sort order: ascending
+       obsdataout:
+         engine:
+           type: H5File
+           obsfile: Data/gnssro_output_2025070400_refnlpep2d.nc4
+       observed variables: []
+       derived variables: [nonLocalPseudoExcessPhase]
+       simulated variables: [nonLocalPseudoExcessPhase]
+     obs operator:
+       name: GnssroRefNLPEP2D 
+       obs options:
+         res: 11.0
+         top_2d: 60.0
+         n_horiz: 19
+         refr_algo: "RuegerBevis"
+         use_compress: true
+     obs filters:
+     - filter: Variable Transforms
+       Transform: NonLocalPseudoExcessPhase
+       UseValidDataOnly: false
+       SkipWhenNoObs: false
+       res: 11.0
+       top_2d: 60.0
+       n_horiz: 19
+       group: "ObsValue"
+       refractivity variable: "atmosphericRefractivity"
+       nlpep variable: "nonLocalPseudoExcessPhase"
+     - filter: ROobserror
+       filter variables:
+       - name: nonLocalPseudoExcessPhase
+       errmodel: NCEP
+     - filter: Variable Transforms
+       Transform: GnssroRefractivityGradient
+       UseValidDataOnly: false
+       SkipWhenNoObs: false
+       group: "ObsValue"
+       min super refraction: 0.5
+       calculate ducting flag: true
+       calculate profile ducting flag: true
+     - filter: Domain Check
+       where:
+       - variable:
+           name: MetaData/height
+         minvalue: 0
+         maxvalue: 62000
+
+References:
+^^^^^^^^^^^
+
+* Sokolovskiy, S., Y. Kuo, and W. Wang, 2005: Assessing the Accuracy of a Linearized Observation Operator for Assimilation of Radio Occultation Data: Case Simulations with a High-Resolution Weather Model. Mon. Wea. Rev., 133, 2200–2212, https://doi.org/10.1175/MWR2948.1.
+
 Ground Based GNSS observation operator (Met Office)
 ---------------------------------------------------
 

docs/inside/jedi-components/ufo/variabletransforms/transforms/gnssrorefractivitygradient.rst

@@ -0,0 +1,73 @@
+
+.. _VT-gnssrorefractivitygradient:
+
+================================================
+Refractivity Vertical Gradient from Refractivity 
+================================================
+
+This transform computes the vertical derivative of refractivity
+with respect to height in a GNSS RO profile. It can optionally
+compute categorical values for ducting and near-ducting conditions.
+
+This is typically used with the GnssroRefNLPEP2D forward operator.
+It is used for diagnostic purposes only.
+
+---------------
+Input Variables
+---------------
+
+- atmosphericRefractivity (N-units)
+- height (m), typically a geometric height
+
+----------------
+Output Variables
+----------------
+
+- atmosphericRefractivityGradient (N-units/meter)
+- sampleDuctingFlag (nondim, integer 0 - 6)
+- profileDuctingFlag (nondim, integer 0 - 6) 
+
+------------------------
+Configuration Parameters
+------------------------
+
+- ``group``: Group name where input refractivity observation is found (default: :code:`ObsValue`)
+- ``refractivity variable``: Name of variable containing input refractivity observation (default: :code:`atmosphericRefractivity`)
+- ``height variable``: Name of variable containing the input height (default: :code:`height`)
+- ``grad variable``: Name of variable containing the output refractivity gradient (default: :code:`atmosphericRefractivityGradient`)
+- ``min super refraction``: Lowest fraction of 1/Rearth consider to be super-refracting (default: :code:`0.5`) 
+- ``calculate ducting flag``: Boolean to control computation and output of a categorical value indicating near-ducting and ducting conditions for each computed gradient. Default is false. Output variable is named :code:`sampleDuctingFlag`.
+- ``calculate profile ducting flag``: Boolean to control computation and output of a categorical value indicating near-ducting and ducting conditions the most extreme sample in a RO profile. All levels in a single profile get the same value. Default is false. Output variable is named :code:`profileDuctingFlag`.
+
+------------------
+Example yaml block
+------------------
+
+.. code-block:: yaml
+
+    obs filters:
+    - filter: Variable Transforms
+      Transform: GnssroRefractivityGradient
+      group: "ObsValue"
+      min super refraction: 0.5
+      calculate ducting flag: true
+      calculate profile ducting flag: true
+
+-------
+Methods
+-------
+
+Iterate over all levels in each RO profile, computing dN/dz using the difference in refractivity
+and height between adjacent levels. dN/dz is set to missing for the topmost level.
+
+If the ducting flag outputs are enabled, bin the absolute value of dN/dz values into the following categories:
+
+- 0 if abs(dN/dz) < :code:`min super refraction` * 1/Re
+- 1 if abs(dN/dz) >= :code:`min super refraction` * 1/Re and < 0.6/Re
+- 2 if abs(dN/dz) >= 0.6/Re and < 0.7/Re
+- 3 if abs(dN/dz) >= 0.7/Re and < 0.8/Re
+- 4 if abs(dN/dz) >= 0.8/Re and < 0.9/Re
+- 5 if abs(dN/dz) >= 0.9/Re and < 1.0/Re
+- 6 if abs(dN/dz) >= 1.0/Re
+
+Here, Re is the radius of the earth in meters.

docs/inside/jedi-components/ufo/variabletransforms/transforms/nonlocalpseudoexcessphase.rst

@@ -0,0 +1,64 @@
+
+.. _VT-nonlocalpseudoexcessphase:
+
+===============================================
+Non-Local Pseudo Excess Phase from Refractivity 
+===============================================
+
+This transform computes non-local pseudo excess phase from refractivity in GNSS RO observations.
+This is typically used with the GnssroRefNLPEP2D forward operator.
+
+---------------
+Input Variables
+---------------
+
+- atmosphericRefractivity (N-units)
+- geometric height (m)
+
+----------------
+Output Variables
+----------------
+
+- nonLocalPseudoExcessPhase (meters) 
+
+------------------------
+Configuration Parameters
+------------------------
+
+- ``group``: Group name where input refractivity observation is found (default: :code:`ObsValue`)
+- ``refractivity variable``: Name of variable containing input refractivity observation (default: :code:`atmosphericRefractivity`)
+- ``nlpep variable``: Name of variable containing the output non-local pseudo excess phase (default: :code:`nonLocalPseudoExcessPhase`)
+- ``ray_path_gen_type``: Method for generating ray path. Must match the GnssroRefNLPEP2D forward operator and has same defaults.
+- ``ray_length``: approximate length of ray path. Must match the GnssroRefNLPEP2D forward operator and has same defaults.
+- ``res``: Horizontal resolution of nodes in ray path. Must match the GnssroRefNLPEP2D forward operator and has same defaults.
+- ``top_2d``: Highest geometric height where 2D NLPEP is computed. Must match the GnssroRefNLPEP2D forward operator and has same defaults.
+- ``n_horiz``: Number of nodes in the ray path. Must match the GnssroRefNLPEP2D forward operator and has same defaults. 
+
+------------------
+Example yaml block
+------------------
+
+.. code-block:: yaml
+
+    obs filters:
+    - filter: Variable Transforms
+      Transform: NonLocalPseudoExcessPhase
+      res: 11.0
+      top_2d: 60.0
+      n_horiz: 19
+      group: "ObsValue"
+      refractivity variable: "atmosphericRefractivity"
+      nlpep variable: "nonLocalPseudoExcessPhase"
+
+-------
+Methods
+-------
+
+Iterates over all nodes in the ray path defined for each RO observation sample (tangent point).
+Determines refractivity at the node of a ray as a function of geometric height by vertically 
+interpolating from the observations at heights above the tangent point in the same Radio Occultation 
+profile. It converts the refractivity :code:`N` to excess index of refraction :code:`n-1` and 
+integrates it over the length of the ray segment associated with that node. 
+
+If the refractivity for any node in a ray is missing, the non-local pseudo excess phase for the 
+entire ray is set to missing.