Skip to content

Commit

Permalink
Merge QZSS into tightly (GPS+SBAS)
Browse files Browse the repository at this point in the history
  • Loading branch information
fenrir-naru committed Mar 28, 2024
2 parents 60a76fa + 75638dd commit e364460
Show file tree
Hide file tree
Showing 7 changed files with 413 additions and 31 deletions.
24 changes: 21 additions & 3 deletions tool/INS_GPS/GNSS_Data.h
Original file line number Diff line number Diff line change
Expand Up @@ -41,6 +41,7 @@

#include "navigation/GPS.h"
#include "navigation/SBAS.h"
#include "navigation/QZSS.h"

template <class FloatT>
struct GNSS_Data {
Expand Down Expand Up @@ -69,17 +70,29 @@ struct GNSS_Data {
struct gps_ephemeris_raw_t : public gps_ephemeris_t::raw_t {
bool set_iodc;
int iode_subframe2, iode_subframe3;
bool is_qzss;
typedef typename gps_ephemeris_t::raw_t super_t;
gps_ephemeris_raw_t()
: super_t(), set_iodc(false), iode_subframe2(-1), iode_subframe3(-1) {}
} gps_ephemeris[32];
: super_t(), set_iodc(false), iode_subframe2(-1), iode_subframe3(-1), is_qzss(false) {}
operator gps_ephemeris_t() const {
return is_qzss
? (gps_ephemeris_t)(
reinterpret_cast<
const typename QZSS_SpaceNode<FloatT>::SatelliteProperties::Ephemeris::raw_t &>(*this))
: super_t::operator gps_ephemeris_t();
}
} gps_ephemeris[32], qzss_ephemeris[10];

typedef typename gps_t::Ionospheric_UTC_Parameters gps_iono_utc_t;

Loader() : gps(NULL), sbas(NULL) {
for(unsigned int i(0); i < sizeof(gps_ephemeris) / sizeof(gps_ephemeris[0]); ++i){
gps_ephemeris[i].svid = i + 1;
}
for(unsigned int i(0); i < sizeof(qzss_ephemeris) / sizeof(qzss_ephemeris[0]); ++i){
qzss_ephemeris[i].svid = i + 193;
qzss_ephemeris[i].is_qzss = true;
}
}

/**
Expand Down Expand Up @@ -128,7 +141,10 @@ struct GNSS_Data {
int subframe_no(parser_t::subframe_id(buf));

if(subframe_no <= 3){
gps_ephemeris_raw_t &eph(gps_ephemeris[data.subframe.sv_number - 1]);
gps_ephemeris_raw_t &eph(
(data.subframe.gnssID == observer_t::gnss_svid_t::QZSS)
? qzss_ephemeris[data.subframe.sv_number - 1]
: gps_ephemeris[data.subframe.sv_number - 1]);

switch(subframe_no){
case 1: eph.template update_subframe1<2, 0>(buf); eph.set_iodc = true; break;
Expand Down Expand Up @@ -184,6 +200,8 @@ struct GNSS_Data {
// TODO: other satellite systems
case observer_t::gnss_svid_t::SBAS:
return load_sbas(data);
case observer_t::gnss_svid_t::QZSS:
if(data.subframe.bytes != sizeof(data.subframe.buffer)){return false;}
case observer_t::gnss_svid_t::GPS:
return load_gps(data);
}
Expand Down
73 changes: 65 additions & 8 deletions tool/INS_GPS/GNSS_Receiver.h
Original file line number Diff line number Diff line change
Expand Up @@ -118,6 +118,7 @@ struct GNSS_Receiver {
if(out_rinex_nav){
typename RINEX_NAV_Writer<FloatT>::space_node_list_t list = {&gps.space_node};
list.sbas = &sbas.space_node;
list.qzss = &gps.space_node;
RINEX_NAV_Writer<FloatT>::write_all(*out_rinex_nav, list);
}
}
Expand All @@ -130,10 +131,32 @@ struct GNSS_Receiver {
struct measurement_items_t : public gps_solver_t::measurement_items_t {
// TODO
};
struct ephemeris_proxy_t {
struct item_t {
const void *impl;
typename base_t::satellite_t (*impl_select)(
const void *,
const typename base_t::prn_t &, const typename base_t::gps_time_t &);
} gps, qzss;
static typename base_t::satellite_t forward(
const void *ptr,
const typename base_t::prn_t &prn, const typename base_t::gps_time_t &t){
const ephemeris_proxy_t *proxy(static_cast<const ephemeris_proxy_t *>(ptr));
const item_t &target(((prn >= 193) && (prn <= 202)) ? proxy->qzss : proxy->gps);
return target.impl_select(target.impl, prn, t);
}
ephemeris_proxy_t(gps_solver_t &solver){
gps.impl = qzss.impl = solver.satellites.impl;
gps.impl_select = qzss.impl_select = solver.satellites.impl_select;
solver.satellites.impl = this;
solver.satellites.impl_select = forward;
}
} ephemeris_proxy;
solver_t(const GNSS_Receiver &rcv)
: base_t(),
gps(rcv.data.gps.space_node),
sbas(rcv.data.sbas.space_node) {
sbas(rcv.data.sbas.space_node),
ephemeris_proxy(gps) {
#if !defined(BUILD_WITHOUT_GNSS_MULTI_CONSTELLATION)
// ionospheric and tropospheric correction methods
typename base_t::range_correction_t ionospheric, tropospheric;
Expand All @@ -153,8 +176,9 @@ struct GNSS_Receiver {
// Proxy functions
const base_t &select(const typename base_t::prn_t &serial) const {
switch(system_t::serial2system(serial)){
case system_t::GPS: return gps;
case system_t::SBAS: return sbas;
case system_t::GPS: return gps;
case system_t::SBAS: return sbas;
case system_t::QZSS: return gps;
}
return *this;
}
Expand Down Expand Up @@ -200,14 +224,16 @@ struct GNSS_Receiver {
void update_ephemeris_source(const data_t &data){
#if !defined(BUILD_WITHOUT_SP3)
typename data_t::sp3_t::satellite_count_t cnt(data.sp3.satellite_count());
if(cnt.gps > 0){data.sp3.push(gps.satellites, SP3_Product<FloatT>::SYSTEM_GPS);}
if(cnt.gps > 0){data.sp3.push(ephemeris_proxy.gps, SP3_Product<FloatT>::SYSTEM_GPS);}
if(cnt.sbas > 0){data.sp3.push(sbas.satellites, SP3_Product<FloatT>::SYSTEM_SBAS);}
if(cnt.qzss > 0){data.sp3.push(ephemeris_proxy.qzss, SP3_Product<FloatT>::SYSTEM_QZSS);}
#endif
#if !defined(BUILD_WITHOUT_RINEX_CLK)
// RINEX clock has higher priority to be applied than SP3
typename data_t::clk_t::count_t cnt2(data.clk.count());
if(cnt2.gps > 0){data.clk.push(gps.satellites, data_t::clk_t::SYSTEM_GPS);}
if(cnt2.gps > 0){data.clk.push(ephemeris_proxy.gps, data_t::clk_t::SYSTEM_GPS);}
if(cnt2.sbas > 0){data.clk.push(sbas.satellites, data_t::clk_t::SYSTEM_SBAS);}
if(cnt2.qzss > 0){data.clk.push(ephemeris_proxy.qzss, data_t::clk_t::SYSTEM_QZSS);}
#endif
}
} solver_GNSS;
Expand Down Expand Up @@ -388,6 +414,7 @@ struct GNSS_Receiver {
std::istream &in(options.spec2istream(value));
typename RINEX_NAV_Reader<FloatT>::space_node_list_t list = {&data.gps.space_node};
list.sbas = &data.sbas.space_node;
list.qzss = &data.gps.space_node;
int ephemeris(RINEX_NAV_Reader<FloatT>::read_all(in, list));
if(ephemeris < 0){
std::cerr << "(error!) Invalid format!" << std::endl;
Expand Down Expand Up @@ -418,6 +445,7 @@ struct GNSS_Receiver {
typename data_t::sp3_t::satellite_count_t cnt(data.sp3.satellite_count());
if(cnt.gps > 0){std::cerr << "SP3 GPS satellites: " << cnt.gps << std::endl;}
if(cnt.sbas > 0){std::cerr << "SP3 SBAS satellites: " << cnt.sbas << std::endl;}
if(cnt.qzss > 0){std::cerr << "SP3 QZSS satellites: " << cnt.qzss << std::endl;}
}
solver_GNSS.update_ephemeris_source(data);
return true;
Expand Down Expand Up @@ -454,6 +482,7 @@ struct GNSS_Receiver {
typename data_t::clk_t::count_t cnt(data.clk.count());
if(cnt.gps > 0){std::cerr << "RINEX clock GPS satellites: " << cnt.gps << std::endl;}
if(cnt.sbas > 0){std::cerr << "RINEX clock SBAS satellites: " << cnt.sbas << std::endl;}
if(cnt.qzss > 0){std::cerr << "RINEX clock QZSS satellites: " << cnt.qzss << std::endl;}
}
solver_GNSS.update_ephemeris_source(data);
return true;
Expand Down Expand Up @@ -549,9 +578,16 @@ data.sbas.solver_options. expr

switch(sys){
case system_t::GPS:
select_all
? data.gps.solver_options.exclude_prn.set(without)
: data.gps.solver_options.exclude_prn.set(sv_id, without);
case system_t::QZSS:
if(select_all){
static const int id_table[][2] = {{1, 32}, {193, 202}};
int idx(sys == system_t::QZSS ? 1 : 0);
for(int i(id_table[idx][0]), i_max(id_table[idx][1]); i <= i_max; ++i){
data.gps.solver_options.exclude_prn.set(i, without);
}
}else{
data.gps.solver_options.exclude_prn.set(sv_id, without);
}
break;
case system_t::SBAS:
select_all
Expand Down Expand Up @@ -611,6 +647,7 @@ data.sbas.solver_options. expr
<< ',' << "GPS_PRN(1-32)"
#if !defined(BUILD_WITHOUT_GNSS_MULTI_CONSTELLATION)
<< ',' << "SBAS_PRN(120-158)"
<< ',' << "QZSS_PRN(193-202)"
#endif
;
}
Expand All @@ -626,6 +663,10 @@ data.sbas.solver_options. expr
out << ',' << "range_residual(SBAS:" << i << ')'
<< ',' << "weight(SBAS:" << i << ')';
}
for(int i(193); i <= 202; ++i){
out << ',' << "range_residual(QZSS:" << i << ')'
<< ',' << "weight(QZSS:" << i << ')';
}
#endif
return out;
}
Expand Down Expand Up @@ -704,6 +745,7 @@ data.sbas.solver_options. expr
<< ',' << mask_printer_t(src.used_satellite_mask, 1, 32)
#if !defined(BUILD_WITHOUT_GNSS_MULTI_CONSTELLATION)
<< ',' << mask_printer_t(src.used_satellite_mask, 120, 158) // SBAS
<< ',' << mask_printer_t(src.used_satellite_mask, 193, 202) // QZSS
#endif
;
}else{
Expand All @@ -727,6 +769,7 @@ data.sbas.solver_options. expr
{1, 32}, // GPS
#if !defined(BUILD_WITHOUT_GNSS_MULTI_CONSTELLATION)
{120, 158}, // SBAS
{193, 202}, // QZSS
#endif
};
unsigned int i_row(0);
Expand Down Expand Up @@ -803,6 +846,12 @@ data.sbas.solver_options. expr
<< ',' << "azimuth(SBAS:" << (*it)->prn << ')'
<< ',' << "elevation(SBAS:" << (*it)->prn << ')';
}
for(int i(193); i <= 202; ++i){
out << ',' << "L1_range(QZSS:" << i << ')'
<< ',' << "L1_rate(QZSS:" << i << ')'
<< ',' << "azimuth(QZSS:" << i << ')'
<< ',' << "elevation(QZSS:" << i << ')';
}
#endif
return out;
}
Expand Down Expand Up @@ -902,6 +951,14 @@ data.sbas.solver_options. expr
it != it_end; ++it){
out << ',' << p((*it)->prn, cmd_sbas) << ',' << p.az_el((*it)->prn);
}
static const cmd_t cmd_qzss[] = {
{items_t::L1_PSEUDORANGE, true}, // range
{items_t::L1_RANGE_RATE, false}, // rate
{items_t::L1_DOPPLER, false, print_t::l1_doppler2rate}, // fallback to using doppler
};
for(int i(193); i <= 202; ++i){
out << ',' << p(i, cmd_qzss) << ',' << p.az_el(i);
}
#endif
return out;
}
Expand Down
16 changes: 10 additions & 6 deletions tool/misc/receiver_debug.rb
Original file line number Diff line number Diff line change
Expand Up @@ -170,8 +170,8 @@ def initialize(options = {})
opt.residual_mask = 1E4 # 10 km (without residual filter, practically)
}
output_options = {
:system => [[:GPS, 1..32]],
:satellites => (1..32).to_a, # [idx, ...] or [[idx, label], ...] is acceptable
:system => [[:GPS, 1..32], [:QZSS, 193..202]],
:satellites => (1..32).to_a + (193..202).to_a, # [idx, ...] or [[idx, label], ...] is acceptable
:FDE => false,
}
options = options.reject{|k, v|
Expand Down Expand Up @@ -280,6 +280,8 @@ def v.to_b; !(self =~ /^(?:false|0|f|off)$/i); end
prns = [svid || (120..158).to_a].flatten
update_output.call(:SBAS, prns)
prns.each{|prn| @solver.sbas_options.send(mode, prn)}
elsif check_sys_svid.call(:QZSS, 193..202) then
[svid || (193..202).to_a].flatten.each{|prn| @solver.gps_options.send(mode, prn)}
else
raise "Unknown satellite: #{spec}"
end
Expand Down Expand Up @@ -407,13 +409,13 @@ def run(meas, t_meas, ref_pos = @base_station)
}

def register_ephemeris(t_meas, sys, prn, bcast_data)
@eph_list ||= Hash[*(1..32).collect{|prn|
@eph_list ||= Hash[*((1..32).to_a + (193..202).to_a).collect{|prn|
eph = GPS::Ephemeris::new
eph.svid = prn
[prn, eph]
}.flatten(1)]
case sys
when :GPS
when :GPS, :QZSS
return unless eph = @eph_list[prn]
sn = @solver.gps_space_node
subframe, iodc_or_iode = eph.parse(bcast_data)
Expand Down Expand Up @@ -521,10 +523,12 @@ def parse_ubx(ubx_fname, &b)
sys, svid = gnss_serial.call(*loader.call(36, 2).reverse)
sigid = (packet[6 + 13] != 0) ? loader.call(38, 1, "C") : 0 # sigID if version(>0); @see UBX-18010854
case sys
when :GPS
when :GPS
sigid = {0 => :L1, 3 => :L2CL, 4 => :L2CM}[sigid]
when :SBAS
sigid = :L1
when :QZSS
sigid = {0 => :L1, 5 => :L2CL, 4 => :L2CM}[sigid]
else; next
end
next unless sigid
Expand Down Expand Up @@ -726,7 +730,7 @@ def attach_rinex_clk(fname)
files.collect!{|fname, ftype|
raise "File not found: #{fname}" unless File::exist?(fname)
ftype ||= case fname
when /\.\d{2}[nh]$/; :rinex_nav
when /\.\d{2}[nhq]$/; :rinex_nav
when /\.\d{2}o$/; :rinex_obs
when /\.ubx$/; :ubx
when /\.sp3$/; :sp3
Expand Down
95 changes: 95 additions & 0 deletions tool/navigation/QZSS.h
Original file line number Diff line number Diff line change
@@ -0,0 +1,95 @@
/*
* Copyright (c) 2021, M.Naruoka (fenrir)
* All rights reserved.
*
* Redistribution and use in source and binary forms, with or without modification,
* are permitted provided that the following conditions are met:
*
* - Redistributions of source code must retain the above copyright notice,
* this list of conditions and the following disclaimer.
* - Redistributions in binary form must reproduce the above copyright notice,
* this list of conditions and the following disclaimer in the documentation
* and/or other materials provided with the distribution.
* - Neither the name of the naruoka.org nor the names of its contributors
* may be used to endorse or promote products derived from this software
* without specific prior written permission.
*
* THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
* AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO,
* THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
* ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS
* BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY,
* OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
* PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
* LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
* HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT,
* STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
* ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE,
* EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
*
*/

#ifndef __QZSS_H__
#define __QZSS_H__

/** @file
* @brief QZSS ICD definitions
*/

#include "GPS.h"
#include <limits>

template <class FloatT>
struct QZSS_SpaceNode {
struct SatelliteProperties {
struct Ephemeris {
typedef typename GPS_SpaceNode<FloatT>::SatelliteProperties::Ephemeris eph_gps_t;
struct raw_t : public eph_gps_t::raw_t {
typedef typename eph_gps_t::raw_t super_t;

operator eph_gps_t() const {
eph_gps_t res(super_t::operator eph_gps_t());
res.fit_interval = super_t::fit_interval_flag
? (std::numeric_limits<FloatT>::max)()
: (2 * 60 * 60); // Fit interval (ICD:4.1.2.4 Subframe 2); should be zero
return res;
};

raw_t &operator=(const eph_gps_t &eph){
super_t::operator=(eph);
super_t::fit_interval_flag = (eph.fit_interval > 2 * 60 * 60);
return *this;
}
};
};
struct Almanac {
typedef typename GPS_SpaceNode<FloatT>::SatelliteProperties::Almanac alm_gps_t;
struct raw_t : public alm_gps_t::raw_t {
typedef typename alm_gps_t::raw_t super_t;
template <int PaddingBits_MSB, int PaddingBits_LSB, class InputT>
void update(const InputT *src){
super_t::template update<PaddingBits_MSB, PaddingBits_LSB, InputT>(src);
if((super_t::svid >= 1) && (super_t::svid <= 10)){
super_t::svid += 192;
}
}
template <int PaddingBits_MSB, int PaddingBits_LSB, class BufferT>
void dump(BufferT *dst, const typename GPS_SpaceNode<FloatT>::u8_t &sf = 4){
super_t::template dump<PaddingBits_MSB, PaddingBits_LSB, BufferT>(dst);
typedef typename GPS_SpaceNode<FloatT>::template BroadcastedMessage<
BufferT, (int)sizeof(BufferT) * CHAR_BIT - PaddingBits_MSB - PaddingBits_LSB, PaddingBits_MSB>
deparse_t;
if((super_t::svid >= 193) && (super_t::svid <= 202)){
deparse_t::subframe_id_set(dst, sf); // sf = 4 or 5
deparse_t::sv_page_id_set(dst, super_t::svid - 192); // [1, 10]
}else{
return;
}
deparse_t::data_id_set(dst, 3); // always 3 @see 4.1.2.6.2. QZS Almanac
}
};
};
};
};

#endif /* __QZSS_H__ */
Loading

0 comments on commit e364460

Please sign in to comment.