From b9ac964654ef4aaf54833966fd0ae67fa10291e7 Mon Sep 17 00:00:00 2001 From: Luis Esteve Date: Thu, 16 May 2013 16:25:51 +0000 Subject: [PATCH] FEATURE FE02: Adding the function that calculates the threshold from false alarm probability. git-svn-id: https://svn.code.sf.net/p/gnss-sdr/code/trunk@357 64b25241-fba3-4117-9849-534c7e92360d --- conf/gnss-sdr.conf | 14 ++--- .../galileo_e1_pcps_ambiguous_acquisition.cc | 39 ++++++++++++-- .../galileo_e1_pcps_ambiguous_acquisition.h | 2 + .../adapters/gps_l1_ca_pcps_acquisition.cc | 54 +++++++++++++++---- .../adapters/gps_l1_ca_pcps_acquisition.h | 3 ++ .../gnuradio_blocks/pcps_acquisition_cc.cc | 2 +- src/algorithms/channel/adapters/channel.cc | 6 ++- 7 files changed, 96 insertions(+), 24 deletions(-) diff --git a/conf/gnss-sdr.conf b/conf/gnss-sdr.conf index da5e83ca8..3668412b7 100644 --- a/conf/gnss-sdr.conf +++ b/conf/gnss-sdr.conf @@ -247,42 +247,36 @@ Channel0.signal=1C ;#satellite: Satellite PRN ID for this channel. Disable this option to random search Channel0.satellite=15 -Channel0.repeat_satellite=false ;######### CHANNEL 1 CONFIG ############ Channel1.system=GPS Channel1.signal=1C Channel1.satellite=18 -Channel1.repeat_satellite=false ;######### CHANNEL 2 CONFIG ############ Channel2.system=GPS Channel2.signal=1C Channel2.satellite=16 -Channel2.repeat_satellite=false ;######### CHANNEL 3 CONFIG ############ Channel3.system=GPS Channel3.signal=1C Channel3.satellite=21 -Channel3.repeat_satellite=false ;######### CHANNEL 4 CONFIG ############ Channel4.system=GPS Channel4.signal=1C Channel4.satellite=3 -Channel4.repeat_satellite=false ;######### CHANNEL 5 CONFIG ############ Channel5.system=GPS Channel5.signal=1C ;Channel5.satellite=21 -;Channel5.repeat_satellite=false ;######### ACQUISITION GLOBAL CONFIG ############ @@ -301,16 +295,18 @@ Acquisition.sampled_ms=1 ;######### ACQUISITION CHANNELS CONFIG ###### ;######### ACQUISITION CH 0 CONFIG ############ -;#implementation: Acquisition algorithm selection for this channel: [GPS_L1_CA_PCPS_Acquisition] +;#implementation: Acquisition algorithm selection for this channel: [GPS_L1_CA_PCPS_Acquisition] or [Galileo_E1_PCPS_Ambiguous_Acquisition] Acquisition0.implementation=GPS_L1_CA_PCPS_Acquisition ;#threshold: Acquisition threshold Acquisition0.threshold=50 +;#pfa: Acquisition false alarm probability. This option overrides the threshold option. Only use with implementations: [GPS_L1_CA_PCPS_Acquisition] or [Galileo_E1_PCPS_Ambiguous_Acquisition] +Acquisition0.pfa=0.001 ;#doppler_max: Maximum expected Doppler shift [Hz] Acquisition0.doppler_max=10000 ;#doppler_max: Doppler step in the grid search [Hz] Acquisition0.doppler_step=250 -;#repeat_satellite: Use only jointly with the satellte PRN ID option. - +;#repeat_satellite: Use only jointly with the satellite PRN ID option. +Acquisition0.repeat_satellite = false ;######### ACQUISITION CH 1 CONFIG ############ Acquisition1.implementation=GPS_L1_CA_PCPS_Acquisition diff --git a/src/algorithms/acquisition/adapters/galileo_e1_pcps_ambiguous_acquisition.cc b/src/algorithms/acquisition/adapters/galileo_e1_pcps_ambiguous_acquisition.cc index 08aa92d28..53a08d877 100644 --- a/src/algorithms/acquisition/adapters/galileo_e1_pcps_ambiguous_acquisition.cc +++ b/src/algorithms/acquisition/adapters/galileo_e1_pcps_ambiguous_acquisition.cc @@ -39,6 +39,8 @@ #include #include #include +#include + using google::LogMessage; @@ -55,7 +57,7 @@ GalileoE1PcpsAmbiguousAcquisition::GalileoE1PcpsAmbiguousAcquisition( DLOG(INFO) << "role " << role; - item_type_ = configuration->property(role + ".item_type", + item_type_ = configuration_->property(role + ".item_type", default_item_type); fs_in_ = configuration_->property("GNSS-SDR.internal_fs_hz", 4000000); @@ -118,8 +120,17 @@ GalileoE1PcpsAmbiguousAcquisition::set_channel(unsigned int channel) void GalileoE1PcpsAmbiguousAcquisition::set_threshold(float threshold) { - threshold_ = threshold; - if (item_type_.compare("gr_complex") == 0) + float pfa = configuration_->property(role_+ boost::lexical_cast(channel_) + ".pfa", 0.0); + if(pfa==0.0) + { + threshold_ = threshold; + } + else + { + threshold_ = calculate_threshold(pfa); + } + + if (item_type_.compare("gr_complex") == 0) { acquisition_cc_->set_threshold(threshold_); } @@ -214,6 +225,28 @@ GalileoE1PcpsAmbiguousAcquisition::reset() } } +float GalileoE1PcpsAmbiguousAcquisition::calculate_threshold(float pfa) +{ + //Calculate the threshold + + unsigned int frequency_bins = 0; + for (int doppler = (int)(-doppler_max_); doppler <= (int)doppler_max_; doppler += doppler_step_) + { + frequency_bins++; + } + + DLOG(INFO) <<"Pfa = "<< pfa; + + unsigned int ncells = vector_length_*frequency_bins; + double exponent = 1/(double)ncells; + double val = pow(1.0-pfa,exponent); + boost::math::exponential_distribution mydist (0.5); + float threshold = (float)quantile(mydist,val); + + DLOG(INFO) << "Threshold = " << threshold; + + return threshold; +} void diff --git a/src/algorithms/acquisition/adapters/galileo_e1_pcps_ambiguous_acquisition.h b/src/algorithms/acquisition/adapters/galileo_e1_pcps_ambiguous_acquisition.h index dac5ffd0e..754682c13 100644 --- a/src/algorithms/acquisition/adapters/galileo_e1_pcps_ambiguous_acquisition.h +++ b/src/algorithms/acquisition/adapters/galileo_e1_pcps_ambiguous_acquisition.h @@ -148,6 +148,8 @@ private: unsigned int out_streams_; gr_msg_queue_sptr queue_; concurrent_queue *channel_internal_queue_; + + float calculate_threshold(float pfa); }; #endif /* GNSS_SDR_GALILEO_E1_PCPS_AMBIGUOUS_ACQUISITION_H_ */ diff --git a/src/algorithms/acquisition/adapters/gps_l1_ca_pcps_acquisition.cc b/src/algorithms/acquisition/adapters/gps_l1_ca_pcps_acquisition.cc index 256a46a4b..af896f1e8 100644 --- a/src/algorithms/acquisition/adapters/gps_l1_ca_pcps_acquisition.cc +++ b/src/algorithms/acquisition/adapters/gps_l1_ca_pcps_acquisition.cc @@ -40,6 +40,7 @@ #include #include #include +#include using google::LogMessage; @@ -47,26 +48,30 @@ GpsL1CaPcpsAcquisition::GpsL1CaPcpsAcquisition( ConfigurationInterface* configuration, std::string role, unsigned int in_streams, unsigned int out_streams, gr_msg_queue_sptr queue) : - role_(role), in_streams_(in_streams), out_streams_(out_streams), queue_( + role_(role), in_streams_(in_streams), out_streams_(out_streams), queue_( queue) { + configuration_=configuration; std::string default_item_type = "gr_complex"; std::string default_dump_filename = "./data/acquisition.dat"; DLOG(INFO) << "role " << role; - item_type_ = configuration->property(role + ".item_type", +// std::cout << "role " << role_ << std::endl; + + item_type_ = configuration_->property(role + ".item_type", default_item_type); - fs_in_ = configuration->property("GNSS-SDR.internal_fs_hz", 2048000); - if_ = configuration->property(role + ".ifreq", 0); - dump_ = configuration->property(role + ".dump", false); - shift_resolution_ = configuration->property(role + ".doppler_max", 15); - sampled_ms_ = configuration->property(role + ".sampled_ms", 1); - dump_filename_ = configuration->property(role + ".dump_filename", + fs_in_ = configuration_->property("GNSS-SDR.internal_fs_hz", 2048000); + if_ = configuration_->property(role + ".ifreq", 0); + dump_ = configuration_->property(role + ".dump", false); + shift_resolution_ = configuration_->property(role + ".doppler_max", 15); + sampled_ms_ = configuration_->property(role + ".sampled_ms", 1); + dump_filename_ = configuration_->property(role + ".dump_filename", default_dump_filename); + //--- Find number of samples per spreading code ------------------------- vector_length_ = round(fs_in_ / (GPS_L1_CA_CODE_RATE_HZ / GPS_L1_CA_CODE_LENGTH_CHIPS)); @@ -113,7 +118,16 @@ void GpsL1CaPcpsAcquisition::set_channel(unsigned int channel) void GpsL1CaPcpsAcquisition::set_threshold(float threshold) { - threshold_ = threshold; + float pfa = configuration_->property(role_+ boost::lexical_cast(channel_) + ".pfa", 0.0); + if(pfa==0.0) + { + threshold_ = threshold; + } + else + { + threshold_ = calculate_threshold(pfa); + } + if (item_type_.compare("gr_complex") == 0) { acquisition_cc_->set_threshold(threshold_); @@ -194,6 +208,28 @@ void GpsL1CaPcpsAcquisition::reset() } } +float GpsL1CaPcpsAcquisition::calculate_threshold(float pfa) +{ + //Calculate the threshold + + unsigned int frequency_bins = 0; + for (int doppler = (int)(-doppler_max_); doppler <= (int)doppler_max_; doppler += doppler_step_) + { + frequency_bins++; + } + + DLOG(INFO) <<"Pfa = "<< pfa; + + unsigned int ncells = vector_length_*frequency_bins; + double exponent = 1/(double)ncells; + double val = pow(1.0-pfa,exponent); + boost::math::exponential_distribution mydist (0.5); + float threshold = (float)quantile(mydist,val); + + DLOG(INFO) << "Threshold = " << threshold; + + return threshold; +} void GpsL1CaPcpsAcquisition::connect(gr_top_block_sptr top_block) { diff --git a/src/algorithms/acquisition/adapters/gps_l1_ca_pcps_acquisition.h b/src/algorithms/acquisition/adapters/gps_l1_ca_pcps_acquisition.h index 966bfb803..8697ad82e 100644 --- a/src/algorithms/acquisition/adapters/gps_l1_ca_pcps_acquisition.h +++ b/src/algorithms/acquisition/adapters/gps_l1_ca_pcps_acquisition.h @@ -127,6 +127,7 @@ public: void reset(); private: + ConfigurationInterface* configuration_; pcps_acquisition_cc_sptr acquisition_cc_; gr::blocks::stream_to_vector::sptr stream_to_vector_; size_t item_size_; @@ -150,6 +151,8 @@ private: unsigned int out_streams_; gr_msg_queue_sptr queue_; concurrent_queue *channel_internal_queue_; + + float calculate_threshold(float pfa); }; #endif /* GNSS_SDR_GPS_L1_CA_PCPS_ACQUISITION_H_ */ diff --git a/src/algorithms/acquisition/gnuradio_blocks/pcps_acquisition_cc.cc b/src/algorithms/acquisition/gnuradio_blocks/pcps_acquisition_cc.cc index f9423d8aa..e1757ad97 100644 --- a/src/algorithms/acquisition/gnuradio_blocks/pcps_acquisition_cc.cc +++ b/src/algorithms/acquisition/gnuradio_blocks/pcps_acquisition_cc.cc @@ -188,7 +188,7 @@ int pcps_acquisition_cc::general_work(int noutput_items, d_input_power = d_input_power / (float)d_fft_size; // 2- Doppler frequency search loop - for (doppler = (int)(-d_doppler_max); doppler < (int)d_doppler_max; doppler += d_doppler_step) + for (doppler = (int)(-d_doppler_max); doppler <= (int)d_doppler_max; doppler += d_doppler_step) { // doppler search steps // Perform the carrier wipe-off diff --git a/src/algorithms/channel/adapters/channel.cc b/src/algorithms/channel/adapters/channel.cc index 749b9abfb..a45ec3132 100644 --- a/src/algorithms/channel/adapters/channel.cc +++ b/src/algorithms/channel/adapters/channel.cc @@ -64,14 +64,16 @@ Channel::Channel(ConfigurationInterface *configuration, unsigned int channel, acq_->set_gnss_synchro(&gnss_synchro_); trk_->set_gnss_synchro(&gnss_synchro_); - acq_->set_threshold(configuration->property("Acquisition" - + boost::lexical_cast(channel_) + ".threshold", 0.0)); +// IMPORTANT: Do not change the order of the following 3 methods + acq_->set_doppler_max(configuration->property("Acquisition" + boost::lexical_cast(channel_) + ".doppler_max", 10000)); acq_->set_doppler_step(configuration->property("Acquisition" + boost::lexical_cast(channel_) + ".doppler_step", 250)); + acq_->set_threshold(configuration->property("Acquisition" + + boost::lexical_cast(channel_) + ".threshold", 0.0)); repeat_ = configuration->property("Acquisition" + boost::lexical_cast< std::string>(channel_) + ".repeat_satellite", false);