Skip to content

Commit

Permalink
mp:
Browse files Browse the repository at this point in the history 
 - implement the mp code biases estimator in the form
   of a trait
 - add documentation

Signed-off-by: Guillaume W. Bres <[email protected]>
  • Loading branch information
@gwbres
gwbres committed Jun 11, 2023
1 parent a55ffe2 commit 9761c87
Show file tree
Hide file tree
Showing 8 changed files with 205 additions and 211 deletions.
Binary file added doc/plots/esbc00dnk_g13_dcb_mp.png
View file
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
45 changes: 27 additions & 18 deletions rinex-cli/doc/processing.md
View file
Original file line number Diff line number Diff line change
Expand Up @@ -21,6 +21,11 @@ This tool supports standard GNSS recombinations, especially for modern RINEX.
Refer to this page for [thorough documentation](gnss-combination.md).
It is important to understand GNSS signal recombinations and what they can represent.

The RINEX tool suites can estimate the following code biases:

- the DCB code biases
- the MP code biases

Differential Code Biases (DBCs)
===============================

Expand Down Expand Up @@ -55,20 +60,9 @@ rinex-cli \

<img align="center" width="650" src="https://github.com/gwbres/rinex/blob/main/doc/plots/esbc00dnk_ph_dcbs.png">

Differential Processing
=======================

When moving to more advanced RINEX processing,
Navigation RINEX (Ephemeris) must be provided with `--nav`.

In this mode, `--fp` is expected to be an Observation RINEX.

Let's remind the user that in this mode, `--sv-epoch` helps
exhibit which vehicles share Ephemeris and Observations for a given epoch
or epoch range. This feature is very important do determine
which vehicle is a good candidate for the operations that follow.

## Code Multipath (MP) analysis
Code Multipath biases
=====================

MP ratios (so called "CMC" for Code Minus Carrier metrics)
are formed by combining Phase and PR observations sampled against different carrier frequencies
Expand All @@ -83,18 +77,33 @@ A very compelling use case of MP code analysis
can be
[found here](https://www.taoglas.com/wp-content/uploads/pdf/Multipath-Analysis-Using-Code-Minus-Carrier-Technique-in-GNSS-Antennas-_WhitePaper_VP__Final-1.pdf).


MP analysis is summoned with `--mp`.

In this example, `ESBC00DNK_R_2020` vehicle GPS#13
has enough data to compute MP ratio for codes "1C", "2W" and "5Q".
This operation requires combining the associated Ephemeris, that we also provide
for demonstration purposes:
has enough data to compute MP ratios:

```bash
rinex-cli \
--fp test_resources/CRNX/V3/ESBC00DNK_R_20201770000_01D_30S_MO.crx.gz \
--nav test_resources/NAV/V3/ESBC00DNK_R_20201770000_01D_MN.rnx.gz \
-P G13 --mp
-P G13 --dcb --mp
```

<img align="center" width="650" src="https://github.com/gwbres/rinex/blob/main/doc/plots/esbc00dnk_g13_dcb_mp.png">

Differential Processing
=======================

When moving to more advanced RINEX processing,
Navigation RINEX (Ephemeris) must be provided with `--nav`.

In this mode, `--fp` is expected to be an Observation RINEX.

Let's remind the user that in this mode, `--sv-epoch` helps
exhibit which vehicles share Ephemeris and Observations for a given epoch
or epoch range. This feature is very important do determine
which vehicle is a good candidate for the operations that follow.

## Cycle slips analysis

Under development
11 changes: 6 additions & 5 deletions rinex-cli/src/main.rs
View file
Original file line number Diff line number Diff line change
Expand Up @@ -88,19 +88,20 @@ pub fn main() -> Result<(), rinex::Error> {
/*
* Code Multipath analysis
*/
/*if cli.multipath() {
if cli.multipath() {
let data = ctx
.primary_rinex
.observation_align_phase_origins()
.observation_code_multipath();
plot::plot_gnss_recombination(
.observation_phase_align_origin()
.observation_phase_carrier_cycles()
.mp();
plot::plot_gnss_dcb(
&mut plot_ctx,
"Code Multipath Biases",
"Meters of delay",
&data,
);
info!("--mp analysis");
}*/
}
/*
* [GF] recombination visualization requested
*/
Expand Down
10 changes: 9 additions & 1 deletion rinex/src/algorithm/dcb.rs
View file
Original file line number Diff line number Diff line change
@@ -1,4 +1,4 @@
//! GNSS code biases estimator
//! GNSS code biases estimators
use crate::prelude::*;
use std::collections::{BTreeMap, HashMap};

Expand All @@ -19,3 +19,11 @@ pub trait Dcb {
/// ```
fn dcb(&self) -> HashMap<String, BTreeMap<Sv, BTreeMap<(Epoch, EpochFlag), f64>>>;
}

/// Code multipath analysis (MP_i), cf.
/// phase data model <https://github.com/gwbres/rinex/blob/main/rinex-cli/doc/gnss-combination.md>.
pub trait Mp {
/// Returns Code Biases estimates due to multipath reflection,
/// sorted per (unique) signal combinations aod for each individual Sv.
fn mp(&self) -> HashMap<String, BTreeMap<Sv, BTreeMap<(Epoch, EpochFlag), f64>>>;
}
2 changes: 1 addition & 1 deletion rinex/src/algorithm/mod.rs
View file
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,7 @@ mod processing;
mod target;

pub use combination::{Combination, Combine};
pub use dcb::Dcb;
pub use dcb::{Dcb, Mp};
pub use ionospheric::IonoDelayDetector;
pub use target::TargetItem;

Expand Down
193 changes: 7 additions & 186 deletions rinex/src/lib.rs
View file
Original file line number Diff line number Diff line change
Expand Up @@ -96,7 +96,7 @@ use prelude::*;
pub use merge::Merge;
pub use split::Split;

use algorithm::{Combination, Combine, Dcb, IonoDelayDetector, Smooth};
use algorithm::{Combination, Combine, Dcb, IonoDelayDetector, Mp, Smooth};

#[macro_use]
extern crate horrorshow;
Expand Down Expand Up @@ -1931,191 +1931,6 @@ impl Rinex {
HashMap::new()
}
}
/// Code multipath analysis (MP_i), cf.
/// phase data model <https://github.com/gwbres/rinex/blob/main/rinex-cli/doc/gnss-combination.md>.
pub fn observation_code_multipath(
&self,
) -> HashMap<String, HashMap<Sv, BTreeMap<(Epoch, EpochFlag), f64>>> {
let mut ret: HashMap<String, HashMap<Sv, BTreeMap<(Epoch, EpochFlag), f64>>> =
HashMap::new();
if let Some(record) = self.record.as_obs() {
/*
* Determine mean value of all datasets
*/
let mut mean: HashMap<Sv, HashMap<Observable, (u32, f64)>> = HashMap::new();
for (_epoch, (_, vehicles)) in record {
for (sv, observations) in vehicles {
if let Some(data) = mean.get_mut(&sv) {
for (observable, obs_data) in observations {
if observable.is_phase_observable()
|| observable.is_pseudorange_observable()
{
if let Some((count, buf)) = data.get_mut(observable) {
*count += 1;
*buf += obs_data.obs;
} else {
data.insert(observable.clone(), (1, obs_data.obs));
}
}
}
} else {
for (observable, obs_data) in observations {
if observable.is_phase_observable()
|| observable.is_pseudorange_observable()
{
let mut map: HashMap<Observable, (u32, f64)> = HashMap::new();
map.insert(observable.clone(), (1, obs_data.obs));
mean.insert(*sv, map);
}
}
}
}
}
mean.iter_mut()
.map(|(_, data)| {
data.iter_mut()
.map(|(_, (n, buf))| {
*buf = *buf / *n as f64;
})
.count()
})
.count();
//println!("MEAN VALUES {:?}", mean); //DEBUG
/*
* Run algorithm
*/
let mut associated: HashMap<String, String> = HashMap::new(); // Ph code to associate to this Mpx
// for operation consistency
for (epoch, (_, vehicles)) in record {
for (sv, observations) in vehicles {
let _mean_sv = mean.get(&sv).unwrap();
for (lhs_observable, lhs_data) in observations {
if lhs_observable.is_pseudorange_observable() {
let pr_i = lhs_data.obs; // - mean_sv.get(lhs_code).unwrap().1;
let lhs_code = lhs_observable.to_string();
let mp_code = &lhs_code[2..]; //TODO will not work on RINEX2
let lhs_carrier = &lhs_code[1..2];
let mut ph_i: Option<f64> = None;
let mut ph_j: Option<f64> = None;
/*
* This will restrict combinations to
* 1 => 2
* and M => 1
*/
let rhs_carrier = match lhs_carrier {
"1" => "2",
_ => "1",
};
/*
* locate related L_i PH code
*/
for (observable, data) in observations {
let ph_code = format!("L{}", mp_code);
let code = observable.to_string();
if code.eq(&ph_code) {
ph_i = Some(data.obs); // - mean_sv.get(code).unwrap().1);
break; // DONE
}
}
/*
* locate another L_j PH code
*/
if let Some(to_locate) = associated.get(mp_code) {
/*
* We already have an association, keep it consistent throughout
* operations
*/
for (observable, data) in observations {
let code = observable.to_string();
let carrier_code = &code[1..2];
if carrier_code == rhs_carrier {
// correct carrier signal
if code.eq(to_locate) {
// match
ph_j = Some(data.obs); // - mean_sv.get(code).unwrap().1);
break; // DONE
}
}
}
} else {
// first: prefer the same code against rhs carrier
let to_locate = format!("L{}{}", rhs_carrier, &mp_code[1..]);
for (observable, data) in observations {
let code = observable.to_string();
let carrier_code = &code[1..2];
if carrier_code == rhs_carrier {
// correct carrier
if code.eq(&to_locate) {
// match
ph_j = Some(data.obs); // - mean_sv.get(code).unwrap().1);
associated.insert(mp_code.to_string(), code.clone());
break; // DONE
}
}
}
if ph_j.is_none() {
/*
* Same code against different carrier does not exist
* try to grab another PH code, against rhs carrier
*/
for (observable, data) in observations {
let code = observable.to_string();
let carrier_code = &code[1..2];
if carrier_code == rhs_carrier {
if observable.is_phase_observable() {
ph_j = Some(data.obs); // - mean_sv.get(code).unwrap().1);
associated
.insert(mp_code.to_string(), code.clone());
break; // DONE
}
}
}
}
}
if ph_i.is_none() || ph_j.is_none() {
break; // incomplete associations, do not proceed further
}
let ph_i = ph_i.unwrap();
let ph_j = ph_j.unwrap();
let lhs_carrier = lhs_observable.carrier(sv.constellation).unwrap();
let rhs_carrier =
lhs_observable //rhs_observable TODO
.carrier(sv.constellation)
.unwrap();
/*let gamma = (lhs_carrier.frequency() / rhs_carrier.frequency()).powf(2.0);
let alpha = (gamma +1.0_f64) / (gamma - 1.0_f64);
let beta = 2.0_f64 / (gamma - 1.0_f64);
let mp = pr_i - alpha * ph_i + beta * ph_j;*/

let alpha = 2.0_f64 * rhs_carrier.frequency().powf(2.0)
/ (lhs_carrier.frequency().powf(2.0)
- rhs_carrier.frequency().powf(2.0));
let mp = pr_i - ph_i - alpha * (ph_i - ph_j);
if let Some(data) = ret.get_mut(mp_code) {
if let Some(data) = data.get_mut(&sv) {
data.insert(*epoch, mp);
} else {
let mut bmap: BTreeMap<(Epoch, EpochFlag), f64> =
BTreeMap::new();
bmap.insert(*epoch, mp);
data.insert(*sv, bmap);
}
} else {
let mut bmap: BTreeMap<(Epoch, EpochFlag), f64> = BTreeMap::new();
let mut map: HashMap<Sv, BTreeMap<(Epoch, EpochFlag), f64>> =
HashMap::new();
bmap.insert(*epoch, mp);
map.insert(*sv, bmap);
ret.insert(mp_code.to_string(), map);
}
}
}
}
}
}
ret
}

/*
/// Single step /stage, in high order phase differencing
/// algorithm, which we use in case of old receiver data / old RINEX
Expand Down Expand Up @@ -2722,6 +2537,12 @@ impl Dcb for Rinex {
}
}

impl Mp for Rinex {
fn mp(&self) -> HashMap<String, BTreeMap<Sv, BTreeMap<(Epoch, EpochFlag), f64>>> {
self.record.dcb()
}
}

impl Combine for Rinex {
fn combine(
&self,
Expand Down
Loading

0 comments on commit 9761c87

Please sign in to comment.