1
0
mirror of https://github.com/gnss-sdr/gnss-sdr synced 2024-12-15 12:40:35 +00:00

Improving and extending GPS L1 CA observables unit test

This commit is contained in:
Javier Arribas 2018-07-26 19:25:10 +02:00
parent ab5f5689f0
commit 83413f2eaf
3 changed files with 360 additions and 258 deletions

View File

@ -0,0 +1,40 @@
/*!
* \file tracking_tests_flags.h
* \brief Helper file for unit testing
* \author Javier Arribas, 2018. jarribas(at)cttc.es
*
* -------------------------------------------------------------------------
*
* Copyright (C) 2010-2018 (see AUTHORS file for a list of contributors)
*
* GNSS-SDR is a software defined Global Navigation
* Satellite Systems receiver
*
* This file is part of GNSS-SDR.
*
* GNSS-SDR is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* GNSS-SDR is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with GNSS-SDR. If not, see <https://www.gnu.org/licenses/>.
*
* -------------------------------------------------------------------------
*/
#ifndef GNSS_SDR_OBSERVABLE_TESTS_FLAGS_H_
#define GNSS_SDR_OBSERVABLE_TESTS_FLAGS_H_
#include <gflags/gflags.h>
#include <limits>
DEFINE_double(skip_obs_transitory_s, 30.0, "Skip the initial observable outputs to avoid transitory results [s]");
#endif

View File

@ -45,6 +45,7 @@ DEFINE_string(filename_raw_data, "signal_out.bin", "Filename of output raw data
DEFINE_int32(fs_gen_sps, 2600000, "Sampling frequency [sps]"); DEFINE_int32(fs_gen_sps, 2600000, "Sampling frequency [sps]");
DEFINE_int32(test_satellite_PRN, 1, "PRN of the satellite under test (must be visible during the observation time)"); DEFINE_int32(test_satellite_PRN, 1, "PRN of the satellite under test (must be visible during the observation time)");
DEFINE_int32(test_satellite_PRN2, 2, "PRN of the satellite under test (must be visible during the observation time)"); DEFINE_int32(test_satellite_PRN2, 2, "PRN of the satellite under test (must be visible during the observation time)");
DEFINE_string(test_satellite_PRN_list, "1,2,3,6,9,10,12,17,20,23,28", "List of PRN of the satellites under test (must be visible during the observation time)");
DEFINE_double(CN0_dBHz, std::numeric_limits<double>::infinity(), "Enable noise generator and set the CN0 [dB-Hz]"); DEFINE_double(CN0_dBHz, std::numeric_limits<double>::infinity(), "Enable noise generator and set the CN0 [dB-Hz]");
#endif #endif

View File

@ -59,6 +59,7 @@
#include "gnss_sdr_sample_counter.h" #include "gnss_sdr_sample_counter.h"
#include <matio.h> #include <matio.h>
#include "test_flags.h" #include "test_flags.h"
#include "observable_tests_flags.h"
#include "gnuplot_i.h" #include "gnuplot_i.h"
@ -188,24 +189,27 @@ public:
bool save_mat_xy(std::vector<double>& x, std::vector<double>& y, std::string filename); bool save_mat_xy(std::vector<double>& x, std::vector<double>& y, std::string filename);
void check_results_carrier_phase( void check_results_carrier_phase(
arma::mat& true_ch0, arma::mat& true_ch0,
arma::mat& true_ch1,
arma::vec& true_tow_s, arma::vec& true_tow_s,
arma::mat& measured_ch0, arma::mat& measured_ch0,
arma::mat& measured_ch1); std::string data_title);
void check_results_code_psudorange( void check_results_carrier_doppler(
arma::mat& true_ch0,
arma::vec& true_tow_s,
arma::mat& measured_ch0,
std::string data_title);
void check_results_code_pseudorange(
arma::mat& true_ch0, arma::mat& true_ch0,
arma::mat& true_ch1, arma::mat& true_ch1,
arma::vec& true_tow_s, arma::vec& true_tow_s,
arma::mat& measured_ch0, arma::mat& measured_ch0,
arma::mat& measured_ch1); arma::mat& measured_ch1,
std::string data_title);
HybridObservablesTest() HybridObservablesTest()
{ {
factory = std::make_shared<GNSSBlockFactory>(); factory = std::make_shared<GNSSBlockFactory>();
config = std::make_shared<InMemoryConfiguration>(); config = std::make_shared<InMemoryConfiguration>();
item_size = sizeof(gr_complex); item_size = sizeof(gr_complex);
gnss_synchro_ch0 = Gnss_Synchro();
gnss_synchro_ch1 = Gnss_Synchro();
} }
~HybridObservablesTest() ~HybridObservablesTest()
@ -217,8 +221,7 @@ public:
gr::top_block_sptr top_block; gr::top_block_sptr top_block;
std::shared_ptr<GNSSBlockFactory> factory; std::shared_ptr<GNSSBlockFactory> factory;
std::shared_ptr<InMemoryConfiguration> config; std::shared_ptr<InMemoryConfiguration> config;
Gnss_Synchro gnss_synchro_ch0; std::vector<Gnss_Synchro> gnss_synchro_vec;
Gnss_Synchro gnss_synchro_ch1;
size_t item_size; size_t item_size;
}; };
@ -268,18 +271,6 @@ int HybridObservablesTest::generate_signal()
void HybridObservablesTest::configure_receiver() void HybridObservablesTest::configure_receiver()
{ {
gnss_synchro_ch0.Channel_ID = 0;
gnss_synchro_ch0.System = 'G';
std::string signal = "1C";
signal.copy(gnss_synchro_ch0.Signal, 2, 0);
gnss_synchro_ch0.PRN = FLAGS_test_satellite_PRN;
gnss_synchro_ch1.Channel_ID = 1;
gnss_synchro_ch1.System = 'G';
signal.copy(gnss_synchro_ch1.Signal, 2, 0);
gnss_synchro_ch1.PRN = FLAGS_test_satellite_PRN2;
config->set_property("GNSS-SDR.internal_fs_sps", std::to_string(baseband_sampling_freq)); config->set_property("GNSS-SDR.internal_fs_sps", std::to_string(baseband_sampling_freq));
// Set Tracking // Set Tracking
@ -290,7 +281,7 @@ void HybridObservablesTest::configure_receiver()
config->set_property("Tracking_1C.dll_bw_hz", "0.20"); config->set_property("Tracking_1C.dll_bw_hz", "0.20");
config->set_property("Tracking_1C.pll_bw_narrow_hz", "1.0"); config->set_property("Tracking_1C.pll_bw_narrow_hz", "1.0");
config->set_property("Tracking_1C.dll_bw_narrow_hz", "0.1"); config->set_property("Tracking_1C.dll_bw_narrow_hz", "0.1");
config->set_property("Tracking_1C.extend_correlation_symbols", "20"); config->set_property("Tracking_1C.extend_correlation_symbols", "1");
config->set_property("Tracking_1C.early_late_space_chips", "0.5"); config->set_property("Tracking_1C.early_late_space_chips", "0.5");
config->set_property("TelemetryDecoder_1C.dump", "true"); config->set_property("TelemetryDecoder_1C.dump", "true");
@ -299,36 +290,31 @@ void HybridObservablesTest::configure_receiver()
void HybridObservablesTest::check_results_carrier_phase( void HybridObservablesTest::check_results_carrier_phase(
arma::mat& true_ch0, arma::mat& true_ch0,
arma::mat& true_ch1,
arma::vec& true_tow_s, arma::vec& true_tow_s,
arma::mat& measured_ch0, arma::mat& measured_ch0,
arma::mat& measured_ch1) std::string data_title)
{ {
//1. True value interpolation to match the measurement times //1. True value interpolation to match the measurement times
double t0 = std::max(measured_ch0(0, 0), measured_ch1(0, 0)); double t0 = measured_ch0(0, 0);
int size1 = measured_ch0.col(0).n_rows; int size1 = measured_ch0.col(0).n_rows;
int size2 = measured_ch1.col(0).n_rows; double t1 = measured_ch0(size1 - 1, 0);
double t1 = std::min(measured_ch0(size1 - 1, 0), measured_ch1(size2 - 1, 0));
arma::vec t = arma::linspace<arma::vec>(t0, t1, floor((t1 - t0) * 1e3)); arma::vec t = arma::linspace<arma::vec>(t0, t1, floor((t1 - t0) * 1e3));
//conversion between arma::vec and std:vector
arma::vec t_from_start = arma::linspace<arma::vec>(0, t1 - t0, floor((t1 - t0) * 1e3));
std::vector<double> time_vector(t_from_start.colptr(0), t_from_start.colptr(0) + t_from_start.n_rows);
arma::vec true_ch0_phase_interp; arma::vec true_ch0_phase_interp;
arma::vec true_ch1_phase_interp;
arma::interp1(true_tow_s, true_ch0.col(3), t, true_ch0_phase_interp); arma::interp1(true_tow_s, true_ch0.col(3), t, true_ch0_phase_interp);
arma::interp1(true_tow_s, true_ch1.col(3), t, true_ch1_phase_interp);
arma::vec meas_ch0_phase_interp; arma::vec meas_ch0_phase_interp;
arma::vec meas_ch1_phase_interp;
arma::interp1(measured_ch0.col(0), measured_ch0.col(3), t, meas_ch0_phase_interp); arma::interp1(measured_ch0.col(0), measured_ch0.col(3), t, meas_ch0_phase_interp);
arma::interp1(measured_ch1.col(0), measured_ch1.col(3), t, meas_ch1_phase_interp);
//2. RMSE //2. RMSE
arma::vec err_ch0_cycles; arma::vec err_ch0_cycles;
arma::vec err_ch1_cycles;
//compute error without the accumulated carrier phase offsets (which depends on the receiver starting time) //compute error without the accumulated carrier phase offsets (which depends on the receiver starting time)
err_ch0_cycles = meas_ch0_phase_interp - true_ch0_phase_interp - meas_ch0_phase_interp(0) + true_ch0_phase_interp(0); err_ch0_cycles = meas_ch0_phase_interp - true_ch0_phase_interp - meas_ch0_phase_interp(0) + true_ch0_phase_interp(0);
err_ch1_cycles = meas_ch1_phase_interp - true_ch1_phase_interp - meas_ch1_phase_interp(0) + true_ch1_phase_interp(0);
arma::vec err2_ch0 = arma::square(err_ch0_cycles); arma::vec err2_ch0 = arma::square(err_ch0_cycles);
double rmse_ch0 = sqrt(arma::mean(err2_ch0)); double rmse_ch0 = sqrt(arma::mean(err2_ch0));
@ -341,21 +327,9 @@ void HybridObservablesTest::check_results_carrier_phase(
double max_error_ch0 = arma::max(err_ch0_cycles); double max_error_ch0 = arma::max(err_ch0_cycles);
double min_error_ch0 = arma::min(err_ch0_cycles); double min_error_ch0 = arma::min(err_ch0_cycles);
arma::vec err2_ch1 = arma::square(err_ch1_cycles);
double rmse_ch1 = sqrt(arma::mean(err2_ch1));
//3. Mean err and variance
double error_mean_ch1 = arma::mean(err_ch1_cycles);
double error_var_ch1 = arma::var(err_ch1_cycles);
// 4. Peaks
double max_error_ch1 = arma::max(err_ch1_cycles);
double min_error_ch1 = arma::min(err_ch1_cycles);
//5. report //5. report
std::streamsize ss = std::cout.precision(); std::streamsize ss = std::cout.precision();
std::cout << std::setprecision(10) << "Channel 0 Carrier phase RMSE = " std::cout << std::setprecision(10) << data_title << " Accumulated Carrier phase RMSE = "
<< rmse_ch0 << ", mean = " << error_mean_ch0 << rmse_ch0 << ", mean = " << error_mean_ch0
<< ", stdev = " << sqrt(error_var_ch0) << ", stdev = " << sqrt(error_var_ch0)
<< " (max,min) = " << max_error_ch0 << " (max,min) = " << max_error_ch0
@ -363,31 +337,115 @@ void HybridObservablesTest::check_results_carrier_phase(
<< " [cycles]" << std::endl; << " [cycles]" << std::endl;
std::cout.precision(ss); std::cout.precision(ss);
//plots
Gnuplot g3("linespoints");
g3.set_title(data_title + "Accumulated Carrier phase error [cycles]");
g3.set_grid();
g3.set_xlabel("Time [s]");
g3.set_ylabel("Carrier Phase error [cycles]");
//conversion between arma::vec and std:vector
std::vector<double> error_vec(err_ch0_cycles.colptr(0), err_ch0_cycles.colptr(0) + err_ch0_cycles.n_rows);
g3.cmd("set key box opaque");
g3.plot_xy(time_vector, error_vec,
"Delta pseudorrange error");
g3.set_legend();
g3.savetops(data_title + "Carrier_phase_error");
if (FLAGS_show_plots)
{
g3.showonscreen(); // window output
}
else
{
g3.disablescreen();
}
ASSERT_LT(rmse_ch0, 5e-2); ASSERT_LT(rmse_ch0, 5e-2);
ASSERT_LT(error_mean_ch0, 5e-2); ASSERT_LT(error_mean_ch0, 5e-2);
ASSERT_GT(error_mean_ch0, -5e-2); ASSERT_GT(error_mean_ch0, -5e-2);
ASSERT_LT(error_var_ch0, 5e-2); ASSERT_LT(error_var_ch0, 5e-2);
ASSERT_LT(max_error_ch0, 5e-2); ASSERT_LT(max_error_ch0, 5e-2);
ASSERT_GT(min_error_ch0, -5e-2); ASSERT_GT(min_error_ch0, -5e-2);
//5. report
ss = std::cout.precision();
std::cout << std::setprecision(10) << "Channel 1 Carrier phase RMSE = "
<< rmse_ch1 << ", mean = " << error_mean_ch1
<< ", stdev = " << sqrt(error_var_ch1)
<< " (max,min) = " << max_error_ch1
<< "," << min_error_ch1
<< " [cycles]" << std::endl;
std::cout.precision(ss);
ASSERT_LT(rmse_ch1, 5e-2);
ASSERT_LT(error_mean_ch1, 5e-2);
ASSERT_GT(error_mean_ch1, -5e-2);
ASSERT_LT(error_var_ch1, 5e-2);
ASSERT_LT(max_error_ch1, 5e-2);
ASSERT_GT(min_error_ch1, -5e-2);
} }
void HybridObservablesTest::check_results_carrier_doppler(
arma::mat& true_ch0,
arma::vec& true_tow_s,
arma::mat& measured_ch0,
std::string data_title)
{
//1. True value interpolation to match the measurement times
double t0 = measured_ch0(0, 0);
int size1 = measured_ch0.col(0).n_rows;
double t1 = measured_ch0(size1 - 1, 0);
arma::vec t = arma::linspace<arma::vec>(t0, t1, floor((t1 - t0) * 1e3));
//conversion between arma::vec and std:vector
arma::vec t_from_start = arma::linspace<arma::vec>(0, t1 - t0, floor((t1 - t0) * 1e3));
std::vector<double> time_vector(t_from_start.colptr(0), t_from_start.colptr(0) + t_from_start.n_rows);
arma::vec true_ch0_doppler_interp;
arma::interp1(true_tow_s, true_ch0.col(2), t, true_ch0_doppler_interp);
arma::vec meas_ch0_doppler_interp;
arma::interp1(measured_ch0.col(0), measured_ch0.col(2), t, meas_ch0_doppler_interp);
//2. RMSE
arma::vec err_ch0_hz;
//compute error
err_ch0_hz = meas_ch0_doppler_interp - true_ch0_doppler_interp;
arma::vec err2_ch0 = arma::square(err_ch0_hz);
double rmse_ch0 = sqrt(arma::mean(err2_ch0));
//3. Mean err and variance
double error_mean_ch0 = arma::mean(err_ch0_hz);
double error_var_ch0 = arma::var(err_ch0_hz);
// 4. Peaks
double max_error_ch0 = arma::max(err_ch0_hz);
double min_error_ch0 = arma::min(err_ch0_hz);
//5. report
std::streamsize ss = std::cout.precision();
std::cout << std::setprecision(10) << data_title << "Carrier Doppler RMSE = "
<< rmse_ch0 << ", mean = " << error_mean_ch0
<< ", stdev = " << sqrt(error_var_ch0)
<< " (max,min) = " << max_error_ch0
<< "," << min_error_ch0
<< " [Hz]" << std::endl;
std::cout.precision(ss);
//plots
Gnuplot g3("linespoints");
g3.set_title(data_title + "Carrier Doppler error [Hz]");
g3.set_grid();
g3.set_xlabel("Time [s]");
g3.set_ylabel("Carrier Doppler error [Hz]");
//conversion between arma::vec and std:vector
std::vector<double> error_vec(err_ch0_hz.colptr(0), err_ch0_hz.colptr(0) + err_ch0_hz.n_rows);
g3.cmd("set key box opaque");
g3.plot_xy(time_vector, error_vec,
"Delta pseudorrange error");
g3.set_legend();
g3.savetops(data_title + "Carrier_doppler_error");
if (FLAGS_show_plots)
{
g3.showonscreen(); // window output
}
else
{
g3.disablescreen();
}
ASSERT_LT(rmse_ch0, 5);
ASSERT_LT(error_mean_ch0, 5);
ASSERT_GT(error_mean_ch0, -5);
ASSERT_LT(error_var_ch0, 10);
ASSERT_LT(max_error_ch0, 10);
ASSERT_GT(min_error_ch0, -5);
}
bool HybridObservablesTest::save_mat_xy(std::vector<double>& x, std::vector<double>& y, std::string filename) bool HybridObservablesTest::save_mat_xy(std::vector<double>& x, std::vector<double>& y, std::string filename)
{ {
@ -424,12 +482,13 @@ bool HybridObservablesTest::save_mat_xy(std::vector<double>& x, std::vector<doub
} }
} }
void HybridObservablesTest::check_results_code_psudorange( void HybridObservablesTest::check_results_code_pseudorange(
arma::mat& true_ch0, arma::mat& true_ch0,
arma::mat& true_ch1, arma::mat& true_ch1,
arma::vec& true_tow_s, arma::vec& true_tow_s,
arma::mat& measured_ch0, arma::mat& measured_ch0,
arma::mat& measured_ch1) arma::mat& measured_ch1,
std::string data_title)
{ {
//1. True value interpolation to match the measurement times //1. True value interpolation to match the measurement times
@ -475,7 +534,7 @@ void HybridObservablesTest::check_results_code_psudorange(
//5. report //5. report
std::streamsize ss = std::cout.precision(); std::streamsize ss = std::cout.precision();
std::cout << std::setprecision(10) << "Delta Observables RMSE = " std::cout << std::setprecision(10) << data_title << "Delta Observables RMSE = "
<< rmse << ", mean = " << error_mean << rmse << ", mean = " << error_mean
<< ", stdev = " << sqrt(error_var) << ", stdev = " << sqrt(error_var)
<< " (max,min) = " << max_error << " (max,min) = " << max_error
@ -484,9 +543,8 @@ void HybridObservablesTest::check_results_code_psudorange(
std::cout.precision(ss); std::cout.precision(ss);
//plots //plots
Gnuplot g3("linespoints"); Gnuplot g3("linespoints");
g3.set_title("Delta Pseudorange error [m]"); g3.set_title(data_title + "Delta Pseudorange error [m]");
g3.set_grid(); g3.set_grid();
g3.set_xlabel("Time [s]"); g3.set_xlabel("Time [s]");
g3.set_ylabel("Pseudorange error [m]"); g3.set_ylabel("Pseudorange error [m]");
@ -496,8 +554,7 @@ void HybridObservablesTest::check_results_code_psudorange(
g3.plot_xy(time_vector, range_error_m, g3.plot_xy(time_vector, range_error_m,
"Delta pseudorrange error"); "Delta pseudorrange error");
g3.set_legend(); g3.set_legend();
g3.savetops("Delta_pseudorrange_error"); g3.savetops(data_title + "Delta_pseudorrange_error");
g3.savetopdf("Delta_pseudorrange_error", 18);
if (FLAGS_show_plots) if (FLAGS_show_plots)
{ {
g3.showonscreen(); // window output g3.showonscreen(); // window output
@ -531,135 +588,125 @@ TEST_F(HybridObservablesTest, ValidationOfResults)
std::chrono::time_point<std::chrono::system_clock> start, end; std::chrono::time_point<std::chrono::system_clock> start, end;
std::chrono::duration<double> elapsed_seconds(0); std::chrono::duration<double> elapsed_seconds(0);
Gnss_Synchro tmp_gnss_synchro;
tmp_gnss_synchro.System = 'G';
std::string signal = "1C";
signal.copy(tmp_gnss_synchro.Signal, 2, 0);
std::istringstream ss(FLAGS_test_satellite_PRN_list);
std::string token;
while (std::getline(ss, token, ','))
{
tmp_gnss_synchro.PRN = boost::lexical_cast<int>(token);
gnss_synchro_vec.push_back(tmp_gnss_synchro);
}
configure_receiver(); configure_receiver();
//open true observables log file written by the simulator //open true observables log file written by the simulator
tracking_true_obs_reader true_obs_data_ch0; std::vector<std::shared_ptr<tracking_true_obs_reader>> true_reader_vec;
tracking_true_obs_reader true_obs_data_ch1; std::vector<std::shared_ptr<TrackingInterface>> tracking_ch_vec;
int test_satellite_PRN = FLAGS_test_satellite_PRN; std::vector<std::shared_ptr<TelemetryDecoderInterface>> tlm_ch_vec;
int test_satellite_PRN2 = FLAGS_test_satellite_PRN2;
std::cout << "Testing satellite PRNs " << test_satellite_PRN << "," << test_satellite_PRN2 << std::endl;
std::string true_obs_file = std::string("./gps_l1_ca_obs_prn");
true_obs_file.append(std::to_string(test_satellite_PRN));
true_obs_file.append(".dat");
ASSERT_NO_THROW({
if (true_obs_data_ch0.open_obs_file(true_obs_file) == false)
{
throw std::exception();
};
}) << "Failure opening true observables file";
true_obs_file = std::string("./gps_l1_ca_obs_prn"); std::vector<gr::blocks::null_sink::sptr> null_sink_vec;
true_obs_file.append(std::to_string(test_satellite_PRN2)); for (unsigned int n = 0; n < gnss_synchro_vec.size(); n++)
true_obs_file.append(".dat"); {
ASSERT_NO_THROW({ //set channels ids
if (true_obs_data_ch1.open_obs_file(true_obs_file) == false) gnss_synchro_vec.at(n).Channel_ID = n;
{ //read true data from the generator logs
throw std::exception(); true_reader_vec.push_back(std::make_shared<tracking_true_obs_reader>());
}; std::cout << "Loading true observable data for PRN " << gnss_synchro_vec.at(n).PRN << std::endl;
}) << "Failure opening true observables file"; std::string true_obs_file = std::string("./gps_l1_ca_obs_prn");
true_obs_file.append(std::to_string(gnss_synchro_vec.at(n).PRN));
true_obs_file.append(".dat");
ASSERT_NO_THROW({
if (true_reader_vec.back()->open_obs_file(true_obs_file) == false)
{
throw std::exception();
};
}) << "Failure opening true observables file";
// load acquisition data based on the first epoch of the true observations
ASSERT_NO_THROW({
if (true_reader_vec.back()->read_binary_obs() == false)
{
throw std::exception();
};
}) << "Failure reading true observables file";
//restart the epoch counter
true_reader_vec.back()->restart();
std::cout << "Initial Doppler [Hz]=" << true_reader_vec.back()->doppler_l1_hz << " Initial code delay [Chips]="
<< true_reader_vec.back()->prn_delay_chips << std::endl;
gnss_synchro_vec.at(n).Acq_delay_samples = (GPS_L1_CA_CODE_LENGTH_CHIPS - true_reader_vec.back()->prn_delay_chips / GPS_L1_CA_CODE_LENGTH_CHIPS) * baseband_sampling_freq * GPS_L1_CA_CODE_PERIOD;
gnss_synchro_vec.at(n).Acq_doppler_hz = true_reader_vec.back()->doppler_l1_hz;
gnss_synchro_vec.at(n).Acq_samplestamp_samples = 0;
//create the tracking channels
tracking_ch_vec.push_back(std::make_shared<GpsL1CaDllPllTracking>(config.get(), "Tracking_1C", 1, 1));
//create the telemetry decoders
tlm_ch_vec.push_back(std::make_shared<GpsL1CaTelemetryDecoder>(config.get(), "TelemetryDecoder_1C", 1, 1));
//create null sinks for observables output
null_sink_vec.push_back(gr::blocks::null_sink::make(sizeof(Gnss_Synchro)));
ASSERT_NO_THROW({
tlm_ch_vec.back()->set_channel(gnss_synchro_vec.at(n).Channel_ID);
tlm_ch_vec.back()->set_satellite(Gnss_Satellite(std::string("GPS"), gnss_synchro_vec.at(n).PRN));
}) << "Failure setting gnss_synchro.";
ASSERT_NO_THROW({
tracking_ch_vec.back()->set_channel(gnss_synchro_vec.at(n).Channel_ID);
}) << "Failure setting channel.";
ASSERT_NO_THROW({
tracking_ch_vec.back()->set_gnss_synchro(&gnss_synchro_vec.at(n));
}) << "Failure setting gnss_synchro.";
}
top_block = gr::make_top_block("Telemetry_Decoder test"); top_block = gr::make_top_block("Telemetry_Decoder test");
std::shared_ptr<TrackingInterface> tracking_ch0 = std::make_shared<GpsL1CaDllPllTracking>(config.get(), "Tracking_1C", 1, 1); boost::shared_ptr<HybridObservablesTest_msg_rx> dummy_msg_rx_trk = HybridObservablesTest_msg_rx_make();
std::shared_ptr<TrackingInterface> tracking_ch1 = std::make_shared<GpsL1CaDllPllTracking>(config.get(), "Tracking_1C", 1, 1); boost::shared_ptr<HybridObservablesTest_tlm_msg_rx> dummy_tlm_msg_rx = HybridObservablesTest_tlm_msg_rx_make();
boost::shared_ptr<HybridObservablesTest_msg_rx> msg_rx_ch0 = HybridObservablesTest_msg_rx_make();
boost::shared_ptr<HybridObservablesTest_msg_rx> msg_rx_ch1 = HybridObservablesTest_msg_rx_make();
// load acquisition data based on the first epoch of the true observations
ASSERT_NO_THROW({
if (true_obs_data_ch0.read_binary_obs() == false)
{
throw std::exception();
};
}) << "Failure reading true observables file";
ASSERT_NO_THROW({
if (true_obs_data_ch1.read_binary_obs() == false)
{
throw std::exception();
};
}) << "Failure reading true observables file";
//restart the epoch counter
true_obs_data_ch0.restart();
true_obs_data_ch1.restart();
std::cout << "Initial Doppler [Hz]=" << true_obs_data_ch0.doppler_l1_hz << " Initial code delay [Chips]=" << true_obs_data_ch0.prn_delay_chips << std::endl;
gnss_synchro_ch0.Acq_delay_samples = (GPS_L1_CA_CODE_LENGTH_CHIPS - true_obs_data_ch0.prn_delay_chips / GPS_L1_CA_CODE_LENGTH_CHIPS) * baseband_sampling_freq * GPS_L1_CA_CODE_PERIOD;
gnss_synchro_ch0.Acq_doppler_hz = true_obs_data_ch0.doppler_l1_hz;
gnss_synchro_ch0.Acq_samplestamp_samples = 0;
std::cout << "Initial Doppler [Hz]=" << true_obs_data_ch1.doppler_l1_hz << " Initial code delay [Chips]=" << true_obs_data_ch1.prn_delay_chips << std::endl;
gnss_synchro_ch1.Acq_delay_samples = (GPS_L1_CA_CODE_LENGTH_CHIPS - true_obs_data_ch1.prn_delay_chips / GPS_L1_CA_CODE_LENGTH_CHIPS) * baseband_sampling_freq * GPS_L1_CA_CODE_PERIOD;
gnss_synchro_ch1.Acq_doppler_hz = true_obs_data_ch1.doppler_l1_hz;
gnss_synchro_ch1.Acq_samplestamp_samples = 0;
//telemetry decoders
std::shared_ptr<TelemetryDecoderInterface> tlm_ch0(new GpsL1CaTelemetryDecoder(config.get(), "TelemetryDecoder_1C", 1, 1));
std::shared_ptr<TelemetryDecoderInterface> tlm_ch1(new GpsL1CaTelemetryDecoder(config.get(), "TelemetryDecoder_1C", 1, 1));
ASSERT_NO_THROW({
tlm_ch0->set_channel(0);
tlm_ch1->set_channel(1);
tlm_ch0->set_satellite(Gnss_Satellite(std::string("GPS"), gnss_synchro_ch0.PRN));
tlm_ch1->set_satellite(Gnss_Satellite(std::string("GPS"), gnss_synchro_ch1.PRN));
}) << "Failure setting gnss_synchro.";
boost::shared_ptr<HybridObservablesTest_tlm_msg_rx> tlm_msg_rx_ch1 = HybridObservablesTest_tlm_msg_rx_make();
boost::shared_ptr<HybridObservablesTest_tlm_msg_rx> tlm_msg_rx_ch2 = HybridObservablesTest_tlm_msg_rx_make();
//Observables //Observables
std::shared_ptr<ObservablesInterface> observables(new HybridObservables(config.get(), "Observables", 3, 2)); std::shared_ptr<ObservablesInterface> observables(new HybridObservables(config.get(), "Observables", tracking_ch_vec.size() + 1, tracking_ch_vec.size()));
ASSERT_NO_THROW({ for (unsigned int n = 0; n < tracking_ch_vec.size(); n++)
tracking_ch0->set_channel(gnss_synchro_ch0.Channel_ID); {
tracking_ch1->set_channel(gnss_synchro_ch1.Channel_ID); ASSERT_NO_THROW({
}) << "Failure setting channel."; tracking_ch_vec.at(n)->connect(top_block);
}) << "Failure connecting tracking to the top_block.";
}
ASSERT_NO_THROW({
tracking_ch0->set_gnss_synchro(&gnss_synchro_ch0);
tracking_ch1->set_gnss_synchro(&gnss_synchro_ch1);
}) << "Failure setting gnss_synchro.";
ASSERT_NO_THROW({
tracking_ch0->connect(top_block);
tracking_ch1->connect(top_block);
}) << "Failure connecting tracking to the top_block.";
ASSERT_NO_THROW({ ASSERT_NO_THROW({
std::string file = "./" + filename_raw_data; std::string file = "./" + filename_raw_data;
const char* file_name = file.c_str(); const char* file_name = file.c_str();
gr::blocks::file_source::sptr file_source = gr::blocks::file_source::make(sizeof(int8_t), file_name, false); gr::blocks::file_source::sptr file_source = gr::blocks::file_source::make(sizeof(int8_t), file_name, false);
gr::blocks::interleaved_char_to_complex::sptr gr_interleaved_char_to_complex = gr::blocks::interleaved_char_to_complex::make(); gr::blocks::interleaved_char_to_complex::sptr gr_interleaved_char_to_complex = gr::blocks::interleaved_char_to_complex::make();
gr::blocks::null_sink::sptr sink_ch0 = gr::blocks::null_sink::make(sizeof(Gnss_Synchro));
gr::blocks::null_sink::sptr sink_ch1 = gr::blocks::null_sink::make(sizeof(Gnss_Synchro));
gnss_sdr_sample_counter_sptr samp_counter = gnss_sdr_make_sample_counter(static_cast<double>(baseband_sampling_freq), sizeof(gr_complex)); gnss_sdr_sample_counter_sptr samp_counter = gnss_sdr_make_sample_counter(static_cast<double>(baseband_sampling_freq), sizeof(gr_complex));
top_block->connect(file_source, 0, gr_interleaved_char_to_complex, 0); top_block->connect(file_source, 0, gr_interleaved_char_to_complex, 0);
top_block->connect(gr_interleaved_char_to_complex, 0, samp_counter, 0); top_block->connect(gr_interleaved_char_to_complex, 0, samp_counter, 0);
//ch0 for (unsigned int n = 0; n < tracking_ch_vec.size(); n++)
top_block->connect(gr_interleaved_char_to_complex, 0, tracking_ch0->get_left_block(), 0); {
top_block->connect(tracking_ch0->get_right_block(), 0, tlm_ch0->get_left_block(), 0); top_block->connect(gr_interleaved_char_to_complex, 0, tracking_ch_vec.at(n)->get_left_block(), 0);
top_block->connect(tlm_ch0->get_right_block(), 0, observables->get_left_block(), 0); top_block->connect(tracking_ch_vec.at(n)->get_right_block(), 0, tlm_ch_vec.at(n)->get_left_block(), 0);
top_block->msg_connect(tracking_ch0->get_right_block(), pmt::mp("events"), msg_rx_ch0, pmt::mp("events")); top_block->connect(tlm_ch_vec.at(n)->get_right_block(), 0, observables->get_left_block(), n);
//ch1 top_block->msg_connect(tracking_ch_vec.at(n)->get_right_block(), pmt::mp("events"), dummy_msg_rx_trk, pmt::mp("events"));
top_block->connect(gr_interleaved_char_to_complex, 0, tracking_ch1->get_left_block(), 0); top_block->connect(observables->get_right_block(), n, null_sink_vec.at(n), 0);
top_block->connect(tracking_ch1->get_right_block(), 0, tlm_ch1->get_left_block(), 0); }
top_block->connect(tlm_ch1->get_right_block(), 0, observables->get_left_block(), 1); //connect sample counter and timmer to the last channel in observables block (extra channel)
top_block->msg_connect(tracking_ch1->get_right_block(), pmt::mp("events"), msg_rx_ch1, pmt::mp("events")); top_block->connect(samp_counter, 0, observables->get_left_block(), tracking_ch_vec.size());
top_block->connect(observables->get_right_block(), 0, sink_ch0, 0);
top_block->connect(observables->get_right_block(), 1, sink_ch1, 0);
top_block->connect(samp_counter, 0, observables->get_left_block(), 2);
}) << "Failure connecting the blocks."; }) << "Failure connecting the blocks.";
tracking_ch0->start_tracking(); for (unsigned int n = 0; n < tracking_ch_vec.size(); n++)
tracking_ch1->start_tracking(); {
tracking_ch_vec.at(n)->start_tracking();
}
EXPECT_NO_THROW({ EXPECT_NO_THROW({
start = std::chrono::system_clock::now(); start = std::chrono::system_clock::now();
@ -684,40 +731,36 @@ TEST_F(HybridObservablesTest, ValidationOfResults)
std::cout << "True observation epochs = " << nepoch << std::endl; std::cout << "True observation epochs = " << nepoch << std::endl;
// Matrices for storing columnwise true GPS time, Range, Doppler and Carrier phase // Matrices for storing columnwise true GPS time, Range, Doppler and Carrier phase
arma::mat true_ch0 = arma::zeros<arma::mat>(nepoch, 4); std::vector<arma::mat> true_obs_vec;
arma::mat true_ch1 = arma::zeros<arma::mat>(nepoch, 4);
true_observables.restart(); true_observables.restart();
long int epoch_counter = 0; long int epoch_counter = 0;
for (unsigned int n = 0; n < tracking_ch_vec.size(); n++)
{
true_obs_vec.push_back(arma::zeros<arma::mat>(nepoch, 4));
}
ASSERT_NO_THROW({ ASSERT_NO_THROW({
while (true_observables.read_binary_obs()) while (true_observables.read_binary_obs())
{ {
if (round(true_observables.prn[0]) != gnss_synchro_ch0.PRN) for (unsigned int n = 0; n < true_obs_vec.size(); n++)
{ {
std::cout << "True observables SV PRN does not match " << round(true_observables.prn[1]) << std::endl; if (round(true_observables.prn[n]) != gnss_synchro_vec.at(n).PRN)
throw std::exception(); {
std::cout << "True observables SV PRN does not match measured ones: "
<< round(true_observables.prn[n]) << " vs. " << gnss_synchro_vec.at(n).PRN << std::endl;
throw std::exception();
}
true_obs_vec.at(n)(epoch_counter, 0) = true_observables.gps_time_sec[n];
true_obs_vec.at(n)(epoch_counter, 1) = true_observables.dist_m[n];
true_obs_vec.at(n)(epoch_counter, 2) = true_observables.doppler_l1_hz[n];
true_obs_vec.at(n)(epoch_counter, 3) = true_observables.acc_carrier_phase_l1_cycles[n];
} }
if (round(true_observables.prn[1]) != gnss_synchro_ch1.PRN)
{
std::cout << "True observables SV PRN does not match " << round(true_observables.prn[1]) << std::endl;
throw std::exception();
}
true_ch0(epoch_counter, 0) = true_observables.gps_time_sec[0];
true_ch0(epoch_counter, 1) = true_observables.dist_m[0];
true_ch0(epoch_counter, 2) = true_observables.doppler_l1_hz[0];
true_ch0(epoch_counter, 3) = true_observables.acc_carrier_phase_l1_cycles[0];
true_ch1(epoch_counter, 0) = true_observables.gps_time_sec[1];
true_ch1(epoch_counter, 1) = true_observables.dist_m[1];
true_ch1(epoch_counter, 2) = true_observables.doppler_l1_hz[1];
true_ch1(epoch_counter, 3) = true_observables.acc_carrier_phase_l1_cycles[1];
epoch_counter++; epoch_counter++;
} }
}); });
//read measured values //read measured values
observables_dump_reader estimated_observables(2); //two channels observables_dump_reader estimated_observables(tracking_ch_vec.size());
ASSERT_NO_THROW({ ASSERT_NO_THROW({
if (estimated_observables.open_obs_file(std::string("./observables.dat")) == false) if (estimated_observables.open_obs_file(std::string("./observables.dat")) == false)
{ {
@ -726,96 +769,114 @@ TEST_F(HybridObservablesTest, ValidationOfResults)
}) << "Failure opening dump observables file"; }) << "Failure opening dump observables file";
nepoch = static_cast<unsigned int>(estimated_observables.num_epochs()); nepoch = static_cast<unsigned int>(estimated_observables.num_epochs());
std::cout << "Measured observation epochs = " << nepoch << std::endl; std::cout << "Measured observations epochs = " << nepoch << std::endl;
// Matrices for storing columnwise measured RX_time, TOW, Doppler, Carrier phase and Pseudorange // Matrices for storing columnwise measured RX_time, TOW, Doppler, Carrier phase and Pseudorange
arma::mat measured_ch0 = arma::zeros<arma::mat>(nepoch, 5); std::vector<arma::mat> measured_obs_vec;
arma::mat measured_ch1 = arma::zeros<arma::mat>(nepoch, 5); std::vector<long int> epoch_counters_vec;
for (unsigned int n = 0; n < tracking_ch_vec.size(); n++)
{
measured_obs_vec.push_back(arma::zeros<arma::mat>(nepoch, 5));
epoch_counters_vec.push_back(0);
}
estimated_observables.restart(); estimated_observables.restart();
epoch_counter = 0;
long int epoch_counter2 = 0;
while (estimated_observables.read_binary_obs()) while (estimated_observables.read_binary_obs())
{ {
if (static_cast<bool>(estimated_observables.valid[0])) for (unsigned int n = 0; n < measured_obs_vec.size(); n++)
{ {
measured_ch0(epoch_counter, 0) = estimated_observables.RX_time[0]; if (static_cast<bool>(estimated_observables.valid[n]))
measured_ch0(epoch_counter, 1) = estimated_observables.TOW_at_current_symbol_s[0]; {
measured_ch0(epoch_counter, 2) = estimated_observables.Carrier_Doppler_hz[0]; measured_obs_vec.at(n)(epoch_counters_vec.at(n), 0) = estimated_observables.RX_time[n];
measured_ch0(epoch_counter, 3) = estimated_observables.Acc_carrier_phase_hz[0]; measured_obs_vec.at(n)(epoch_counters_vec.at(n), 1) = estimated_observables.TOW_at_current_symbol_s[n];
measured_ch0(epoch_counter, 4) = estimated_observables.Pseudorange_m[0]; measured_obs_vec.at(n)(epoch_counters_vec.at(n), 2) = estimated_observables.Carrier_Doppler_hz[n];
epoch_counter++; measured_obs_vec.at(n)(epoch_counters_vec.at(n), 3) = estimated_observables.Acc_carrier_phase_hz[n];
} measured_obs_vec.at(n)(epoch_counters_vec.at(n), 4) = estimated_observables.Pseudorange_m[n];
if (static_cast<bool>(estimated_observables.valid[1])) epoch_counters_vec.at(n)++;
{ }
measured_ch1(epoch_counter2, 0) = estimated_observables.RX_time[1];
measured_ch1(epoch_counter2, 1) = estimated_observables.TOW_at_current_symbol_s[1];
measured_ch1(epoch_counter2, 2) = estimated_observables.Carrier_Doppler_hz[1];
measured_ch1(epoch_counter2, 3) = estimated_observables.Acc_carrier_phase_hz[1];
measured_ch1(epoch_counter2, 4) = estimated_observables.Pseudorange_m[1];
epoch_counter2++;
} }
} }
//Cut measurement tail zeros //Cut measurement tail zeros
arma::uvec index = arma::find(measured_ch0.col(0) > 0.0, 1, "last"); arma::uvec index;
if ((index.size() > 0) and index(0) < (nepoch - 1)) for (unsigned int n = 0; n < measured_obs_vec.size(); n++)
{ {
measured_ch0.shed_rows(index(0) + 1, nepoch - 1); index = arma::find(measured_obs_vec.at(n).col(0) > 0.0, 1, "last");
} if ((index.size() > 0) and index(0) < (nepoch - 1))
index = arma::find(measured_ch1.col(0) > 0.0, 1, "last"); {
if ((index.size() > 0) and index(0) < (nepoch - 1)) measured_obs_vec.at(n).shed_rows(index(0) + 1, nepoch - 1);
{ }
measured_ch1.shed_rows(index(0) + 1, nepoch - 1);
} }
//Cut measurement initial transitory of the measurements //Cut measurement initial transitory of the measurements
double initial_transitory_s = 30.0; double initial_transitory_s = FLAGS_skip_obs_transitory_s;
index = arma::find(measured_ch0.col(0) >= (measured_ch0(0, 0) + initial_transitory_s), 1, "first"); for (unsigned int n = 0; n < measured_obs_vec.size(); n++)
if ((index.size() > 0) and (index(0) > 0))
{ {
measured_ch0.shed_rows(0, index(0)); index = arma::find(measured_obs_vec.at(n).col(0) >= (measured_obs_vec.at(n)(0, 0) + initial_transitory_s), 1, "first");
} if ((index.size() > 0) and (index(0) > 0))
index = arma::find(measured_ch1.col(0) >= (measured_ch1(0, 0) + initial_transitory_s), 1, "first"); {
if ((index.size() > 0) and (index(0) > 0)) measured_obs_vec.at(n).shed_rows(0, index(0));
{ }
measured_ch1.shed_rows(0, index(0));
index = arma::find(measured_obs_vec.at(n).col(0) >= true_obs_vec.at(n)(0, 0), 1, "first");
if ((index.size() > 0) and (index(0) > 0))
{
measured_obs_vec.at(n).shed_rows(0, index(0));
}
} }
index = arma::find(measured_ch0.col(0) >= true_ch0(0, 0), 1, "first");
if ((index.size() > 0) and (index(0) > 0))
{
measured_ch0.shed_rows(0, index(0));
}
index = arma::find(measured_ch1.col(0) >= true_ch1(0, 0), 1, "first");
if ((index.size() > 0) and (index(0) > 0))
{
measured_ch1.shed_rows(0, index(0));
}
//Correct the clock error using true values (it is not possible for a receiver to correct //Correct the clock error using true values (it is not possible for a receiver to correct
//the receiver clock offset error at the observables level because it is required the //the receiver clock offset error at the observables level because it is required the
//decoding of the ephemeris data and solve the PVT equations) //decoding of the ephemeris data and solve the PVT equations)
//Find the reference satellite (the nearest) and compute the receiver time offset at observable level //Find the reference satellite (the nearest) and compute the receiver time offset at observable level
double min_pr = std::numeric_limits<double>::max();
unsigned int min_pr_ch_id = 0;
arma::vec receiver_time_offset_s; arma::vec receiver_time_offset_s;
if (measured_ch0(0, 4) < measured_ch1(0, 4)) for (unsigned int n = 0; n < measured_obs_vec.size(); n++)
{ {
receiver_time_offset_s = true_ch0.col(1) / GPS_C_m_s - GPS_STARTOFFSET_ms / 1000.0; if (measured_obs_vec.at(n)(0, 4) < min_pr)
{
min_pr = measured_obs_vec.at(n)(0, 4);
min_pr_ch_id = n;
}
} }
else
receiver_time_offset_s = true_obs_vec.at(min_pr_ch_id).col(1) / GPS_C_m_s - GPS_STARTOFFSET_ms / 1000.0;
for (unsigned int n = 0; n < measured_obs_vec.size(); n++)
{ {
receiver_time_offset_s = true_ch1.col(1) / GPS_C_m_s - GPS_STARTOFFSET_ms / 1000.0; arma::vec corrected_reference_TOW_s = true_obs_vec.at(min_pr_ch_id).col(0) - receiver_time_offset_s;
std::cout << "[CH " << n << "] Receiver time offset " << receiver_time_offset_s(0) * 1e3 << " [ms]" << std::endl;
//Compare measured observables
if (min_pr_ch_id != n)
{
check_results_code_pseudorange(true_obs_vec.at(n),
true_obs_vec.at(min_pr_ch_id),
corrected_reference_TOW_s,
measured_obs_vec.at(n),
measured_obs_vec.at(min_pr_ch_id),
"[CH " + std::to_string(n) + "] PRN " + std::to_string(gnss_synchro_vec.at(n).PRN) + " ");
}
else
{
std::cout << "[CH " << std::to_string(n) << "] PRN " << std::to_string(gnss_synchro_vec.at(n).PRN) << " is the reference satellite" << std::endl;
}
std::cout << "true_obs_vec.at(n): " << true_obs_vec.at(n)(0, 2) << std::endl;
check_results_carrier_phase(true_obs_vec.at(n),
corrected_reference_TOW_s,
measured_obs_vec.at(n),
"[CH " + std::to_string(n) + "] PRN " + std::to_string(gnss_synchro_vec.at(n).PRN) + " ");
check_results_carrier_doppler(true_obs_vec.at(n),
corrected_reference_TOW_s,
measured_obs_vec.at(n),
"[CH " + std::to_string(n) + "] PRN " + std::to_string(gnss_synchro_vec.at(n).PRN) + " ");
} }
arma::vec corrected_reference_TOW_s = true_ch0.col(0) - receiver_time_offset_s;
std::cout << "Receiver time offset: " << receiver_time_offset_s(0) * 1e3 << " [ms]" << std::endl;
//Compare measured observables
check_results_code_psudorange(true_ch0, true_ch1, corrected_reference_TOW_s, measured_ch0, measured_ch1);
check_results_carrier_phase(true_ch0, true_ch1, corrected_reference_TOW_s, measured_ch0, measured_ch1);
std::cout << "Test completed in " << elapsed_seconds.count() << " [s]" << std::endl; std::cout << "Test completed in " << elapsed_seconds.count() << " [s]" << std::endl;
} }