diff --git a/src/algorithms/libs/gnss_sdr_flags.cc b/src/algorithms/libs/gnss_sdr_flags.cc index 96ce9a705..0e8e1cf59 100644 --- a/src/algorithms/libs/gnss_sdr_flags.cc +++ b/src/algorithms/libs/gnss_sdr_flags.cc @@ -30,11 +30,16 @@ #include "gnss_sdr_flags.h" +#include // for exists #include #include +#include DEFINE_string(c, "-", "Path to the configuration file (if set, overrides --config_file)."); +DEFINE_string(config_file, std::string(GNSSSDR_INSTALL_DIR "/share/gnss-sdr/conf/default.conf"), + "Path to the configuration file."); + DEFINE_string(s, "-", "If defined, path to the file containing the signal samples (overrides the configuration file and --signal_source)."); @@ -62,71 +67,126 @@ DEFINE_double(pll_bw_hz, 0.0, "If defined, bandwidth of the PLL low pass filter, #if GFLAGS_GREATER_2_0 +static bool ValidateC(const char* flagname, const std::string & value) +{ + if (boost::filesystem::exists( value ) or value.compare("-") == 0) // value is ok + return true; + std::cout << "Invalid value for flag -" << flagname << ". The file '" << value << "' does not exist." << std::endl; + std::cout << "GNSS-SDR program ended." << std::endl; + return false; +} + +static bool ValidateConfigFile(const char* flagname, const std::string & value) +{ + if (boost::filesystem::exists( value ) or value.compare(std::string(GNSSSDR_INSTALL_DIR "/share/gnss-sdr/conf/default.conf")) == 0) // value is ok + return true; + std::cout << "Invalid value for flag -" << flagname << ". The file '" << value << "' does not exist." << std::endl; + std::cout << "GNSS-SDR program ended." << std::endl; + return false; +} + +static bool ValidateS(const char* flagname, const std::string & value) +{ + if (boost::filesystem::exists( value ) or value.compare("-") == 0) // value is ok + return true; + std::cout << "Invalid value for flag -" << flagname << ". The file '" << value << "' does not exist." << std::endl; + std::cout << "GNSS-SDR program ended." << std::endl; + return false; +} + +static bool ValidateSignalSource(const char* flagname, const std::string & value) +{ + if (boost::filesystem::exists( value ) or value.compare("-") == 0) // value is ok + return true; + std::cout << "Invalid value for flag -" << flagname << ". The file '" << value << "' does not exist." << std::endl; + std::cout << "GNSS-SDR program ended." << std::endl; + return false; +} + static bool ValidateDopplerMax(const char* flagname, int32_t value) { - if (value >= 0 && value < 1000000) // value is ok + const int32_t max_value = 1000000; + if (value >= 0 && value < max_value) // value is ok return true; - std::cout << "Invalid value for " << flagname << ": " << value << std::endl; + std::cout << "Invalid value for flag -" << flagname << ": " << value << ". Allowed range is 0 < " << flagname << " < " << max_value << " Hz." << std::endl; + std::cout << "GNSS-SDR program ended." << std::endl; return false; } static bool ValidateDopplerStep(const char* flagname, int32_t value) { - if (value >= 0 && value < 10000) // value is ok + const int32_t max_value = 10000; + if (value >= 0 && value < max_value) // value is ok return true; - std::cout << "Invalid value for " << flagname << ": " << value << std::endl; + std::cout << "Invalid value for flag -" << flagname << ": " << value << ". Allowed range is 0 < " << flagname << " < " << max_value << " Hz." << std::endl; + std::cout << "GNSS-SDR program ended." << std::endl; return false; } static bool ValidateCn0Samples(const char* flagname, int32_t value) { - if (value > 0 && value < 10000) // value is ok + const int32_t max_value = 10000; + if (value > 0 && value < max_value) // value is ok return true; - std::cout << "Invalid value for " << flagname << ": " << value << std::endl; + std::cout << "Invalid value for flag -" << flagname << ": " << value << ". Allowed range is 0 < " << flagname << " < " << max_value << " samples." << std::endl; + std::cout << "GNSS-SDR program ended." << std::endl; return false; } static bool ValidateCn0Min(const char* flagname, int32_t value) { - if (value > 0 && value < 100) // value is ok + const int32_t max_value = 100; + if (value > 0 && value < max_value) // value is ok return true; - std::cout << "Invalid value for " << flagname << ": " << value << std::endl; + std::cout << "Invalid value for flag -" << flagname << ": " << value << ". Allowed range is 0 < " << flagname << " < " << max_value << " dB-Hz." << std::endl; + std::cout << "GNSS-SDR program ended." << std::endl; return false; } static bool ValidateMaxLockFail(const char* flagname, int32_t value) { - if (value > 0 && value < 10000) // value is ok + const int32_t max_value = 10000; + if (value > 0 && value < max_value) // value is ok return true; - std::cout << "Invalid value for " << flagname << ": " << value << std::endl; + std::cout << "Invalid value for flag -" << flagname << ": " << value << ". Allowed range is 0 < " << flagname << " < " << max_value << " fails." << std::endl; + std::cout << "GNSS-SDR program ended." << std::endl; return false; } static bool ValidateCarrierLockTh(const char* flagname, double value) { - if (value > 0.0 && value < 1.508) // value is ok + const double max_value = 1.508; + if (value > 0.0 && value < max_value) // value is ok return true; - std::cout << "Invalid value for " << flagname << ": " << value << std::endl; + std::cout << "Invalid value for flag -" << flagname << ": " << value << ". Allowed range is 0 < " << flagname << " < " << max_value << " rad." << std::endl; + std::cout << "GNSS-SDR program ended." << std::endl; return false; } static bool ValidateDllBw(const char* flagname, double value) { - if (value >= 0.0 && value < 10000.0) // value is ok + const double max_value = 10000.0; + if (value >= 0.0 && value < max_value) // value is ok return true; - std::cout << "Invalid value for " << flagname << ": " << value << std::endl; + std::cout << "Invalid value for flag -" << flagname << ": " << value << ". Allowed range is 0 < " << flagname << " < " << max_value << " Hz." << std::endl; + std::cout << "GNSS-SDR program ended." << std::endl; return false; } static bool ValidatePllBw(const char* flagname, double value) { - if (value >= 0.0 && value < 10000.0) // value is ok + const double max_value = 10000.0; + if (value >= 0.0 && value < max_value) // value is ok return true; - std::cout << "Invalid value for " << flagname << ": " << value << std::endl; + std::cout << "Invalid value for flag -" << flagname << ": " << value << ". Allowed range is 0 < " << flagname << " < " << max_value << " Hz." << std::endl; + std::cout << "GNSS-SDR program ended." << std::endl; return false; } - +DEFINE_validator(c, &ValidateC); +DEFINE_validator(config_file, &ValidateConfigFile); +DEFINE_validator(s, &ValidateS); +DEFINE_validator(signal_source, &ValidateSignalSource); DEFINE_validator(doppler_max, &ValidateDopplerMax); DEFINE_validator(doppler_step, &ValidateDopplerStep); DEFINE_validator(cn0_samples, &ValidateCn0Samples); diff --git a/src/core/receiver/CMakeLists.txt b/src/core/receiver/CMakeLists.txt index aa96cb4dd..b87b2bd1a 100644 --- a/src/core/receiver/CMakeLists.txt +++ b/src/core/receiver/CMakeLists.txt @@ -90,8 +90,6 @@ if(ENABLE_FMCOMMS2) set(OPT_RECEIVER_INCLUDE_DIRS ${OPT_RECEIVER_INCLUDE_DIRS} ${IIO_INCLUDE_DIRS}) endif(ENABLE_FMCOMMS2) -add_definitions(-DGNSSSDR_INSTALL_DIR="${CMAKE_INSTALL_PREFIX}") - if(${PC_GNURADIO_RUNTIME_VERSION} VERSION_GREATER "3.7.13" ) add_definitions( -DGR_GREATER_38=1 ) endif(${PC_GNURADIO_RUNTIME_VERSION} VERSION_GREATER "3.7.13" ) diff --git a/src/core/receiver/control_thread.cc b/src/core/receiver/control_thread.cc index 073807b43..04a76b971 100644 --- a/src/core/receiver/control_thread.cc +++ b/src/core/receiver/control_thread.cc @@ -66,8 +66,6 @@ extern concurrent_queue global_gps_acq_assist_queue; using google::LogMessage; -DEFINE_string(config_file, std::string(GNSSSDR_INSTALL_DIR "/share/gnss-sdr/conf/default.conf"), - "Path to the configuration file"); ControlThread::ControlThread() { diff --git a/src/tests/CMakeLists.txt b/src/tests/CMakeLists.txt index 2994eb0c4..aac3aa98f 100644 --- a/src/tests/CMakeLists.txt +++ b/src/tests/CMakeLists.txt @@ -605,6 +605,7 @@ if(NOT ${ENABLE_PACKAGING}) target_link_libraries(control_thread_test ${Boost_LIBRARIES} ${GFlags_LIBS} ${GLOG_LIBRARIES} + ${GNURADIO_RUNTIME_LIBRARIES} ${GTEST_LIBRARIES} gnss_sp_libs gnss_system_parameters @@ -628,6 +629,7 @@ if(NOT ${ENABLE_PACKAGING}) target_link_libraries(flowgraph_test ${Boost_LIBRARIES} ${GFlags_LIBS} ${GLOG_LIBRARIES} + ${GNURADIO_RUNTIME_LIBRARIES} ${GTEST_LIBRARIES} gnss_sp_libs gnss_rx diff --git a/src/utils/front-end-cal/CMakeLists.txt b/src/utils/front-end-cal/CMakeLists.txt index 0402cfa45..73f1c1f2a 100644 --- a/src/utils/front-end-cal/CMakeLists.txt +++ b/src/utils/front-end-cal/CMakeLists.txt @@ -33,6 +33,7 @@ include_directories( ${CMAKE_SOURCE_DIR}/src/core/libs/supl/asn-supl ${CMAKE_SOURCE_DIR}/src/algorithms/acquisition/adapters ${CMAKE_SOURCE_DIR}/src/algorithms/acquisition/gnuradio_blocks + ${CMAKE_SOURCE_DIR}/src/algorithms/libs ${GLOG_INCLUDE_DIRS} ${GFlags_INCLUDE_DIRS} ${GNURADIO_RUNTIME_INCLUDE_DIRS} @@ -59,6 +60,7 @@ target_link_libraries(front_end_cal_lib ${MAC_LIBRARIES} ${VOLK_GNSSSDR_LIBRARIES} ${ORC_LIBRARIES} ${GNSS_SDR_OPTIONAL_LIBS} rx_core_lib + gnss_sdr_flags gnss_rx channel_fsm ) diff --git a/src/utils/front-end-cal/main.cc b/src/utils/front-end-cal/main.cc index 50af2fdd2..a482e29ce 100644 --- a/src/utils/front-end-cal/main.cc +++ b/src/utils/front-end-cal/main.cc @@ -32,26 +32,7 @@ #define FRONT_END_CAL_VERSION "0.0.1" #endif -#include -#include -#include // for ctime -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include +#include "front_end_cal.h" #include "concurrent_map.h" #include "concurrent_queue.h" #include "file_configuration.h" @@ -72,21 +53,32 @@ #include "galileo_utc_model.h" #include "sbas_ephemeris.h" #include "gnss_sdr_supl_client.h" +#include "gnss_sdr_flags.h" +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include // for ctime +#include +#include +#include +#include -#include "front_end_cal.h" - using google::LogMessage; DECLARE_string(log_dir); -std::string s1_(GNSSSDR_INSTALL_DIR); -std::string s2_("/share/gnss-sdr/conf/front-end-cal.conf"); -std::string s3_ = s1_ + s2_; - -DEFINE_string(config_file, s3_, - "Path to the file containing the configuration parameters"); - concurrent_map global_gps_ephemeris_map; concurrent_map global_gps_iono_map; concurrent_map global_gps_utc_model_map; diff --git a/src/utils/reproducibility/ieee-access/L2-access18.conf b/src/utils/reproducibility/ieee-access/L2-access18.conf new file mode 100644 index 000000000..fff74d193 --- /dev/null +++ b/src/utils/reproducibility/ieee-access/L2-access18.conf @@ -0,0 +1,154 @@ + +[GNSS-SDR] + +;######### GLOBAL OPTIONS ################## +;internal_fs_sps: Internal signal sampling frequency after the signal conditioning stage [samples per second]. +GNSS-SDR.internal_fs_sps=3000000 + + +;######### SIGNAL_SOURCE CONFIG ############ +;#implementation: Use [File_Signal_Source] or [UHD_Signal_Source] or [GN3S_Signal_Source] (experimental) +SignalSource.implementation=File_Signal_Source + +;#filename: path to file with the captured GNSS signal samples to be processed +SignalSource.filename=./data/L125_III1b_210s_L2_3Msps.bin ; <- Available at https://zenodo.org/record/1184601 + +;#item_type: Type and resolution for each of the signal samples. Use only gr_complex in this version. +SignalSource.item_type=ibyte + +;#sampling_frequency: Original Signal sampling frequency in [Hz] +SignalSource.sampling_frequency=3000000 + +;#samples: Number of samples to be processed. Notice that 0 indicates the entire file. +SignalSource.samples=0 + +;#repeat: Repeat the processing file. Disable this option in this version +SignalSource.repeat=false + +;#dump: Dump the Signal source data to a file. Disable this option in this version +SignalSource.dump=false + +;#enable_throttle_control: Enabling this option tells the signal source to keep the delay between samples in post processing. +; it helps to not overload the CPU, but the processing time will be longer. +SignalSource.enable_throttle_control=false + + +;######### SIGNAL_CONDITIONER CONFIG ############ + +SignalConditioner.implementation=Signal_Conditioner + +;######### DATA_TYPE_ADAPTER CONFIG ############ + +DataTypeAdapter.implementation=Ibyte_To_Complex + +;######### INPUT_FILTER CONFIG ############ + +InputFilter.implementation=Pass_Through + + +;######### RESAMPLER CONFIG ############ + +Resampler.implementation=Pass_Through +Resampler.item_type=gr_complex + + +;######### CHANNELS GLOBAL CONFIG ############ + +Channels_2S.count=10 + +Channels.in_acquisition=1 + +Channel0.signal=2S +Channel1.signal=2S +Channel2.signal=2S +Channel3.signal=2S +Channel4.signal=2S +Channel5.signal=2S +Channel6.signal=2S +Channel7.signal=2S +Channel8.signal=2S +Channel9.signal=2S + + +;######### ACQUISITION GLOBAL CONFIG ############ + +Acquisition_2S.implementation=GPS_L2_M_PCPS_Acquisition + +Acquisition_2S.item_type=gr_complex + +Acquisition_2S.doppler_max=4500 + +Acquisition_2S.doppler_step=125 + +Acquisition_2S.use_CFAR_algorithm=false + +Acquisition_2S.threshold=19.5 + +Acquisition_2S.blocking=true + + +;######### TRACKING GLOBAL CONFIG ############ + + +Tracking_2S.implementation=GPS_L2_M_DLL_PLL_Tracking + +Tracking_2S.item_type=gr_complex + +Tracking_2S.pll_bw_hz=4.0; + +Tracking_2S.dll_bw_hz=0.75; + +Tracking_2S.early_late_space_chips=0.5; + +Tracking_2S.dump=true + +Tracking_2S.dump_filename=./data/track_ch_ + +;######### TELEMETRY DECODER CONFIG ############ + +TelemetryDecoder_2S.implementation=GPS_L2C_Telemetry_Decoder + +;######### OBSERVABLES CONFIG ############ + +Observables.implementation=Hybrid_Observables + +;######### PVT CONFIG ############ + +PVT.implementation=RTKLIB_PVT + +PVT.positioning_mode=Single; options: Single, Static, Kinematic, PPP_Static, PPP_Kinematic +PVT.iono_model=OFF; options: OFF, Broadcast, SBAS, Iono-Free-LC, Estimate_STEC, IONEX +PVT.trop_model=OFF; options: OFF, Saastamoinen, SBAS, Estimate_ZTD, Estimate_ZTD_Grad + +;#output_rate_ms: Period between two PVT outputs. Notice that the minimum period is equal to the tracking integration time [ms] +PVT.output_rate_ms=100 + +;#display_rate_ms: Position console print (std::out) interval [ms]. Notice that output_rate_ms<=display_rate_ms. +PVT.display_rate_ms=500 + +;# KML, GeoJSON, NMEA and RTCM output configuration +;#dump_filename: Log path and filename without extension. Notice that PVT will add ".dat" to the binary dump and ".kml" to GoogleEarth dump. +PVT.dump_filename=./data/PVT + +;#nmea_dump_filename: NMEA log path and filename +PVT.nmea_dump_filename=./gnss_sdr_pvt.nmea + +;#flag_nmea_tty_port: Enable or disable the NMEA log to a serial TTY port (Can be used with real hardware or virtual one) +PVT.flag_nmea_tty_port=false + +;#nmea_dump_devname: serial device descriptor for NMEA logging +PVT.nmea_dump_devname=/dev/pts/4 + +PVT.flag_rtcm_server=false + +PVT.rtcm_tcp_port=2101 + +PVT.rtcm_station_id=1234 + +PVT.flag_rtcm_tty_port=false + +PVT.rtcm_dump_devname=/dev/pts/1 + +PVT.dump=true + +PVT.elevation_mask=5 diff --git a/src/utils/reproducibility/ieee-access/README.md b/src/utils/reproducibility/ieee-access/README.md new file mode 100644 index 000000000..3194502c5 --- /dev/null +++ b/src/utils/reproducibility/ieee-access/README.md @@ -0,0 +1,42 @@ +Continuous Reproducibility in GNSS Signal Processing +---------------------------------------------------- + + +This folder contains files required for the reproduction of the experiment proposed in: + +C. Fernández-Prades, J. Vilà-Valls, J. Arribas and A. Ramos, *Continuous Reproducibility in GNSS Signal Processing*, submitted to IEEE Access, Feb. 2018. + +The dataset used in this paper is available at + +The sample format is `ibyte`: Interleaved (I&Q) stream of samples of type signed integer, 8-bit two’s complement number ranging from -128 to 127.  + +The figure appearing in that paper can be automatically generated with the pipeline available at https://gitlab.com/gnss-sdr/gnss-sdr/pipelines + +After the **Build** stage, which compiles the source code in several versions of the most popular GNU/Linux distributions, and the **Test** stage, which executes GNSS-SDR’s QA code, the **Deploy** stage creates and publishes an image of a software container ready to execute the experiment. This container is available by doing: + +``` +$ docker pull carlesfernandez/docker-gnsssdr:access18 +``` + +Then, in the **Experiment** stage, a job installs the image created in the previous step, grabs the data file, executes the experiment and produces a figure with the obtained results. + +The steps to reproduce the experiment in your own machine (with [Docker](https://www.docker.com) already installed and running) are: + +``` +$ docker pull carlesfernandez/docker-gnsssdr:access18 +$ git clone https://github.com/gnss-sdr/gnss-sdr +$ cd gnss-sdr +$ git checkout next +$ mkdir -p exp-access18/data +$ cd ex-access18/data +$ curl https://zenodo.org/record/1184601/files/L2_signal_samples.tar.xz --output L2_signal_samples.tar.xz +$ tar xvfJ L2_signal_samples.tar.xz +$ echo "3a04c1eeb970776bb77f5e3b7eaff2df L2_signal_samples.tar.xz" > data.md5 +$ md5sum -c data.md5 +$ cd .. +$ cp ../src/utils/reproducibility/ieee-access18/L2-access18.conf . +$ cp ../src/utils/reproducibility/ieee-access18/plot_dump.m . +$ cp -r ../src/utils/matlab/libs/geoFunctions . +$ octave --no-gui plot_dump.m +$ epspdf Figure2.eps Figure2.pdf +``` diff --git a/src/utils/reproducibility/ieee-access/plot_dump.m b/src/utils/reproducibility/ieee-access/plot_dump.m new file mode 100644 index 000000000..6bbb388c5 --- /dev/null +++ b/src/utils/reproducibility/ieee-access/plot_dump.m @@ -0,0 +1,239 @@ +% /*! +% * \file plot_dump.m +% * \brief Read GNSS-SDR Tracking dump binary file and plot some internal +% variables +% * \author Antonio Ramos, 2018. antonio.ramos(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 . +% * +% * ------------------------------------------------------------------------- +% */ + +clear all; +clc; + +n_channel = 0; +symbol_period = 20e-3; +filename = 'track_ch_'; + +fontsize = 12; + +addpath('./data') % Path to gnss-sdr dump files (Tracking and PVT) +addpath('./geoFunctions') + +load([filename int2str(n_channel) '.mat']); +t = (0 : length(abs_P) - 1) * symbol_period; +hf = figure('visible', 'off'); +set(hf, 'paperorientation', 'landscape'); +subplot(3, 3, [1,3]) +plot(t, abs_E, t, abs_P, t, abs_L) +xlabel('Time [s]','fontname','Times','fontsize', fontsize) +ylabel('Correlation result','fontname','Times','fontsize', fontsize) +legend('Early', 'Prompt', 'Late') +grid on + + +subplot(3, 3, 7) +plot(Prompt_I./1000, Prompt_Q./1000, 'linestyle', 'none', 'marker', '.') +xlabel('I','fontname','Times','fontsize', fontsize) +ylabel('Q','fontname','Times','fontsize', fontsize) +axis equal +grid on + +subplot(3, 3, [4,6]) +plot(t, Prompt_I) +xlabel('Time [s]','fontname','Times','fontsize', fontsize) +ylabel('Navigation data bits','fontname','Times','fontsize', fontsize) +grid on + + +fileID = fopen('data/PVT_ls_pvt.dat', 'r'); +dinfo = dir('data/PVT_ls_pvt.dat'); +filesize = dinfo.bytes; +aux = 1; +while ne(ftell(fileID), filesize) + navsol.RX_time(aux) = fread(fileID, 1, 'double'); + navsol.X(aux) = fread(fileID, 1, 'double'); + navsol.Y(aux) = fread(fileID, 1, 'double'); + navsol.Z(aux) = fread(fileID, 1, 'double'); + navsol.user_clock(aux) = fread(fileID, 1, 'double'); + navsol.lat(aux) = fread(fileID, 1, 'double'); + navsol.long(aux) = fread(fileID, 1, 'double'); + navsol.height(aux) = fread(fileID, 1, 'double'); + aux = aux + 1; +end +fclose(fileID); + + +mean_Latitude=mean(navsol.lat); +mean_Longitude=mean(navsol.long); +mean_h=mean(navsol.height); +utmZone = findUtmZone(mean_Latitude,mean_Longitude); +[ref_X_cart,ref_Y_cart,ref_Z_cart]=geo2cart(dms2mat(deg2dms(mean_Latitude)), dms2mat(deg2dms(mean_Longitude)), mean_h, 5); +[mean_utm_X,mean_utm_Y,mean_utm_Z]=cart2utm(ref_X_cart,ref_Y_cart,ref_Z_cart,utmZone); + + +numPoints=length(navsol.X); +aux=0; +for n=1:numPoints + aux=aux+1; + [E(aux),N(aux),U(aux)]=cart2utm(navsol.X(n), navsol.Y(n), navsol.Z(n), utmZone); +end + +v_2d=[E;N].'; %2D East Nort position vectors +v_3d=[E;N;U].'; %2D East Nort position vectors + + +%% ACCURACY + +% 2D ------------------- + +sigma_E_accuracy=sqrt((1/(numPoints-1))*sum((v_2d(:,1)-mean_utm_X).^2)); +sigma_N_accuracy=sqrt((1/(numPoints-1))*sum((v_2d(:,2)-mean_utm_Y).^2)); + +sigma_ratio_2d_accuracy=sigma_N_accuracy/sigma_E_accuracy + +% if sigma_ratio=1 -> Prob in circle with r=DRMS -> 65% +DRMS_accuracy=sqrt(sigma_E_accuracy^2+sigma_N_accuracy^2) +% if sigma_ratio=1 -> Prob in circle with r=2DRMS -> 95% +TWO_DRMS_accuracy=2*DRMS_accuracy +% if sigma_ratio>0.3 -> Prob in circle with r=CEP -> 50% +CEP_accuracy=0.62*sigma_E_accuracy+0.56*sigma_N_accuracy + +% 3D ------------------- + +sigma_U_accuracy=sqrt((1/(numPoints-1))*sum((v_3d(:,3)-mean_utm_Z).^2)); + +% if sigma_ratio=1 -> Prob in circle with r=DRMS -> 50% +SEP_accuracy=0.51*sqrt(sigma_E_accuracy^2+sigma_N_accuracy^2+sigma_U_accuracy^2) + +% if sigma_ratio=1 -> Prob in circle with r=DRMS -> 61% +MRSE_accuracy=sqrt(sigma_E_accuracy^2+sigma_N_accuracy^2+sigma_U_accuracy^2) +% if sigma_ratio=1 -> Prob in circle with r=2DRMS -> 95% +TWO_MRSE_accuracy=2*MRSE_accuracy + + + +%% PRECISION + +% 2D analysis +% Simulated X,Y measurements +%v1=randn(1000,2); + +% 2D Mean and Variance +mean_2d = [mean(v_2d(:,1)) ; mean(v_2d(:,2))]; +sigma_2d = [sqrt(var(v_2d(:,1))) ; sqrt(var(v_2d(:,2)))]; + +sigma_ratio_2d=sigma_2d(2)/sigma_2d(1) + +% if sigma_ratio=1 -> Prob in circle with r=DRMS -> 65% +DRMS=sqrt(sigma_2d(1)^2+sigma_2d(2)^2) +% if sigma_ratio=1 -> Prob in circle with r=2DRMS -> 95% +TWO_DRMS=2*DRMS +% if sigma_ratio>0.3 -> Prob in circle with r=CEP -> 50% +CEP=0.62*sigma_2d(1)+0.56*sigma_2d(2) + + +% Mean and Variance +mean_3d=[mean(v_3d(:,1)) ; mean(v_3d(:,2)) ; mean(v_3d(:,3))]; +sigma_3d=[sqrt(var(v_3d(:,1))) ; sqrt(var(v_3d(:,2))) ; sqrt(var(v_3d(:,3)))]; + +% absolute mean error +% 2D + +error_2D_vec=[mean_utm_X-mean_2d(1) mean_utm_Y-mean_2d(2)]; +error_2D_m=norm(error_2D_vec) + +error_3D_vec=[mean_utm_X-mean_3d(1) mean_utm_Y-mean_3d(2) mean_utm_Z-mean_3d(3)]; +error_3D_m=norm(error_3D_vec) + +% RMSE 2D + +RMSE_X=sqrt(mean((v_3d(:,1)-mean_utm_X).^2)) +RMSE_Y=sqrt(mean((v_3d(:,2)-mean_utm_Y).^2)) +RMSE_Z=sqrt(mean((v_3d(:,3)-mean_utm_Z).^2)) + + +RMSE_2D=sqrt(mean((v_2d(:,1)-mean_utm_X).^2+(v_2d(:,2)-mean_utm_Y).^2)) + +RMSE_3D=sqrt(mean((v_3d(:,1)-mean_utm_X).^2+(v_3d(:,2)-mean_utm_Y).^2+(v_3d(:,3)-mean_utm_Z).^2)) + +% SCATTER PLOT +subplot(3,3,8) +scatter(v_2d(:,1)-mean_2d(1),v_2d(:,2)-mean_2d(2)); +hold on; + +plot(0,0,'k*'); + + +[x,y,z] = cylinder([TWO_DRMS TWO_DRMS],200); +plot(x(1,:),y(1,:),[0 0.6 0],'Color',[0 0.6 0]); +str = strcat('2DRMS=',num2str(TWO_DRMS), ' m'); +text(cosd(65)*TWO_DRMS,sind(65)*TWO_DRMS,str,'Color',[0 0.6 0]); + + +[x,y,z] = cylinder([CEP CEP],200); + +plot(x(1,:),y(1,:),'r--'); +str = strcat('CEP=',num2str(CEP), ' m'); +text(cosd(80)*CEP,sind(80)*CEP,str,'Color','r'); + +grid on +axis equal; +xlabel('North [m]','fontname','Times','fontsize', fontsize) +ylabel('East [m]','fontname','Times','fontsize', fontsize) + +% 3D analysis +% Simulated X,Y,Z measurements + +% if sigma_ratio=1 -> Prob in circle with r=DRMS -> 50% +SEP=0.51*sqrt(sigma_3d(1)^2+sigma_3d(2)^2+sigma_3d(3)^2) + +% if sigma_ratio=1 -> Prob in circle with r=DRMS -> 61% +MRSE=sqrt(sigma_3d(1)^2+sigma_3d(2)^2+sigma_3d(3)^2) +% if sigma_ratio=1 -> Prob in circle with r=2DRMS -> 95% +TWO_MRSE=2*MRSE + + + +% SCATTER PLOT +subplot(3,3,9) +scatter3(v_3d(:,1)-mean_3d(1),v_3d(:,2)-mean_3d(2), v_3d(:,3)-mean_3d(3)); + +hold on; + +[x,y,z] = sphere(); +hSurface=surf(MRSE*x,MRSE*y,MRSE*z); % sphere centered at origin + +set(hSurface,'facecolor','none','edgecolor',[0 0.6 0],'edgealpha',1,'facealpha',1); + +%axis equal; +xlabel('North [m]','fontname','Times','fontsize', fontsize) +ylabel('East [m]','fontname','Times','fontsize', fontsize) +zlabel('Up [m]','fontname','Times','fontsize', fontsize) +str = strcat('MRSE=',num2str(MRSE), ' m'); +text(cosd(45)*MRSE,sind(45)*MRSE,20,str,'Color',[0 0.6 0]); + +hh=findall(hf,'-property','FontName'); +set(hh,'FontName','Times'); +print(hf, 'Figure2.eps', '-depsc') +close(hf);