webrtc常用算法之指数平滑和卡尔曼滤波

本文作者汇总了一下webrtc中的常用算法。关于算法的实现原理，列了一些网上资料，如果想深入研究的同学，可以自己看看。

主要算法有：

指数平滑、双线性插值、卡尔曼滤波、贝叶斯预测

今天先分享一下指数平滑和卡尔曼滤波的学习资料。

1、指数平滑，一次指数、二次指数平滑

主要用来平滑一些抖动数据，如rtt、丢包率等。参考链接：https://blog.51cto.com/u_15054042/4745278

公式：

webrtc中的实现：

float ExpFilter::Apply(float exp, float sample) {
  if (filtered_ == kValueUndefined) {
    // Initialize filtered value.
    filtered_ = sample;
  } else if (exp == 1.0) {
    filtered_ = alpha_ * filtered_ + (1 - alpha_) * sample;
  } else {
    float alpha = std::pow(alpha_, exp);
    filtered_ = alpha * filtered_ + (1 - alpha) * sample;
  }
  if (max_ != kValueUndefined && filtered_ > max_) {
    filtered_ = max_;
  }
  return filtered_;
}

基本思想：

       预测值是以前观测值的加权和，且对不同的数据给予不同的权，新数据给较大的权，旧数据给较小的权。

基本原理：

    指数平滑法是移动平均法中的一种，其特点在于给过去的观测值不一样的权重。较近期观测值的权重数比较远期观测值的权数要大。

根据平滑次数不同，指数平滑法分为一次指数平滑法、二次指数平滑法和三次指数平滑法等。

但它们的基本思想都是：预测值是以前观测值的加权和，且对不同的数据给予不同的权数，新数据给予较大的权数，旧数据给予较小的权数。

指数平滑法对实际序列具有平滑作用，权系数（平滑系数）a 越小，平滑作用越强，但对实际数据的变动反应较迟缓。 

指数平滑值出现滞后偏差的程度随着权系数（平滑系数）a 的增大而减少，但当时间序列的变动出现直线趋势时，用一次指数平滑法来进行预测仍将存在着明显的滞后偏差。

平滑之后的数据存在偏差怎么办？需要修正，修正的方法是在一次指数平滑的基础上再进行二次指数平滑。

2、卡尔曼滤波首推B站DR_CAN的教程

链接：https://www.bilibili.com/video/BV1ez4y1X7eR/

webrtc中卡尔曼算法的实现，主要用来平滑帧的抖动延时。

jitter_estimator.cc：

// Updates Kalman estimate of the channel.
// The caller is expected to sanity check the inputs.
void VCMJitterEstimator::KalmanEstimateChannel(int64_t frameDelayMS,
                                               int32_t deltaFSBytes) {
  double Mh[2];
  double hMh_sigma;
  double kalmanGain[2];
  double measureRes;
  double t00, t01;

  // Kalman filtering

  // Prediction
  // M = M + Q
  _thetaCov[0][0] += _Qcov[0][0];
  _thetaCov[0][1] += _Qcov[0][1];
  _thetaCov[1][0] += _Qcov[1][0];
  _thetaCov[1][1] += _Qcov[1][1];

  // Kalman gain
  // 这里 M 是PK-
  // h'是H的转置
  // K = M*h'/(sigma2n + h*M*h') = M*h'/(1 + h*M*h')
  // h = [dFS 1] 
  // Mh = M*h' 分子部分
  // hMh_sigma = h*M*h' + R 分母部分
  Mh[0] = _thetaCov[0][0] * deltaFSBytes + _thetaCov[0][1];
  Mh[1] = _thetaCov[1][0] * deltaFSBytes + _thetaCov[1][1];
  // sigma weights measurements with a small deltaFS as noisy and
  // measurements with large deltaFS as good
  if (_maxFrameSize < 1.0) {
    return;
  }
  double sigma = (300.0 * exp(-fabs(static_cast<double>(deltaFSBytes)) /
                              (1e0 * _maxFrameSize)) +
                  1) *
                 sqrt(_varNoise);
  if (sigma < 1.0) {
    sigma = 1.0;
  }
  hMh_sigma = deltaFSBytes * Mh[0] + Mh[1] + sigma;
  if ((hMh_sigma < 1e-9 && hMh_sigma >= 0) ||
      (hMh_sigma > -1e-9 && hMh_sigma <= 0)) {
    assert(false);
    return;
  }
  kalmanGain[0] = Mh[0] / hMh_sigma;
  kalmanGain[1] = Mh[1] / hMh_sigma;

  // Correction
  // theta = theta + K*(dT - h*theta)
  measureRes = frameDelayMS - (deltaFSBytes * _theta[0] + _theta[1]);
  _theta[0] += kalmanGain[0] * measureRes;
  _theta[1] += kalmanGain[1] * measureRes;

  if (_theta[0] < _thetaLow) {
    _theta[0] = _thetaLow;
  }

  // M = (I - K*h)*M
  t00 = _thetaCov[0][0];
  t01 = _thetaCov[0][1];
  _thetaCov[0][0] = (1 - kalmanGain[0] * deltaFSBytes) * t00 -
                    kalmanGain[0] * _thetaCov[1][0];
  _thetaCov[0][1] = (1 - kalmanGain[0] * deltaFSBytes) * t01 -
                    kalmanGain[0] * _thetaCov[1][1];
  _thetaCov[1][0] = _thetaCov[1][0] * (1 - kalmanGain[1]) -
                    kalmanGain[1] * deltaFSBytes * t00;
  _thetaCov[1][1] = _thetaCov[1][1] * (1 - kalmanGain[1]) -
                    kalmanGain[1] * deltaFSBytes * t01;

  // Covariance matrix, must be positive semi-definite.
  assert(_thetaCov[0][0] + _thetaCov[1][1] >= 0 &&
         _thetaCov[0][0] * _thetaCov[1][1] -
                 _thetaCov[0][1] * _thetaCov[1][0] >=
             0 &&
         _thetaCov[0][0] >= 0);
}

作者：音视频之路