Implement Theil-Sen's method for fitting a line to noisy data (used in bandwidth estimation).

webrtc/modules/congestion_controller/delay_based_bwe.cc

 /*
  *  Copyright (c) 2016 The WebRTC project authors. All Rights Reserved.
  *
  *  Use of this source code is governed by a BSD-style license
  *  that can be found in the LICENSE file in the root of the source
  *  tree. An additional intellectual property rights grant can be found
  *  in the file PATENTS.  All contributing project authors may
  *  be found in the AUTHORS file in the root of the source tree.
  */
 
 #include "webrtc/modules/congestion_controller/delay_based_bwe.h"
 
 #include <algorithm>
 #include <cmath>
 #include <string>
 
 #include "webrtc/base/checks.h"
 #include "webrtc/base/constructormagic.h"
 #include "webrtc/base/logging.h"
 #include "webrtc/base/thread_annotations.h"
 #include "webrtc/modules/congestion_controller/include/congestion_controller.h"
 #include "webrtc/modules/pacing/paced_sender.h"
 #include "webrtc/modules/remote_bitrate_estimator/include/remote_bitrate_estimator.h"
+#include "webrtc/modules/remote_bitrate_estimator/test/bwe_test_logging.h"
 #include "webrtc/system_wrappers/include/field_trial.h"
 #include "webrtc/system_wrappers/include/metrics.h"
 #include "webrtc/typedefs.h"
 
 namespace {
 constexpr int kTimestampGroupLengthMs = 5;
 constexpr int kAbsSendTimeFraction = 18;
 constexpr int kAbsSendTimeInterArrivalUpshift = 8;
 constexpr int kInterArrivalShift =
     kAbsSendTimeFraction + kAbsSendTimeInterArrivalUpshift;
 constexpr double kTimestampToMs =
     1000.0 / static_cast<double>(1 << kInterArrivalShift);
 // This ssrc is used to fulfill the current API but will be removed
 // after the API has been changed.
 constexpr uint32_t kFixedSsrc = 0;
 constexpr int kInitialRateWindowMs = 500;
 constexpr int kRateWindowMs = 150;
 
+// Parameters for linear least squares fit of regression line to noisy data.
 constexpr size_t kDefaultTrendlineWindowSize = 15;
 constexpr double kDefaultTrendlineSmoothingCoeff = 0.9;
 constexpr double kDefaultTrendlineThresholdGain = 4.0;
 
+// Parameters for Theil-Sen robust fitting of line to noisy data.
+constexpr size_t kDefaultTheilSenWindowSize = 20;
+constexpr double kDefaultTheilSenThresholdGain = 4.0;
+
 const char kBitrateEstimateExperiment[] = "WebRTC-ImprovedBitrateEstimate";
 const char kBweTrendlineFilterExperiment[] = "WebRTC-BweTrendlineFilter";
+const char kBweTheilSenFilterExperiment[] = "WebRTC-BweTheilSenFilter";
 
 bool BitrateEstimateExperimentIsEnabled() {
   return webrtc::field_trial::FindFullName(kBitrateEstimateExperiment) ==
          "Enabled";
 }
 
 bool TrendlineFilterExperimentIsEnabled() {
   std::string experiment_string =
       webrtc::field_trial::FindFullName(kBweTrendlineFilterExperiment);
   // The experiment is enabled iff the field trial string begins with "Enabled".
   return experiment_string.find("Enabled") == 0;
 }
 
+bool TheilSenFilterExperimentIsEnabled() {
+  std::string experiment_string =
+      webrtc::field_trial::FindFullName(kBweTheilSenFilterExperiment);
+  // The experiment is enabled iff the field trial string begins with "Enabled".
+  return experiment_string.find("Enabled") == 0;
+}
+
 bool ReadTrendlineFilterExperimentParameters(size_t* window_points,

[brandtr, 2016/11/28 16:27:19]

window_size? Here and in function below.

[terelius, 2016/12/02 16:45:52]

Done.

                                              double* smoothing_coef,
                                              double* threshold_gain) {
   RTC_DCHECK(TrendlineFilterExperimentIsEnabled());
+  RTC_DCHECK(!TheilSenFilterExperimentIsEnabled());
+  RTC_DCHECK(window_points != nullptr);
+  RTC_DCHECK(smoothing_coef != nullptr);
+  RTC_DCHECK(threshold_gain != nullptr);
   std::string experiment_string =
       webrtc::field_trial::FindFullName(kBweTrendlineFilterExperiment);
   int parsed_values = sscanf(experiment_string.c_str(), "Enabled-%zu,%lf,%lf",
                              window_points, smoothing_coef, threshold_gain);
   if (parsed_values == 3) {
     RTC_CHECK_GT(*window_points, 1) << "Need at least 2 points to fit a line.";
     RTC_CHECK(0 <= *smoothing_coef && *smoothing_coef <= 1)
         << "Coefficient needs to be between 0 and 1 for weighted average.";
     RTC_CHECK_GT(*threshold_gain, 0) << "Threshold gain needs to be positive.";
     return true;
   }
   LOG(LS_WARNING) << "Failed to parse parameters for BweTrendlineFilter "
                      "experiment from field trial string. Using default.";
   *window_points = kDefaultTrendlineWindowSize;
   *smoothing_coef = kDefaultTrendlineSmoothingCoeff;
   *threshold_gain = kDefaultTrendlineThresholdGain;
   return false;
 }
 
+bool ReadTheilSenFilterExperimentParameters(size_t* window_points,
+                                            double* threshold_gain) {
+  RTC_DCHECK(!TrendlineFilterExperimentIsEnabled());
+  RTC_DCHECK(TheilSenFilterExperimentIsEnabled());
+  RTC_DCHECK(window_points != nullptr);
+  RTC_DCHECK(threshold_gain != nullptr);
+  std::string experiment_string =
+      webrtc::field_trial::FindFullName(kBweTheilSenFilterExperiment);
+  int parsed_values = sscanf(experiment_string.c_str(), "Enabled-%zu,%lf",
+                             window_points, threshold_gain);
+  if (parsed_values == 2) {
+    RTC_CHECK_GT(*window_points, 1) << "Need at least 2 points to fit a line.";
+    RTC_CHECK_GT(*threshold_gain, 0) << "Threshold gain needs to be positive.";
+    return true;
+  }
+  LOG(LS_WARNING) << "Failed to parse parameters for BweTheilSenFilter "
+                     "experiment from field trial string. Using default.";
+  *window_points = kDefaultTheilSenWindowSize;
+  *threshold_gain = kDefaultTheilSenThresholdGain;
+  return false;
+}
 }  // namespace
 
 namespace webrtc {
 
 DelayBasedBwe::BitrateEstimator::BitrateEstimator()
     : sum_(0),
       current_win_ms_(0),
       prev_time_ms_(-1),
       bitrate_estimate_(-1.0f),
       bitrate_estimate_var_(50.0f),
@@ skipping 67 matching lines @@
   return bitrate_sample;
 }
 
 rtc::Optional<uint32_t> DelayBasedBwe::BitrateEstimator::bitrate_bps() const {
   if (bitrate_estimate_ < 0.f)
     return rtc::Optional<uint32_t>();
   return rtc::Optional<uint32_t>(bitrate_estimate_ * 1000);
 }
 
 DelayBasedBwe::DelayBasedBwe(Clock* clock)
-    : clock_(clock),
+    : in_trendline_experiment_(TrendlineFilterExperimentIsEnabled()),
+      in_theil_sen_experiment_(TheilSenFilterExperimentIsEnabled()),
+      clock_(clock),
       inter_arrival_(),
       kalman_estimator_(),
       trendline_estimator_(),
       detector_(OverUseDetectorOptions()),
       receiver_incoming_bitrate_(),
       last_update_ms_(-1),
       last_seen_packet_ms_(-1),
       uma_recorded_(false),
       trendline_window_size_(kDefaultTrendlineWindowSize),
       trendline_smoothing_coeff_(kDefaultTrendlineSmoothingCoeff),
       trendline_threshold_gain_(kDefaultTrendlineThresholdGain),
-      in_trendline_experiment_(TrendlineFilterExperimentIsEnabled()) {
+      theil_sen_window_size_(kDefaultTheilSenWindowSize),
+      theil_sen_threshold_gain_(kDefaultTheilSenThresholdGain) {
   if (in_trendline_experiment_) {
     ReadTrendlineFilterExperimentParameters(&trendline_window_size_,
                                             &trendline_smoothing_coeff_,
                                             &trendline_threshold_gain_);
   }
+  if (in_theil_sen_experiment_) {
+    ReadTheilSenFilterExperimentParameters(&trendline_window_size_,
+                                           &trendline_threshold_gain_);
+  }
+
   network_thread_.DetachFromThread();
 }
 
 DelayBasedBwe::Result DelayBasedBwe::IncomingPacketFeedbackVector(
     const std::vector<PacketInfo>& packet_feedback_vector) {
   RTC_DCHECK(network_thread_.CalledOnValidThread());
   if (!uma_recorded_) {
     RTC_HISTOGRAM_ENUMERATION(kBweTypeHistogram,
                               BweNames::kSendSideTransportSeqNum,
                               BweNames::kBweNamesMax);
@@ skipping 17 matching lines @@
   // Reset if the stream has timed out.
   if (last_seen_packet_ms_ == -1 ||
       now_ms - last_seen_packet_ms_ > kStreamTimeOutMs) {
     inter_arrival_.reset(
         new InterArrival((kTimestampGroupLengthMs << kInterArrivalShift) / 1000,
                          kTimestampToMs, true));
     kalman_estimator_.reset(new OveruseEstimator(OverUseDetectorOptions()));
     trendline_estimator_.reset(new TrendlineEstimator(
         trendline_window_size_, trendline_smoothing_coeff_,
         trendline_threshold_gain_));
+    theil_sen_estimator_.reset(new TheilSenEstimator(
+        theil_sen_window_size_, theil_sen_threshold_gain_));
   }
   last_seen_packet_ms_ = now_ms;
 
   uint32_t send_time_24bits =
       static_cast<uint32_t>(
           ((static_cast<uint64_t>(info.send_time_ms) << kAbsSendTimeFraction) +
            500) /
           1000) &
       0x00FFFFFF;
   // Shift up send time to use the full 32 bits that inter_arrival works with,
   // so wrapping works properly.
   uint32_t timestamp = send_time_24bits << kAbsSendTimeInterArrivalUpshift;
 
   uint32_t ts_delta = 0;
   int64_t t_delta = 0;
   int size_delta = 0;
   if (inter_arrival_->ComputeDeltas(timestamp, info.arrival_time_ms, now_ms,
                                     info.payload_size, &ts_delta, &t_delta,
                                     &size_delta)) {
     double ts_delta_ms = (1000.0 * ts_delta) / (1 << kInterArrivalShift);
     if (in_trendline_experiment_) {
       trendline_estimator_->Update(t_delta, ts_delta_ms, info.arrival_time_ms);
       detector_.Detect(trendline_estimator_->trendline_slope(), ts_delta_ms,
                        trendline_estimator_->num_of_deltas(),
                        info.arrival_time_ms);
-
+    } else if (in_theil_sen_experiment_) {
+      theil_sen_estimator_->Update(t_delta, ts_delta_ms, info.arrival_time_ms);
+      detector_.Detect(theil_sen_estimator_->trendline_slope(), ts_delta_ms,
+                       theil_sen_estimator_->num_of_deltas(),
+                       info.arrival_time_ms);
     } else {
       kalman_estimator_->Update(t_delta, ts_delta_ms, size_delta,
                                 detector_.State(), info.arrival_time_ms);
       detector_.Detect(kalman_estimator_->offset(), ts_delta_ms,
                        kalman_estimator_->num_of_deltas(),
                        info.arrival_time_ms);
     }
   }
 
   int probing_bps = 0;
@@ skipping 18 matching lines @@
         UpdateEstimate(info.arrival_time_ms, now_ms, acked_bitrate_bps,
                        &result.target_bitrate_bps);
   }
   if (!result.updated &&
       (last_update_ms_ == -1 ||
        now_ms - last_update_ms_ > rate_control_.GetFeedbackInterval())) {
     result.updated =
         UpdateEstimate(info.arrival_time_ms, now_ms, acked_bitrate_bps,
                        &result.target_bitrate_bps);
   }
-  if (result.updated)
+  if (result.updated) {
     last_update_ms_ = now_ms;
+    BWE_TEST_LOGGING_PLOT(1, "target_bitrate_bps", now_ms,
+                          result.target_bitrate_bps);
+  }
 
   return result;
 }
 
 bool DelayBasedBwe::UpdateEstimate(int64_t arrival_time_ms,
                                    int64_t now_ms,
                                    rtc::Optional<uint32_t> acked_bitrate_bps,
                                    uint32_t* target_bitrate_bps) {
   // TODO(terelius): RateControlInput::noise_var is deprecated and will be
   // removed. In the meantime, we set it to zero.
@@ skipping 22 matching lines @@
   *bitrate_bps = rate_control_.LatestEstimate();
   return true;
 }
 
 void DelayBasedBwe::SetMinBitrate(int min_bitrate_bps) {
   // Called from both the configuration thread and the network thread. Shouldn't
   // be called from the network thread in the future.
   rate_control_.SetMinBitrate(min_bitrate_bps);
 }
 }  // namespace webrtc