From 5e056ccd6ecb881325bc8b82c3dec92daa197ef1 Mon Sep 17 00:00:00 2001 From: "Guillaume W. Bres" Date: Mon, 4 Sep 2023 08:14:06 +0200 Subject: [PATCH 01/19] move deprecated file All qc analysis moved into a module. Signed-off-by: Guillaume W. Bres --- rinex/src/qc/sampling.rs | 128 --------------------------------------- 1 file changed, 128 deletions(-) delete mode 100644 rinex/src/qc/sampling.rs diff --git a/rinex/src/qc/sampling.rs b/rinex/src/qc/sampling.rs deleted file mode 100644 index 7ad496aff..000000000 --- a/rinex/src/qc/sampling.rs +++ /dev/null @@ -1,128 +0,0 @@ -use crate::prelude::*; -use horrorshow::RenderBox; - -/// Sampling QC report -#[derive(Debug, Clone, PartialEq)] -#[cfg_attr(feature = "serde", derive(Serialize, Deserialize))] -pub struct QcReport { - /// First epoch identified - pub first_epoch: Epoch, - /// Last epoch identified - pub last_epoch: Epoch, - /// Time line - pub time_line: Duration, - /// Dominant sample rate - pub sample_rate: Duration, - /// Unusual data gaps - pub gaps: Vec<(Epoch, Duration)>, -} - -impl QcReport { - pub fn new(rnx: &Rinex) -> Self { - let first_epoch = rnx - .first_epoch() - .expect("Sampling QC expects a RINEX indexed by epochs"); - let last_epoch = rnx - .last_epoch() - .expect("Sampling QC expects a RINEX indexed by epochs"); - let sample_rate = rnx - .sampling_interval() - .expect("failed to determine sample rate"); - Self { - first_epoch, - last_epoch, - sample_rate, - time_line: last_epoch - first_epoch, - gaps: rnx.data_gaps(), - } - } - fn gap_analysis(gaps: &Vec<(Epoch, Duration)>) -> Box { - box_html! { - @ if gaps.len() == 0 { - tr { - th { - b { - : "Gap analysis: " - } - } - td { - : "Data missing" - } - } - } else { - tr { - td { - b { - : "Gap Aanalysis:" - } - } - } - tr { - td { - : "" - } - td { - : "Epoch (start)" - } - td { - : "Duration" - } - } - @ for (epoch, duration) in gaps { - tr { - td { - : "" - } - td { - : epoch.to_string() - } - td { - : duration.to_string() - } - } - } - } - } - } - pub fn to_inline_html(&self) -> Box { - box_html! { - table(id="sampling") { - tr { - td { - : "First Epoch:" - } - td { - : self.first_epoch.to_string() - } - } - tr { - td { - : "Last Epoch:" - } - td { - : self.last_epoch.to_string() - } - } - tr { - td { - : "Time line:" - } - td { - : self.time_line.to_string() - } - } - tr { - td { - : "Sample rate:" - } - td { - : self.sample_rate.to_string() - } - } - table(id="gap-analysis") { - : Self::gap_analysis(&self.gaps) - } - } - } - } -} From dedd15b53b39feb96c02c24dae4da844a8cf1c88 Mon Sep 17 00:00:00 2001 From: "Guillaume W. Bres" Date: Mon, 4 Sep 2023 08:14:35 +0200 Subject: [PATCH 02/19] working on OBS Data QC Signed-off-by: Guillaume W. Bres --- rinex/src/lib.rs | 46 ++++- rinex/src/qc/analysis/obs.rs | 282 ++++++++++++++++-------------- rinex/src/qc/analysis/sampling.rs | 18 +- rinex/src/qc/mod.rs | 50 +++++- rinex/src/qc/opts.rs | 3 + 5 files changed, 258 insertions(+), 141 deletions(-) diff --git a/rinex/src/lib.rs b/rinex/src/lib.rs index 59867e605..2969149ed 100644 --- a/rinex/src/lib.rs +++ b/rinex/src/lib.rs @@ -1688,7 +1688,7 @@ impl Rinex { svnn.iter() .flat_map(|(_sv, observables)| observables.keys()) }) - .fold(vec![], |mut list, items| { + /*.fold(vec![], |mut list, items| { // create a unique list for item in items { if !list.contains(&item) { @@ -1697,7 +1697,8 @@ impl Rinex { } list }) - .into_iter(), + .into_iter(),*/ + .unique() ) } else if let Some(_) = self.record.as_meteo() { Box::new( @@ -1810,6 +1811,47 @@ impl Rinex { #[cfg(feature = "obs")] #[cfg_attr(docrs, doc(cfg(feature = "obs")))] impl Rinex { + /// Returns a Unique Iterator over identified [`Carrier`]s + pub fn carrier(&self) -> Box> { + Box::new( + self.observation() + .flat_map(|(_, (_, sv))| { + sv.iter() + .flat_map(|(_, obsvervations)| { + observations + .keys() + .map(|observable| { + if let Ok(carrier) = observable.carrier(sv.constellation) { + Some(carrier) + } else { + None + } + }) + .unique() + }) + }) + ) + } + pub fn code(&self) -> Box> { + Box::new( + self.observation() + .flat_map(|(_, (_, sv))| { + sv.iter() + .flat_map(|(_, obsvervations)| { + observations + .keys() + .map(|observable| { + if let Ok(code) = observable.code() { + Some(code) + } else { + None + } + }) + .unique() + }) + }) + ) + } /// Returns ([`Epoch`] [`EpochFlag`]) iterator, where each {`EpochFlag`] /// validates or invalidates related [`Epoch`] /// ``` diff --git a/rinex/src/qc/analysis/obs.rs b/rinex/src/qc/analysis/obs.rs index b4255842b..908334467 100644 --- a/rinex/src/qc/analysis/obs.rs +++ b/rinex/src/qc/analysis/obs.rs @@ -33,104 +33,94 @@ fn report_signals(list: &Vec) -> String { /* * Epoch anomalies formatter */ -fn report_anomalies(anomalies: &Vec<(Epoch, EpochFlag)>) -> Box { - if anomalies.len() == 0 { - box_html! { - table(class="table is-bordered") { - th { - : "Anomalies" - } +fn report_anomalies( + cs: &Vec, + power: &Vec, + other: &Vec<(Epoch, EpochFlag)>, +) -> Box { + box_html! { + table(class="table is-bordered") { + tr { td { - : "None" + : "Power Failures" } } - } - } else { - box_html! { - table(class="table is-bordered") { - thead { + @ if power.is_empty() { + tr { th { - : "Anomalies" + : "None" } - th { - : "Power failure" + } + } else { + tr { + td { + : format!("{:?}", power) + } + } + tr { + td { + : "Longest" + } + td { + : power_failures.max_by(|(_, d1), (_, d2)| d1.cmp(d2)).unwrap().to_string() + } + td { + : "Average Power failure" } + td { + : "TODO" + } + } + } + tr { + td { + : "Cycle slip(s)" + } + } + tr { + @ if cs.is_empty() { th { - : "Antenna movement detected" + : "None" + } + } else { + td { + : format!("{:?}", cs) } + } + } + tr { + @ if other.is_empty() { th { - : "Kinematic" + : "Other anomalies" } th { - : "External event" + : "None" } + } else { th { - : "Cycle Slips" + : "Other anomalies" } - } - tbody { - @ for (epoch, _flag) in anomalies { - tr { - td { - : epoch.to_string() - } - /*@match flag { - EpochFlag::PowerFailure => { - td { - : "x" - } - }, - EpochFlag::AntennaBeingMoved => { - td { - : "" - } - td { - : "x" - } - }, - EpochFlag::NewSiteOccupation => { - td { - : "" - } - td { - : "" - } - td { - : "x" - } - }, - EpochFlag::ExternalEvent => { - td { - : "" - } - td { - : "" - } - td { - : "" - } - td { - : "x" - } - }, - EpochFlag::CycleSlip => { - td { - : "" - } - td { - : "" - } - td { - : "" - } - td { - : "" - } - td { - : "x" - } - }, - }*/ + td { + : "Epoch" + } + td { + : "Event" + } + @ for (e, event) in other { + td { + : "" + } + td { + : e.to_string() + } + @ if event == EpochFlag::AntennaBeingMoved { + : "Antenna being moved" + } else if event == EpochFlag::Kinematic { + : "Kinematic" + } else if event == EpochFlag::NewSiteOccupation { + : "New Site Occupation" + } else if event == EpochFlag::ExternalEvent { + : "External Event" } } } @@ -249,6 +239,15 @@ fn report_ssi_statistics( } } +fn derivative_dt(data: Vec<(Epoch, f64>), window: Duration) -> Vec<(Epoch, f64)> { + let mut acc = 0.0_f64; + let mut prev_epoch : Option = None + let mut ret: Vec<(Epoch, f64)> = Vec::new(); + for (epoch, value) in data { + let dt = + } +} + #[derive(Debug, Clone)] /// OBS RINEX specific QC analysis. /// Full OBS RINEX analysis requires both the "obs" and "processing" features. @@ -261,8 +260,12 @@ pub struct QcObsAnalysis { codes: Vec, /// true if doppler observations were identified has_doppler: bool, - /// Abnormal events, by chronological epochs - anomalies: Vec<(Epoch, EpochFlag)>, + /// CS anomalies + cs_anomalies: Vec, + /// Epochs where power failures took place, and their duration + power_failures: Vec, + /// Other abnormal events, by chronological epochs + other_anomalies: Vec<(Epoch, EpochFlag)>, /// Total number of epochs identified total_epochs: usize, /// Epochs with at least 1 observation @@ -275,31 +278,31 @@ pub struct QcObsAnalysis { #[cfg(feature = "obs")] /// SSi statistical analysis (mean, stddev, skew) ssi_stats: HashMap, - #[cfg(feature = "processing")] - /// Receiver clock drift analysis - clock_drift: Option, + /// RX clock drift + clock_drift: Vec, } impl QcObsAnalysis { pub fn new(rnx: &Rinex, opts: &QcOpts) -> Self { - let sv: Vec<_> = rnx.sv().collect(); - let obs = rnx.header.obs.as_ref().unwrap(); - let mut observables = obs.codes.clone(); - // TODO: improve this .get(0) - let observables = observables.get_mut(&sv[0].constellation).unwrap(); - let mut signals: Vec = Vec::new(); - let mut codes: Vec = Vec::new(); - let mut anomalies: Vec<(Epoch, EpochFlag)> = Vec::new(); - let mut total_epochs: usize = 0; + let mut observables : Vec<_> = rnx.observable().collect(); + let has_doppler = observables.iter() + .fold(|boolean, obs| { + boolean |= obs.is_doppler_observable() + }); + let mut signals : Vec<_> = rnx.signal().collect(); + let mut codes : Vec<_> = rnx.code().collect(); + let anomalies = rnx.epoch_anomalies(); + + let mut total_epochs = self.epoch().count(); let mut epoch_with_obs: Vec = Vec::new(); + let mut complete_epochs: HashMap = HashMap::new(); let mut min_max_snr = ( (Sv::default(), Epoch::default(), Snr::DbHz54), (Sv::default(), Epoch::default(), Snr::DbHz0), ); - let mut ssi_stats: HashMap = HashMap::new(); - let mut clock_drift: Option = None; + let mut ssi_stats: HashMap = HashMap::new(); if let Some(rec) = rnx.record.as_obs() { let mask: Filter = Filter::from(MaskFilter { @@ -307,34 +310,18 @@ impl QcObsAnalysis { item: TargetItem::ClockItem, }); let clk_data = rec.filter(mask); - //TODO: run analysis on smoothed derivative - //let der = clk_data.derivative(); - //let mov = der.moving_average(Duration::from_seconds(600.0), None); } if let Some(r) = rnx.record.as_obs() { total_epochs = r.len(); for ((epoch, flag), (clk, svs)) in r { - if let Some(_clk) = clk {} - if !flag.is_ok() { - anomalies.push((*epoch, *flag)); - } for (sv, observables) in svs { - if observables.len() > 0 && !epoch_with_obs.contains(&epoch) { + if !observables.is_empty() { epoch_with_obs.push(*epoch); } for (observable, observation) in observables { - let code = observable.code().unwrap(); - let carrier = observable.carrier(sv.constellation).unwrap(); - if !signals.contains(&carrier) { - signals.push(carrier); - } - if !codes.contains(&code) { - codes.push(code); - } - if let Some(snr) = observation.snr { if snr < min_max_snr.0 .2 { min_max_snr.0 .0 = *sv; @@ -350,7 +337,6 @@ impl QcObsAnalysis { } } } - /* * Now that signals have been determined, * determine observation completion @@ -449,28 +435,49 @@ impl QcObsAnalysis { } } } - + /* + * sort, prior reporting + */ codes.sort(); - signals.sort(); observables.sort(); + signals.sort(); Self { - observables: { observables.iter().map(|v| v.to_string()).collect() }, - has_doppler: { - let mut ret = false; - for obs in observables.iter() { - if obs.is_doppler_observable() { - ret = true; - break; - } - } - ret - }, codes, signals, - anomalies, + observables, + has_doppler, + cs_anomalies: { + anomalies.filter(|(e, flag) { + if flag == EpochFlag::CycleSlip { + Some(e) + } else { + None + } + }).collect() + }, + power_failures: { + anomalies.filter(|(e, flag) { + if flag == EpochFlag::PowerFailure { + Some(e) + } else { + None + } + }) + .collect() + }, + other_anomalies: { + anomalies.filter(|(e, flag) { + if flag != EpochFlag::PowerFailure && flag != EpochFlag::CycleSlip { + Some(e) + } else { + None + } + }) + .collect() + }, total_epochs, - total_with_obs: epoch_with_obs.len(), + total_with_obs: epoch_with_obs.unique().len(), complete_epochs: { let mut ret: Vec<(Carrier, usize)> = complete_epochs.iter().map(|(k, v)| (*k, *v)).collect(); @@ -479,7 +486,12 @@ impl QcObsAnalysis { }, min_max_snr, ssi_stats, - clock_drift, + clock_drift: { + let rx_clock : Vec<_> = rnx.recvr_clock() + .collect(); + let rx_clock_drift: Vec<(Epoch, f64)> = derivative_dt(rx_clock); + moving_average(rx_clock_drift, opts.clock_avg_window) + } } } } @@ -533,7 +545,7 @@ impl HtmlReport for QcObsAnalysis { } tr { table { - : report_anomalies(&self.anomalies) + : report_anomalies(&self.cs_anomalies, &self.power_failures, &self.other_anomalies) } } tr { diff --git a/rinex/src/qc/analysis/sampling.rs b/rinex/src/qc/analysis/sampling.rs index 3266f6962..9862341d9 100644 --- a/rinex/src/qc/analysis/sampling.rs +++ b/rinex/src/qc/analysis/sampling.rs @@ -92,13 +92,27 @@ impl HtmlReport for QcSamplingAnalysis { th { : "Sample rate" } + th { + : "Dominan Sample rate" + } + } + tr { @ if let Some(rate) = self.sample_rate { td { : format!("{} ({:.3} Hz)", rate, 1.0 / rate.to_unit(Unit::Second)) } } else { td { - : "Unknown" + : "Not Specified" + } + } + @ if let Some(rate) = self.dominant_sample_rate { + td { + : format!("{} ({:.3} Hz)", rate, 1.0 / rate.to_unit(Unit::Second)) + } + } else { + td { + : "Undetermined" } } } @@ -107,7 +121,7 @@ impl HtmlReport for QcSamplingAnalysis { : "Gap analysis" } } - @ if self.gaps.len() == 0 { + @ if self.gaps.is_empty() { tr { th { : "No Data Gaps detected" diff --git a/rinex/src/qc/mod.rs b/rinex/src/qc/mod.rs index bbe5c9c0f..2c094bb4c 100644 --- a/rinex/src/qc/mod.rs +++ b/rinex/src/qc/mod.rs @@ -186,10 +186,56 @@ impl QcReport { * Report all analysis that were performed */ div(id="analysis") { - @ for analysis in &analysis { + @ if analysis.len() == 1 { + /* + * Single analysis case + */ table(class="table is-bordered") { + thead { + tr { + td { + : "Multi constellation analysis" + } + } + } tbody { - : analysis.to_inline_html() + : self.analysis[0].to_inline_html() + } + } + } else { + /* + * Report all analysis + * and emphasize how they were sorted (self.opts.classfication) + */ + @ for i in 0..self.analysis.len() { + table(class="table is-bordered") { + @ match self.opts.classification { + QcClassification::GNSS => { + tr { + td { + : format!("{} analysis", self.constellation().nth(i).unwrap()) + } + } + }, + QcClassification::Sv => { + tr { + td { + : format!("{} analysis", self.sv().nth(i).unwrap()) + } + } + }, + QcClassification::Physics => { + tr { + td { + : format!("{} analysis", self.observable().nth(i).unwrap()) + } + } + }, + } + + tbody { + : self.analysis[0].to_inline_html() + } } } } diff --git a/rinex/src/qc/opts.rs b/rinex/src/qc/opts.rs index 61ce82b90..14faa3e01 100644 --- a/rinex/src/qc/opts.rs +++ b/rinex/src/qc/opts.rs @@ -161,6 +161,8 @@ pub struct QcOpts { pub gap_tolerance: Option, /// Manually defined Ground position (ECEF) pub ground_position: Option, + /// Window duration to be used, during RX clock drift analysis + pub clock_drift_window: Duration, } impl QcOpts { @@ -198,6 +200,7 @@ impl Default for QcOpts { min_snr_db: 20.0, // dB elev_mask: None, classification: QcClassification::default(), + clock_window: Duration::from_seconds(10.0 * 60.0), } } } From 92a354c7c381accb4002bb5468da00f81aa7a675 Mon Sep 17 00:00:00 2001 From: "Guillaume W. Bres" Date: Mon, 4 Sep 2023 09:36:21 +0200 Subject: [PATCH 03/19] fix a couple of syntax issues Signed-off-by: Guillaume W. Bres --- rinex/src/lib.rs | 64 ++++++++++------------ rinex/src/qc/analysis/obs.rs | 102 +++++++++++++++++++---------------- rinex/src/qc/mod.rs | 62 ++++++++++----------- 3 files changed, 114 insertions(+), 114 deletions(-) diff --git a/rinex/src/lib.rs b/rinex/src/lib.rs index 2969149ed..74ec74ac3 100644 --- a/rinex/src/lib.rs +++ b/rinex/src/lib.rs @@ -1688,7 +1688,7 @@ impl Rinex { svnn.iter() .flat_map(|(_sv, observables)| observables.keys()) }) - /*.fold(vec![], |mut list, items| { + .fold(vec![], |mut list, items| { // create a unique list for item in items { if !list.contains(&item) { @@ -1697,8 +1697,7 @@ impl Rinex { } list }) - .into_iter(),*/ - .unique() + .into_iter(), ) } else if let Some(_) = self.record.as_meteo() { Box::new( @@ -1813,43 +1812,38 @@ impl Rinex { impl Rinex { /// Returns a Unique Iterator over identified [`Carrier`]s pub fn carrier(&self) -> Box> { - Box::new( - self.observation() - .flat_map(|(_, (_, sv))| { - sv.iter() - .flat_map(|(_, obsvervations)| { - observations - .keys() - .map(|observable| { - if let Ok(carrier) = observable.carrier(sv.constellation) { - Some(carrier) - } else { - None - } - }) - .unique() - }) - }) - ) + Box::new(self.observation().flat_map(|(_, (_, sv))| { + sv.iter().flat_map(|(_, observations)| { + observations + .keys() + .filter_map(|observable| { + if let Ok(carrier) = observable.carrier(sv.constellation) { + Some(carrier) + } else { + None + } + }) + .fold(vec![], |mut list, item| { + if !list.contains(&item) { + list.push(item); + } + list + }) + .into_iter() + }) + })) } - pub fn code(&self) -> Box> { + pub fn code(&self) -> Box> { Box::new( self.observation() .flat_map(|(_, (_, sv))| { - sv.iter() - .flat_map(|(_, obsvervations)| { - observations - .keys() - .map(|observable| { - if let Ok(code) = observable.code() { - Some(code) - } else { - None - } - }) - .unique() - }) + sv.iter().flat_map(|(_, observations)| { + observations + .keys() + .filter_map(|observable| observable.code()) }) + }) + .unique(), ) } /// Returns ([`Epoch`] [`EpochFlag`]) iterator, where each {`EpochFlag`] diff --git a/rinex/src/qc/analysis/obs.rs b/rinex/src/qc/analysis/obs.rs index 908334467..2f096bc1c 100644 --- a/rinex/src/qc/analysis/obs.rs +++ b/rinex/src/qc/analysis/obs.rs @@ -35,7 +35,7 @@ fn report_signals(list: &Vec) -> String { */ fn report_anomalies( cs: &Vec, - power: &Vec, + power: &Vec, other: &Vec<(Epoch, EpochFlag)>, ) -> Box { box_html! { @@ -50,7 +50,7 @@ fn report_anomalies( th { : "None" } - } + } } else { tr { td { @@ -62,7 +62,7 @@ fn report_anomalies( : "Longest" } td { - : power_failures.max_by(|(_, d1), (_, d2)| d1.cmp(d2)).unwrap().to_string() + : power.iter().max_by(|(_, d1), (_, d2)| d1.cmp(d2)).unwrap().to_string() } td { : "Average Power failure" @@ -239,13 +239,20 @@ fn report_ssi_statistics( } } -fn derivative_dt(data: Vec<(Epoch, f64>), window: Duration) -> Vec<(Epoch, f64)> { +fn derivative_dt(data: Vec<(Epoch, f64)>) -> Vec<(Epoch, f64)> { let mut acc = 0.0_f64; - let mut prev_epoch : Option = None + let mut prev_epoch: Option = None; let mut ret: Vec<(Epoch, f64)> = Vec::new(); - for (epoch, value) in data { - let dt = - } + for (epoch, value) in data {} + ret +} + +fn moving_average(data: Vec<(Epoch, f64)>, window: Duration) -> Vec<(Epoch, f64)> { + let mut acc = 0.0_f64; + let mut prev_epoch: Option = None; + let mut ret: Vec<(Epoch, f64)> = Vec::new(); + for (epoch, value) in data {} + ret } #[derive(Debug, Clone)] @@ -279,23 +286,24 @@ pub struct QcObsAnalysis { /// SSi statistical analysis (mean, stddev, skew) ssi_stats: HashMap, /// RX clock drift - clock_drift: Vec, + clock_drift: Vec<(Epoch, f64)>, } impl QcObsAnalysis { pub fn new(rnx: &Rinex, opts: &QcOpts) -> Self { - let mut observables : Vec<_> = rnx.observable().collect(); - let has_doppler = observables.iter() - .fold(|boolean, obs| { - boolean |= obs.is_doppler_observable() - }); - let mut signals : Vec<_> = rnx.signal().collect(); - let mut codes : Vec<_> = rnx.code().collect(); + let mut observables: Vec<_> = rnx.observable().collect(); + let has_doppler = observables + .iter() + .fold(|boolean, obs| boolean |= obs.is_doppler_observable()); + + let mut signals: Vec<_> = rnx.signal().collect(); + let mut codes: Vec<_> = rnx.code().map(|c| c.to_string()).collect(); + let anomalies = rnx.epoch_anomalies(); - let mut total_epochs = self.epoch().count(); + let mut total_epochs = rnx.epoch().count(); let mut epoch_with_obs: Vec = Vec::new(); - + let mut complete_epochs: HashMap = HashMap::new(); let mut min_max_snr = ( (Sv::default(), Epoch::default(), Snr::DbHz54), @@ -315,7 +323,6 @@ impl QcObsAnalysis { if let Some(r) = rnx.record.as_obs() { total_epochs = r.len(); for ((epoch, flag), (clk, svs)) in r { - for (sv, observables) in svs { if !observables.is_empty() { epoch_with_obs.push(*epoch); @@ -448,33 +455,37 @@ impl QcObsAnalysis { observables, has_doppler, cs_anomalies: { - anomalies.filter(|(e, flag) { - if flag == EpochFlag::CycleSlip { - Some(e) - } else { - None - } - }).collect() + anomalies + .filter(|(e, flag)| { + if *flag == EpochFlag::CycleSlip { + Some(e) + } else { + None + } + }) + .collect() }, power_failures: { - anomalies.filter(|(e, flag) { - if flag == EpochFlag::PowerFailure { - Some(e) - } else { - None - } - }) - .collect() + anomalies + .filter(|(e, flag)| { + if *flag == EpochFlag::PowerFailure { + Some(e) + } else { + None + } + }) + .collect() }, other_anomalies: { - anomalies.filter(|(e, flag) { - if flag != EpochFlag::PowerFailure && flag != EpochFlag::CycleSlip { - Some(e) - } else { - None - } - }) - .collect() + anomalies + .filter(|(e, flag)| { + if *flag != EpochFlag::PowerFailure && *flag != EpochFlag::CycleSlip { + Some(e) + } else { + None + } + }) + .collect() }, total_epochs, total_with_obs: epoch_with_obs.unique().len(), @@ -487,11 +498,10 @@ impl QcObsAnalysis { min_max_snr, ssi_stats, clock_drift: { - let rx_clock : Vec<_> = rnx.recvr_clock() - .collect(); + let rx_clock: Vec<_> = rnx.recvr_clock().collect(); let rx_clock_drift: Vec<(Epoch, f64)> = derivative_dt(rx_clock); - moving_average(rx_clock_drift, opts.clock_avg_window) - } + moving_average(rx_clock_drift, opts.clock_drift_window) + }, } } } diff --git a/rinex/src/qc/mod.rs b/rinex/src/qc/mod.rs index 2c094bb4c..854b998e3 100644 --- a/rinex/src/qc/mod.rs +++ b/rinex/src/qc/mod.rs @@ -199,45 +199,41 @@ impl QcReport { } } tbody { - : self.analysis[0].to_inline_html() + : analysis[0].to_inline_html() } - } + } } else { /* * Report all analysis * and emphasize how they were sorted (self.opts.classfication) */ - @ for i in 0..self.analysis.len() { - table(class="table is-bordered") { - @ match self.opts.classification { - QcClassification::GNSS => { - tr { - td { - : format!("{} analysis", self.constellation().nth(i).unwrap()) - } - } - }, - QcClassification::Sv => { - tr { - td { - : format!("{} analysis", self.sv().nth(i).unwrap()) - } - } - }, - QcClassification::Physics => { - tr { - td { - : format!("{} analysis", self.observable().nth(i).unwrap()) - } - } - }, - } - - tbody { - : self.analysis[0].to_inline_html() - } - } - } + //@ for i in 0..self.analysis.len() { + // table(class="table is-bordered") { + // @ match self.opts.classification { + // QcClassification::GNSS => { + // tr { + // td { + // : format!("{} analysis", self.constellation().nth(i).unwrap()) + // } + // } + // }, + // QcClassification::Sv => { + // tr { + // td { + // : format!("{} analysis", self.sv().nth(i).unwrap()) + // } + // } + // }, + // QcClassification::Physics => { + // tr { + // td { + // : format!("{} analysis", self.observable().nth(i).unwrap()) + // } + // } + // }, + // } + // } + //} } }//div=analysis } From a67ad8fc7d2361e345c590e3aa2889c3a3443753 Mon Sep 17 00:00:00 2001 From: "Guillaume W. Bres" Date: Mon, 4 Sep 2023 09:57:34 +0200 Subject: [PATCH 04/19] fix a couple of syntax issues Signed-off-by: Guillaume W. Bres --- rinex/src/lib.rs | 2 +- rinex/src/qc/analysis/obs.rs | 50 +++++++++++++++++++++--------------- rinex/src/qc/opts.rs | 10 +++++++- 3 files changed, 40 insertions(+), 22 deletions(-) diff --git a/rinex/src/lib.rs b/rinex/src/lib.rs index 74ec74ac3..c36c745ec 100644 --- a/rinex/src/lib.rs +++ b/rinex/src/lib.rs @@ -1813,7 +1813,7 @@ impl Rinex { /// Returns a Unique Iterator over identified [`Carrier`]s pub fn carrier(&self) -> Box> { Box::new(self.observation().flat_map(|(_, (_, sv))| { - sv.iter().flat_map(|(_, observations)| { + sv.iter().flat_map(|(sv, observations)| { observations .keys() .filter_map(|observable| { diff --git a/rinex/src/qc/analysis/obs.rs b/rinex/src/qc/analysis/obs.rs index 2f096bc1c..92fd9b18c 100644 --- a/rinex/src/qc/analysis/obs.rs +++ b/rinex/src/qc/analysis/obs.rs @@ -62,7 +62,8 @@ fn report_anomalies( : "Longest" } td { - : power.iter().max_by(|(_, d1), (_, d2)| d1.cmp(d2)).unwrap().to_string() + //: power.iter().max_by(|(_, d1), (_, d2)| d1.cmp(d2)).unwrap().to_string() + : "TODO" } td { : "Average Power failure" @@ -113,13 +114,13 @@ fn report_anomalies( td { : e.to_string() } - @ if event == EpochFlag::AntennaBeingMoved { + @ if *event == EpochFlag::AntennaBeingMoved { : "Antenna being moved" - } else if event == EpochFlag::Kinematic { - : "Kinematic" - } else if event == EpochFlag::NewSiteOccupation { + //} else if *event == EpochFlag::Kinematic { + // : "Kinematic" + } else if *event == EpochFlag::NewSiteOccupation { : "New Site Occupation" - } else if event == EpochFlag::ExternalEvent { + } else if *event == EpochFlag::ExternalEvent { : "External Event" } } @@ -291,10 +292,14 @@ pub struct QcObsAnalysis { impl QcObsAnalysis { pub fn new(rnx: &Rinex, opts: &QcOpts) -> Self { - let mut observables: Vec<_> = rnx.observable().collect(); - let has_doppler = observables - .iter() - .fold(|boolean, obs| boolean |= obs.is_doppler_observable()); + let has_doppler = rnx.observable().fold(|mut is_doppler, observables| { + for obs in observables { + is_doppler |= obs.is_doppler_observable(); + } + is_doppler + }); + + let mut observables: Vec = rnx.observable().map(|obs| obs.to_string()).collect(); let mut signals: Vec<_> = rnx.signal().collect(); let mut codes: Vec<_> = rnx.code().map(|c| c.to_string()).collect(); @@ -325,7 +330,9 @@ impl QcObsAnalysis { for ((epoch, flag), (clk, svs)) in r { for (sv, observables) in svs { if !observables.is_empty() { - epoch_with_obs.push(*epoch); + if !epoch_with_obs.contains(&epoch) { + epoch_with_obs.push(*epoch); + } } for (observable, observation) in observables { @@ -456,8 +463,8 @@ impl QcObsAnalysis { has_doppler, cs_anomalies: { anomalies - .filter(|(e, flag)| { - if *flag == EpochFlag::CycleSlip { + .filter_map(|(e, flag)| { + if flag == EpochFlag::CycleSlip { Some(e) } else { None @@ -467,8 +474,8 @@ impl QcObsAnalysis { }, power_failures: { anomalies - .filter(|(e, flag)| { - if *flag == EpochFlag::PowerFailure { + .filter_map(|(e, flag)| { + if flag == EpochFlag::PowerFailure { Some(e) } else { None @@ -478,9 +485,9 @@ impl QcObsAnalysis { }, other_anomalies: { anomalies - .filter(|(e, flag)| { - if *flag != EpochFlag::PowerFailure && *flag != EpochFlag::CycleSlip { - Some(e) + .filter_map(|(e, flag)| { + if flag != EpochFlag::PowerFailure && flag != EpochFlag::CycleSlip { + Some((e, flag)) } else { None } @@ -488,7 +495,7 @@ impl QcObsAnalysis { .collect() }, total_epochs, - total_with_obs: epoch_with_obs.unique().len(), + total_with_obs: epoch_with_obs.len(), complete_epochs: { let mut ret: Vec<(Carrier, usize)> = complete_epochs.iter().map(|(k, v)| (*k, *v)).collect(); @@ -498,7 +505,10 @@ impl QcObsAnalysis { min_max_snr, ssi_stats, clock_drift: { - let rx_clock: Vec<_> = rnx.recvr_clock().collect(); + let rx_clock: Vec<_> = rnx + .recvr_clock() + .map(|((e, flag), value)| (e, value)) + .collect(); let rx_clock_drift: Vec<(Epoch, f64)> = derivative_dt(rx_clock); moving_average(rx_clock_drift, opts.clock_drift_window) }, diff --git a/rinex/src/qc/opts.rs b/rinex/src/qc/opts.rs index 14faa3e01..afc3a08ab 100644 --- a/rinex/src/qc/opts.rs +++ b/rinex/src/qc/opts.rs @@ -200,7 +200,7 @@ impl Default for QcOpts { min_snr_db: 20.0, // dB elev_mask: None, classification: QcClassification::default(), - clock_window: Duration::from_seconds(10.0 * 60.0), + clock_drift_window: Duration::from_seconds(3600.0), } } } @@ -255,6 +255,14 @@ impl HtmlReport for QcOpts { } } } + tr { + th { + : "Clock Drift Window" + } + td { + : self.clock_drift_window.to_string() + } + } } } } From 69d6589876897ebd296b801e447d873effbbf0f6 Mon Sep 17 00:00:00 2001 From: "Guillaume W. Bres" Date: Mon, 4 Sep 2023 10:15:52 +0200 Subject: [PATCH 05/19] fix a couple of syntax issues Signed-off-by: Guillaume W. Bres --- rinex/src/qc/analysis/obs.rs | 11 +++-------- rinex/src/tests/obs.rs | 2 +- 2 files changed, 4 insertions(+), 9 deletions(-) diff --git a/rinex/src/qc/analysis/obs.rs b/rinex/src/qc/analysis/obs.rs index 92fd9b18c..4e87442ea 100644 --- a/rinex/src/qc/analysis/obs.rs +++ b/rinex/src/qc/analysis/obs.rs @@ -292,16 +292,11 @@ pub struct QcObsAnalysis { impl QcObsAnalysis { pub fn new(rnx: &Rinex, opts: &QcOpts) -> Self { - let has_doppler = rnx.observable().fold(|mut is_doppler, observables| { - for obs in observables { - is_doppler |= obs.is_doppler_observable(); - } - is_doppler - }); + let doppler_obs = rnx.observable().filter(|obs| obs.is_doppler_observable()); let mut observables: Vec = rnx.observable().map(|obs| obs.to_string()).collect(); - let mut signals: Vec<_> = rnx.signal().collect(); + let mut signals: Vec<_> = rnx.carrier().collect(); let mut codes: Vec<_> = rnx.code().map(|c| c.to_string()).collect(); let anomalies = rnx.epoch_anomalies(); @@ -460,7 +455,7 @@ impl QcObsAnalysis { codes, signals, observables, - has_doppler, + has_doppler: doppler_obs.count() > 0, cs_anomalies: { anomalies .filter_map(|(e, flag)| { diff --git a/rinex/src/tests/obs.rs b/rinex/src/tests/obs.rs index 4ab469df9..6af66d117 100644 --- a/rinex/src/tests/obs.rs +++ b/rinex/src/tests/obs.rs @@ -1184,7 +1184,7 @@ mod test { let mut irnss_sv: Vec = rnx .sv() .filter_map(|sv| { - if (sv.constellation == Constellation::IRNSS) { + if sv.constellation == Constellation::IRNSS { Some(sv) } else { None From b81468b967d1aac4795b354d69fa63d6e37f41a2 Mon Sep 17 00:00:00 2001 From: "Guillaume W. Bres" Date: Mon, 4 Sep 2023 11:53:10 +0200 Subject: [PATCH 06/19] working on averaging and derivation methods Signed-off-by: Guillaume W. Bres --- rinex/src/algorithm/averaging.rs | 67 +++++++--------- rinex/src/algorithm/derivative.rs | 43 +++++++++++ rinex/src/algorithm/mod.rs | 4 + rinex/src/lib.rs | 4 +- rinex/src/qc/analysis/obs.rs | 122 ++++++++++++++++-------------- 5 files changed, 140 insertions(+), 100 deletions(-) create mode 100644 rinex/src/algorithm/derivative.rs diff --git a/rinex/src/algorithm/averaging.rs b/rinex/src/algorithm/averaging.rs index eb7c163cb..e59d0ece0 100644 --- a/rinex/src/algorithm/averaging.rs +++ b/rinex/src/algorithm/averaging.rs @@ -1,58 +1,43 @@ -use crate::prelude::*; +use hifitime::{Duration, Epoch}; -pub enum AverageType { - MovingAverage, - Exponential, - Cummulative, +fn moving_average( + data: Vec<(Epoch, T)>, + window: Duration, +) -> Vec<(Epoch, T)> { + let mut acc = T::default(); + let mut prev_epoch: Option = None; + let mut ret: Vec<(Epoch, T)> = Vec::new(); + for (epoch, value) in data {} + ret } -pub struct Averager { - buffer: Vec, - next_epoch: Option, - window: Duration, - avgtype: AverageType, +#[derive(Debug, Clone, Copy)] +pub enum Averager { + MovingAverage(Duration), } -impl Averager { - pub fn new(window: Duration) -> Self { - Self { - buffer: Vec::new(), - next_epoch: None, - window, - } +impl Default for Averager { + fn default() -> Self { + Self::MovingAverage(Duration::from_seconds(600.0_f64)) } - pub fn average(&mut self, data: (f64, Epoch)) -> Option { - match self.avgtype { - AverageType::MovingAverage => self.moving_average(data), - AverageType::Exponential => self.moving_average(data), - AverageType::Cummulative => self.moving_average(data), - } +} + +impl Averager { + pub fn mov(window: Duration) -> Self { + Self::MovingAverage(window) } - pub fn moving_average(&mut self, data: (f64, Epoch)) -> Option { - self.buffer.push(data.0); - if self.next_epoch.is_none() { - self.next_epoch = Some(data.1 + self.window); - } - if let Some(next_epoch) = self.next_epoch { - if data.1 >= next_epoch { - self.next_epoch = Some(data.1 + self.window); - let mut avg = 0.0_f64; - for b in &self.buffer { - avg += b; - } - let ret = avg / self.buffer.len() as f64; - self.buffer.clear(); - return Some(ret); - } + pub fn eval(&self, input: Vec<(Epoch, T)>) -> Vec<(Epoch, T)> { + match self { + Self::MovingAverage(dt) => moving_average(input, *dt), } - None } } #[cfg(test)] mod test { + use super::*; #[test] fn test_moving_average() { - + let mov = Averager::mov(Duration::from_seconds(10.0_f64)); } } diff --git a/rinex/src/algorithm/derivative.rs b/rinex/src/algorithm/derivative.rs new file mode 100644 index 000000000..453b850da --- /dev/null +++ b/rinex/src/algorithm/derivative.rs @@ -0,0 +1,43 @@ +use hifitime::Epoch; + +pub struct Derivative { + order: usize, +} + +/* + * Derivative of an number of array sorted by chronological Epoch + */ +pub(crate) fn derivative + std::marker::Copy>( + input: Vec<(Epoch, T)>, + buf: &mut Vec<(Epoch, T)>, +) { + let mut prev: Option<(Epoch, T)> = None; + for (e, value) in input { + if let Some((prev_e, prev_v)) = prev { + let dt = e - prev_e; + let dy = value - prev_v; + buf.push((e, dy)); + } + prev = Some((e, value)); + } +} + +/* + * Derivative^2 of an number of array sorted by chronological Epoch + */ +impl Derivative { + pub fn new(order: usize) -> Self { + Self { order } + } + pub fn eval + std::marker::Copy>( + input: Vec<(Epoch, T)>, + ) -> Vec<(Epoch, T)> { + let mut buf: Vec<(Epoch, T)> = Vec::with_capacity(input.len()); + derivative(input, &mut buf); + //for i in 1..self.order { + // derivative(&ret, &mut ret); + //} + //ret + buf + } +} diff --git a/rinex/src/algorithm/mod.rs b/rinex/src/algorithm/mod.rs index 1b77a046e..96e068401 100644 --- a/rinex/src/algorithm/mod.rs +++ b/rinex/src/algorithm/mod.rs @@ -1,3 +1,5 @@ +mod averaging; +mod derivative; mod filters; mod target; @@ -7,3 +9,5 @@ pub use filters::{ Decimate, DecimationFilter, DecimationType, Filter, InterpFilter, InterpMethod, Interpolate, Mask, MaskFilter, MaskOperand, Preprocessing, Smooth, SmoothingFilter, SmoothingType, }; + +pub use averaging::Averager; diff --git a/rinex/src/lib.rs b/rinex/src/lib.rs index c36c745ec..4de29c234 100644 --- a/rinex/src/lib.rs +++ b/rinex/src/lib.rs @@ -1811,7 +1811,7 @@ impl Rinex { #[cfg_attr(docrs, doc(cfg(feature = "obs")))] impl Rinex { /// Returns a Unique Iterator over identified [`Carrier`]s - pub fn carrier(&self) -> Box> { + pub fn carrier(&self) -> Box + '_> { Box::new(self.observation().flat_map(|(_, (_, sv))| { sv.iter().flat_map(|(sv, observations)| { observations @@ -1833,7 +1833,7 @@ impl Rinex { }) })) } - pub fn code(&self) -> Box> { + pub fn code(&self) -> Box + '_> { Box::new( self.observation() .flat_map(|(_, (_, sv))| { diff --git a/rinex/src/qc/analysis/obs.rs b/rinex/src/qc/analysis/obs.rs index 4e87442ea..abf7fd3d4 100644 --- a/rinex/src/qc/analysis/obs.rs +++ b/rinex/src/qc/analysis/obs.rs @@ -33,11 +33,11 @@ fn report_signals(list: &Vec) -> String { /* * Epoch anomalies formatter */ -fn report_anomalies( - cs: &Vec, - power: &Vec, - other: &Vec<(Epoch, EpochFlag)>, -) -> Box { +fn report_anomalies<'a>( + cs: &'a Vec, + power: &'a Vec, + other: &'a Vec<(Epoch, EpochFlag)>, +) -> Box { box_html! { table(class="table is-bordered") { tr { @@ -114,15 +114,28 @@ fn report_anomalies( td { : e.to_string() } - @ if *event == EpochFlag::AntennaBeingMoved { - : "Antenna being moved" - //} else if *event == EpochFlag::Kinematic { - // : "Kinematic" - } else if *event == EpochFlag::NewSiteOccupation { - : "New Site Occupation" - } else if *event == EpochFlag::ExternalEvent { - : "External Event" - } + //@ match event { + // EpochFlag::AntennaBeingMoved => { + // td { + // : "Antenna being moved" + // } + // }, + // _ => {}, + //} + // } + //td { + // @ match event { + // EpochFlag::NewSiteOccupation => { + // : "New Site Occupation" + // } + // EpochFlag::ExternalEvent => { + // : "External Event" + // } + // _ => { + // : "Other" + // } + // } + //} } } } @@ -248,14 +261,6 @@ fn derivative_dt(data: Vec<(Epoch, f64)>) -> Vec<(Epoch, f64)> { ret } -fn moving_average(data: Vec<(Epoch, f64)>, window: Duration) -> Vec<(Epoch, f64)> { - let mut acc = 0.0_f64; - let mut prev_epoch: Option = None; - let mut ret: Vec<(Epoch, f64)> = Vec::new(); - for (epoch, value) in data {} - ret -} - #[derive(Debug, Clone)] /// OBS RINEX specific QC analysis. /// Full OBS RINEX analysis requires both the "obs" and "processing" features. @@ -286,6 +291,7 @@ pub struct QcObsAnalysis { #[cfg(feature = "obs")] /// SSi statistical analysis (mean, stddev, skew) ssi_stats: HashMap, + #[cfg(feature = "obs")] /// RX clock drift clock_drift: Vec<(Epoch, f64)>, } @@ -299,7 +305,38 @@ impl QcObsAnalysis { let mut signals: Vec<_> = rnx.carrier().collect(); let mut codes: Vec<_> = rnx.code().map(|c| c.to_string()).collect(); - let anomalies = rnx.epoch_anomalies(); + let cs_anomalies: Vec<_> = rnx + .epoch_anomalies() + .filter_map(|(e, flag)| { + if flag == EpochFlag::CycleSlip { + Some(e) + } else { + None + } + }) + .collect(); + + let power_failures: Vec<_> = rnx + .epoch_anomalies() + .filter_map(|(e, flag)| { + if flag == EpochFlag::PowerFailure { + Some(e) + } else { + None + } + }) + .collect(); + + let other_anomalies: Vec<_> = rnx + .epoch_anomalies() + .filter_map(|(e, flag)| { + if flag != EpochFlag::PowerFailure && flag != EpochFlag::CycleSlip { + Some((e, flag)) + } else { + None + } + }) + .collect(); let mut total_epochs = rnx.epoch().count(); let mut epoch_with_obs: Vec = Vec::new(); @@ -456,39 +493,9 @@ impl QcObsAnalysis { signals, observables, has_doppler: doppler_obs.count() > 0, - cs_anomalies: { - anomalies - .filter_map(|(e, flag)| { - if flag == EpochFlag::CycleSlip { - Some(e) - } else { - None - } - }) - .collect() - }, - power_failures: { - anomalies - .filter_map(|(e, flag)| { - if flag == EpochFlag::PowerFailure { - Some(e) - } else { - None - } - }) - .collect() - }, - other_anomalies: { - anomalies - .filter_map(|(e, flag)| { - if flag != EpochFlag::PowerFailure && flag != EpochFlag::CycleSlip { - Some((e, flag)) - } else { - None - } - }) - .collect() - }, + cs_anomalies, + power_failures, + other_anomalies, total_epochs, total_with_obs: epoch_with_obs.len(), complete_epochs: { @@ -505,7 +512,8 @@ impl QcObsAnalysis { .map(|((e, flag), value)| (e, value)) .collect(); let rx_clock_drift: Vec<(Epoch, f64)> = derivative_dt(rx_clock); - moving_average(rx_clock_drift, opts.clock_drift_window) + let mov = Averager::mov(opts.clock_drift_window); + mov.eval(rx_clock_drift) }, } } From fbb6f50cfe19044694e979c390e681223dfadaa7 Mon Sep 17 00:00:00 2001 From: "Guillaume W. Bres" Date: Mon, 4 Sep 2023 18:17:13 +0200 Subject: [PATCH 07/19] default clock window Signed-off-by: Guillaume W. Bres --- rinex/src/qc/opts.rs | 7 ++++++- 1 file changed, 6 insertions(+), 1 deletion(-) diff --git a/rinex/src/qc/opts.rs b/rinex/src/qc/opts.rs index afc3a08ab..f9792ddc6 100644 --- a/rinex/src/qc/opts.rs +++ b/rinex/src/qc/opts.rs @@ -162,6 +162,7 @@ pub struct QcOpts { /// Manually defined Ground position (ECEF) pub ground_position: Option, /// Window duration to be used, during RX clock drift analysis + #[serde(default = "default_drift_window")] pub clock_drift_window: Duration, } @@ -192,6 +193,10 @@ impl QcOpts { } } +fn default_drift_window() -> Duration { + Duration::from_seconds(3600.0) +} + impl Default for QcOpts { fn default() -> Self { Self { @@ -200,7 +205,7 @@ impl Default for QcOpts { min_snr_db: 20.0, // dB elev_mask: None, classification: QcClassification::default(), - clock_drift_window: Duration::from_seconds(3600.0), + clock_drift_window: default_drift_window(), } } } From e84992c556527ceec05ebe9b5bff3df8c7ac658d Mon Sep 17 00:00:00 2001 From: "Guillaume W. Bres" Date: Mon, 4 Sep 2023 19:03:24 +0200 Subject: [PATCH 08/19] improving qc Signed-off-by: Guillaume W. Bres --- rinex/src/qc/analysis/obs.rs | 3 +- rinex/src/qc/analysis/sampling.rs | 38 +++++++------- rinex/src/qc/mod.rs | 82 +++++++++++++------------------ rinex/src/qc/opts.rs | 2 +- 4 files changed, 54 insertions(+), 71 deletions(-) diff --git a/rinex/src/qc/analysis/obs.rs b/rinex/src/qc/analysis/obs.rs index abf7fd3d4..6fc8aa0eb 100644 --- a/rinex/src/qc/analysis/obs.rs +++ b/rinex/src/qc/analysis/obs.rs @@ -1,4 +1,5 @@ use crate::{carrier, observation::Snr, prelude::*, Carrier}; +use itertools::Itertools; use super::{pretty_array, QcOpts}; use std::collections::HashMap; @@ -302,7 +303,7 @@ impl QcObsAnalysis { let mut observables: Vec = rnx.observable().map(|obs| obs.to_string()).collect(); - let mut signals: Vec<_> = rnx.carrier().collect(); + let mut signals: Vec<_> = rnx.carrier().unique().collect(); let mut codes: Vec<_> = rnx.code().map(|c| c.to_string()).collect(); let cs_anomalies: Vec<_> = rnx diff --git a/rinex/src/qc/analysis/sampling.rs b/rinex/src/qc/analysis/sampling.rs index 9862341d9..4c175113e 100644 --- a/rinex/src/qc/analysis/sampling.rs +++ b/rinex/src/qc/analysis/sampling.rs @@ -90,28 +90,28 @@ impl HtmlReport for QcSamplingAnalysis { } tr { th { - : "Sample rate" + : "Sample rate (Header)" } - th { - : "Dominan Sample rate" - } - } - tr { @ if let Some(rate) = self.sample_rate { td { : format!("{} ({:.3} Hz)", rate, 1.0 / rate.to_unit(Unit::Second)) } } else { - td { - : "Not Specified" + th { + : "Unspecified" } } + } + tr { + th { + : "Dominant Sample rate" + } @ if let Some(rate) = self.dominant_sample_rate { td { : format!("{} ({:.3} Hz)", rate, 1.0 / rate.to_unit(Unit::Second)) } } else { - td { + th { : "Undetermined" } } @@ -120,21 +120,19 @@ impl HtmlReport for QcSamplingAnalysis { th { : "Gap analysis" } - } - @ if self.gaps.is_empty() { - tr { + + @ if self.gaps.is_empty() { th { - : "No Data Gaps detected" + : "No gaps detected" } - } - } else { - @ for (epoch, dt) in &self.gaps { + } else { tr { td { - : format!("Start : {}", epoch) - } - td { - : format!("Duration: {}", dt) + @ for (epoch, dt) in &self.gaps { + p { + : format!("Start : {}, Duration: {}", epoch, dt) + } + } } } } diff --git a/rinex/src/qc/mod.rs b/rinex/src/qc/mod.rs index 854b998e3..3d2b99b50 100644 --- a/rinex/src/qc/mod.rs +++ b/rinex/src/qc/mod.rs @@ -153,30 +153,36 @@ impl QcReport { } }//div=header div(id="context") { - h3(class="title") { - : "Context" - } table(class="table is-bordered") { + thead { + th { + : "Context" + } + } tbody { : context.to_inline_html() } } }//div=context div(id="parameters") { - h3(class="title") { - : "Parameters" - } table(class="table is-bordered") { + thead { + th { + : "Parameters" + } + } tbody { : opts.to_inline_html() } } } //div=parameters div(id="header") { - h3(class="title") { - : "File Header" - } table(class="table is-bordered") { + thead { + th { + : "File Header" + } + } tbody { : context.primary_data().header.to_inline_html() } @@ -186,54 +192,32 @@ impl QcReport { * Report all analysis that were performed */ div(id="analysis") { - @ if analysis.len() == 1 { - /* - * Single analysis case - */ + /* + * Report all analysis + * and emphasize how they were sorted (self.opts.classfication) + */ + @ for i in 0..analysis.len() { table(class="table is-bordered") { thead { - tr { - td { - : "Multi constellation analysis" + @ if opts.classification == QcClassification::GNSS { + th { + : format!("{} analysis", context.primary_data().constellation().nth(i).unwrap()) + } + } else if opts.classification == QcClassification::Sv { + th { + : format!("{} analysis", context.primary_data().sv().nth(i).unwrap()) + } + + } else if opts.classification == QcClassification::Physics { + th { + : format!("{} analysis", context.primary_data().nth(i).unwrap()) } } } tbody { - : analysis[0].to_inline_html() + : analysis[i].to_inline_html() } } - } else { - /* - * Report all analysis - * and emphasize how they were sorted (self.opts.classfication) - */ - //@ for i in 0..self.analysis.len() { - // table(class="table is-bordered") { - // @ match self.opts.classification { - // QcClassification::GNSS => { - // tr { - // td { - // : format!("{} analysis", self.constellation().nth(i).unwrap()) - // } - // } - // }, - // QcClassification::Sv => { - // tr { - // td { - // : format!("{} analysis", self.sv().nth(i).unwrap()) - // } - // } - // }, - // QcClassification::Physics => { - // tr { - // td { - // : format!("{} analysis", self.observable().nth(i).unwrap()) - // } - // } - // }, - // } - // } - //} } }//div=analysis } diff --git a/rinex/src/qc/opts.rs b/rinex/src/qc/opts.rs index f9792ddc6..bf56c300d 100644 --- a/rinex/src/qc/opts.rs +++ b/rinex/src/qc/opts.rs @@ -119,7 +119,7 @@ impl Default for ProcessingOpts { } /// Qc Report classification method -#[derive(Default, Debug, Clone)] +#[derive(Default, Debug, Clone, PartialEq)] #[cfg_attr(feature = "serde", derive(Serialize, Deserialize))] /// Classify the QC report by desired data set pub enum QcClassification { From 169c8a243704dfe91679807c4a9feb2a22d3e630 Mon Sep 17 00:00:00 2001 From: "Guillaume W. Bres" Date: Tue, 5 Sep 2023 08:31:45 +0200 Subject: [PATCH 09/19] improving qc Signed-off-by: Guillaume W. Bres --- rinex/src/hardware.rs | 4 ++-- rinex/src/header.rs | 30 ++++++++++++++++++++++++------ rinex/src/qc/analysis/mod.rs | 32 +++++++++++++++++++------------- rinex/src/qc/mod.rs | 15 +++++---------- 4 files changed, 50 insertions(+), 31 deletions(-) diff --git a/rinex/src/hardware.rs b/rinex/src/hardware.rs index 78ed15509..cc57709ab 100644 --- a/rinex/src/hardware.rs +++ b/rinex/src/hardware.rs @@ -106,7 +106,7 @@ impl HtmlReport for Antenna { table(class="table is-bordered") { tr { th { - : "Antenna model" + : "Model" } th { : "SN#" @@ -167,7 +167,7 @@ impl HtmlReport for Rcvr { } fn to_inline_html(&self) -> Box { box_html! { - table(class="table is-bordered") { + table(class="table is-bordered; style=\"margin-bottom: 20px\"") { tr { th { : "Model" diff --git a/rinex/src/header.rs b/rinex/src/header.rs index afc574841..7a046b3a6 100644 --- a/rinex/src/header.rs +++ b/rinex/src/header.rs @@ -1833,14 +1833,32 @@ impl HtmlReport for Header { } fn to_inline_html(&self) -> Box { box_html! { - @ if let Some(antenna) = &self.rcvr_antenna { - div(id="antenna") { - : antenna.to_inline_html() + tr { + th { + : "Antenna" + } + @ if let Some(antenna) = &self.rcvr_antenna { + td { + : antenna.to_inline_html() + } + } else { + td { + : "No information" + } } } - @ if let Some(rcvr) = &self.rcvr { - div(id="receiver") { - : rcvr.to_inline_html() + tr { + th { + : "Receiver" + } + @ if let Some(rcvr) = &self.rcvr { + td { + : rcvr.to_inline_html() + } + } else { + td { + : "No information" + } } } } diff --git a/rinex/src/qc/analysis/mod.rs b/rinex/src/qc/analysis/mod.rs index 504629ad8..8d4859a36 100644 --- a/rinex/src/qc/analysis/mod.rs +++ b/rinex/src/qc/analysis/mod.rs @@ -66,32 +66,38 @@ impl HtmlReport for QcAnalysis { } fn to_inline_html(&self) -> Box { box_html! { - div(id="analysis") { + div(id="analysis; style=\"display: flex; flex-direction: column; gap: 30px\"") { div(id="sampling") { - h3(class="title") { - : "Sampling" - } - table(class="table is-bordered") { + table(class="table is-bordered; style=\"margin-bottom: 30px\"") { + thead { + th { + : "Sampling" + } + } tbody { : self.sampling.to_inline_html() } } } div(id="sv") { - h3(class="title") { - : "Sv" - } - table(class="table is-bordered") { + table(class="table is-bordered; style=\"margin-bottom: 30px\"") { + thead { + th { + : "Sv" + } + } tbody { : self.sv.to_inline_html() } } } div(id="observations") { - h3(class="title") { - : "Observations" - } - table(class="table is-bordered") { + table(class="table is-bordered; style=\"margin-bottom: 30px\"") { + thead { + th { + : "Observations" + } + } tbody { : self.observ.to_inline_html() } diff --git a/rinex/src/qc/mod.rs b/rinex/src/qc/mod.rs index 3d2b99b50..8cff1a53c 100644 --- a/rinex/src/qc/mod.rs +++ b/rinex/src/qc/mod.rs @@ -153,7 +153,7 @@ impl QcReport { } }//div=header div(id="context") { - table(class="table is-bordered") { + table(class="table is-bordered; style=\"margin-bottom: 20px\"") { thead { th { : "Context" @@ -165,7 +165,7 @@ impl QcReport { } }//div=context div(id="parameters") { - table(class="table is-bordered") { + table(class="table is-bordered; style=\"margin-bottom: 20px\"") { thead { th { : "Parameters" @@ -177,12 +177,7 @@ impl QcReport { } } //div=parameters div(id="header") { - table(class="table is-bordered") { - thead { - th { - : "File Header" - } - } + table(class="table is-bordered; style=\"margin-bottom: 20px\"") { tbody { : context.primary_data().header.to_inline_html() } @@ -197,7 +192,7 @@ impl QcReport { * and emphasize how they were sorted (self.opts.classfication) */ @ for i in 0..analysis.len() { - table(class="table is-bordered") { + table(class="table is-bordered; style=\"margin-bottom: 20px\"") { thead { @ if opts.classification == QcClassification::GNSS { th { @@ -210,7 +205,7 @@ impl QcReport { } else if opts.classification == QcClassification::Physics { th { - : format!("{} analysis", context.primary_data().nth(i).unwrap()) + : format!("{} analysis", context.primary_data().observable().nth(i).unwrap()) } } } From 402044e9c430e23d3cf6697aa2d66034227006ab Mon Sep 17 00:00:00 2001 From: "Guillaume W. Bres" Date: Tue, 5 Sep 2023 09:08:26 +0200 Subject: [PATCH 10/19] improving qc Signed-off-by: Guillaume W. Bres --- rinex/src/qc/analysis/mod.rs | 7 +- rinex/src/qc/analysis/obs.rs | 287 +++++++++++++++++++---------------- rinex/src/qc/analysis/sv.rs | 2 +- rinex/src/qc/mod.rs | 31 +++- 4 files changed, 189 insertions(+), 138 deletions(-) diff --git a/rinex/src/qc/analysis/mod.rs b/rinex/src/qc/analysis/mod.rs index 8d4859a36..2c330e917 100644 --- a/rinex/src/qc/analysis/mod.rs +++ b/rinex/src/qc/analysis/mod.rs @@ -1,4 +1,4 @@ -use super::{pretty_array, HtmlReport, QcOpts}; +use super::{pretty_array, table_lengthy_td, HtmlReport, QcOpts}; use crate::prelude::*; use horrorshow::{helper::doctype, RenderBox}; @@ -93,11 +93,6 @@ impl HtmlReport for QcAnalysis { } div(id="observations") { table(class="table is-bordered; style=\"margin-bottom: 30px\"") { - thead { - th { - : "Observations" - } - } tbody { : self.observ.to_inline_html() } diff --git a/rinex/src/qc/analysis/obs.rs b/rinex/src/qc/analysis/obs.rs index 6fc8aa0eb..5357986d9 100644 --- a/rinex/src/qc/analysis/obs.rs +++ b/rinex/src/qc/analysis/obs.rs @@ -41,102 +41,100 @@ fn report_anomalies<'a>( ) -> Box { box_html! { table(class="table is-bordered") { - tr { - td { - : "Power Failures" + thead { + th { + : "Anomalies" } } - @ if power.is_empty() { + tbody { tr { th { - : "None" + : "Power Failures" } - } - } else { - tr { - td { - : format!("{:?}", power) + @ if power.is_empty() { + td { + : "None" + } + } else { + td { + : format!("{:?}", power) + } + tr { + th { + : "Longest" + } + td { + //: power.iter().max_by(|(_, d1), (_, d2)| d1.cmp(d2)).unwrap().to_string() + : "TODO" + } + td { + : "Average Power failure" + } + td { + : "TODO" + } + } } } tr { - td { - : "Longest" - } - td { - //: power.iter().max_by(|(_, d1), (_, d2)| d1.cmp(d2)).unwrap().to_string() - : "TODO" - } - td { - : "Average Power failure" - } - td { - : "TODO" - } - } - } - tr { - td { - : "Cycle slip(s)" - } - } - tr { - @ if cs.is_empty() { th { - : "None" + : "Cycle slip(s)" } - } else { - td { - : format!("{:?}", cs) + @ if cs.is_empty() { + td { + : "None" + } + } else { + td { + : format!("{:?}", cs) + } } } - } - tr { - @ if other.is_empty() { - th { - : "Other anomalies" - } - th { - : "None" - } - } else { + tr { th { : "Other anomalies" } - td { - : "Epoch" - } - td { - : "Event" - } - @ for (e, event) in other { + @ if other.is_empty() { + td { + : "None" + } + } else { td { - : "" + : "Epoch" } td { - : e.to_string() + : "Event" + } + @ for (e, event) in other { + td { + : "" + } + td { + : e.to_string() + } + //@ match event { + // EpochFlag::AntennaBeingMoved => { + // td { + // : "Antenna being moved" + // } + // }, + // _ => {}, + //} + // } + //td { + // @ match event { + // EpochFlag::NewSiteOccupation => { + // : "New Site Occupation" + // } + // EpochFlag::ExternalEvent => { + // : "External Event" + // } + // _ => { + // : "Other" + // } + // } + //} } - //@ match event { - // EpochFlag::AntennaBeingMoved => { - // td { - // : "Antenna being moved" - // } - // }, - // _ => {}, - //} - // } - //td { - // @ match event { - // EpochFlag::NewSiteOccupation => { - // : "New Site Occupation" - // } - // EpochFlag::ExternalEvent => { - // : "External Event" - // } - // _ => { - // : "Other" - // } - // } - //} } } } @@ -157,10 +155,15 @@ fn report_epoch_completion( ) -> Box { box_html! { table(class="table is-bordered") { + thead { + th { + : "Epochs" + } + } tbody { tr { th { - : "Total Epochs" + : "Total#" } td { : total.to_string() @@ -198,6 +201,58 @@ fn report_epoch_completion( } } +/* + * SNR analysis report + */ +fn report_snr_analysis( + min: (Sv, Epoch, Snr), + max: (Sv, Epoch, Snr), + ssi_stats: &HashMap, +) -> Box { + box_html! { + tr { + th { + : "Minimum" + } + td { + p { + : min.0.to_string() + } + p { + : format!("{:e}", min.2) + } + p { + : format!("@{}", min.1) + } + } + } + tr { + th { + : "Best" + } + td { + p { + : max.0.to_string() + } + p { + : format!("{:e}", max.2) + } + p { + : format!("@{}", max.1) + } + } + } + tr { + th { + : "SSI Statistics" + } + td { + : report_ssi_statistics(ssi_stats) + } + } + } +} + /* * Reports statistical analysis results for SSx observations */ @@ -287,8 +342,11 @@ pub struct QcObsAnalysis { /// Complete epochs, with respect to given signal complete_epochs: Vec<(Carrier, usize)>, #[cfg(feature = "obs")] - /// Min. Max. SNR (sv @ epoch) - min_max_snr: ((Sv, Epoch, Snr), (Sv, Epoch, Snr)), + /// Min. SNR (sv @ epoch) + min_snr: (Sv, Epoch, Snr), + #[cfg(feature = "obs")] + /// Max. SNR (sv @ epoch) + max_snr: (Sv, Epoch, Snr), #[cfg(feature = "obs")] /// SSi statistical analysis (mean, stddev, skew) ssi_stats: HashMap, @@ -343,10 +401,9 @@ impl QcObsAnalysis { let mut epoch_with_obs: Vec = Vec::new(); let mut complete_epochs: HashMap = HashMap::new(); - let mut min_max_snr = ( - (Sv::default(), Epoch::default(), Snr::DbHz54), - (Sv::default(), Epoch::default(), Snr::DbHz0), - ); + + let mut min_snr = (Sv::default(), Epoch::default(), Snr::DbHz54); + let mut max_snr = (Sv::default(), Epoch::default(), Snr::DbHz0); let mut ssi_stats: HashMap = HashMap::new(); @@ -370,15 +427,15 @@ impl QcObsAnalysis { for (observable, observation) in observables { if let Some(snr) = observation.snr { - if snr < min_max_snr.0 .2 { - min_max_snr.0 .0 = *sv; - min_max_snr.0 .1 = *epoch; - min_max_snr.0 .2 = snr; + if snr < min_snr.2 { + min_snr.0 = *sv; + min_snr.1 = *epoch; + min_snr.2 = snr; } - if snr > min_max_snr.1 .2 { - min_max_snr.1 .0 = *sv; - min_max_snr.1 .1 = *epoch; - min_max_snr.1 .2 = snr; + if snr > max_snr.2 { + max_snr.0 = *sv; + max_snr.1 = *epoch; + max_snr.2 = snr; } } } @@ -505,7 +562,8 @@ impl QcObsAnalysis { ret.sort(); ret }, - min_max_snr, + min_snr, + max_snr, ssi_stats, clock_drift: { let rx_clock: Vec<_> = rnx @@ -577,49 +635,18 @@ impl HtmlReport for QcObsAnalysis { : report_epoch_completion(self.total_epochs, self.total_with_obs, &self.complete_epochs) } } - table(class="table is-bordered") { - thead { - tr { - th { - : "Worst SNR" - } + tr { + table { + thead { th { - : "Best SNR" + : "SNR" } } - } - tbody { - tr { - td { - p { - : self.min_max_snr.0.0.to_string() - } - b { - : format!("{:e}", self.min_max_snr.0.2) - } - p { - : format!("@{}", self.min_max_snr.0.1) - } - } - td { - p { - : self.min_max_snr.1.0.to_string() - } - b { - : format!("{:e}", self.min_max_snr.1.2) - } - p { - : format!("@{}", self.min_max_snr.1.1) - } - } + tbody { + : report_snr_analysis(self.min_snr, self.max_snr, &self.ssi_stats) } } } - tr { - table { - : report_ssi_statistics(&self.ssi_stats) - } - } } } } diff --git a/rinex/src/qc/analysis/sv.rs b/rinex/src/qc/analysis/sv.rs index 01996cbc5..a482656ea 100644 --- a/rinex/src/qc/analysis/sv.rs +++ b/rinex/src/qc/analysis/sv.rs @@ -1,4 +1,4 @@ -use super::{pretty_array, QcOpts}; +use super::{pretty_array, table_lengthy_td, QcOpts}; use crate::prelude::*; #[derive(Debug, Clone)] diff --git a/rinex/src/qc/mod.rs b/rinex/src/qc/mod.rs index 8cff1a53c..5d143c99e 100644 --- a/rinex/src/qc/mod.rs +++ b/rinex/src/qc/mod.rs @@ -13,6 +13,35 @@ use analysis::QcAnalysis; #[cfg(feature = "processing")] use crate::preprocessing::*; +/* + * Methods used when reporting lenghty vectors or data subsets in a table. + * Makes tables cleaner and nicer by wrapping string content, into several paragraphs. + */ +pub(crate) fn table_lengthy_td( + list: &Vec, + max_items: usize, +) -> Box { + let mut content = String::with_capacity(64 * max_items); + let mut paragraphs: Vec = Vec::new(); + + for i in 0..list.len() { + content.push_str(&format!("{}, ", list[i])); + if i.rem_euclid(max_items) == 0 { + paragraphs.push(content.clone()); + content.clear(); + } else if i == list.len() - 1 { + paragraphs.push(content.clone()); + } + } + box_html! { + @ for paragraph in paragraphs { + p { + : paragraph.to_string() + } + } + } +} + /* * Array (CSV) pretty formatter */ @@ -139,7 +168,7 @@ impl QcReport { h2(class="title") { : "RINEX Quality Check summary" } - table(class="table is-bordered") { + table(class="table is-bordered; style=\"margin-bottom: 20px\"") { tbody { tr { th { From 82330a84c6fd22e7c1977443234570b64206e124 Mon Sep 17 00:00:00 2001 From: "Guillaume W. Bres" Date: Tue, 5 Sep 2023 09:16:35 +0200 Subject: [PATCH 11/19] cli configuration Signed-off-by: Guillaume W. Bres --- rinex-cli/config/{basic.json => basic_gnss.json} | 0 rinex-cli/config/basic_manual_gap.json | 7 ------- rinex-cli/config/sv_manual_gap.json | 4 ++++ rinex-cli/src/cli.rs | 1 + 4 files changed, 5 insertions(+), 7 deletions(-) rename rinex-cli/config/{basic.json => basic_gnss.json} (100%) delete mode 100644 rinex-cli/config/basic_manual_gap.json create mode 100644 rinex-cli/config/sv_manual_gap.json diff --git a/rinex-cli/config/basic.json b/rinex-cli/config/basic_gnss.json similarity index 100% rename from rinex-cli/config/basic.json rename to rinex-cli/config/basic_gnss.json diff --git a/rinex-cli/config/basic_manual_gap.json b/rinex-cli/config/basic_manual_gap.json deleted file mode 100644 index 82618e5d0..000000000 --- a/rinex-cli/config/basic_manual_gap.json +++ /dev/null @@ -1,7 +0,0 @@ -{ - "reporting": "PerConstellation", - "statistical": { - "window": "10%" - }, - "manual_gap": "1O seconds" -} diff --git a/rinex-cli/config/sv_manual_gap.json b/rinex-cli/config/sv_manual_gap.json new file mode 100644 index 000000000..764e59183 --- /dev/null +++ b/rinex-cli/config/sv_manual_gap.json @@ -0,0 +1,4 @@ +{ + "classification": "Sv", + "min_snr_db": 25.0 +} diff --git a/rinex-cli/src/cli.rs b/rinex-cli/src/cli.rs index f6e557d4c..97fe0de78 100644 --- a/rinex-cli/src/cli.rs +++ b/rinex-cli/src/cli.rs @@ -298,6 +298,7 @@ Refer to README")) if let Ok(content) = std::fs::read_to_string(path) { let opts = serde_json::from_str(&content); if let Ok(opts) = opts { + info!("qc parameter file \"{}\"", path); opts } else { error!("failed to parse parameter file \"{}\"", path); From 51ace4360ea2a7ae320de7dc48e1a0170cd992cb Mon Sep 17 00:00:00 2001 From: "Guillaume W. Bres" Date: Tue, 5 Sep 2023 09:48:51 +0200 Subject: [PATCH 12/19] improving qc Signed-off-by: Guillaume W. Bres --- rinex-cli/src/main.rs | 2 -- rinex/src/algorithm/derivative.rs | 1 + rinex/src/algorithm/mod.rs | 1 + rinex/src/qc/analysis/obs.rs | 51 +++++++++++++++++++++++++------ 4 files changed, 44 insertions(+), 11 deletions(-) diff --git a/rinex-cli/src/main.rs b/rinex-cli/src/main.rs index 5c2bea3fd..772947d13 100644 --- a/rinex-cli/src/main.rs +++ b/rinex-cli/src/main.rs @@ -471,8 +471,6 @@ pub fn main() -> Result<(), rinex::Error> { */ if qc { info!("entering qc mode"); - let mut qc_opts = cli.qc_config(); - /* * QC Config / versus command line * let the possibility to define some parameters diff --git a/rinex/src/algorithm/derivative.rs b/rinex/src/algorithm/derivative.rs index 453b850da..3367fd257 100644 --- a/rinex/src/algorithm/derivative.rs +++ b/rinex/src/algorithm/derivative.rs @@ -30,6 +30,7 @@ impl Derivative { Self { order } } pub fn eval + std::marker::Copy>( + &self, input: Vec<(Epoch, T)>, ) -> Vec<(Epoch, T)> { let mut buf: Vec<(Epoch, T)> = Vec::with_capacity(input.len()); diff --git a/rinex/src/algorithm/mod.rs b/rinex/src/algorithm/mod.rs index 96e068401..b7b60c5f7 100644 --- a/rinex/src/algorithm/mod.rs +++ b/rinex/src/algorithm/mod.rs @@ -11,3 +11,4 @@ pub use filters::{ }; pub use averaging::Averager; +pub use derivative::Derivative; diff --git a/rinex/src/qc/analysis/obs.rs b/rinex/src/qc/analysis/obs.rs index 5357986d9..eeaadc893 100644 --- a/rinex/src/qc/analysis/obs.rs +++ b/rinex/src/qc/analysis/obs.rs @@ -31,6 +31,34 @@ fn report_signals(list: &Vec) -> String { s } +/* + * Report RX Clock drift analysis + */ +fn report_clock_drift(data: &Vec<(Epoch, f64)>) -> Box { + box_html! { + table(class="table is-bordered") { + tr { + th { + : "Epoch" + } + th { + : "Mean Clock drift [s/s]" + } + } + @ for (epoch, drift) in data { + tr { + td { + : epoch.to_string() + } + td { + : format!("{:e}", drift) + } + } + } + } + } +} + /* * Epoch anomalies formatter */ @@ -309,14 +337,6 @@ fn report_ssi_statistics( } } -fn derivative_dt(data: Vec<(Epoch, f64)>) -> Vec<(Epoch, f64)> { - let mut acc = 0.0_f64; - let mut prev_epoch: Option = None; - let mut ret: Vec<(Epoch, f64)> = Vec::new(); - for (epoch, value) in data {} - ret -} - #[derive(Debug, Clone)] /// OBS RINEX specific QC analysis. /// Full OBS RINEX analysis requires both the "obs" and "processing" features. @@ -570,7 +590,8 @@ impl QcObsAnalysis { .recvr_clock() .map(|((e, flag), value)| (e, value)) .collect(); - let rx_clock_drift: Vec<(Epoch, f64)> = derivative_dt(rx_clock); + let der = Derivative::new(1); + let rx_clock_drift: Vec<(Epoch, f64)> = der.eval(rx_clock); let mov = Averager::mov(opts.clock_drift_window); mov.eval(rx_clock_drift) }, @@ -647,6 +668,18 @@ impl HtmlReport for QcObsAnalysis { } } } + tr { + table { + thead { + th { + : "(RX) Clock Drift" + } + } + tbody { + : report_clock_drift(&self.clock_drift) + } + } + } } } } From dc8e2c21015f305b19a10a580dfd5773966494d9 Mon Sep 17 00:00:00 2001 From: "Guillaume W. Bres" Date: Tue, 5 Sep 2023 18:55:47 +0200 Subject: [PATCH 13/19] qc doc Signed-off-by: Guillaume W. Bres --- rinex-cli/config/advanced_repair.json | 12 -- rinex-cli/config/advanced_study.json | 12 -- .../{basic_gnss.json => gnss_snr30db.json} | 0 rinex-cli/config/sv_manual_gap.json | 2 +- rinex-cli/doc/qc.md | 121 +++++++++--------- 5 files changed, 63 insertions(+), 84 deletions(-) delete mode 100644 rinex-cli/config/advanced_repair.json delete mode 100644 rinex-cli/config/advanced_study.json rename rinex-cli/config/{basic_gnss.json => gnss_snr30db.json} (100%) diff --git a/rinex-cli/config/advanced_repair.json b/rinex-cli/config/advanced_repair.json deleted file mode 100644 index d5d4027be..000000000 --- a/rinex-cli/config/advanced_repair.json +++ /dev/null @@ -1,12 +0,0 @@ -{ - "reporting": "PerConstellation", - "statistical": { - "window": "2 hours" - }, - "processing": { - "cs": "study", - "iono_rate_tolerance": 400E-2, - "iono_rate_tolerance_dt": "1 min", - "clk_drift_window": "10 mins" - } -} diff --git a/rinex-cli/config/advanced_study.json b/rinex-cli/config/advanced_study.json deleted file mode 100644 index d5d4027be..000000000 --- a/rinex-cli/config/advanced_study.json +++ /dev/null @@ -1,12 +0,0 @@ -{ - "reporting": "PerConstellation", - "statistical": { - "window": "2 hours" - }, - "processing": { - "cs": "study", - "iono_rate_tolerance": 400E-2, - "iono_rate_tolerance_dt": "1 min", - "clk_drift_window": "10 mins" - } -} diff --git a/rinex-cli/config/basic_gnss.json b/rinex-cli/config/gnss_snr30db.json similarity index 100% rename from rinex-cli/config/basic_gnss.json rename to rinex-cli/config/gnss_snr30db.json diff --git a/rinex-cli/config/sv_manual_gap.json b/rinex-cli/config/sv_manual_gap.json index 764e59183..22aa92e90 100644 --- a/rinex-cli/config/sv_manual_gap.json +++ b/rinex-cli/config/sv_manual_gap.json @@ -1,4 +1,4 @@ { "classification": "Sv", - "min_snr_db": 25.0 + "gap_tolerance": "1 h" } diff --git a/rinex-cli/doc/qc.md b/rinex-cli/doc/qc.md index 241c9a22b..dd560ca9f 100644 --- a/rinex-cli/doc/qc.md +++ b/rinex-cli/doc/qc.md @@ -1,10 +1,10 @@ Quality Check (QC) ================== -RINEX quality check is a special mode of this tool, activated with the `--qc` option. +RINEX quality check is a special mode, activated with `--qc`. QC is first developed for Observation files analysis, but this tool -will accept other RINEX files, for which it will compute basic statistical anlysis. +will accept other RINEX files, for which it will compute basic statistical analysis. QC is a well known procedure in data preprocessing. There are a few differences people who are very familiar with `teqc` must @@ -17,7 +17,7 @@ to advanced "teqc" users. Among them: * Quick GNSS filters (`-G`, `-R`, ...) still exist * In augmented mode, `-no_orbit X` is feasible if you know how to operate the [preprocessor](doc/preprocessing.md) -* Similar position anlysis and reporting +* Similar position analysis and reporting * Similar signals analysis and reporting ## Differences with `teqc` @@ -31,13 +31,14 @@ Although we support as many Navigation files to be specified, when the _augmente This allows for example Glonass and other NAV context to be correctly defined. 1. `--fp [FILE]` for the Observation file -2. `--nav [FILE1] [FILE2]..` : pass NAV contexts +2. `--nav [FILE1] --nav [FILE2]..` : define NAV context +3. `--sp3 [FILE1] --sp3 [FILE2]..` : pass SP3 data We will generate products the [workspace](https://github.com/gwbres/rinex/tree/rinex-cli/workspace) folder, that includes QC reports. -Unlike teqc, we do not support BINEX nor SP3 input data/files as of today. +Unlike teqc, we do not support BINEX input data/files as of today. Unlike teqc we do not limit ourselves to the analysis of GPS and Glonass constellations. @@ -57,10 +58,10 @@ For example, this would apply to ## QC specific command line options * `--qc-only`: ensures the tool will only perform the QCs, -other graphs or record analysis is turned off. This is the most efficient -QC mode. +other graphs or record analysis is turned off. +This command line option makes the quality checks and reporting faster. -* `--qc-config`: pass a configuration file for QC reporting and calculations management (see down below) +* `--qc-cfg`: QC analysis and reporting fine tuning (see down below). ## Basic QC (No NAV) @@ -82,83 +83,85 @@ Run this configuration for the most basic QC: ```bash rinex-cli \ - -P gps,glo \ + -P GPS,GLO \ --qc-only \ - --qc-conf rinex-cli/config/basic.json \ + --qc-cfg rinex-cli/config/gnss_snr30db.json \ --fp test_resources/CRNX/V3/ESBC00DNK_R_20201770000_01D_30S_MO.crx.gz ``` -`--qc-conf` is independent to `--qc` activation. If you activate QC without any configuration, the default configuration is used. -basic.json specifies that we want to report on a constellation basis. -If you compare this report to the previous one +## QC Configuration file + +The QC configuration file allows precise control of how analysis are performed. +We provide example configurations in `rinex-cli/config`. + +All parameters are optional, in a configuration file description. That means, +the default value is used if one parameter is omitted. For example, the classification +method is "mandatory" to a QC analysis. If not defined, the GNSS classification method is used. + + +### QC Classification method + +- "GNSS": the provided context is splitted into all constellations that remain +after pre processing. You can one sub report per constellation. +It is most convenient when one or two constellations are left out, for example with `-P GLO,GPS`. +The report is designed in alphabetical order. + +- "Sv" : split by Satellite vehicle. It is most convenient when a few vehicles are isolated, +for example either by retaining a single constellation : `-P GPS`, or by retaining +specific vehicles: `-P G08,G15,C19,C31,R08`. -1. you get one table per constellation. We only retained GPS and Glonass in this example, -therefore the report is made of two tables. -2. All statistical analysis are made on constellations separately and independently -3. A "25%" statistical window is specified. "window" accepts either a Duration -or a percentage of the file. Here we use the latter, and 25% means we will have 4 statistical -analysis performed over the course of 6 hours, because this file is 24h long. -4. You also have less information in this basic configuration, because most calculations are turned off +- "Physics": split by Observable. -Try this configuration now: +### Sv classification case + +The Sv classification case is interesting, because the Sampling/Gap anlalysis emphasizes +when isolated vehicles go out of sight. + +For example, compare the report generated by this "default" run : ```bash -rinex-cli \ - -F gps,glo \ - --qc-only \ - --qc-conf rinex-cli/config/basic_manual_gap.json \ - --fp test_resources/CRNX/V3/ESBC00DNK_R_20201770000_01D_30S_MO.crx.gz +./target/release/rinex-cli \ + --fp test_resources/CRNX/V3/MOJN00DNK_R_20201770000_01D_30S_MO.crx.gz \ + -P G08,G15,G16,R23,C19,C09 \ + --qc ``` -A "10%" statistical slot is specified, you get more statistical analysis because the time slot -now spans approximately 2 hours. -Also "manual_gap" is specified and set to "10 seconds". That means that a duration -of 10 seconds is now considered as an abnormal gap, in the data gap analysis. -In default configuration, there is no manual gap. That means an abnormal gap -is any abnormal duration above the dominant epoch interval (sample rate). +To this one : + +```bash +./target/release/rinex-cli \ + --fp test_resources/CRNX/V3/MOJN00DNK_R_20201770000_01D_30S_MO.crx.gz \ + -P G08,G15,G16,R23,C19,C09 \ + --qc --qc-cfg rinex-cli/config/sv_manual_gap.json +``` -When no configuration files are given, the default configuration is used +### SNR parametrization -```json TODO -``` -## Advanced configurations +### CS analysis -Now let's move on to more "advanced" configuration, in which basically all -calculations are active and customized +TODO -```bash -rinex-cli \ - -P gps,glo \ - --qc-separate \ - --qc-conf rinex-cli/config/basic_manual_gap.json \ - --fp test_resources/CRNX/V3/ESBC00DNK_R_20201770000_01D_30S_MO.crx.gz -``` +## QC with Navigation Data -## Basic QC (with NAV) +Providing Navigation Data is the only way for complete QCs. +With such context, we can now use the "elevation mask" in the configuration file: -Let's go back to our basic demo and provide Navigation context: +Compare these two runs. In the second one, Sv loss of sight is larger because +the go out of sight more rapidly, due to the stringent elevation criteria : ```bash rinex-cli \ -P gps,glo \ - --qc-separate \ - --qc-conf rinex-cli/config/basic.json \ + -P G08,G15,G16,R23,C19,C09 \ + --qc --qc-cfg rinex-cli/config/sv_manual_gap_ev35.json --fp test_resources/CRNX/V3/ESBC00DNK_R_20201770000_01D_30S_MO.crx.gz \ --nav test_resources/OBS/V3/ESBC00DNK_R_20201770000_01D_MN.rnx.gz ``` - -Navigation context is fully taken into account in advanced calculations -```bash -rinex-cli \ - -F gps,glo \ - --qc-separate \ - --qc-conf rinex-cli/config/advanced_study.json \ - --fp test_resources/CRNX/V3/ESBC00DNK_R_20201770000_01D_30S_MO.crx.gz \ - --nav test_resources/OBS/V3/ESBC00DNK_R_20201770000_01D_MN.rnx.gz -`` +## QC on other RINEX +TODO From e30993d08c8da1a7e8443e362eb3e5b673af9c94 Mon Sep 17 00:00:00 2001 From: "Guillaume W. Bres" Date: Tue, 5 Sep 2023 19:04:43 +0200 Subject: [PATCH 14/19] improving qc Signed-off-by: Guillaume W. Bres --- rinex/src/qc/analysis/obs.rs | 60 ++++++++++++++++++------------------ 1 file changed, 30 insertions(+), 30 deletions(-) diff --git a/rinex/src/qc/analysis/obs.rs b/rinex/src/qc/analysis/obs.rs index eeaadc893..7cec2d2bf 100644 --- a/rinex/src/qc/analysis/obs.rs +++ b/rinex/src/qc/analysis/obs.rs @@ -36,22 +36,28 @@ fn report_signals(list: &Vec) -> String { */ fn report_clock_drift(data: &Vec<(Epoch, f64)>) -> Box { box_html! { - table(class="table is-bordered") { - tr { - th { - : "Epoch" - } - th { - : "Mean Clock drift [s/s]" - } + @ if data.is_empty() { + td { + : "Unfeasible (Data Missing)" } - @ for (epoch, drift) in data { + } else { + table(class="table is-bordered") { tr { - td { - : epoch.to_string() + th { + : "Epoch" } - td { - : format!("{:e}", drift) + th { + : "Mean Clock drift [s/s]" + } + } + @ for (epoch, drift) in data { + tr { + td { + : epoch.to_string() + } + td { + : format!("{:e}", drift) + } } } } @@ -486,7 +492,6 @@ impl QcObsAnalysis { continue; // phase should have SNR information attached to it } } - /* * Signal condition */ @@ -533,7 +538,6 @@ impl QcObsAnalysis { } } } - /* * SSI statistical analysis */ @@ -657,25 +661,21 @@ impl HtmlReport for QcObsAnalysis { } } tr { - table { - thead { - th { - : "SNR" - } - } - tbody { - : report_snr_analysis(self.min_snr, self.max_snr, &self.ssi_stats) - } + th { + : "SNR" } } tr { table { - thead { - th { - : "(RX) Clock Drift" - } - } - tbody { + : report_snr_analysis(self.min_snr, self.max_snr, &self.ssi_stats) + } + } + tr { + th { + : "(RX) Clock Drift" + } + td { + table { : report_clock_drift(&self.clock_drift) } } From f5d95054b67c8e475384333c35e380a337bf75b3 Mon Sep 17 00:00:00 2001 From: "Guillaume W. Bres" Date: Wed, 6 Sep 2023 09:10:45 +0200 Subject: [PATCH 15/19] improving qc Signed-off-by: Guillaume W. Bres --- rinex-cli/config/sv_manual_gap.json | 5 +- rinex/src/qc/analysis/obs.rs | 77 ++++++++++++++++------------- rinex/src/qc/mod.rs | 5 ++ rinex/src/qc/opts.rs | 6 +-- 4 files changed, 55 insertions(+), 38 deletions(-) diff --git a/rinex-cli/config/sv_manual_gap.json b/rinex-cli/config/sv_manual_gap.json index 22aa92e90..d7b4f68f6 100644 --- a/rinex-cli/config/sv_manual_gap.json +++ b/rinex-cli/config/sv_manual_gap.json @@ -1,4 +1,7 @@ { "classification": "Sv", - "gap_tolerance": "1 h" + "gap_tolerance": { + "centuries": 0, + "nanoseconds": 7200000000000 + } } diff --git a/rinex/src/qc/analysis/obs.rs b/rinex/src/qc/analysis/obs.rs index 7cec2d2bf..d1527c740 100644 --- a/rinex/src/qc/analysis/obs.rs +++ b/rinex/src/qc/analysis/obs.rs @@ -37,8 +37,15 @@ fn report_signals(list: &Vec) -> String { fn report_clock_drift(data: &Vec<(Epoch, f64)>) -> Box { box_html! { @ if data.is_empty() { - td { - : "Unfeasible (Data Missing)" + table(class="table is-bordered") { + tr { + th { + : "Unfeasible" + } + td { + : "Missing Data" + } + } } } else { table(class="table is-bordered") { @@ -74,13 +81,6 @@ fn report_anomalies<'a>( other: &'a Vec<(Epoch, EpochFlag)>, ) -> Box { box_html! { - table(class="table is-bordered") { - thead { - th { - : "Anomalies" - } - } - tbody { tr { th { : "Power Failures" @@ -171,8 +171,6 @@ fn report_anomalies<'a>( } } } - } - } } } @@ -189,12 +187,6 @@ fn report_epoch_completion( ) -> Box { box_html! { table(class="table is-bordered") { - thead { - th { - : "Epochs" - } - } - tbody { tr { th { : "Total#" @@ -230,7 +222,6 @@ fn report_epoch_completion( } } } - } } } } @@ -651,31 +642,49 @@ impl HtmlReport for QcObsAnalysis { } } tr { - table { - : report_anomalies(&self.cs_anomalies, &self.power_failures, &self.other_anomalies) - } - } - tr { - table { - : report_epoch_completion(self.total_epochs, self.total_with_obs, &self.complete_epochs) + table(class="table is-bordered") { + thead { + th { + : "Anomalies" + } + } + tbody { + : report_anomalies(&self.cs_anomalies, &self.power_failures, &self.other_anomalies) + } } } tr { - th { - : "SNR" + table(class="table is-bordered") { + thead { + th { + : "Epochs" + } + } + tbody { + : report_epoch_completion(self.total_epochs, self.total_with_obs, &self.complete_epochs) + } } } tr { - table { - : report_snr_analysis(self.min_snr, self.max_snr, &self.ssi_stats) + table(class="table is-bordered") { + thead { + th { + : "SNR" + } + } + tbody { + : report_snr_analysis(self.min_snr, self.max_snr, &self.ssi_stats) + } } } tr { - th { - : "(RX) Clock Drift" - } - td { - table { + table(class="table is-bordered") { + thead { + th { + : "(RX) Clock Drift" + } + } + tbody { : report_clock_drift(&self.clock_drift) } } diff --git a/rinex/src/qc/mod.rs b/rinex/src/qc/mod.rs index 5d143c99e..8474a61d2 100644 --- a/rinex/src/qc/mod.rs +++ b/rinex/src/qc/mod.rs @@ -207,6 +207,11 @@ impl QcReport { } //div=parameters div(id="header") { table(class="table is-bordered; style=\"margin-bottom: 20px\"") { + thead { + th { + : "File Header" + } + } tbody { : context.primary_data().header.to_inline_html() } diff --git a/rinex/src/qc/opts.rs b/rinex/src/qc/opts.rs index bf56c300d..5bb3acb2e 100644 --- a/rinex/src/qc/opts.rs +++ b/rinex/src/qc/opts.rs @@ -251,11 +251,11 @@ impl HtmlReport for QcOpts { : "Data gap" } @ if let Some(tol) = self.gap_tolerance { - th { - : format!("Tolerance: ({})", tol) + td { + : format!("{} tolerance", tol) } } else { - th { + td { : "No tolerance" } } From 59718845404d4cc75b42f311cf793e061cbe2c94 Mon Sep 17 00:00:00 2001 From: "Guillaume W. Bres" Date: Wed, 6 Sep 2023 10:35:33 +0200 Subject: [PATCH 16/19] improving qc Signed-off-by: Guillaume W. Bres --- rinex/Cargo.toml | 4 +- rinex/src/lib.rs | 95 ++++-------------- rinex/src/qc/analysis/obs.rs | 186 ++++++++++++++--------------------- 3 files changed, 94 insertions(+), 191 deletions(-) diff --git a/rinex/Cargo.toml b/rinex/Cargo.toml index c38dbce5d..5ffb344bb 100644 --- a/rinex/Cargo.toml +++ b/rinex/Cargo.toml @@ -15,11 +15,11 @@ rust-version = "1.64" [features] default = [] # no features by default sbas = ["geo", "wkt"] -obs = ["statrs"] +obs = [] meteo = [] nav = [] processing = [] -qc = ["horrorshow", "processing", "sp3"] # rinex Quality Check (mainly OBS RINEX) +qc = ["horrorshow", "processing", "sp3", "statrs"] # rinex Quality Check (mainly OBS RINEX) [package.metadata.docs.rs] all-features = true diff --git a/rinex/src/lib.rs b/rinex/src/lib.rs index 4de29c234..d64987fe2 100644 --- a/rinex/src/lib.rs +++ b/rinex/src/lib.rs @@ -1802,11 +1802,13 @@ impl Rinex { } } +#[cfg(feature = "obs")] +use crate::observation::Snr; + /* * OBS RINEX specific methods: only available on crate feature. * Either specific Iterators, or meaningful data we can extract. */ - #[cfg(feature = "obs")] #[cfg_attr(docrs, doc(cfg(feature = "obs")))] impl Rinex { @@ -2059,6 +2061,21 @@ impl Rinex { }) })) } + /// Returns an Iterator over signal SNR indications. + /// All observation that did not come with such indication are filtered out. + pub fn snr(&self) -> Box + '_> { + Box::new(self.observation().flat_map(|(e, (_, vehicles))| { + vehicles.iter().flat_map(|(sv, observations)| { + observations.iter().filter_map(|(obs, obsdata)| { + if let Some(snr) = obsdata.snr { + Some((*e, *sv, obs, snr)) + } else { + None + } + }) + }) + })) + } } #[cfg(feature = "nav")] @@ -2851,82 +2868,6 @@ impl Decimate for Rinex { } } -#[cfg(feature = "obs")] -use crate::observation::Observation; - -#[cfg(feature = "obs")] -#[cfg_attr(docrs, doc(cfg(feature = "obs")))] -impl Observation for Rinex { - fn min(&self) -> (Option, HashMap>) { - if let Some(r) = self.record.as_obs() { - r.min() - } else { - (None, HashMap::new()) - } - } - fn min_observable(&self) -> HashMap { - if let Some(r) = self.record.as_obs() { - r.min_observable() - } else if let Some(r) = self.record.as_meteo() { - r.min_observable() - } else { - HashMap::new() - } - } - fn max(&self) -> (Option, HashMap>) { - if let Some(r) = self.record.as_obs() { - r.max() - } else { - (None, HashMap::new()) - } - } - fn max_observable(&self) -> HashMap { - if let Some(r) = self.record.as_obs() { - r.max_observable() - } else if let Some(r) = self.record.as_meteo() { - r.max_observable() - } else { - HashMap::new() - } - } - fn mean(&self) -> (Option, HashMap>) { - if let Some(r) = self.record.as_obs() { - r.mean() - } else { - (None, HashMap::new()) - } - } - fn mean_observable(&self) -> HashMap { - if let Some(r) = self.record.as_obs() { - r.mean_observable() - } else if let Some(r) = self.record.as_meteo() { - r.mean_observable() - } else { - HashMap::new() - } - } - fn std_dev(&self) -> (Option, HashMap>) { - if let Some(r) = self.record.as_obs() { - r.std_dev() - } else { - (None, HashMap::new()) - } - } - fn std_dev_observable(&self) -> HashMap { - HashMap::new() - } - fn std_var(&self) -> (Option, HashMap>) { - if let Some(r) = self.record.as_obs() { - r.std_var() - } else { - (None, HashMap::new()) - } - } - fn std_var_observable(&self) -> HashMap { - HashMap::new() - } -} - #[cfg(feature = "obs")] use observation::Dcb; diff --git a/rinex/src/qc/analysis/obs.rs b/rinex/src/qc/analysis/obs.rs index d1527c740..20f46d59b 100644 --- a/rinex/src/qc/analysis/obs.rs +++ b/rinex/src/qc/analysis/obs.rs @@ -4,6 +4,8 @@ use itertools::Itertools; use super::{pretty_array, QcOpts}; use std::collections::HashMap; +use statrs::statistics::Statistics; + #[cfg(feature = "processing")] use crate::preprocessing::*; // include preprocessing toolkit when feasible, // for complete analysis @@ -229,61 +231,34 @@ fn report_epoch_completion( /* * SNR analysis report */ -fn report_snr_analysis( - min: (Sv, Epoch, Snr), - max: (Sv, Epoch, Snr), - ssi_stats: &HashMap, +fn report_snr_statistics( + snr_stats: &HashMap, ) -> Box { box_html! { - tr { - th { - : "Minimum" - } - td { - p { - : min.0.to_string() - } - p { - : format!("{:e}", min.2) - } - p { - : format!("@{}", min.1) + @ for (observable, _) in snr_stats { + tr { + th { + : observable.to_string() } } } tr { th { - : "Best" + : "Worst" } - td { - p { - : max.0.to_string() - } - p { - : format!("{:e}", max.2) - } - p { - : format!("@{}", max.1) + @ for (_, (min, _)) in snr_stats { + td { + : format!("{:e} @{}", min.0, min.1) } } } - tr { - th { - : "SSI Statistics" - } - td { - : report_ssi_statistics(ssi_stats) - } - } } } /* * Reports statistical analysis results for SSx observations */ -fn report_ssi_statistics( - ssi_stats: &HashMap, -) -> Box { +fn report_ssi_statistics(ssi_stats: &HashMap) -> Box { box_html! { table(class="table is-bordered") { thead { @@ -303,7 +278,7 @@ fn report_ssi_statistics( th { : "Mean" } - @ for (_, (mean, _, _)) in ssi_stats { + @ for (_, (mean, _)) in ssi_stats { td { : format!("{:.3} dB", mean) } @@ -313,19 +288,9 @@ fn report_ssi_statistics( th { : "Deviation" // (σ)" } - @ for (_, (_, sigma, _)) in ssi_stats { - td { - : format!("{:.3} dB", sigma) - } - } - } - tr { - th { - : "Skewness" //(μ / σγ)" - } - @ for (_, (_, _, sk)) in ssi_stats { + @ for (_, (_, std)) in ssi_stats { td { - : format!("{:.3}", sk) + : format!("{:.3} dB", std) } } } @@ -359,14 +324,11 @@ pub struct QcObsAnalysis { /// Complete epochs, with respect to given signal complete_epochs: Vec<(Carrier, usize)>, #[cfg(feature = "obs")] - /// Min. SNR (sv @ epoch) - min_snr: (Sv, Epoch, Snr), + /// Min. & Max. SNR (signal @ epoch) + snr_stats: HashMap, #[cfg(feature = "obs")] - /// Max. SNR (sv @ epoch) - max_snr: (Sv, Epoch, Snr), - #[cfg(feature = "obs")] - /// SSi statistical analysis (mean, stddev, skew) - ssi_stats: HashMap, + /// SSI statistical analysis (mean, stddev) + ssi_stats: HashMap, #[cfg(feature = "obs")] /// RX clock drift clock_drift: Vec<(Epoch, f64)>, @@ -416,22 +378,8 @@ impl QcObsAnalysis { let mut total_epochs = rnx.epoch().count(); let mut epoch_with_obs: Vec = Vec::new(); - let mut complete_epochs: HashMap = HashMap::new(); - let mut min_snr = (Sv::default(), Epoch::default(), Snr::DbHz54); - let mut max_snr = (Sv::default(), Epoch::default(), Snr::DbHz0); - - let mut ssi_stats: HashMap = HashMap::new(); - - if let Some(rec) = rnx.record.as_obs() { - let mask: Filter = Filter::from(MaskFilter { - operand: MaskOperand::Equals, - item: TargetItem::ClockItem, - }); - let clk_data = rec.filter(mask); - } - if let Some(r) = rnx.record.as_obs() { total_epochs = r.len(); for ((epoch, flag), (clk, svs)) in r { @@ -441,21 +389,6 @@ impl QcObsAnalysis { epoch_with_obs.push(*epoch); } } - - for (observable, observation) in observables { - if let Some(snr) = observation.snr { - if snr < min_snr.2 { - min_snr.0 = *sv; - min_snr.1 = *epoch; - min_snr.2 = snr; - } - if snr > max_snr.2 { - max_snr.0 = *sv; - max_snr.1 = *epoch; - max_snr.2 = snr; - } - } - } } } /* @@ -529,31 +462,49 @@ impl QcObsAnalysis { } } } - /* - * SSI statistical analysis - */ - let mut mean_ssi: HashMap<_, _> = r.mean_observable(); - mean_ssi.retain(|obs, _| obs.is_ssi_observable()); - for (obs, mean) in mean_ssi { - ssi_stats.insert(obs.clone(), (mean, 0.0_f64, 0.0_f64)); - } - - let mut stddev_ssi: HashMap<_, _> = r.mean_observable(); // TODO r.stddev_observable(); - stddev_ssi.retain(|obs, _| obs.is_ssi_observable()); - for (obs, stddev) in stddev_ssi { - if let Some((_, dev, _)) = ssi_stats.get_mut(&obs) { - *dev = stddev; - } + } + // append ssi: drop vehicle differentiation + let mut ssi: HashMap> = HashMap::new(); + for (_, _, obs, value) in rnx.ssi() { + if let Some(values) = ssi.get_mut(&obs) { + values.push(value); + } else { + ssi.insert(obs.clone(), vec![value]); } - - let mut skew_ssi: HashMap<_, _> = r.mean_observable(); // TODO r.skewness_observable(); - skew_ssi.retain(|obs, _| obs.is_ssi_observable()); - for (obs, skew) in skew_ssi { - if let Some((_, _, sk)) = ssi_stats.get_mut(&obs) { - *sk = skew; - } + } + /* + * SSI statistical analysis: {mean, stddev,} + * per signal: we do not differentiate vehicles + */ + let ssi_stats: HashMap = ssi + .iter() + .map(|(obs, values)| (obs.clone(), (values.mean(), values.std_dev()))) + .collect(); + // append snr: drop vehicle differentiation + let mut snr: HashMap> = HashMap::new(); + for (e, _, obs, value) in rnx.snr() { + let value_f64: f64 = (value as u8).into(); + if let Some(values) = snr.get_mut(&obs) { + values.push((e.0, value_f64)); + } else { + snr.insert(obs.clone(), vec![(e.0, value_f64)]); } } + /* + * SNR analysis: {min, max} + * per signal: we do not differentiate vehicles + */ + let snr_stats: HashMap = snr + .iter() + .map(|(obs, values)| { + let min = Statistics::min(values.1); + let max = Statistics::max(values.1); + ( + obs.clone(), + ((Epoch::default(), min), (Epoch::default(), max)), + ) + }) + .collect(); /* * sort, prior reporting */ @@ -577,8 +528,7 @@ impl QcObsAnalysis { ret.sort(); ret }, - min_snr, - max_snr, + snr_stats, ssi_stats, clock_drift: { let rx_clock: Vec<_> = rnx @@ -673,7 +623,19 @@ impl HtmlReport for QcObsAnalysis { } } tbody { - : report_snr_analysis(self.min_snr, self.max_snr, &self.ssi_stats) + : report_snr_statistics(&self.snr_stats) + } + } + } + tr { + table(class="table is-bordered") { + thead { + th { + : "SSI" + } + } + tbody { + : report_ssi_statistics(&self.ssi_stats) } } } From 6bcc8426928b925a735b2f183e5cfa602016dafe Mon Sep 17 00:00:00 2001 From: "Guillaume W. Bres" Date: Thu, 7 Sep 2023 10:05:13 +0200 Subject: [PATCH 17/19] fix snr reporting Signed-off-by: Guillaume W. Bres --- rinex/src/qc/analysis/obs.rs | 95 +++++++++++++++++++++++++----------- 1 file changed, 66 insertions(+), 29 deletions(-) diff --git a/rinex/src/qc/analysis/obs.rs b/rinex/src/qc/analysis/obs.rs index 20f46d59b..7e15e4c74 100644 --- a/rinex/src/qc/analysis/obs.rs +++ b/rinex/src/qc/analysis/obs.rs @@ -1,5 +1,6 @@ use crate::{carrier, observation::Snr, prelude::*, Carrier}; use itertools::Itertools; +use std::str::FromStr; use super::{pretty_array, QcOpts}; use std::collections::HashMap; @@ -206,8 +207,19 @@ fn report_epoch_completion( } } tr { - th { - : "Complete" + td { + b { + : "Complete" + } + p { + : "Epochs with at least Phase + PR" + } + p { + : "in dual frequency, with" + } + p { + : "both SNR and elev above masks" + } } @ for (signal, count) in complete { td { @@ -215,12 +227,8 @@ fn report_epoch_completion( : format!("L1/{}", signal) } p { - : count.to_string() - } - b { - : format!("{}%", count * 100/total) + : format!("{} ({}%)", count, count * 100 / total) } - //: format!("L1/{}: {} ({}%)", signal, count, count *100/total) } } } @@ -235,10 +243,32 @@ fn report_snr_statistics( snr_stats: &HashMap, ) -> Box { box_html! { - @ for (observable, _) in snr_stats { - tr { - th { - : observable.to_string() + tr { + td { + : "" + } + @ for (observable, _) in snr_stats { + @ if observable.is_phase_observable() || observable.is_pseudorange_observable() { + td { + : observable.to_string() + } + } + } + } + tr { + th { + : "Best" + } + @ for (observable, (_, max)) in snr_stats { + @ if observable.is_phase_observable() || observable.is_pseudorange_observable() { + td { + b { + : format!("{:e}", Snr::from_str(&format!("{}", max.1)).unwrap()) + } + p { + : format!("@{}", max.0) + } + } } } } @@ -246,9 +276,16 @@ fn report_snr_statistics( th { : "Worst" } - @ for (_, (min, _)) in snr_stats { - td { - : format!("{:e} @{}", min.0, min.1) + @ for (observable, (min, _)) in snr_stats { + @ if observable.is_phase_observable() || observable.is_pseudorange_observable() { + td { + b { + : format!("{:e}", Snr::from_str(&format!("{}", min.1)).unwrap()) + } + p { + : format!("@{}", min.0) + } + } } } } @@ -482,29 +519,29 @@ impl QcObsAnalysis { .collect(); // append snr: drop vehicle differentiation let mut snr: HashMap> = HashMap::new(); - for (e, _, obs, value) in rnx.snr() { - let value_f64: f64 = (value as u8).into(); + for ((e, _), _, obs, snr_value) in rnx.snr() { + let snr_f64: f64 = (snr_value as u8).into(); if let Some(values) = snr.get_mut(&obs) { - values.push((e.0, value_f64)); + values.push((e, snr_f64)); } else { - snr.insert(obs.clone(), vec![(e.0, value_f64)]); + snr.insert(obs.clone(), vec![(e, snr_f64)]); } } /* * SNR analysis: {min, max} * per signal: we do not differentiate vehicles */ - let snr_stats: HashMap = snr - .iter() - .map(|(obs, values)| { - let min = Statistics::min(values.1); - let max = Statistics::max(values.1); - ( - obs.clone(), - ((Epoch::default(), min), (Epoch::default(), max)), - ) - }) - .collect(); + let mut snr_stats: HashMap = HashMap::new(); + for (obs, data) in snr { + let values: Vec = data.iter().map(|(e, value)| *value).collect(); + let min = values.clone().min(); + println!("MIN: {}", min); + let epoch_min = data.iter().find(|(e, value)| *value == min).unwrap().0; + let max = values.clone().max(); + println!("MAX: {}", max); + let epoch_max = data.iter().find(|(e, value)| *value == max).unwrap().0; + snr_stats.insert(obs, ((epoch_min, min), (epoch_max, max))); + } /* * sort, prior reporting */ From 86f642daeaef57cc32a3b5777eb2da6554f7f70f Mon Sep 17 00:00:00 2001 From: "Guillaume W. Bres" Date: Thu, 7 Sep 2023 10:29:21 +0200 Subject: [PATCH 18/19] fix warnings Signed-off-by: Guillaume W. Bres --- rinex/src/algorithm/derivative.rs | 16 +-- rinex/src/algorithm/mod.rs | 4 +- rinex/src/qc/analysis/mod.rs | 4 +- rinex/src/qc/analysis/obs.rs | 162 ++++++++++++++---------------- rinex/src/qc/analysis/sv.rs | 4 +- rinex/src/qc/mod.rs | 3 +- 6 files changed, 91 insertions(+), 102 deletions(-) diff --git a/rinex/src/algorithm/derivative.rs b/rinex/src/algorithm/derivative.rs index 3367fd257..5de98965d 100644 --- a/rinex/src/algorithm/derivative.rs +++ b/rinex/src/algorithm/derivative.rs @@ -7,15 +7,12 @@ pub struct Derivative { /* * Derivative of an number of array sorted by chronological Epoch */ -pub(crate) fn derivative + std::marker::Copy>( - input: Vec<(Epoch, T)>, - buf: &mut Vec<(Epoch, T)>, -) { - let mut prev: Option<(Epoch, T)> = None; +pub(crate) fn derivative(input: Vec<(Epoch, f64)>, buf: &mut Vec<(Epoch, f64)>) { + let mut prev: Option<(Epoch, f64)> = None; for (e, value) in input { if let Some((prev_e, prev_v)) = prev { let dt = e - prev_e; - let dy = value - prev_v; + let dy = (value - prev_v) / dt.to_seconds(); buf.push((e, dy)); } prev = Some((e, value)); @@ -29,11 +26,8 @@ impl Derivative { pub fn new(order: usize) -> Self { Self { order } } - pub fn eval + std::marker::Copy>( - &self, - input: Vec<(Epoch, T)>, - ) -> Vec<(Epoch, T)> { - let mut buf: Vec<(Epoch, T)> = Vec::with_capacity(input.len()); + pub fn eval(&self, input: Vec<(Epoch, f64)>) -> Vec<(Epoch, f64)> { + let mut buf: Vec<(Epoch, f64)> = Vec::with_capacity(input.len()); derivative(input, &mut buf); //for i in 1..self.order { // derivative(&ret, &mut ret); diff --git a/rinex/src/algorithm/mod.rs b/rinex/src/algorithm/mod.rs index b7b60c5f7..5ee754691 100644 --- a/rinex/src/algorithm/mod.rs +++ b/rinex/src/algorithm/mod.rs @@ -1,4 +1,4 @@ -mod averaging; +//mod averaging; mod derivative; mod filters; mod target; @@ -10,5 +10,5 @@ pub use filters::{ Mask, MaskFilter, MaskOperand, Preprocessing, Smooth, SmoothingFilter, SmoothingType, }; -pub use averaging::Averager; +//pub use averaging::Averager; pub use derivative::Derivative; diff --git a/rinex/src/qc/analysis/mod.rs b/rinex/src/qc/analysis/mod.rs index 2c330e917..3d1b862be 100644 --- a/rinex/src/qc/analysis/mod.rs +++ b/rinex/src/qc/analysis/mod.rs @@ -1,6 +1,6 @@ -use super::{pretty_array, table_lengthy_td, HtmlReport, QcOpts}; +use super::{pretty_array, HtmlReport, QcOpts}; use crate::prelude::*; -use horrorshow::{helper::doctype, RenderBox}; +use horrorshow::{helper::doctype, RenderBox}; //table_lengthy_td mod sv; diff --git a/rinex/src/qc/analysis/obs.rs b/rinex/src/qc/analysis/obs.rs index 7e15e4c74..c5d5f3cbd 100644 --- a/rinex/src/qc/analysis/obs.rs +++ b/rinex/src/qc/analysis/obs.rs @@ -11,9 +11,6 @@ use statrs::statistics::Statistics; use crate::preprocessing::*; // include preprocessing toolkit when feasible, // for complete analysis -#[cfg(feature = "obs")] -use crate::observation::Observation; // having this feature unlocks full OBS RINEX analysis - /* * GNSS signal special formatting */ @@ -84,96 +81,91 @@ fn report_anomalies<'a>( other: &'a Vec<(Epoch, EpochFlag)>, ) -> Box { box_html! { + tr { + th { + : "Power Failures" + } + @ if power.is_empty() { + td { + : "None" + } + } else { + td { + : format!("{:?}", power) + } tr { th { - : "Power Failures" + : "Longest" } - @ if power.is_empty() { - td { - : "None" - } - } else { - td { - : format!("{:?}", power) - } - tr { - th { - : "Longest" - } - td { - //: power.iter().max_by(|(_, d1), (_, d2)| d1.cmp(d2)).unwrap().to_string() - : "TODO" - } - td { - : "Average Power failure" - } - td { - : "TODO" - } - } + td { + //: power.iter().max_by(|(_, d1), (_, d2)| d1.cmp(d2)).unwrap().to_string() + : "TODO" } - } - tr { - th { - : "Cycle slip(s)" + td { + : "Average Duration" } - @ if cs.is_empty() { - td { - : "None" - } - } else { - td { - : format!("{:?}", cs) - } + td { + : "TODO" } } - tr { - th { - : "Other anomalies" + } + } + tr { + th { + : "Cycle slip(s)" + } + @ if cs.is_empty() { + td { + : "None" + } + } else { + td { + : format!("{:?}", cs) + } + } + } + tr { + th { + : "Other anomalies" + } + @ if other.is_empty() { + td { + : "None" + } + } else { + td { + : "Epoch" + } + td { + : "Event" + } + @ for (e, event) in other { + td { + : "" + } + td { + : format!("{}", e) } - @ if other.is_empty() { + @ if *event == EpochFlag::AntennaBeingMoved { td { - : "None" + : "Antenna Being Moved" } - } else { + } else if *event == EpochFlag::NewSiteOccupation { td { - : "Epoch" + : "New Site Occupation" } + } else if *event == EpochFlag::ExternalEvent { td { - : "Event" + : "External Event" } - @ for (e, event) in other { - td { - : "" - } - td { - : e.to_string() - } - //@ match event { - // EpochFlag::AntennaBeingMoved => { - // td { - // : "Antenna being moved" - // } - // }, - // _ => {}, - //} - // } - //td { - // @ match event { - // EpochFlag::NewSiteOccupation => { - // : "New Site Occupation" - // } - // EpochFlag::ExternalEvent => { - // : "External Event" - // } - // _ => { - // : "Other" - // } - // } - //} + } else { + td { + : "Other" } } } + } + } } } @@ -419,8 +411,8 @@ impl QcObsAnalysis { if let Some(r) = rnx.record.as_obs() { total_epochs = r.len(); - for ((epoch, flag), (clk, svs)) in r { - for (sv, observables) in svs { + for ((epoch, _flag), (_clk, svs)) in r { + for (_sv, observables) in svs { if !observables.is_empty() { if !epoch_with_obs.contains(&epoch) { epoch_with_obs.push(*epoch); @@ -533,13 +525,13 @@ impl QcObsAnalysis { */ let mut snr_stats: HashMap = HashMap::new(); for (obs, data) in snr { - let values: Vec = data.iter().map(|(e, value)| *value).collect(); + let values: Vec = data.iter().map(|(_e, value)| *value).collect(); let min = values.clone().min(); println!("MIN: {}", min); - let epoch_min = data.iter().find(|(e, value)| *value == min).unwrap().0; + let epoch_min = data.iter().find(|(_e, value)| *value == min).unwrap().0; let max = values.clone().max(); println!("MAX: {}", max); - let epoch_max = data.iter().find(|(e, value)| *value == max).unwrap().0; + let epoch_max = data.iter().find(|(_e, value)| *value == max).unwrap().0; snr_stats.insert(obs, ((epoch_min, min), (epoch_max, max))); } /* @@ -570,12 +562,14 @@ impl QcObsAnalysis { clock_drift: { let rx_clock: Vec<_> = rnx .recvr_clock() - .map(|((e, flag), value)| (e, value)) + .map(|((e, _flag), value)| (e, value)) .collect(); let der = Derivative::new(1); let rx_clock_drift: Vec<(Epoch, f64)> = der.eval(rx_clock); - let mov = Averager::mov(opts.clock_drift_window); - mov.eval(rx_clock_drift) + //TODO + //let mov = Averager::mov(opts.clock_drift_window); + //mov.eval(rx_clock_drift) + rx_clock_drift }, } } diff --git a/rinex/src/qc/analysis/sv.rs b/rinex/src/qc/analysis/sv.rs index a482656ea..dcce4fa05 100644 --- a/rinex/src/qc/analysis/sv.rs +++ b/rinex/src/qc/analysis/sv.rs @@ -1,5 +1,5 @@ -use super::{pretty_array, table_lengthy_td, QcOpts}; -use crate::prelude::*; +use super::{pretty_array, QcOpts}; +use crate::prelude::*; //table_lengthy_td #[derive(Debug, Clone)] pub struct QcSvAnalysis { diff --git a/rinex/src/qc/mod.rs b/rinex/src/qc/mod.rs index 8474a61d2..3648c6595 100644 --- a/rinex/src/qc/mod.rs +++ b/rinex/src/qc/mod.rs @@ -16,7 +16,7 @@ use crate::preprocessing::*; /* * Methods used when reporting lenghty vectors or data subsets in a table. * Makes tables cleaner and nicer by wrapping string content, into several paragraphs. - */ + * pub(crate) fn table_lengthy_td( list: &Vec, max_items: usize, @@ -41,6 +41,7 @@ pub(crate) fn table_lengthy_td( } } } +*/ /* * Array (CSV) pretty formatter From 6c679d6ae9a2784a2c58499cd044f73b79960fc1 Mon Sep 17 00:00:00 2001 From: "Guillaume W. Bres" Date: Thu, 7 Sep 2023 11:37:26 +0200 Subject: [PATCH 19/19] remove "observation" trait, rely on statrs and mark Mp as work in progress Signed-off-by: Guillaume W. Bres --- rinex-cli/src/main.rs | 24 +-- rinex/src/lib.rs | 125 ++++++++---- rinex/src/meteo/record.rs | 66 ------ rinex/src/observation/mod.rs | 113 ----------- rinex/src/observation/record.rs | 305 +--------------------------- rinex/src/tests/mod.rs | 1 - rinex/src/tests/stats.rs | 347 -------------------------------- 7 files changed, 95 insertions(+), 886 deletions(-) delete mode 100644 rinex/src/tests/stats.rs diff --git a/rinex-cli/src/main.rs b/rinex-cli/src/main.rs index 772947d13..671546bea 100644 --- a/rinex-cli/src/main.rs +++ b/rinex-cli/src/main.rs @@ -270,18 +270,18 @@ pub fn main() -> Result<(), rinex::Error> { * Code Multipath analysis */ if cli.multipath() { - let data = ctx - .primary_data() - .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"); + //let data = ctx + // .primary_data() + // .observation_phase_align_origin() + // .observation_phase_carrier_cycles() + // .mp(); + //plot::plot_gnss_dcb( + // &mut plot_ctx, + // "Code Multipath Biases", + // "Meters of delay", + // &data, + //); + warn!("--mp analysis not available yet"); } /* * [GF] recombination visualization requested diff --git a/rinex/src/lib.rs b/rinex/src/lib.rs index d64987fe2..701a10d99 100644 --- a/rinex/src/lib.rs +++ b/rinex/src/lib.rs @@ -984,36 +984,6 @@ impl Rinex { results } */ - /* - /// Extracts Raw Carrier Phase observations, - /// from this Observation record, on an epoch basis an per space vehicle. - /// Does not produce anything if self is not an Observation RINEX. - pub fn observation_phase(&self) -> BTreeMap<(Epoch, EpochFlag), HashMap>> { - let mut ret: BTreeMap>> = BTreeMap::new(); - if let Some(r) = self.record.as_obs() { - for ((e, _), (_, sv)) in record.iter() { - let mut map: BTreeMap> = BTreeMap::new(); - for (sv, obs) in sv.iter() { - let mut v: Vec<(String, f64)> = Vec::new(); - for (observable, data) in obs.iter() { - if observable.is_phase_observation() { - v.push((code.clone(), data.obs)); - } - } - if v.len() > 0 { - // did come with at least 1 Phase obs - map.insert(*sv, v); - } - } - if map.len() > 0 { - // did produce something - results.insert(*e, map); - } - } - } - ret - } - */ /* /// Extracts Carrier phases without Ionospheric path delay contributions, /// by extracting [carrier_phases] and using the differential (dual frequency) compensation. @@ -2882,19 +2852,88 @@ impl Dcb for Rinex { } } -#[cfg(feature = "obs")] -use observation::Mp; - -#[cfg(feature = "obs")] -impl Mp for Rinex { - fn mp(&self) -> HashMap>> { - if let Some(r) = self.record.as_obs() { - r.dcb() - } else { - panic!("wrong rinex type"); - } - } -} +//#[cfg(feature = "obs")] +//use observation::Mp; +// +//#[cfg(feature = "obs")] +//impl Mp for Rinex { +// fn mp(&self) -> HashMap>> { +// /* +// * Determine mean value of all observed Phase and Pseudo Range +// * observations, for all Sv +// */ +// let mut mean: HashMap> = HashMap::new(); +// /* +// * Associate a Phase code to all PR codes +// */ +// let mut associations: HashMap = HashMap::new(); +// +// let pr_codes: Vec = self.observable() +// .filter_map(|obs| +// if obs.is_pseudorange_observable() { +// Some(obs.clone()) +// } else { +// None +// }) +// .collect(); +// +// for observable in self.observable() { +// if !observable.is_phase_observable() { +// if !observable.is_pseudorange_observable() { +// continue; +// } +// } +// // code associations (for future combianations) +// if observable.is_phase_observable() { +// for pr_code in &pr_codes { +// +// } +// } +// +// for sv in self.sv() { +// /* +// * average phase values +// */ +// let phases: Vec = self.carrier_phase() +// .filter_map(|(_, svnn, obs, ph)| { +// if sv == svnn && observable == obs { +// Some(ph) +// } else { +// None +// } +// }) +// .collect(); +// if let Some(data) = mean.get_mut(&sv) { +// data.insert(observable.clone(), phases.mean()); +// } else { +// let mut map: HashMap = HashMap::new(); +// map.insert(observable.clone(), phases.mean()); +// mean.insert(sv, map); +// } +// /* +// * average PR values +// */ +// let pr: Vec = self.pseudo_range() +// .filter_map(|(_, svnn, obs, pr)| { +// if sv == svnn && observable == obs { +// Some(pr) +// } else { +// None +// } +// }) +// .collect(); +// if let Some(data) = mean.get_mut(&sv) { +// data.insert(observable.clone(), pr.mean()); +// } else { +// let mut map: HashMap = HashMap::new(); +// map.insert(observable.clone(), pr.mean()); +// mean.insert(sv, map); +// } +// } +// } +// HashMap::new() +// } +//} #[cfg(feature = "obs")] use observation::Combine; diff --git a/rinex/src/meteo/record.rs b/rinex/src/meteo/record.rs index 2df3cd840..7e15a611f 100644 --- a/rinex/src/meteo/record.rs +++ b/rinex/src/meteo/record.rs @@ -436,69 +436,3 @@ impl Interpolate for Record { unimplemented!("meteo:record:interpolate_mut()") } } - -#[cfg(feature = "obs")] -use crate::observation::{Observation, StatisticalOps}; - -#[cfg(feature = "obs")] -use statrs::statistics::Statistics; - -#[cfg(feature = "obs")] -fn statistical_estimate(rec: &Record, ops: StatisticalOps) -> HashMap { - let mut ret: HashMap = HashMap::new(); - let mut dataset: HashMap> = HashMap::new(); - for (_t, observables) in rec { - for (observable, point) in observables { - if let Some(data) = dataset.get_mut(&observable) { - data.push(*point); - } else { - dataset.insert(observable.clone(), vec![*point]); - } - } - } - for (observable, data) in dataset { - let stats = match ops { - StatisticalOps::Min => data.min(), - StatisticalOps::Max => data.max(), - StatisticalOps::Mean => data.mean(), - StatisticalOps::StdDev => data.std_dev(), - StatisticalOps::StdVar => data.variance(), - }; - ret.insert(observable.clone(), stats); - } - ret -} - -#[cfg(feature = "obs")] -impl Observation for Record { - fn min(&self) -> (Option, HashMap>) { - unimplemented!("only on OBS RINEX!"); - } - fn max(&self) -> (Option, HashMap>) { - unimplemented!("only on OBS RINEX!"); - } - fn mean(&self) -> (Option, HashMap>) { - unimplemented!("only on OBS RINEX!"); - } - fn std_var(&self) -> (Option, HashMap>) { - unimplemented!("only on OBS RINEX!"); - } - fn std_dev(&self) -> (Option, HashMap>) { - unimplemented!("only on OBS RINEX!"); - } - fn min_observable(&self) -> HashMap { - statistical_estimate(self, StatisticalOps::Min) - } - fn max_observable(&self) -> HashMap { - statistical_estimate(self, StatisticalOps::Max) - } - fn std_dev_observable(&self) -> HashMap { - statistical_estimate(self, StatisticalOps::StdDev) - } - fn std_var_observable(&self) -> HashMap { - statistical_estimate(self, StatisticalOps::StdVar) - } - fn mean_observable(&self) -> HashMap { - statistical_estimate(self, StatisticalOps::Mean) - } -} diff --git a/rinex/src/observation/mod.rs b/rinex/src/observation/mod.rs index d311885a1..3725d2acb 100644 --- a/rinex/src/observation/mod.rs +++ b/rinex/src/observation/mod.rs @@ -174,122 +174,9 @@ impl HeaderFields { } } -#[cfg(feature = "obs")] -#[derive(Debug, Clone, Copy)] -pub(crate) enum StatisticalOps { - Min, - Max, - Mean, - StdDev, - StdVar, -} - #[cfg(feature = "obs")] use std::collections::BTreeMap; -/// OBS RINEX specific analysis trait. -/// Include this trait to unlock Observation analysis, mainly statistical analysis. -#[cfg(feature = "obs")] -#[cfg_attr(docrs, doc(cfg(feature = "obs")))] -pub trait Observation { - /// Returns minimum value observed, throughout all epochs, sorted by Observable. - /// This also applies to clock receiver estimate, - /// when requested on OBS RINEX files, not METEO files. - /// ``` - /// use rinex::*; // prelude + macros - /// use rinex::prelude::*; - /// use std::str::FromStr; // observable! - /// use rinex::observation::Observation; // .min_observable() - /// - /// // OBS RINEX example - /// let rinex = Rinex::from_file("../test_resources/OBS/V3/DUTH0630.22O") - /// .unwrap(); - /// let min_values = rinex.min_observable(); - /// for (observable, min_value) in min_values { - /// if observable == observable!("S1C") { - /// // minimum signal strength for carrier 1 - /// assert_eq!(min_value, 37.75); // L1 carrier min (worst) RSSI - /// } - /// } - /// - /// // METEO RINEX example - /// let rinex = Rinex::from_file("../test_resources/MET/V2/clar0020.00m") - /// .unwrap(); - /// let min_values = rinex.min_observable(); - /// for (observable, min_value) in min_values { - /// if observable == Observable::Temperature { - /// assert_eq!(min_value, 8.4); // min value encountered on that day - /// } - /// } - /// ``` - fn min_observable(&self) -> HashMap; - - /// Returns maximal value observed, throughout all epochs, sorted by Observable. - /// See [Self::min_observable()] for API example. - fn max_observable(&self) -> HashMap; - - /// Returns mean observation, throughout all epochs, sorted by Observable. - /// See [Self::min_observable()] for API example. - fn mean_observable(&self) -> HashMap; - - /// Returns standard deviation for all Observables. - /// See [Self::min_observable()] for API example. - fn std_dev_observable(&self) -> HashMap; - - /// Returns standard variance for all Observables. - /// See [Self::min_observable()] for API example. - fn std_var_observable(&self) -> HashMap; - - /// Returns minimum value observed throughout all epochs sorted by - /// Satellite vehicle and Observable. This does not apply to METEO - /// RINEX files. - /// ``` - /// use rinex::*; - /// use rinex::prelude::*; // basics - /// use std::str::FromStr; // sv!, observable! - /// use rinex::observation::Observation; // .min() - /// - /// // OBS RINEX example - /// let rinex = Rinex::from_file("../test_resources/OBS/V3/DUTH0630.22O") - /// .unwrap(); - /// - /// let (min_clock, min_values) = rinex.min(); - /// assert!(min_clock.is_none()); // we don't have an example file with such information yet - /// - /// for (sv, observables) in min_values { - /// if sv == sv!("G08") { - /// for (observable, min_value) in observables { - /// if observable == observable!("S1C") { - /// // minimum signal strength for carrier 1 - /// // for that particular vehicle - /// } - /// } - /// } - /// } - /// ``` - fn min(&self) -> (Option, HashMap>); - - /// Returns maximal value observed throughout all epochs sorted by - /// Satellite vehicle and Observable. This does not apply to METEO - /// RINEX files. See [Self::min()] for API example. - fn max(&self) -> (Option, HashMap>); - - /// Returns mean value observed throughout all epochs sorted by - /// Satellite vehicle and Observable. This does not apply to METEO - /// RINEX files. See [Self::min()] for API example. - fn mean(&self) -> (Option, HashMap>); - - /// Returns observations deviation throughout all epochs sorted by - /// Satellite vehicle and Observable. This does not apply to METEO - /// RINEX files. See [Self::min()] for API example. - fn std_dev(&self) -> (Option, HashMap>); - - /// Returns observations variance throughout all epochs sorted by - /// Satellite vehicle and Observable. This does not apply to METEO - /// RINEX files. See [Self::min()] for API example. - fn std_var(&self) -> (Option, HashMap>); -} - /// GNSS signal recombination trait. /// Import this to recombine OBS RINEX with usual recombination methods. /// This only applies to OBS RINEX records. diff --git a/rinex/src/observation/record.rs b/rinex/src/observation/record.rs index 09c200484..4e7c88aab 100644 --- a/rinex/src/observation/record.rs +++ b/rinex/src/observation/record.rs @@ -1282,167 +1282,6 @@ impl Decimate for Record { } } -#[cfg(feature = "obs")] -use statrs::statistics::Statistics; - -#[cfg(feature = "obs")] -use crate::observation::{Observation, StatisticalOps}; - -#[cfg(feature = "obs")] -/* - * Evaluate a specific statistical estimate on this record data set - */ -fn statistical_estimate( - rec: &Record, - ops: StatisticalOps, -) -> (Option, HashMap>) { - let mut ret: (Option, HashMap>) = (None, HashMap::new()); - - // Vectorize clock offset (in time), so we can invoke statrs() on it - let data: Vec = rec - .iter() - .filter_map(|(_, (clk, _))| { - if let Some(clk) = clk { - Some(*clk) - } else { - None - } - }) - .collect(); - - if data.len() > 0 { - // data exists - match ops { - StatisticalOps::Min => ret.0 = Some(data.min()), - StatisticalOps::Max => ret.0 = Some(data.max()), - StatisticalOps::Mean => ret.0 = Some(data.mean()), - StatisticalOps::StdDev => ret.0 = Some(data.std_dev()), - StatisticalOps::StdVar => ret.0 = Some(data.variance()), - } - } - - // Vectorize data (in time) so we can invoke statrs() on it - let mut data: HashMap>> = HashMap::new(); - - for (_epoch, (_clk, svnn)) in rec { - for (sv, observations) in svnn { - for (observable, obs_data) in observations { - if let Some(data) = data.get_mut(&sv) { - if let Some(data) = data.get_mut(&observable) { - data.push(obs_data.obs); // append - } else { - data.insert(observable.clone(), vec![obs_data.obs]); - } - } else { - let mut map: HashMap> = HashMap::new(); - map.insert(observable.clone(), vec![obs_data.obs]); - data.insert(sv.clone(), map); - } - } - } - } - - // build returned data - for (sv, observables) in data { - for (observable, data) in observables { - // invoke statrs - let stats = match ops { - StatisticalOps::Min => data.min(), - StatisticalOps::Max => data.max(), - StatisticalOps::Mean => data.mean(), - StatisticalOps::StdDev => data.std_dev(), - StatisticalOps::StdVar => data.variance(), - }; - if let Some(data) = ret.1.get_mut(&sv) { - data.insert(observable.clone(), stats); - } else { - let mut map: HashMap = HashMap::new(); - map.insert(observable.clone(), stats); - ret.1.insert(sv.clone(), map); - } - } - } - - ret -} - -#[cfg(feature = "obs")] -/* - * Evaluate a specific statistical estimate on this record data set - */ -fn statistical_observable_estimate(rec: &Record, ops: StatisticalOps) -> HashMap { - let mut ret: HashMap = HashMap::new(); - // For all statistical estimates but StdDev/StdVar - // We can compute across all vehicles then shrink 1D - match ops { - StatisticalOps::StdVar | StatisticalOps::StdDev => {}, - _ => { - // compute across all observables, then per Sv - let (_, data) = statistical_estimate(rec, ops); - let mut data_set: HashMap> = HashMap::new(); - for (_sv, observables) in data { - for (observable, value) in observables { - if let Some(data) = data_set.get_mut(&observable) { - data.push(value); - } else { - data_set.insert(observable.clone(), vec![value]); - } - } - } - for (observable, values) in data_set { - // invoke statrs - match ops { - StatisticalOps::Min => { - ret.insert(observable.clone(), values.min()); - }, - StatisticalOps::Max => { - ret.insert(observable.clone(), values.max()); - }, - StatisticalOps::Mean => { - ret.insert(observable.clone(), values.mean()); - }, - _ => {}, // does not apply - }; - } - }, - } - ret -} - -#[cfg(feature = "obs")] -impl Observation for Record { - fn min(&self) -> (Option, HashMap>) { - statistical_estimate(self, StatisticalOps::Min) - } - fn max(&self) -> (Option, HashMap>) { - statistical_estimate(self, StatisticalOps::Max) - } - fn mean(&self) -> (Option, HashMap>) { - statistical_estimate(self, StatisticalOps::Mean) - } - fn std_var(&self) -> (Option, HashMap>) { - statistical_estimate(self, StatisticalOps::StdVar) - } - fn std_dev(&self) -> (Option, HashMap>) { - statistical_estimate(self, StatisticalOps::StdDev) - } - fn min_observable(&self) -> HashMap { - statistical_observable_estimate(self, StatisticalOps::Min) - } - fn max_observable(&self) -> HashMap { - statistical_observable_estimate(self, StatisticalOps::Max) - } - fn std_dev_observable(&self) -> HashMap { - statistical_observable_estimate(self, StatisticalOps::StdDev) - } - fn std_var_observable(&self) -> HashMap { - statistical_observable_estimate(self, StatisticalOps::StdVar) - } - fn mean_observable(&self) -> HashMap { - statistical_observable_estimate(self, StatisticalOps::Mean) - } -} - #[cfg(feature = "obs")] use crate::observation::Combine; @@ -1834,7 +1673,7 @@ impl Combine for Record { #[cfg(feature = "obs")] use crate::{ carrier, - observation::{Dcb, Mp}, + observation::Dcb, //Mp}, }; #[cfg(feature = "obs")] @@ -2000,148 +1839,6 @@ impl IonoDelay for Record { } } -#[cfg(feature = "obs")] -impl Mp for Record { - fn mp(&self) -> HashMap>> { - let mut ret: HashMap>> = - HashMap::new(); - /* - * Determine mean value of all datasets - */ - let (_, mean) = self.mean(); // drop mean{clk} - //println!("MEAN VALUES {:?}", mean); //DEBUG - /* - * Run algorithm - */ - let mut associated: HashMap = HashMap::new(); // Ph code to associate to this Mpx - // for operation consistency - for (epoch, (_, vehicles)) in self { - 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 = None; - let mut ph_j: Option = 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: BTreeMap> = - BTreeMap::new(); - bmap.insert(*epoch, mp); - map.insert(*sv, bmap); - ret.insert(mp_code.to_string(), map); - } - } - } - } - } - ret - } -} - #[cfg(test)] mod test { use super::*; diff --git a/rinex/src/tests/mod.rs b/rinex/src/tests/mod.rs index a92ab8e64..dacb79326 100644 --- a/rinex/src/tests/mod.rs +++ b/rinex/src/tests/mod.rs @@ -17,4 +17,3 @@ mod parsing; mod production; mod sampling; mod smoothing; -mod stats; diff --git a/rinex/src/tests/stats.rs b/rinex/src/tests/stats.rs deleted file mode 100644 index 5de0dfdbc..000000000 --- a/rinex/src/tests/stats.rs +++ /dev/null @@ -1,347 +0,0 @@ -#[cfg(test)] -mod test { - use crate::*; - use crate::{observation::*, prelude::*}; - use std::str::FromStr; - fn run_test(rinex: &Rinex, ops: &str, expected: Vec<(Sv, Vec<(Observable, f64)>)>) { - let (clk_stats, stats) = match ops { - "min" => rinex.min(), - "max" => rinex.max(), - "mean" => rinex.mean(), - _ => panic!("unknown test"), - }; - assert!(clk_stats.is_none()); // always none on our test pool - - for (sv, fields) in expected { - let sv_data = stats.get(&sv); - assert!( - sv_data.is_some(), - "{}", - format!("\"{}\" stats missing for vehicle \"{}\"", ops, sv) - ); - let sv_data = sv_data.unwrap(); - for (observable, expected) in fields { - let data = sv_data.get(&observable); - assert!( - data.is_some(), - "{}", - format!( - "\"{}\" stats missing for vehicle \"{}\" - \"{}\"", - ops, sv, observable - ) - ); - let data = data.unwrap(); - let e = *data - expected; - assert!( - e.abs() < 1E-5, - "{}", - format!( - "\"{}\" test failed for \"{}\" - \"{}\", expecting {} got {}", - ops, sv, observable, expected, *data - ) - ); - } - } - } - fn run_obsv_test(rinex: &Rinex, ops: &str, expected: Vec<(Observable, f64)>) { - let stats = match ops { - "min" => rinex.min_observable(), - "max" => rinex.max_observable(), - "mean" => rinex.mean_observable(), - _ => panic!("unknown test"), - }; - - for (observable, expected) in expected { - let results = stats.get(&observable); - assert!( - results.is_some(), - "{}", - format!( - "\"{}\" observable stats results missing for observable \"{}\"", - ops, observable - ) - ); - let results = results.unwrap(); - let err = *results - expected; - assert!( - err.abs() < 1E-5, - "{}", - format!( - "\"{}\" observable test failed for - \"{}\", expecting {} - got {}", - ops, observable, expected, results - ) - ); - } - } - #[test] - fn stats() { - let rinex = Rinex::from_file("../test_resources/OBS/V3/DUTH0630.22O").unwrap(); - - let expected: Vec<(Sv, Vec<(Observable, f64)>)> = vec![ - ( - sv!("g01"), - vec![ - ( - observable!("c1c"), - (20243517.560 + 20805393.080 + 21653418.260) / 3.0, - ), - (observable!("d1c"), (-1242.766 - 2193.055 - 2985.516) / 3.0), - (observable!("s1c"), (51.25 + 50.75 + 49.5) / 3.0), - ], - ), - ( - sv!("g03"), - vec![ - ( - observable!("c1c"), - (20619020.680 + 20425456.580 + 20410261.460) / 3.0, - ), - (observable!("d1c"), (852.785 + 328.797 - 244.031) / 3.0), - (observable!("s1c"), (50.75 + 51.0 + 51.0) / 3.0), - ], - ), - ]; - run_test(&rinex, "mean", expected); - - let expected: Vec<(Sv, Vec<(Observable, f64)>)> = vec![ - ( - sv!("g01"), - vec![ - (observable!("c1c"), 20243517.560), - (observable!("d1c"), -2985.516), - (observable!("s1c"), 49.5), - ], - ), - ( - sv!("g03"), - vec![ - (observable!("c1c"), 20410261.460), - (observable!("d1c"), -244.031), - (observable!("s1c"), 50.75), - ], - ), - ]; - run_test(&rinex, "min", expected); - - let expected: Vec<(Sv, Vec<(Observable, f64)>)> = vec![ - ( - sv!("g01"), - vec![ - (observable!("c1c"), 21653418.260), - (observable!("d1c"), -1242.766), - (observable!("s1c"), 51.25), - ], - ), - ( - sv!("g03"), - vec![ - (observable!("c1c"), 20619020.680), - (observable!("d1c"), 852.785), - (observable!("s1c"), 51.0), - ], - ), - ]; - run_test(&rinex, "max", expected); - - /* ********************************************* - * STATS FOR ALL SV / ALL OBS : meteo compatible - **********************************************/ - /*let expected : Vec<(Observable, f64)> = vec![ - (observable!("c1c"), (20243517.560 + 20619020.680 + 21542633.500 + 24438727.980 + 22978068.560 + 23460759.840 + 21923317.180 + 23434790.440 + 22401985.340 + 24991723.280 + 19727826.340 + 23171275.620 + 20662538.580 + 23450513.820 + 23044984.180 + 22909354.040 + 20116780.920 + 19708379.260 + 20805393.080 + 20425456.580 + 20887001.400 + 23371156.300 + 23031543.660 + 23117350.280 + 22726604.680 + 24425563.640 + 22689941.780 + 19677287.000 + 22265147.080 + 21462395.740 + 23740237.340 + 22432243.520 + 21750541.080 + 21199384.320 + 19680274.400 + 21653418.260 + 20410261.460 + 20488105.720 + 23647940.540 + 22436978.380 + 23392660.200 + 23154069.760 + 23689895.760 + 25161827.280 + 23333751.720 + 19831816.600 + 21490078.220 + 22328018.960 + 22235350.560 + 20915624.780 + 22543866.020 + 20147683.700) / 52.0), - ]; - run_obsv_test(&rinex, "mean", expected);*/ - - /* - * Test on METEO RINEX - */ - let rinex = Rinex::from_file("../test_resources/MET/V2/clar0020.00m").unwrap(); - let expected: Vec<(Observable, f64)> = vec![( - observable!("PR"), - (970.5 - + 970.4 - + 970.4 - + 970.3 - + 970.2 - + 970.3 - + 970.1 - + 970.3 - + 970.2 - + 970.2 - + 972.1 - + 972.3 - + 972.4 - + 972.6 - + 972.8 - + 973.1 - + 973.3 - + 973.4 - + 973.6 - + 973.6 - + 973.7 - + 973.7 - + 973.5 - + 973.6 - + 973.6 - + 973.6 - + 973.6 - + 973.5 - + 973.5 - + 973.4 - + 973.3 - + 973.2 - + 972.9 - + 972.8 - + 972.7 - + 972.5 - + 972.4 - + 972.4 - + 972.2 - + 972.2 - + 972.2 - + 972.3 - + 972.2 - + 972.2 - + 972.1 - + 972.2 - + 972.1 - + 972.2 - + 972.2 - + 972.2 - + 972.1 - + 972.2 - + 972.2 - + 972.3 - + 972.3 - + 972.5 - + 972.5) - / 57.0, - )]; - run_obsv_test(&rinex, "mean", expected); - - let expected: Vec<(Observable, f64)> = vec![( - observable!("TD"), - (10.7 - + 10.6 - + 10.4 - + 10.3 - + 10.1 - + 9.9 - + 9.8 - + 9.5 - + 9.6 - + 9.6 - + 8.4 - + 10.1 - + 10.1 - + 10.3 - + 10.2 - + 10.3 - + 10.6 - + 10.8 - + 11.1 - + 11.2 - + 11.4 - + 11.5 - + 11.8 - + 11.8 - + 11.8 - + 12.7 - + 12.7 - + 13.2 - + 13.5 - + 13.0 - + 13.5 - + 13.9 - + 14.0 - + 14.7 - + 15.0 - + 14.8 - + 15.6 - + 15.1 - + 14.7 - + 14.8 - + 14.9 - + 14.9 - + 15.9 - + 15.6 - + 16.2 - + 15.6 - + 15.5 - + 15.9 - + 15.8 - + 15.7 - + 15.7 - + 15.6 - + 15.4 - + 14.9 - + 14.6 - + 14.3 - + 14.2) - / 57.0, - )]; - run_obsv_test(&rinex, "mean", expected); - - let expected: Vec<(Observable, f64)> = vec![( - observable!("HR"), - (71.4 - + 72.2 - + 72.9 - + 72.9 - + 75.0 - + 76.7 - + 78.1 - + 79.7 - + 80.5 - + 79.1 - + 70.7 - + 65.3 - + 66.4 - + 68.0 - + 68.4 - + 68.6 - + 66.8 - + 65.4 - + 63.3 - + 62.7 - + 60.0 - + 57.2 - + 56.0 - + 54.9 - + 52.2 - + 49.3 - + 46.3 - + 43.8 - + 44.8 - + 42.7 - + 41.2 - + 40.2 - + 37.9 - + 36.4 - + 37.3 - + 36.1 - + 32.7 - + 31.6 - + 36.1 - + 35.7 - + 34.0 - + 32.6 - + 29.5 - + 28.0 - + 26.7 - + 25.5 - + 25.0 - + 30.7 - + 25.0 - + 23.5 - + 27.7 - + 28.2 - + 28.4 - + 30.9 - + 33.3 - + 33.3 - + 33.2) - / 57.0, - )]; - run_obsv_test(&rinex, "mean", expected); - } -}