mirrors/gnss-sdr
1
0
Fork 0
mirror of https://github.com/gnss-sdr/gnss-sdr synced 2024-12-14 12:10:34 +00:00
Code Issues Releases Wiki Activity

Working towards a RINEX printer for multiconstellation data. Code

Browse Source 
cleaning
This commit is contained in:
Carles Fernandez 2014-09-05 12:44:51 +02:00
parent a57c5ccf8a
commit 0d0c603420
7 changed files with 468 additions and 667 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
src/algorithms/PVT/gnuradio_blocks/galileo_e1_pvt_cc.cc
View File

@ -169,6 +169,7 @@ int galileo_e1_pvt_cc::general_work (int noutput_items, gr_vector_int &ninput_it
// Almanac data is shared for all the Galileo satellites. Read always at ID=0 // Almanac data is shared for all the Galileo satellites. Read always at ID=0
global_galileo_almanac_map.read(0, d_ls_pvt->galileo_almanac); global_galileo_almanac_map.read(0, d_ls_pvt->galileo_almanac);
} }
// ############ 2 COMPUTE THE PVT ################################ // ############ 2 COMPUTE THE PVT ################################
if (gnss_pseudoranges_map.size() > 0 and d_ls_pvt->galileo_ephemeris_map.size() > 0) if (gnss_pseudoranges_map.size() > 0 and d_ls_pvt->galileo_ephemeris_map.size() > 0)
{ {
@ -178,7 +179,6 @@ int galileo_e1_pvt_cc::general_work (int noutput_items, gr_vector_int &ninput_it
bool pvt_result; bool pvt_result;
pvt_result = d_ls_pvt->get_PVT(gnss_pseudoranges_map, d_rx_time, d_flag_averaging); pvt_result = d_ls_pvt->get_PVT(gnss_pseudoranges_map, d_rx_time, d_flag_averaging);
if (pvt_result == true) if (pvt_result == true)
{ {
d_kml_dump->print_position_galileo(d_ls_pvt, d_flag_averaging); d_kml_dump->print_position_galileo(d_ls_pvt, d_flag_averaging);

318
src/algorithms/PVT/gnuradio_blocks/hybrid_pvt_cc.cc
View File

@ -48,6 +48,7 @@ using google::LogMessage;
extern concurrent_map<Galileo_Ephemeris> global_galileo_ephemeris_map; extern concurrent_map<Galileo_Ephemeris> global_galileo_ephemeris_map;
extern concurrent_map<Galileo_Iono> global_galileo_iono_map; extern concurrent_map<Galileo_Iono> global_galileo_iono_map;
extern concurrent_map<Galileo_Utc_Model> global_galileo_utc_model_map; extern concurrent_map<Galileo_Utc_Model> global_galileo_utc_model_map;
extern concurrent_map<Galileo_Almanac> global_galileo_almanac_map;
extern concurrent_map<Gps_Ephemeris> global_gps_ephemeris_map; extern concurrent_map<Gps_Ephemeris> global_gps_ephemeris_map;
extern concurrent_map<Gps_Iono> global_gps_iono_map; extern concurrent_map<Gps_Iono> global_gps_iono_map;
@ -136,6 +137,7 @@ int hybrid_pvt_cc::general_work (int noutput_items, gr_vector_int &ninput_items,
gr_vector_const_void_star &input_items, gr_vector_void_star &output_items) gr_vector_const_void_star &input_items, gr_vector_void_star &output_items)
{ {
d_sample_counter++; d_sample_counter++;
bool arrived_galileo_almanac = false;
std::map<int,Gnss_Synchro> gnss_pseudoranges_map; std::map<int,Gnss_Synchro> gnss_pseudoranges_map;
@ -173,6 +175,12 @@ int hybrid_pvt_cc::general_work (int noutput_items, gr_vector_int &ninput_items,
global_galileo_iono_map.read(0, d_ls_pvt->galileo_iono); global_galileo_iono_map.read(0, d_ls_pvt->galileo_iono);
} }
if (global_galileo_almanac_map.size() > 0)
{
// Almanac data is shared for all the Galileo satellites. Read always at ID=0
arrived_galileo_almanac = global_galileo_almanac_map.read(0, d_ls_pvt->galileo_almanac);
}
// ############ 1. READ GPS EPHEMERIS/UTC_MODE/IONO FROM GLOBAL MAPS #### // ############ 1. READ GPS EPHEMERIS/UTC_MODE/IONO FROM GLOBAL MAPS ####
if (global_gps_ephemeris_map.size() > 0) if (global_gps_ephemeris_map.size() > 0)
@ -202,299 +210,49 @@ int hybrid_pvt_cc::general_work (int noutput_items, gr_vector_int &ninput_items,
// compute on the fly PVT solution // compute on the fly PVT solution
if ((d_sample_counter % d_output_rate_ms) == 0) if ((d_sample_counter % d_output_rate_ms) == 0)
{ {
//bool pvt_result;
d_ls_pvt->get_PVT(gnss_pseudoranges_map, d_rx_time, d_flag_averaging); d_ls_pvt->get_PVT(gnss_pseudoranges_map, d_rx_time, d_flag_averaging);
//std::cout << "pvt_result = " << pvt_result << std::endl;
if (d_ls_pvt->b_valid_position == true) bool pvt_result;
pvt_result = d_ls_pvt->get_PVT(gnss_pseudoranges_map, d_rx_time, d_flag_averaging);
if (pvt_result == true)
{ {
valid_solution_counter ++;
d_kml_dump->print_position_hybrid(d_ls_pvt, d_flag_averaging); d_kml_dump->print_position_hybrid(d_ls_pvt, d_flag_averaging);
/* *********************************** COMPUTE STATISTICS OVER THE FIRST 1000 SOLUTOINS **********************************/
if (valid_solution_counter<=1000)
//if (d_ls_pvt->d_valid_observations >= 14)
//if (d_ls_pvt->d_valid_GPS_obs >= 7 or d_ls_pvt->d_valid_GAL_obs >=7)
{
GDOP_vect.push_back (d_ls_pvt->d_GDOP);
//PDOP_vect.push_back (d_ls_pvt->d_PDOP);
X_vect.push_back (d_ls_pvt->X); // add the last X value to X_vect
Y_vect.push_back (d_ls_pvt->Y); // add the last X value to X_vect
Z_vect.push_back (d_ls_pvt->Z);
longitude_vect.push_back (d_ls_pvt->d_longitude_d);
latitude_vect.push_back (d_ls_pvt->d_latitude_d);
h_vect.push_back (d_ls_pvt->d_height_m);
tot_obs_vect.push_back (d_ls_pvt->d_valid_observations);
Gal_obs_vect.push_back (d_ls_pvt->d_valid_GAL_obs);
GPS_obs_vect.push_back (d_ls_pvt->d_valid_GPS_obs);
time_vect.push_back (d_ls_pvt->d_position_UTC_time);
//valid_solution_16_sat_counter ++;
}
//if(valid_solution_16_sat_counter ==1000)
/* if (valid_solution_counter==1000)
{
//compute ECEF average and standard deviation over the first 1000 solutions
GDOP_sum = std::accumulate( GDOP_vect.begin(), GDOP_vect.end(), 0.0);
GDOP_mean = GDOP_sum / GDOP_vect.size();
longitude_vect_sum = std::accumulate( longitude_vect.begin(), longitude_vect.end(), 0.0);
longitude_mean = longitude_vect_sum / longitude_vect.size();
//IFEN true solutions
double ref_longitude = 11.80800563;
double ref_latitude = 48.17149767;
double ref_X = 4171691.011;
double ref_Y = 872120.003;
double ref_Z = 4730005.761;
// REAL CAPTURE reference solutions (obtained with 4 Galileo)
// double ref_longitude= 1.987686994;
// double ref_latitude= 41.274786935;
// double ref_X=4797680.560;
// double ref_Y= 166506.414;
// double ref_Z=4185453.947;
//compute mean value for precision
latitude_vect_sum = std::accumulate( latitude_vect.begin(), latitude_vect.end(), 0.0);
latitude_mean = latitude_vect_sum / latitude_vect.size();
X_vect_sum = std::accumulate(X_vect.begin(), X_vect.end(), 0.0);
X_vect_mean = X_vect_sum / X_vect.size();
X_vect_sq_sum = std::inner_product(X_vect.begin(), X_vect.end(), X_vect.begin(), 0.0);
X_vect_stdev = std::sqrt(X_vect_sq_sum / X_vect.size() - X_vect_mean * X_vect_mean);
Y_vect_sum = std::accumulate(Y_vect.begin(), Y_vect.end(), 0.0);
Y_vect_mean = Y_vect_sum / Y_vect.size();
Y_vect_sq_sum = std::inner_product(Y_vect.begin(), Y_vect.end(), Y_vect.begin(), 0.0);
Y_vect_stdev = std::sqrt(Y_vect_sq_sum / Y_vect.size() - Y_vect_mean * Y_vect_mean);
Z_vect_sum = std::accumulate(Z_vect.begin(), Z_vect.end(), 0.0);
Z_vect_mean = Z_vect_sum / Z_vect.size();
Z_vect_sq_sum = std::inner_product(Z_vect.begin(), Z_vect.end(), Z_vect.begin(), 0.0);
Z_vect_stdev = std::sqrt(Z_vect_sq_sum / Z_vect.size() - Z_vect_mean * Z_vect_mean);
std::fstream file_solutions;
std::ofstream file_statitics;
//std::fstream file_RES_X;
file_solutions.open ("GNSS_SDR_solutions.txt", std::fstream::out);
//file_statitics.open ("GNSS_SDR_statitics.txt", std::ofstream::out | std::ofstream::app);
file_statitics.open ("GNSS_SDR_statitics.txt", std::ofstream::out);
file_solutions.setf(file_solutions.fixed, file_solutions.floatfield);
file_statitics.setf(file_statitics.fixed, file_statitics.floatfield);
file_statitics << std::setprecision(9);
file_solutions << std::setw(13) << "time" << std::setw(26) <<"X [m]" << std::setw(20)<< "Y [m]"<< std::setw(20)<< "Z [m]" << // X Y Z solutions
std::setw(20) << "Long [deg]" << std::setw(15) << "Lat [deg]" << std::setw(15) << "h [m]" << // long, lat, h solutions
std::setw(18) << "E(Acc) [m]" << std::setw(18) << "N(Acc) [m]" << std::setw(18) <<"Up(Acc) [m]" << //ENU residual for accuracy columns 7:8:9
std::setw(18) << "E(Pre) [m]" << std::setw(18) << "N(Pre) [m]" << std::setw(18) <<"Up(Pre) [m]" << //ENU residual for precision columns 10:11:12
std::setw(18) << "Tot Sat" << std::setw(18) << "Gal"<< std::setw(18) << "GPS" << //number of satellites
std::setw(18) << "GDOP" << std::endl; //GDOP
std::cout << "X mean over 1000 solutions [m]= " << X_vect_mean << ", X dev st over 1000 solutins= " << X_vect_stdev << std::endl;
std::cout << "Y mean over 1000 solutions [m]= " << Y_vect_mean << ", Y dev st over 1000 solutins= " << Y_vect_stdev << std::endl;
std::cout << "Z mean over 1000 solutions [m]= " << Z_vect_mean << ", Z dev st over 1000 solutins= " << Z_vect_stdev << std::endl;
//file_statitics << "X mean over 1000 solutions [m]= " << X_vect_mean << ", X dev st over 1000 solutins= " << X_vect_stdev << std::endl;
//file_statitics << "Y mean over 1000 solutions [m]= " << Y_vect_mean << ", Y dev st over 1000 solutins= " << Y_vect_stdev << std::endl;
//file_statitics << "Z mean over 1000 solutions [m]= " << Z_vect_mean << ", Z dev st over 1000 solutins= " << Z_vect_stdev << std::endl << std::endl;
file_statitics << "Num of GPS observation " << d_ls_pvt->d_valid_GPS_obs << std::endl;
file_statitics << "Num of GALILEO observation " << d_ls_pvt->d_valid_GAL_obs << std::endl;
file_statitics << "GDOP mean= " << GDOP_mean << std::endl << std::endl;
file_statitics << "ENU computed at (IFEN true coordinates): ref Longitude = " << ref_longitude << ", Ref Latitude = " << ref_latitude << " for Accuracy"<<std::endl;
file_statitics << "ENU computed at (average coordinates) mean Longitude = " << longitude_mean << ", mean Latitude = " << latitude_mean << " for Precision"<<std::endl;
//file_statitics << "Residual computed with respect (IFEN true coordinates) X Y Z= " << ref_X << '\t'<< ref_Y << '\t'<< ref_Z << std::endl;
// 1- Rotation matrix from ECEF coordinates to ENU coordinates
// ref: http://www.navipedia.net/index.php/Transformations_between_ECEF_and_ENU_coordinates
arma::mat R = arma::zeros(3,3);//matrix with REFERENCE position
arma::mat M = arma::zeros(3,3);//matrix with MEAN position
//matrix with MEAN position
M(0,0) = -sin(GPS_TWO_PI*(longitude_mean/360.0));
M(0,1) = -sin(GPS_TWO_PI*(latitude_mean/360.0))*cos(GPS_TWO_PI*(longitude_mean/360.0));
M(0,2) = cos(GPS_TWO_PI*(latitude_mean/360.0))*cos(GPS_TWO_PI*(longitude_mean/360.0));
M(1,0) = cos((GPS_TWO_PI*longitude_mean)/360.0);
M(1,1) = -sin((GPS_TWO_PI*latitude_mean)/360.0)*sin((GPS_TWO_PI*longitude_mean)/360.0);
M(1,2) = cos((GPS_TWO_PI*latitude_mean/360.0))*sin((GPS_TWO_PI*longitude_mean)/360.0);
M(2,0) = 0;
M(2,1) = cos((GPS_TWO_PI*latitude_mean)/360.0);
M(2,2) = sin((GPS_TWO_PI*latitude_mean/360.0));
//matrix with REFERENCE position
R(0,0) = -sin(GPS_TWO_PI*(ref_longitude/360.0));
R(0,1) = -sin(GPS_TWO_PI*(ref_latitude/360.0))*cos(GPS_TWO_PI*(ref_longitude/360.0));
R(0,2) = cos(GPS_TWO_PI*(ref_latitude/360.0))*cos(GPS_TWO_PI*(ref_longitude/360.0));
R(1,0) = cos((GPS_TWO_PI*ref_longitude)/360.0);
R(1,1) = -sin((GPS_TWO_PI*ref_latitude)/360.0)*sin((GPS_TWO_PI*ref_longitude)/360.0);
R(1,2) = cos((GPS_TWO_PI*ref_latitude/360.0))*sin((GPS_TWO_PI*ref_longitude)/360.0);
R(2,0) = 0;
R(2,1) = cos((GPS_TWO_PI*ref_latitude)/360.0);
R(2,2) = sin((GPS_TWO_PI*ref_latitude/360.0));
for(unsigned i=0; i<X_vect.size(); i++ )
{
//compute residual with respect to the REFERENCE SOLUTIONS (true IFEN coordinates)
X_vect_res.push_back (X_vect[i] - ref_X); //accuracy
Y_vect_res.push_back (Y_vect[i] - ref_Y); //accuracy
Z_vect_res.push_back (Z_vect[i] - ref_Z); //accuracy
//compute residual with respect to the MEAN SOLUTIONS
X_vect_res_precision.push_back (X_vect[i] - X_vect_mean); //precision
Y_vect_res_precision.push_back (Y_vect[i] - Y_vect_mean); //precision
Z_vect_res_precision.push_back (Z_vect[i] - Z_vect_mean); //precision
//apply to residual value the rotation from ECEF to ENU IN REFERENCE POSITION (matrix R)
E_res.push_back (R(0,0)*X_vect_res[i] + R(1,0)*Y_vect_res[i] + R(2,0) * Z_vect_res[i]);
N_res.push_back (R(0,1)*X_vect_res[i] + R(1,1)*Y_vect_res[i] + R(2,1) * Z_vect_res[i]);
Up_res.push_back (R(0,2)*X_vect_res[i] + R(1,2)*Y_vect_res[i] + R(2,2) * Z_vect_res[i]);
//apply to residual value the rotation from ECEF to ENU MEAN POSITOIN (matrix R)
E_res_precision.push_back (M(0,0)*X_vect_res_precision[i] + M(1,0)*Y_vect_res_precision[i] + M(2,0) * Z_vect_res_precision[i]);
N_res_precision.push_back (M(0,1)*X_vect_res_precision[i] + M(1,1)*Y_vect_res_precision[i] + M(2,1) * Z_vect_res_precision[i]);
Up_res_precision.push_back (M(0,2)*X_vect_res_precision[i] + M(1,2)*Y_vect_res_precision[i] + M(2,2) * Z_vect_res_precision[i]);
file_solutions << std::setprecision(5);
file_solutions << std::setw(13)<< time_vect[i]<<
std::setw(22)<< X_vect[i]<< std::setw(20)<< Y_vect[i]<< std::setw(20)<< Z_vect[i] << // X Y Z solutions
std::setw(15)<< longitude_vect [i]<<std::setw(15)<< latitude_vect [i]<< std::setw(15)<< h_vect [i]<< // long, lat, h solutions
std::setw(18)<< E_res[i]<< std::setw(18)<< N_res[i]<< std::setw(18)<< Up_res[i]<< //ENU residual for accuracy columns 7:8:9
std::setw(18)<< E_res_precision[i]<< std::setw(18)<< N_res_precision[i]<< std::setw(18)<<Up_res_precision[i]<< //ENU residual for precision columns 10:11:12
std::setw(18)<< std::setprecision(0)<<tot_obs_vect[i] << std::setw(18)<< std::setprecision(0)<<Gal_obs_vect[i] << std::setw(18)<<std::setprecision(0)<< GPS_obs_vect[i] << // number of satellites
std::setw(18)<<std::setprecision(2)<<GDOP_vect[i] <<std::endl;
}
//compute ENU average and standard deviation over the first 1000 residual solutions ACCURACY
E_res_sum = std::accumulate(E_res.begin(), E_res.end(), 0.0);
E_res_mean = E_res_sum / E_res.size();
E_res_sq_sum = std::inner_product(E_res.begin(), E_res.end(), E_res.begin(), 0.0);
E_res_stdev = std::sqrt(E_res_sq_sum / E_res.size() - E_res_mean * E_res_mean);
N_res_sum = std::accumulate(N_res.begin(), N_res.end(), 0.0);
N_res_mean = N_res_sum / N_res.size();
N_res_sq_sum = std::inner_product(N_res.begin(), N_res.end(), N_res.begin(), 0.0);
N_res_stdev = std::sqrt(N_res_sq_sum / N_res.size() - N_res_mean * N_res_mean);
Up_res_sum = std::accumulate(Up_res.begin(), Up_res.end(), 0.0);
Up_res_mean = Up_res_sum / Up_res.size();
Up_res_sq_sum = std::inner_product(Up_res.begin(), Up_res.end(), Up_res.begin(), 0.0);
Up_res_stdev = std::sqrt(Up_res_sq_sum / Up_res.size() - Up_res_mean * Up_res_mean);
//compute ENU average and standard deviation over the first 1000 residual solutions PRECISION
E_res_sum_precision = std::accumulate(E_res_precision.begin(), E_res_precision.end(), 0.0);
E_res_mean_precision = E_res_sum_precision / E_res_precision.size();
E_res_sq_sum_precision = std::inner_product(E_res_precision.begin(), E_res_precision.end(), E_res_precision.begin(), 0.0);
E_res_stdev_precision = std::sqrt(E_res_sq_sum_precision / E_res_precision.size() - E_res_mean_precision * E_res_mean_precision);
N_res_sum_precision = std::accumulate(N_res_precision.begin(), N_res_precision.end(), 0.0);
N_res_mean_precision = N_res_sum_precision / N_res_precision.size();
N_res_sq_sum_precision = std::inner_product(N_res_precision.begin(), N_res_precision.end(), N_res_precision.begin(), 0.0);
N_res_stdev_precision = std::sqrt(N_res_sq_sum_precision / N_res_precision.size() - N_res_mean_precision * N_res_mean_precision);
Up_res_sum_precision = std::accumulate(Up_res_precision.begin(), Up_res_precision.end(), 0.0);
Up_res_mean_precision = Up_res_sum_precision / Up_res_precision.size();
Up_res_sq_sum_precision = std::inner_product(Up_res_precision.begin(), Up_res_precision.end(), Up_res_precision.begin(), 0.0);
Up_res_stdev_precision = std::sqrt(Up_res_sq_sum_precision / Up_res_precision.size() - Up_res_mean_precision * Up_res_mean_precision);
DRMS = std::sqrt( E_res_stdev * E_res_stdev + N_res_stdev*N_res_stdev);
DUE_DRMS = 2*DRMS;
CEP = 0.62 * N_res_stdev + 0.56 * E_res_stdev;
MRSE = std::sqrt( E_res_stdev * E_res_stdev + N_res_stdev*N_res_stdev + Up_res_stdev*Up_res_stdev);
SEP = 0.51* (E_res_stdev * E_res_stdev + N_res_stdev*N_res_stdev + Up_res_stdev*Up_res_stdev);
file_statitics << std::endl;
file_statitics << "ACCURACY (respect true position)" << std::endl;
file_statitics << "East offset [m] = " << E_res_mean << ", East st. dev = " << E_res_stdev << std::endl;
file_statitics << "Nord offset [m] = " << N_res_mean << ",Noth st. dev = " << N_res_stdev << std::endl;
file_statitics << "Up offset [m] = " << Up_res_mean << ", Up st. dev = " << Up_res_stdev << std::endl;
file_statitics << std::endl;
file_statitics << "DRMS= " << DRMS << std::endl;
file_statitics << "DUE_DRMS= " << DUE_DRMS << std::endl;
file_statitics << "CEP= " << CEP << std::endl;
file_statitics << "MRSE= " << MRSE << std::endl;
file_statitics << "SEP= " << SEP << std::endl;
file_statitics << std::endl;
file_statitics << "PRECISION (respect average solution)" << std::endl;
file_statitics << "East offset [m] = " << E_res_mean_precision << ", East st. dev = " << E_res_stdev_precision << std::endl;
file_statitics << "Nord offset [m] = " << N_res_mean_precision << ", ,Noth st. dev = " << N_res_stdev_precision << std::endl;
file_statitics << "Up offset [m]= " << Up_res_mean_precision << ", Up st. dev = " << Up_res_stdev_precision << std::endl;
file_statitics << "----------------------------------------------------------------------------------------------" << std::endl;
file_solutions.close();
file_statitics.close();
} / * END *********************************** COMPUTE STATISTICS OVER THE FIRST 1000 SOLUTOINS ********************************** /
//ToDo: Implement Galileo RINEX and Galileo NMEA outputs //ToDo: Implement Galileo RINEX and Galileo NMEA outputs
// d_nmea_printer->Print_Nmea_Line(d_ls_pvt, d_flag_averaging); // d_nmea_printer->Print_Nmea_Line(d_ls_pvt, d_flag_averaging);
// //
// if (!b_rinex_header_writen) // & we have utc data in nav message! if (!b_rinex_header_writen) // & we have utc data in nav message!
// { {
// std::map<int,Gps_Ephemeris>::iterator gps_ephemeris_iter; std::map<int, Galileo_Ephemeris>::iterator galileo_ephemeris_iter;
// gps_ephemeris_iter = d_ls_pvt->gps_ephemeris_map.begin(); galileo_ephemeris_iter = d_ls_pvt->galileo_ephemeris_map.begin();
// if (gps_ephemeris_iter != d_ls_pvt->gps_ephemeris_map.end()) std::map<int, Gps_Ephemeris>::iterator gps_ephemeris_iter;
// { gps_ephemeris_iter = d_ls_pvt->gps_ephemeris_map.begin();
// rp->rinex_obs_header(rp->obsFile, gps_ephemeris_iter->second,d_rx_time); if ((galileo_ephemeris_iter != d_ls_pvt->galileo_ephemeris_map.end()) && (gps_ephemeris_iter != d_ls_pvt->gps_ephemeris_map.end()) )
// rp->rinex_nav_header(rp->navFile,d_ls_pvt->gps_iono, d_ls_pvt->gps_utc_model); {
// b_rinex_header_writen = true; // do not write header anymore if (arrived_galileo_almanac)
// } {
// } //rp->rinex_obs_header(rp->obsFile, galileo_ephemeris_iter->second, d_rx_time);
// if(b_rinex_header_writen) // Put here another condition to separate annotations (e.g 30 s) rp->rinex_nav_header(rp->navMixFile, d_ls_pvt->gps_iono, d_ls_pvt->gps_utc_model, d_ls_pvt->galileo_iono, d_ls_pvt->galileo_utc_model, d_ls_pvt->galileo_almanac);
b_rinex_header_writen = true; // do not write header anymore
}
}
}
//if(b_rinex_header_writen) // Put here another condition to separate annotations (e.g 30 s)
// { // {
// // Limit the RINEX navigation output rate to 1/6 seg // // Limit the RINEX navigation output rate to 1/6 seg
// // Notice that d_sample_counter period is 1ms (for GPS correlators) // // Notice that d_sample_counter period is 4ms (for Galileo correlators)
// if ((d_sample_counter-d_last_sample_nav_output)>=6000) // if ((d_sample_counter - d_last_sample_nav_output) >= 6000)
// { // {
// rp->log_rinex_nav(rp->navFile, d_ls_pvt->gps_ephemeris_map); // rp->log_rinex_nav(rp->navGalFile, d_ls_pvt->galileo_ephemeris_map);
// d_last_sample_nav_output=d_sample_counter; // d_last_sample_nav_output = d_sample_counter;
// } // }
// std::map<int,Gps_Ephemeris>::iterator gps_ephemeris_iter; // std::map<int, Galileo_Ephemeris>::iterator galileo_ephemeris_iter;
// gps_ephemeris_iter = d_ls_pvt->gps_ephemeris_map.begin(); // galileo_ephemeris_iter = d_ls_pvt->galileo_ephemeris_map.begin();
// if (gps_ephemeris_iter != d_ls_pvt->gps_ephemeris_map.end()) // if (galileo_ephemeris_iter != d_ls_pvt->galileo_ephemeris_map.end())
// { // {
// rp->log_rinex_obs(rp->obsFile, gps_ephemeris_iter->second, d_rx_time, gnss_pseudoranges_map); // rp->log_rinex_obs(rp->obsFile, galileo_ephemeris_iter->second, d_rx_time, gnss_pseudoranges_map);
// }
// } // }
// } */
} }
} }

93
src/algorithms/PVT/gnuradio_blocks/hybrid_pvt_cc.h
View File

@ -116,99 +116,6 @@ private:
std::shared_ptr<Nmea_Printer> d_nmea_printer; std::shared_ptr<Nmea_Printer> d_nmea_printer;
double d_rx_time; double d_rx_time;
double d_TOW_at_curr_symbol_constellation; double d_TOW_at_curr_symbol_constellation;
/* *********variable used for final statistics****/
std::vector<double> GDOP_vect;
std::vector<double> PDOP_vect;
double PDOP_sum;
double PDOP_mean;
double GDOP_sum;
double GDOP_mean;
std::vector<double> tot_obs_vect;
std::vector<double> Gal_obs_vect;
std::vector<double> GPS_obs_vect;
std::vector<boost::posix_time::ptime> time_vect;
std::vector<double> longitude_vect;
std::vector<double> latitude_vect;
std::vector<double> h_vect;
std::vector<double> X_vect;
std::vector<double> X_vect_res;
std::vector<double> X_vect_res_precision;
std::vector<double> Y_vect;
std::vector<double> Y_vect_res;
std::vector<double> Y_vect_res_precision;
std::vector<double> Z_vect;
std::vector<double> Z_vect_res;
std::vector<double> Z_vect_res_precision;
std::vector<double> E_res;
std::vector<double> N_res;
std::vector<double> Up_res;
std::vector<double> E_res_precision;
std::vector<double> N_res_precision;
std::vector<double> Up_res_precision;
double longitude_vect_sum;
double longitude_mean;
double latitude_vect_sum;
double latitude_mean;
double X_vect_sum;
double X_vect_mean;
double X_vect_sq_sum;
double X_vect_stdev;
double Y_vect_sum;
double Y_vect_mean;
double Y_vect_sq_sum;
double Y_vect_stdev;
double Z_vect_sum;
double Z_vect_mean;
double Z_vect_sq_sum;
double Z_vect_stdev;
//precision
double E_res_sum_precision;
double E_res_mean_precision;
double E_res_sq_sum_precision;
double E_res_stdev_precision;
double N_res_sum_precision;
double N_res_mean_precision;
double N_res_sq_sum_precision;
double N_res_stdev_precision;
double Up_res_sum_precision;
double Up_res_mean_precision;
double Up_res_sq_sum_precision;
double Up_res_stdev_precision;
//accuracy
double E_res_sum;
double E_res_mean;
double E_res_sq_sum;
double E_res_stdev;
double N_res_sum;
double N_res_mean;
double N_res_sq_sum;
double N_res_stdev;
double Up_res_sum;
double Up_res_mean;
double Up_res_sq_sum;
double Up_res_stdev;
double DRMS; //2D-65%
double DUE_DRMS; //2D-95%
double CEP; //2D-50%
double MRSE; //3D-61%
double SEP; //3D-50%
/* END variable used for final statistics */
std::shared_ptr<hybrid_ls_pvt> d_ls_pvt; std::shared_ptr<hybrid_ls_pvt> d_ls_pvt;
bool pseudoranges_pairCompare_min(std::pair<int,Gnss_Synchro> a, std::pair<int,Gnss_Synchro> b); bool pseudoranges_pairCompare_min(std::pair<int,Gnss_Synchro> a, std::pair<int,Gnss_Synchro> b);

77
src/algorithms/PVT/libs/hybrid_ls_pvt.cc
View File

@ -236,16 +236,14 @@ bool hybrid_ls_pvt::get_PVT(std::map<int,Gnss_Synchro> gnss_pseudoranges_map, do
std::map<int,Gps_Ephemeris>::iterator gps_ephemeris_iter; std::map<int,Gps_Ephemeris>::iterator gps_ephemeris_iter;
int valid_pseudoranges = gnss_pseudoranges_map.size(); int valid_pseudoranges = gnss_pseudoranges_map.size();
arma::mat W = arma::eye(valid_pseudoranges, valid_pseudoranges); //channels weights matrix arma::mat W = arma::eye(valid_pseudoranges, valid_pseudoranges); // channels weights matrix
arma::vec obs = arma::zeros(valid_pseudoranges); // pseudoranges observation vector arma::vec obs = arma::zeros(valid_pseudoranges); // pseudoranges observation vector
arma::mat satpos = arma::zeros(3, valid_pseudoranges); //satellite positions matrix arma::mat satpos = arma::zeros(3, valid_pseudoranges); // satellite positions matrix
int Galileo_week_number = 0; int Galileo_week_number = 0;
int GPS_week; int GPS_week;
double utc = 0; double utc = 0;
double utc_tx_corrected = 0; //utc computed at tx_time_corrected, added for Galileo constellation (in GPS utc is directly computed at TX_time_corrected_s) double utc_tx_corrected = 0; //utc computed at tx_time_corrected, added for Galileo constellation (in GPS utc is directly computed at TX_time_corrected_s)
//double SV_clock_drift_s = 0;
// double SV_relativistic_clock_corr_s = 0;
double TX_time_corrected_s; double TX_time_corrected_s;
double SV_clock_bias_s = 0; double SV_clock_bias_s = 0;
@ -282,26 +280,14 @@ bool hybrid_ls_pvt::get_PVT(std::map<int,Gnss_Synchro> gnss_pseudoranges_map, do
//double day_in_week = 7; //double day_in_week = 7;
// t = WN*sec_in_day*day_in_week + TOW; // t is Galileo System Time to use to compute satellite positions // t = WN*sec_in_day*day_in_week + TOW; // t is Galileo System Time to use to compute satellite positions
//JAVIER VERSION:
double Rx_time = hybrid_current_time; double Rx_time = hybrid_current_time;
//std::cout<<"Gal_Rx_time = "<< Rx_time << std::endl;
//to compute satellite position we need GST = WN+TOW (everything expressed in seconds)
//double Rx_time = galileo_current_time + Galileo_week_number*sec_in_day*day_in_week;
double Tx_time = Rx_time - gnss_pseudoranges_iter->second.Pseudorange_m/GALILEO_C_m_s; double Tx_time = Rx_time - gnss_pseudoranges_iter->second.Pseudorange_m/GALILEO_C_m_s;
//std::cout<<"Gal_Tx_time = "<< Tx_time << std::endl;
// 2- compute the clock drift using the clock model (broadcast) for this SV // 2- compute the clock drift using the clock model (broadcast) for this SV
// SV_clock_drift_s = galileo_ephemeris_iter->second.sv_clock_drift(Tx_time);
SV_clock_bias_s = galileo_ephemeris_iter->second.sv_clock_drift(Tx_time); SV_clock_bias_s = galileo_ephemeris_iter->second.sv_clock_drift(Tx_time);
// 3- compute the relativistic clock drift using the clock model (broadcast) for this SV // 3- compute the current ECEF position for this SV using corrected TX time
//SV_relativistic_clock_corr_s = galileo_ephemeris_iter->second.sv_clock_relativistic_term(Tx_time);
// 4- compute the current ECEF position for this SV using corrected TX time
//SV_clock_bias_s = SV_clock_drift_s + SV_relativistic_clock_corr_s;
TX_time_corrected_s = Tx_time - SV_clock_bias_s; TX_time_corrected_s = Tx_time - SV_clock_bias_s;
//std::cout<<"Gal_TX_time_corrected_s = "<< TX_time_corrected_s << std::endl;
galileo_ephemeris_iter->second.satellitePosition(TX_time_corrected_s); galileo_ephemeris_iter->second.satellitePosition(TX_time_corrected_s);
satpos(0,obs_counter) = galileo_ephemeris_iter->second.d_satpos_X; satpos(0,obs_counter) = galileo_ephemeris_iter->second.d_satpos_X;
@ -337,7 +323,6 @@ bool hybrid_ls_pvt::get_PVT(std::map<int,Gnss_Synchro> gnss_pseudoranges_map, do
<< " [m] Y=" << galileo_ephemeris_iter->second.d_satpos_Y << " [m] Y=" << galileo_ephemeris_iter->second.d_satpos_Y
<< " [m] Z=" << galileo_ephemeris_iter->second.d_satpos_Z << " [m] Z=" << galileo_ephemeris_iter->second.d_satpos_Z
<< " [m] PR_obs=" << obs(obs_counter) << " [m]"; << " [m] PR_obs=" << obs(obs_counter) << " [m]";
//std::cout<<"E"<<galileo_ephemeris_iter->second.i_satellite_PRN<< " satellite position computed"<< std::endl;
} }
else // the ephemeris are not available for this SV else // the ephemeris are not available for this SV
@ -347,11 +332,8 @@ bool hybrid_ls_pvt::get_PVT(std::map<int,Gnss_Synchro> gnss_pseudoranges_map, do
obs(obs_counter) = 1; // to avoid algorithm problems (divide by zero) obs(obs_counter) = 1; // to avoid algorithm problems (divide by zero)
DLOG(INFO) << "No ephemeris data for SV " << gnss_pseudoranges_iter->first; DLOG(INFO) << "No ephemeris data for SV " << gnss_pseudoranges_iter->first;
} }
} }
else if (gnss_pseudoranges_iter->second.System == 'G') else if (gnss_pseudoranges_iter->second.System == 'G')
{ {
//std::cout << "Satellite System: " << gnss_pseudoranges_iter->second.System <<std::endl; //std::cout << "Satellite System: " << gnss_pseudoranges_iter->second.System <<std::endl;
@ -367,17 +349,12 @@ bool hybrid_ls_pvt::get_PVT(std::map<int,Gnss_Synchro> gnss_pseudoranges_map, do
// COMMON RX TIME PVT ALGORITHM MODIFICATION (Like RINEX files) // COMMON RX TIME PVT ALGORITHM MODIFICATION (Like RINEX files)
// first estimate of transmit time // first estimate of transmit time
double Rx_time = hybrid_current_time; double Rx_time = hybrid_current_time;
//std::cout<<"Gps_Rx_time = "<< Rx_time << std::endl;
double Tx_time = Rx_time - gnss_pseudoranges_iter->second.Pseudorange_m/GPS_C_m_s; double Tx_time = Rx_time - gnss_pseudoranges_iter->second.Pseudorange_m/GPS_C_m_s;
//std::cout<<"Gps_Tx_time = "<< Tx_time << std::endl;
// 2- compute the clock drift using the clock model (broadcast) for this SV // 2- compute the clock drift using the clock model (broadcast) for this SV
SV_clock_bias_s = gps_ephemeris_iter->second.sv_clock_drift(Tx_time); SV_clock_bias_s = gps_ephemeris_iter->second.sv_clock_drift(Tx_time);
// 3- compute the relativistic clock drift using the clock model (broadcast) for this SV // 3- compute the current ECEF position for this SV using corrected TX time
//SV_relativistic_clock_corr_s = gps_ephemeris_iter->second.sv_clock_relativistic_term(Tx_time);
// 4- compute the current ECEF position for this SV using corrected TX time
//SV_clock_bias_s = SV_clock_drift_s + SV_relativistic_clock_corr_s - gps_ephemeris_iter->second.d_TGD;
TX_time_corrected_s = Tx_time - SV_clock_bias_s; TX_time_corrected_s = Tx_time - SV_clock_bias_s;
gps_ephemeris_iter->second.satellitePosition(TX_time_corrected_s); gps_ephemeris_iter->second.satellitePosition(TX_time_corrected_s);
@ -397,15 +374,10 @@ bool hybrid_ls_pvt::get_PVT(std::map<int,Gnss_Synchro> gnss_pseudoranges_map, do
<< " [m] Y=" << gps_ephemeris_iter->second.d_satpos_Y << " [m] Y=" << gps_ephemeris_iter->second.d_satpos_Y
<< " [m] Z=" << gps_ephemeris_iter->second.d_satpos_Z << " [m] Z=" << gps_ephemeris_iter->second.d_satpos_Z
<< " [m] PR_obs=" << obs(obs_counter) << " [m]"; << " [m] PR_obs=" << obs(obs_counter) << " [m]";
//std::cout<<"G"<<gps_ephemeris_iter->second.i_satellite_PRN<< " satellite position computed"<< std::endl;
// compute the UTC time for this SV (just to print the asociated UTC timestamp) // compute the UTC time for this SV (just to print the asociated UTC timestamp)
GPS_week = gps_ephemeris_iter->second.i_GPS_week; GPS_week = gps_ephemeris_iter->second.i_GPS_week;
//std::cout<<"GPS_week = "<< GPS_week << std::endl;
utc = gps_utc_model.utc_time(TX_time_corrected_s, GPS_week); utc = gps_utc_model.utc_time(TX_time_corrected_s, GPS_week);
//std::cout << "GPS UTC at TX_time_corrected_s= " << utc << std::endl;
} }
else // the ephemeris are not available for this SV else // the ephemeris are not available for this SV
{ {
@ -414,18 +386,10 @@ bool hybrid_ls_pvt::get_PVT(std::map<int,Gnss_Synchro> gnss_pseudoranges_map, do
obs(obs_counter) = 1; // to avoid algorithm problems (divide by zero) obs(obs_counter) = 1; // to avoid algorithm problems (divide by zero)
DLOG(INFO) << "No ephemeris data for SV " << gnss_pseudoranges_iter->first; DLOG(INFO) << "No ephemeris data for SV " << gnss_pseudoranges_iter->first;
} }
} }
obs_counter++; obs_counter++;
} }
//std::cout<<"Number of observations in PVT = " << obs_counter << std::endl;
// ******************************************************************************** // ********************************************************************************
// ****** SOLVE LEAST SQUARES****************************************************** // ****** SOLVE LEAST SQUARES******************************************************
// ******************************************************************************** // ********************************************************************************
@ -436,11 +400,7 @@ bool hybrid_ls_pvt::get_PVT(std::map<int,Gnss_Synchro> gnss_pseudoranges_map, do
if (valid_obs >= 4) if (valid_obs >= 4)
{ {
//std::cout<<"LS SOLVER CALLED " << obs_counter << std::endl;
arma::vec mypos; arma::vec mypos;
//std::cout << "satpos=" << satpos;
//std::cout << "obs="<< obs;
//std::cout << "W=" << W;
DLOG(INFO) << "satpos=" << satpos; DLOG(INFO) << "satpos=" << satpos;
DLOG(INFO) << "obs="<< obs; DLOG(INFO) << "obs="<< obs;
DLOG(INFO) << "W=" << W; DLOG(INFO) << "W=" << W;
@ -474,12 +434,6 @@ bool hybrid_ls_pvt::get_PVT(std::map<int,Gnss_Synchro> gnss_pseudoranges_map, do
LOG(INFO) << "Hybrid Position at " << boost::posix_time::to_simple_string(p_time) LOG(INFO) << "Hybrid Position at " << boost::posix_time::to_simple_string(p_time)
<< " is Lat = " << d_latitude_d << " [deg], Long = " << d_longitude_d << " is Lat = " << d_latitude_d << " [deg], Long = " << d_longitude_d
<< " [deg], Height= " << d_height_m << " [m]"; << " [deg], Height= " << d_height_m << " [m]";
//std::cout <<"X = " << mypos(0) << "Y = " << mypos(1) << "Z = " << mypos(2) << std::endl;
X = (double)mypos(0);
Y = (double)mypos(1);
Z = (double)mypos(2);
int count_solutions;
// ###### Compute DOPs ######## // ###### Compute DOPs ########
@ -498,20 +452,7 @@ bool hybrid_ls_pvt::get_PVT(std::map<int,Gnss_Synchro> gnss_pseudoranges_map, do
F(2,1) = cos((GPS_TWO_PI*d_latitude_d)/360.0); F(2,1) = cos((GPS_TWO_PI*d_latitude_d)/360.0);
F(2,2) = sin((GPS_TWO_PI*d_latitude_d/360.0)); F(2,2) = sin((GPS_TWO_PI*d_latitude_d/360.0));
// std::fstream fs;
// fs.open ("append_all_solutions.txt", std::fstream::out | std::fstream::app); //file printed for debug purposes
// fs.setf(fs.fixed, fs.floatfield);
// fs << std::setprecision(5);
// fs << boost::posix_time::to_simple_string(p_time) <<'\t'<< '\t'<< X << '\t'<< '\t'<< Y << '\t'<< '\t'<< Z << '\t'<< '\t' <<
// d_latitude_d << '\t'<< '\t'<< d_longitude_d << '\t'<< '\t'<< d_height_m << '\t'<< '\t' <<
// valid_obs << '\t'<< '\t'<< valid_obs_GALILEO_counter <<'\t'<< '\t'<< valid_obs_GPS_counter << std::endl;
// fs.close();
// 2- Apply the rotation to the latest covariance matrix (available in ECEF from LS) // 2- Apply the rotation to the latest covariance matrix (available in ECEF from LS)
arma::mat Q_ECEF = d_Q.submat(0, 0, 2, 2); arma::mat Q_ECEF = d_Q.submat(0, 0, 2, 2);
arma::mat DOP_ENU = arma::zeros(3, 3); arma::mat DOP_ENU = arma::zeros(3, 3);
@ -866,6 +807,8 @@ void hybrid_ls_pvt::topocent(double *Az, double *El, double *D, arma::vec x, arm
*D = sqrt(dx(0)*dx(0) + dx(1)*dx(1) + dx(2)*dx(2)); *D = sqrt(dx(0)*dx(0) + dx(1)*dx(1) + dx(2)*dx(2));
} }
void hybrid_ls_pvt::tropo(double *ddr_m, double sinel, double hsta_km, double p_mb, double t_kel, double hum, double hp_km, double htkel_km, double hhum_km) void hybrid_ls_pvt::tropo(double *ddr_m, double sinel, double hsta_km, double p_mb, double t_kel, double hum, double hp_km, double htkel_km, double hhum_km)
{ {
/* Inputs: /* Inputs:

4
src/algorithms/PVT/libs/hybrid_ls_pvt.h
View File

@ -85,6 +85,7 @@ public:
std::map<int,Gps_Ephemeris> gps_ephemeris_map; //!< Map storing new Galileo_Ephemeris std::map<int,Gps_Ephemeris> gps_ephemeris_map; //!< Map storing new Galileo_Ephemeris
Galileo_Utc_Model galileo_utc_model; Galileo_Utc_Model galileo_utc_model;
Galileo_Iono galileo_iono; Galileo_Iono galileo_iono;
Galileo_Almanac galileo_almanac;
Gps_Utc_Model gps_utc_model; Gps_Utc_Model gps_utc_model;
Gps_Iono gps_iono; Gps_Iono gps_iono;
@ -96,9 +97,6 @@ public:
double d_latitude_d; //!< Latitude in degrees double d_latitude_d; //!< Latitude in degrees
double d_longitude_d; //!< Longitude in degrees double d_longitude_d; //!< Longitude in degrees
double d_height_m; //!< Height [m] double d_height_m; //!< Height [m]
double X;
double Y;
double Z;
//averaging //averaging
std::deque<double> d_hist_latitude_d; std::deque<double> d_hist_latitude_d;
std::deque<double> d_hist_longitude_d; std::deque<double> d_hist_longitude_d;

200
src/algorithms/PVT/libs/rinex_printer.cc
View File

@ -29,6 +29,7 @@
*/ */
#include "rinex_printer.h" #include "rinex_printer.h"
#include <unistd.h>
#include <algorithm> // for min and max #include <algorithm> // for min and max
#include <cmath> // for floor #include <cmath> // for floor
#include <cstdlib> // for getenv() #include <cstdlib> // for getenv()
@ -54,18 +55,20 @@ Rinex_Printer::Rinex_Printer()
obsfilename = Rinex_Printer::createFilename("RINEX_FILE_TYPE_OBS"); obsfilename = Rinex_Printer::createFilename("RINEX_FILE_TYPE_OBS");
sbsfilename = Rinex_Printer::createFilename("RINEX_FILE_TYPE_SBAS"); sbsfilename = Rinex_Printer::createFilename("RINEX_FILE_TYPE_SBAS");
navGalfilename = Rinex_Printer::createFilename("RINEX_FILE_TYPE_GAL_NAV"); navGalfilename = Rinex_Printer::createFilename("RINEX_FILE_TYPE_GAL_NAV");
navMixfilename = Rinex_Printer::createFilename("RINEX_FILE_TYPE_MIXED_NAV");
Rinex_Printer::navFile.open(navfilename, std::ios::out | std::ios::app); Rinex_Printer::navFile.open(navfilename, std::ios::out | std::ios::app);
Rinex_Printer::obsFile.open(obsfilename, std::ios::out | std::ios::app); Rinex_Printer::obsFile.open(obsfilename, std::ios::out | std::ios::app);
Rinex_Printer::sbsFile.open(sbsfilename, std::ios::out | std::ios::app); Rinex_Printer::sbsFile.open(sbsfilename, std::ios::out | std::ios::app);
Rinex_Printer::navGalFile.open(navGalfilename, std::ios::out | std::ios::app); Rinex_Printer::navGalFile.open(navGalfilename, std::ios::out | std::ios::app);
Rinex_Printer::navMixFile.open(navMixfilename, std::ios::out | std::ios::app);
// RINEX v3.00 codes // RINEX v3.00 codes
satelliteSystem["GPS"] = "G"; satelliteSystem["GPS"] = "G";
satelliteSystem["GLONASS"] = "R"; satelliteSystem["GLONASS"] = "R";
satelliteSystem["SBAS payload"] = "S"; satelliteSystem["SBAS payload"] = "S";
satelliteSystem["Galileo"] = "E"; satelliteSystem["Galileo"] = "E";
satelliteSystem["Compass"] = "C"; satelliteSystem["Beidou"] = "C";
observationCode["GPS_L1_CA"] = "1C"; // "1C" GPS L1 C/A observationCode["GPS_L1_CA"] = "1C"; // "1C" GPS L1 C/A
observationCode["GPS_L1_P"] = "1P"; // "1P" GPS L1 P observationCode["GPS_L1_P"] = "1P"; // "1P" GPS L1 P
@ -141,6 +144,11 @@ Rinex_Printer::Rinex_Printer()
version = 3; version = 3;
stringVersion = "3.01"; stringVersion = "3.01";
} }
if ( FLAGS_RINEX_version.compare("3.02") == 0 )
{
version = 3;
stringVersion = "3.02";
}
else if ( FLAGS_RINEX_version.compare("2.11") == 0 ) else if ( FLAGS_RINEX_version.compare("2.11") == 0 )
{ {
version = 2; version = 2;
@ -164,11 +172,12 @@ Rinex_Printer::Rinex_Printer()
Rinex_Printer::~Rinex_Printer() Rinex_Printer::~Rinex_Printer()
{ {
// close RINEX files // close RINEX files
long posn, poso, poss, posng; long posn, poso, poss, posng, posmn;
posn = navFile.tellp(); posn = navFile.tellp();
poso = obsFile.tellp(); poso = obsFile.tellp();
poss = sbsFile.tellp(); poss = sbsFile.tellp();
posng = navGalFile.tellp(); posng = navGalFile.tellp();
posmn = navMixFile.tellp();
Rinex_Printer::navFile.close(); Rinex_Printer::navFile.close();
Rinex_Printer::obsFile.close(); Rinex_Printer::obsFile.close();
Rinex_Printer::sbsFile.close(); Rinex_Printer::sbsFile.close();
@ -190,6 +199,10 @@ Rinex_Printer::~Rinex_Printer()
{ {
remove(navGalfilename.c_str()); remove(navGalfilename.c_str());
} }
if (posmn == 0)
{
remove(navMixfilename.c_str());
}
} }
@ -289,7 +302,18 @@ std::string Rinex_Printer::getLocalTime()
std::string line; std::string line;
line += std::string("GNSS-SDR"); line += std::string("GNSS-SDR");
line += std::string(12, ' '); line += std::string(12, ' ');
line += Rinex_Printer::leftJustify("CTTC", 20); //put a flag to let the user change this std::string username;
char c_username[20] = {0};
int nGet = getlogin_r(c_username, sizeof(c_username) - 1);
if (nGet == 0)
{
username = c_username;
}
else
{
username = "UNKNOWN USER";
}
line += Rinex_Printer::leftJustify(username, 20);
boost::gregorian::date today = boost::gregorian::day_clock::local_day(); boost::gregorian::date today = boost::gregorian::day_clock::local_day();
boost::local_time::time_zone_ptr zone(new boost::local_time::posix_time_zone("UTC")); boost::local_time::time_zone_ptr zone(new boost::local_time::posix_time_zone("UTC"));
@ -645,6 +669,139 @@ void Rinex_Printer::rinex_nav_header(std::ofstream& out, Gps_Iono iono, Gps_Utc_
out << line << std::endl; out << line << std::endl;
} }
void Rinex_Printer::rinex_nav_header(std::ofstream& out, Gps_Iono gps_iono, Gps_Utc_Model gps_utc_model, Galileo_Iono galileo_iono, Galileo_Utc_Model galileo_utc_model, Galileo_Almanac galileo_almanac)
{
std::string line;
stringVersion = "3.02";
version = 3;
// -------- Line 1
line = std::string(5, ' ');
line += stringVersion;
line += std::string(11, ' ');
line += std::string("N: GNSS NAV DATA");
line += std::string(4, ' ');
line += std::string("M: MIXED");
line += std::string(12, ' ');
line += std::string("RINEX VERSION / TYPE");
Rinex_Printer::lengthCheck(line);
out << line << std::endl;
// -------- Line 2
line.clear();
line += Rinex_Printer::getLocalTime();
line += std::string("PGM / RUN BY / DATE");
line += std::string(1, ' ');
Rinex_Printer::lengthCheck(line);
out << line << std::endl;
// -------- Line COMMENT
line.clear();
line += Rinex_Printer::leftJustify("GNSS NAVIGATION MESSAGE FILE GENERATED BY GNSS-SDR", 60);
line += Rinex_Printer::leftJustify("COMMENT", 20);
Rinex_Printer::lengthCheck(line);
out << line << std::endl;
// -------- Line COMMENT
line.clear();
std::string gnss_sdr_version(GNSS_SDR_VERSION);
line += "GNSS-SDR VERSION ";
line += Rinex_Printer::leftJustify(gnss_sdr_version, 43);
line += Rinex_Printer::leftJustify("COMMENT", 20);
Rinex_Printer::lengthCheck(line);
out << line << std::endl;
// -------- Line COMMENT
line.clear();
line += Rinex_Printer::leftJustify("See http://gnss-sdr.org", 60);
line += Rinex_Printer::leftJustify("COMMENT", 20);
Rinex_Printer::lengthCheck(line);
out << line << std::endl;
// -------- Line ionospheric info 1
line.clear();
line += std::string("GAL ");
line += std::string(1, ' ');
line += Rinex_Printer::rightJustify(Rinex_Printer::doub2for(galileo_iono.ai0_5, 10, 2), 12);
line += Rinex_Printer::rightJustify(Rinex_Printer::doub2for(galileo_iono.ai1_5, 10, 2), 12);
line += Rinex_Printer::rightJustify(Rinex_Printer::doub2for(galileo_iono.ai2_5, 10, 2), 12);
double zero = 0.0;
line += Rinex_Printer::rightJustify(Rinex_Printer::doub2for(zero, 10, 2), 12);
line += std::string(7, ' ');
line += Rinex_Printer::leftJustify("IONOSPHERIC CORR", 20);
Rinex_Printer::lengthCheck(line);
out << line << std::endl;
// -------- Line ionospheric info 2
line.clear();
line += std::string("GPSA");
line += std::string(1, ' ');
line += Rinex_Printer::rightJustify(Rinex_Printer::doub2for(gps_iono.d_alpha0, 10, 2), 12);
line += Rinex_Printer::rightJustify(Rinex_Printer::doub2for(gps_iono.d_alpha1, 10, 2), 12);
line += Rinex_Printer::rightJustify(Rinex_Printer::doub2for(gps_iono.d_alpha2, 10, 2), 12);
line += Rinex_Printer::rightJustify(Rinex_Printer::doub2for(gps_iono.d_alpha3, 10, 2), 12);
line += std::string(7, ' ');
line += Rinex_Printer::leftJustify("IONOSPHERIC CORR", 20);
Rinex_Printer::lengthCheck(line);
out << line << std::endl;
// -------- Line system time correction
line.clear();
line += std::string("GAUT");
line += Rinex_Printer::rightJustify(Rinex_Printer::doub2for(galileo_utc_model.A0_6, 16, 2), 18);
line += Rinex_Printer::rightJustify(Rinex_Printer::doub2for(galileo_utc_model.A1_6, 15, 2), 16);
line += Rinex_Printer::rightJustify(boost::lexical_cast<std::string>(galileo_utc_model.t0t_6), 7);
line += Rinex_Printer::rightJustify(boost::lexical_cast<std::string>(galileo_utc_model.WNot_6), 5);
line += std::string(10, ' ');
line += Rinex_Printer::leftJustify("TIME SYSTEM CORR", 20);
Rinex_Printer::lengthCheck(line);
out << line << std::endl;
// -------- Line system time correction 2
line.clear();
line += std::string("GPGA");
line += Rinex_Printer::rightJustify(Rinex_Printer::doub2for(galileo_almanac.A_0G_10, 16, 2), 18);
line += Rinex_Printer::rightJustify(Rinex_Printer::doub2for(galileo_almanac.A_1G_10, 15, 2), 16);
line += Rinex_Printer::rightJustify(boost::lexical_cast<std::string>(galileo_almanac.t_0G_10), 7);
line += Rinex_Printer::rightJustify(boost::lexical_cast<std::string>(galileo_almanac.WN_0G_10), 5);
line += std::string(10, ' ');
line += Rinex_Printer::leftJustify("TIME SYSTEM CORR", 20);
Rinex_Printer::lengthCheck(line);
out << line << std::endl;
// -------- Line system time correction 3
line.clear();
line += std::string("GPUT");
line += Rinex_Printer::rightJustify(Rinex_Printer::doub2for(gps_utc_model.d_A0, 16, 2), 18);
line += Rinex_Printer::rightJustify(Rinex_Printer::doub2for(gps_utc_model.d_A1, 15, 2), 16);
line += Rinex_Printer::rightJustify(boost::lexical_cast<std::string>(gps_utc_model.d_t_OT), 7);
line += Rinex_Printer::rightJustify(boost::lexical_cast<std::string>(gps_utc_model.i_WN_T + 1024), 5); // valid until 2019
line += std::string(10, ' ');
line += Rinex_Printer::leftJustify("TIME SYSTEM CORR", 20);
Rinex_Printer::lengthCheck(line);
out << line << std::endl;
// -------- Line 6 leap seconds
// For leap second information, see http://www.endruntechnologies.com/leap.htm
line.clear();
line += Rinex_Printer::rightJustify(boost::lexical_cast<std::string>(gps_utc_model.d_DeltaT_LS), 6);
line += Rinex_Printer::rightJustify(boost::lexical_cast<std::string>(gps_utc_model.d_DeltaT_LSF), 6);
line += Rinex_Printer::rightJustify(boost::lexical_cast<std::string>(gps_utc_model.i_WN_LSF), 6);
line += Rinex_Printer::rightJustify(boost::lexical_cast<std::string>(gps_utc_model.i_DN), 6);
line += std::string(36, ' ');
line += Rinex_Printer::leftJustify("LEAP SECONDS", 20);
Rinex_Printer::lengthCheck(line);
out << line << std::endl;
// -------- End of Header
line.clear();
line += std::string(60, ' ');
line += Rinex_Printer::leftJustify("END OF HEADER", 20);
Rinex_Printer::lengthCheck(line);
out << line << std::endl;
}
void Rinex_Printer::rinex_sbs_header(std::ofstream& out) void Rinex_Printer::rinex_sbs_header(std::ofstream& out)
{ {
@ -665,7 +822,18 @@ void Rinex_Printer::rinex_sbs_header(std::ofstream& out)
// -------- Line 2 // -------- Line 2
line.clear(); line.clear();
line += Rinex_Printer::leftJustify("GNSS-SDR", 20); line += Rinex_Printer::leftJustify("GNSS-SDR", 20);
line += Rinex_Printer::leftJustify("CTTC", 20); std::string username;
char c_username[20] = {0};
int nGet = getlogin_r(c_username, sizeof(c_username) - 1);
if (nGet == 0)
{
username = c_username;
}
else
{
username = "UNKNOWN USER";
}
line += Rinex_Printer::leftJustify(username, 20);
// Date of file creation (dd-mmm-yy hhmm) // Date of file creation (dd-mmm-yy hhmm)
boost::local_time::time_zone_ptr zone(new boost::local_time::posix_time_zone("UTC")); boost::local_time::time_zone_ptr zone(new boost::local_time::posix_time_zone("UTC"));
boost::local_time::local_date_time pt = boost::local_time::local_sec_clock::local_time(zone); boost::local_time::local_date_time pt = boost::local_time::local_sec_clock::local_time(zone);
@ -1299,7 +1467,17 @@ void Rinex_Printer::rinex_obs_header(std::ofstream& out, Gps_Ephemeris eph, doub
// -------- Line OBSERVER / AGENCY // -------- Line OBSERVER / AGENCY
line.clear(); line.clear();
std::string username=getenv("USER"); std::string username;
char c_username[20] = {0};
int nGet = getlogin_r(c_username, sizeof(c_username) - 1);
if (nGet == 0)
{
username = c_username;
}
else
{
username = "UNKNOWN USER";
}
line += leftJustify(username, 20); line += leftJustify(username, 20);
line += Rinex_Printer::leftJustify("CTTC", 40); // add flag and property line += Rinex_Printer::leftJustify("CTTC", 40); // add flag and property
line += Rinex_Printer::leftJustify("OBSERVER / AGENCY", 20); line += Rinex_Printer::leftJustify("OBSERVER / AGENCY", 20);
@ -1546,7 +1724,17 @@ void Rinex_Printer::rinex_obs_header(std::ofstream& out, Galileo_Ephemeris eph,
// -------- Line OBSERVER / AGENCY // -------- Line OBSERVER / AGENCY
line.clear(); line.clear();
std::string username=getenv("USER"); std::string username;
char c_username[20] = {0};
int nGet = getlogin_r(c_username, sizeof(c_username) - 1);
if (nGet == 0)
{
username = c_username;
}
else
{
username = "UNKNOWN USER";
}
line += leftJustify(username, 20); line += leftJustify(username, 20);
line += Rinex_Printer::leftJustify("CTTC", 40); // add flag and property line += Rinex_Printer::leftJustify("CTTC", 40); // add flag and property
line += Rinex_Printer::leftJustify("OBSERVER / AGENCY", 20); line += Rinex_Printer::leftJustify("OBSERVER / AGENCY", 20);

7
src/algorithms/PVT/libs/rinex_printer.h
View File

@ -88,6 +88,7 @@ public:
std::ofstream navFile ; //<! Output file stream for RINEX navigation data file std::ofstream navFile ; //<! Output file stream for RINEX navigation data file
std::ofstream sbsFile ; //<! Output file stream for RINEX SBAS raw data file std::ofstream sbsFile ; //<! Output file stream for RINEX SBAS raw data file
std::ofstream navGalFile ; //<! Output file stream for RINEX Galileo navigation data file std::ofstream navGalFile ; //<! Output file stream for RINEX Galileo navigation data file
std::ofstream navMixFile ; //<! Output file stream for RINEX Mixed navigation data file
/*! /*!
* \brief Generates the GPS Navigation Data header * \brief Generates the GPS Navigation Data header
@ -99,6 +100,11 @@ public:
*/ */
void rinex_nav_header(std::ofstream& out, Galileo_Iono iono, Galileo_Utc_Model utc_model, Galileo_Almanac galileo_almanac); void rinex_nav_header(std::ofstream& out, Galileo_Iono iono, Galileo_Utc_Model utc_model, Galileo_Almanac galileo_almanac);
/*!
* \brief Generates the Mixed (GPS/Galileo) Navigation Data header
*/
void rinex_nav_header(std::ofstream& out, Gps_Iono gps_iono, Gps_Utc_Model gps_utc_model, Galileo_Iono galileo_iono, Galileo_Utc_Model galileo_utc_model, Galileo_Almanac galileo_almanac);
/*! /*!
* \brief Generates the GPS Observation data header * \brief Generates the GPS Observation data header
*/ */
@ -194,6 +200,7 @@ private:
std::string obsfilename; std::string obsfilename;
std::string sbsfilename; std::string sbsfilename;
std::string navGalfilename; std::string navGalfilename;
std::string navMixfilename;
/* /*
* Generates the data for the PGM / RUN BY / DATE line * Generates the data for the PGM / RUN BY / DATE line