1）前言
- 前一篇文章分析了
FrameBuffer模塊對視頻幀的插入原理扳肛，以及出隊(送到解碼隊列)的機制走敌。 - 在出隊的過程中涉及到了很多和延遲相關的信息庭惜，沒有分析，諸如渲染時間的計算涝开、幀延遲的計算循帐、抖動的計算等都未進行相應的分析。
- 其中
FrameBuffer延遲JitterDleay對視頻流的單向延遲有重要的影響舀武，很大程度上決定了應用的實時性,本文結合參考文獻對JitterDleay的計算進行深入原理分析拄养。 - 本文首先介紹gcc算法的Arrival-time model模型
- 其次介紹Arrival-time filter（卡爾曼濾波）的公式推導
- 最后分析結合代碼分析卡爾曼濾波在
JitterDleay中的應用，以及了解jitter delay究竟描述的是個啥東西
2）Arrival-time model模型介紹
image-20210610215218805.png
- 上圖為網(wǎng)絡傳輸過程中相鄰兩幀數(shù)據(jù)的到達時間模型
- 已知T(i-1)為第(i-1)幀數(shù)據(jù)發(fā)送時間银舱，T(i)為第(i)幀數(shù)據(jù)發(fā)送時間
- t(i)為第i幀數(shù)據(jù)的接收時間瘪匿，t(i-1)為第(i-1)幀數(shù)據(jù)的接收時間
- 理論上若傳輸過程中無噪聲（網(wǎng)絡噪聲）和其他因素（傳輸大幀引起的延遲）的影響，第i幀傳輸曲線應該如上述紅色線條（這樣jitter delay應該等于0）寻馏，但正由于網(wǎng)絡噪聲和傳輸大幀導致的影響使得第I幀傳輸?shù)臅r延遲可能大于第i-1幀棋弥，也或者小于第i-1幀，所以就有了jitter delay
- webrtc的
FrameBuffer模塊中通過獲取當前的幀間延遲.來告知解碼隊列何時解碼該幀比較合適诚欠，從而來讓接收端所接收到的視頻幀能夠較為平滑的進入到解碼器顽染，從而緩解卡頓現(xiàn)象 - 通過上圖可以理論建模如下：
d(i) = t(i) - t(i-1) - (T(i) - T(i-1)) (2.1.1)
= L(i)/C(i) - L(i-1)/C(i-1) + w(i)
L(i)-L(i-1)
= -------------- + w(i)
C(i)
= dL(i)/C(i) + w(i)
- d(i) 就是最終的幀間延遲
- 其中
C(i)為信道傳輸速率、L(i)為每幀數(shù)據(jù)的數(shù)據(jù)量轰绵，dL(i)為兩幀數(shù)據(jù)的數(shù)據(jù)量的差值粉寞，以上假定傳輸速率恒定 - 而
w(i)為是隨機過程w的一個樣本，它是C(i)左腔、當前發(fā)送數(shù)據(jù)量和當前發(fā)送速率的函數(shù)唧垦，并假定w為高斯白噪聲 - 當傳輸信道過載發(fā)送數(shù)據(jù)時，w(i)的方差將變大液样，傳輸信道空載時w(i)將隨之變小振亮，其他情況w(i)為0
- 若
t(i) - t(i-1) > T(i) - T(i-1)表示幀I相對幀i-1的延遲大了 - 若將網(wǎng)絡排隊延遲
m(i)從w(i)中提取出來使得過程噪聲的均值為0巧还，則會得出如下：
w(i) = m(i) + v(i)
d(i) = dL(i)/C(i) + m(i) + v(i)
-
m(i)表示網(wǎng)絡排隊延遲（如經(jīng)歷路由器的時候的排隊延遲），v(i)表示測量噪聲(如最大幀數(shù)據(jù)量和平均每幀數(shù)據(jù)量的計算誤差坊秸，時間同步等)狞悲。 - 以上就是幀間延遲的基本模型，在該模型中我們實際要求的是利用卡爾曼濾波來求得
C(i)和m(i) - 并通過迭代使得
C(i)和m(i)的誤差最小妇斤，從而得到最優(yōu)的d(i)幀間延遲
3）Kalman filter 建模及理論推導
3.1 卡爾曼濾波-建立狀態(tài)空間方程
theta(i) = A * theta(i-1) + u(i-1) P(u) ~ (0,Q)
= theta(i-1) + u(i-1) (3.1.1)
Q(i) = E{u_bar(i) * u_bar(i)^T}
diag(Q(i)) = [10^-13 10^-3]^T
theta(i)為狀態(tài)空間變量摇锋，其中A為狀態(tài)轉移矩陣，該案例取值為二維單位矩陣站超，u(i-1)為過程噪聲(無法測量)荸恕，其概率分布服從正太分布,數(shù)學期望為0，協(xié)方差矩陣為Q,推薦為對角矩陣
theta(i)包含兩個變量如下
theta_bar(i) = [1/C(i) m(i)]^T
-
1/C(i)為信道傳輸速率的倒數(shù) -
m(i)為幀i的網(wǎng)絡排隊時延(路由器或者交換機處理數(shù)據(jù)包排隊所消耗的時間) -
C(i)和m(i)也是我們實際要求得的目標值
3.2 卡爾曼濾波-建立觀測方程
- 觀測方程如下：
d(i) = H * theta(i) + v(i) P(V) ~ (0,R) (3.2.2)
- 上述觀測方程測量的是時間死相，為一維方程
- H為觀測系數(shù)矩陣融求，其定義如下：
h_bar(i) = [dL(i) 1]^T
H = h_bar(i)^T = [dL(i) 1]
dL(i)為第i幀和第i-1幀的數(shù)據(jù)量之差（detal L(i)）將其變形為google 官方公式如下：
d(i) = h_bar(i)^T * theta_bar(i) + v(i) (3.2.3)
- v(i)為觀測噪聲，同樣其概率分布服從正太分布算撮，數(shù)學期望為0生宛，協(xié)方差矩陣為R
variance var_v = sigma(v,i)^2
R(i) = E{v_bar(i) * v_bar(i)^T}
= var_v
- 將
h_bar向量代入到(3.2.3)式中
d(i) = h_bar(i)^T * theta_bar(i) + v(i)
= [dL(i) 1] * [1/C(i) m(i)]^T + v(i) (3.2.4)
- 有了狀態(tài)方程和測量方程則可以根據(jù)卡爾曼濾波的預測和校正模型進行估計
3.3 卡爾曼濾波-預測計算先驗估計
- 綜合上述公式可求得先驗估計值如下：
- 預測其實就是使用上一次的最優(yōu)結果預測當前的值
- 首先是計算先驗估計的誤差協(xié)方差
theta_hat^-(i) = theta_hat(i-1) + u(i-1) (3.3.1)
-
theta_hat^-(i)第i幀的先驗估計值 - 其次根據(jù)先驗估計值計算先驗估計的誤差協(xié)方差
3.4 卡爾曼濾波-預測計算先驗估計誤差協(xié)方差
e^-(i) = theta(i) - theta_hat^-(i) P(E(i)) ~ (0 , P) (3.4.1)
P^-(i) = {e^-(i) * e^-(i)^T}
= E{(theta(i) - theta_hat^-(i)) * (theta(i) - theta_hat^-(i))^T} (3.4.2)
= A * P(i-1) * A^T + Q (3.4.3)
= P(i-1) + Q
= E(i-1) + Q (3.4.4)
-
e^-(i)表示當前幀的實際值和估計值之間的誤差，也就是先驗估計的誤差 -
P^-(i)表示當前幀的先驗誤差協(xié)方差 -
P(i-1)為上一次的誤差協(xié)方差 - 也就是說先驗估計的誤差協(xié)方差等于上一次的后驗估計的誤差協(xié)方差 + 過程噪聲的的誤差協(xié)方差
- 有了先驗估計的誤差協(xié)方差之后就是計算當前幀的最優(yōu)卡爾曼增益為迭代做準備
3.5 卡爾曼濾波-校正計算卡爾曼增益
- 此處不做數(shù)學理論推導肮柜，直接給出如下公式
P^-(i) * H^T
k_bar(i) = ------------------------------------------------------
H * P^-(i) * H^T + R
P^-(i) * h_bar(i)
= ------------------------------------------------------ (3.5.1)
h_bar(i)^T * P^-(i) * h_bar(i) + R
( E(i-1) + Q(i) ) * h_bar(i)
= ------------------------------------------------------ (3.5.2)
var_v_hat(i) + h_bar(i)^T * (E(i-1) + Q(i)) * h_bar(i)
- 其中
R = var_v_hat(i) - 在webrtc gcc算法中R表示測量誤差協(xié)方差陷舅，它使用指數(shù)平均濾波器，通過如下方法計算R
The variance var_v(i) = sigma_v(i)^2 is estimated using an exponential averaging filter, modified for variable sampling rate
var_v_hat(i) = max(beta * var_v_hat(i-1) + (1-beta) * z(i)^2, 1) (3.5.3)
beta = (1-chi)^(30/(1000 * f_max)) (3.5.4)
- 測量誤差協(xié)方差最小取值為1
- z(i)為網(wǎng)絡殘差(當前測量值 - 上一次估計值）
3.6 卡爾曼濾波-校正計算后驗估計值
theta_hat(i) = theta_hat^-(i) + k_bar(i) * (d(i) - H * theta_hat^-(i))
= theta_hat(i-1) + k_bar(i) * (d(i) - H * theta_hat(i-1)) (3.6.1)
= theta_hat(i-1) + k_bar(i) * d(i) - k_bar(i) * H * theta_hat(i-1)
= (1 - k_bar(i) * H) * theta_hat(i-1) + k_bar(i) * d(i) (3.6.2)
k_bar(i) ~ [0 ~ 1/H]
-
theta_hat(i)為第i幀的估計值审洞，叫做后驗估計值 -
k_bar(i)為當前幀(第i幀)的卡爾曼增益 - 卡爾曼濾波的目標就是去尋找適當?shù)?code>k_bar(i)莱睁，使得
theta_hat(i)的值更加趨向theta(i)（也就是實際值），如何去選擇值肯定是和誤差v(i)以及u(i-1)息息相關的 - 當
k_bar(i)趨近0的時候芒澜，更相信算出來的結果（估計值） - 當
k_bar(i)趨近1/H的時候仰剿，更相信測量出來的結果 - 將上述公式進行替換如下
z(i) = d(i) - h_bar(i)^T * theta_hat(i-1) (3.6.3)
theta_hat(i) = theta_hat(i-1) + z(i) * k_bar(i) (3.6.4)
-
z(i)也叫網(wǎng)絡殘差(當前測量值 - 上一次估計值） - 上述公式則為google給出的官方公式
- 結合
(3.6.1)和(3.5.2)可知當var_v_hat(i)越大時，說明測量誤差較大痴晦，此時卡爾曼增益將越小南吮，最終的估計值將更加趨近與一次的估計值
3.7 卡爾曼濾波-更新誤差協(xié)方差
- 計算后驗估計誤差協(xié)方差
e(i) = theta(i) - theta_hat(i) P(E(i)) ~ (0 , P) (3.7.1)
P(i) = E{e(i) * e(i)^T} (3.7.2)
= E{(theta(i) - theta_hat(i)) * (theta(i) - theta_hat(i))^T} (3.7.3)
= E(i)
e(i)表示當前幀的實際值和估計值之間的誤差P(i)表示當前幀的誤差協(xié)方差，此處叫做后驗誤差協(xié)方差綜合(3.6.1)誊酌、(3.5.1)最后得出
P(i) = (I - k_bar(i) * H) * P^-(i) (3.7.4)
= (I - k_bar(i) * h_bar(i)^T) * (E(i-1) + Q(i)) (3.7.5)
= E(i)
- 其中I為2*2的單位矩陣
3.8 卡爾曼濾波-系統(tǒng)模型圖
-
結合公式3.6.2我們可以得出如下模型圖
image-20210523222453018.png
m(i) = (1 ? K(i)) * m(i?1) + K(i) * (dm(i))
4）WebRTC中JitterDelay的計算和迭代過程
image-20210612222555068.png
- 步驟1部凑、在
FrameBuffer模塊中通過Decode任務隊列，通過重復任務調用FindNextFrame和GetNextFrame獲取待解碼的視頻幀 - 步驟2术辐、然后通過
VCMInterFrameDelay模塊的CalculateDelay()函數(shù)依據(jù)當前找到的幀的rtp時間戳以及接收時間計算幀間延遲得到frameDelayMS - 步驟3砚尽、以
frameDelayMS和當前幀的frameSizeBytes(當前幀的數(shù)據(jù)大小字節(jié))為參數(shù)，調用VCMJitterEstimator::UpdateEstimate()函數(shù)對卡爾曼濾波器進行迭代辉词，并計算出當前幀的最優(yōu)jitterDelay - 步驟4、對卡爾曼濾波進行校準供下一次迭代使用
- 最后通過
VCMJitterEstimator::GetJitterEstimate()函數(shù)返回最優(yōu)估計jitterdelay并將其作用到VCMTiming模塊供期望RenderDelay的計算 - 其中步驟2中的
frameDelayMS的計算對應了第二節(jié)的Arrival-time model示意圖猾骡，并將T(i)代表第i幀的發(fā)送的rtp時間戳瑞躺，T(i-1)為第i-1幀的rtp時間戳敷搪，t(i)為第i幀的本地接收時間，t(i-1)為第i-1幀的本地接收時間 - 步驟2的實現(xiàn)代碼如下：
bool VCMInterFrameDelay::CalculateDelay(uint32_t timestamp,
int64_t* delay,
int64_t currentWallClock) {
.....
// Compute the compensated timestamp difference and convert it to ms and round
// it to closest integer.
_dTS = static_cast<int64_t>(
(timestamp + wrapAroundsSincePrev * (static_cast<int64_t>(1) << 32) -
_prevTimestamp) /
90.0 +
0.5);\
// frameDelay is the difference of dT and dTS -- i.e. the difference of the
// wall clock time difference and the timestamp difference between two
// following frames.
*delay = static_cast<int64_t>(currentWallClock - _prevWallClock - _dTS);
_prevTimestamp = timestamp;
_prevWallClock = currentWallClock;
return true;
}
-
_dTS表示相鄰兩幀的rtp時間戳 -
_prevWallClock為上一幀的本地接收時間 -
currentWallClock為當前幀的本地接收時間 - 接下來首先介紹
JitterDelay的計算流程幢哨，最后分析卡爾曼濾波預測和校正過程
4.1）計算JitterDelay
- 在
FrameBuffer模塊中使用如下函數(shù)獲取
int VCMJitterEstimator::GetJitterEstimate(
double rttMultiplier,
absl::optional<double> rttMultAddCapMs) {
//調用CalculateEstimate()計算當前的jitterDelay,OPERATING_SYSTEM_JITTER默認為10ms
//這就意味著默認最小的jittterDelay至少是10ms?
double jitterMS = CalculateEstimate() + OPERATING_SYSTEM_JITTER;
uint64_t now = clock_->TimeInMicroseconds();
//kNackCountTimeoutMs = 60000
// FrameNacked會更新_latestNackTimestamp單位為微秒
//1分鐘內若所有幀都未丟包則清除
if (now - _latestNackTimestamp > kNackCountTimeoutMs * 1000)
_nackCount = 0;
if (_filterJitterEstimate > jitterMS)
jitterMS = _filterJitterEstimate;
if (_nackCount >= _nackLimit) {//_nackLimit
if (rttMultAddCapMs.has_value()) {
jitterMS +=
std::min(_rttFilter.RttMs() * rttMultiplier, rttMultAddCapMs.value());
} else {
jitterMS += _rttFilter.RttMs() * rttMultiplier;
}
}
....
return rtc::checked_cast<int>(std::max(0.0, jitterMS) + 0.5);
}
-
FrameBuffer只針對未重傳過包的幀進行jitterDelay的迭代和計算赡勘，也就是說假設第i幀數(shù)據(jù)有丟包，那么該幀是不會計算JitterDelay的捞镰，該幀的期望渲染時間計算過程中所用到的jitterDelay值是上一幀的jitterDelay - webrtc jitter delay對丟過包的數(shù)據(jù)幀會在
FrameBuffer::GetNextFrame()函數(shù)中通過jitter_estimator_.FrameNacked()函數(shù)告知VCMJitterEstimator模塊闸与，使得_nackCount變量自增(前提通過WebRTC-AddRttToPlayoutDelay/Enable使能)，超過3會清除 - 當
_nackCount大于等于3的時候在時候jitterDelay會在后續(xù)未丟包的幀中加上一個RTT乘以某個系數(shù) - 通過參閱網(wǎng)上的一些文摘表示岸售，如果某幀有丟包如果不處理該RTT會很容易造成卡頓
- 以上函數(shù)的核心主要是調用
CalculateEstimate()函數(shù)計算jitterDelay
// Calculates the current jitter estimate from the filtered estimates.
double VCMJitterEstimator::CalculateEstimate() {
double ret = _theta[0] * (_maxFrameSize - _avgFrameSize) + NoiseThreshold();
.......
_prevEstimate = ret;
return ret;
}
-
_theta[0]記錄的是信道傳輸速率的倒數(shù)践樱，也就是上述卡爾曼狀態(tài)方程中的1/c(i) -
_maxFrameSize表示自會話開始以來所收到的最大幀大小 -
_avgFrameSize表示當前平均幀大小 -
NoiseThreshold()計算噪聲補償閥值
double VCMJitterEstimator::NoiseThreshold() const {
double noiseThreshold = _noiseStdDevs * sqrt(_varNoise) - _noiseStdDevOffset;
if (noiseThreshold < 1.0) {
noiseThreshold = 1.0;
}
return noiseThreshold;
}
-
_noiseStdDevs表示噪聲標準差系數(shù)取值2.33 -
_varNoise表示測量噪聲方差對應公式(3.2.3)中的v(i)，默認初始值為4.0 -
_noiseStdDevOffset噪聲標準差偏移（噪聲扣除常數(shù)）取值30.0ms
4.2）JitterDelay迭代更新機制
- 通過第4大節(jié)中的截圖看出凸丸，卡爾曼的更新是通過
Framebuffer模塊調用其UpdateEstimate函數(shù)來實現(xiàn)的
// Updates the estimates with the new measurements.
void VCMJitterEstimator::UpdateEstimate(int64_t frameDelayMS,
uint32_t frameSizeBytes,
bool incompleteFrame /* = false */) {
//1)計算當前幀和上一幀的數(shù)據(jù)量之差
int deltaFS = frameSizeBytes - _prevFrameSize;
//2)計算_avgFrameSize平均每幀數(shù)據(jù)的大小
if (_fsCount < kFsAccuStartupSamples) {
_fsSum += frameSizeBytes;
_fsCount++;
} else if (_fsCount == kFsAccuStartupSamples) {//kFsAccuStartupSamples取值為5超過5幀開始計算平均幀大小
// Give the frame size filter.
_avgFrameSize = static_cast<double>(_fsSum) / static_cast<double>(_fsCount);
_fsCount++;
}
/*3）若當前輸入幀的大小比平均每幀數(shù)據(jù)的數(shù)據(jù)量要大,則對其進行滑動平均處理拷邢，比如說如果當前是一個I幀，數(shù)據(jù)量顯然會比較大,
默認incompleteFrame為false屎慢，所以每幀都會計算平均值*/
if (!incompleteFrame || frameSizeBytes > _avgFrameSize) {
//滑動平均算法瞭稼，_phi的取值為0.97，取接近前30幀數(shù)據(jù)大小的平均值腻惠，求得的avgFrameSize值為接近近30幀數(shù)據(jù)的平均大小
double avgFrameSize = _phi * _avgFrameSize + (1 - _phi) * frameSizeBytes;
//如果I幀數(shù)據(jù)量會比較大环肘，如下的判斷會不成立，偏移太大不計賦值_avgFrameSize
if (frameSizeBytes < _avgFrameSize + 2 * sqrt(_varFrameSize)) {
// Only update the average frame size if this sample wasn't a key frame.
_avgFrameSize = avgFrameSize;
}
// Update the variance anyway since we want to capture cases where we only
// get key frames.
//3.1）此處更新平均幀大下的方差集灌，默認方差為100,取其最大值,根據(jù)_varFrameSize可以得出在傳輸過程中每幀數(shù)據(jù)大小的均勻性
// 若方差較大則說明幀的大小偏離平均幀大小的程度越大廷臼，則均勻性也越差
_varFrameSize = VCM_MAX(
_phi * _varFrameSize + (1 - _phi) * (frameSizeBytes - avgFrameSize) *
(frameSizeBytes - avgFrameSize),
1.0);
}
// Update max frameSize estimate.
//4)計算最大幀數(shù)據(jù)量
_maxFrameSize =
VCM_MAX(_psi * _maxFrameSize, static_cast<double>(frameSizeBytes));
if (_prevFrameSize == 0) {
_prevFrameSize = frameSizeBytes;
return;
}
//賦值上一幀數(shù)據(jù)大小
_prevFrameSize = frameSizeBytes;
// Cap frameDelayMS based on the current time deviation noise.
/*5) 根據(jù)當前時間偏移噪聲求frameDelayMS,_varNoise為噪聲方差，默認4.0很顯然該值在傳輸過程中會變化绝页，
time_deviation_upper_bound_為時間偏移上限值荠商，默認為3.5，所以默認初始值計算出來max_time_deviation_ms
為7续誉，對于幀率越高莱没，默認輸入的frameDelayMS會越小，這里和max_time_deviation_ms去最小值酷鸦，當噪聲的方差越大
max_time_deviation_ms的值月越大饰躲，其取值就會越接近取向傳入的frameDelayMS*/
int64_t max_time_deviation_ms =
static_cast<int64_t>(time_deviation_upper_bound_ * sqrt(_varNoise) + 0.5);
frameDelayMS = std::max(std::min(frameDelayMS, max_time_deviation_ms),
-max_time_deviation_ms);
/*6）根據(jù)得出的延遲時間計算樣本與卡爾曼濾波器估計的期望延遲之間的延遲差(反映網(wǎng)絡噪聲的大小),計算公式為
frameDelayMS - (_theta[0] * deltaFSBytes + _theta[1])
當前測量值 - 上一次卡爾曼濾波后的估計值 對應公式3.6.3和3.6.4*/
double deviation = DeviationFromExpectedDelay(frameDelayMS, deltaFS);
// Only update the Kalman filter if the sample is not considered an extreme
// outlier. Even if it is an extreme outlier from a delay point of view, if
// the frame size also is large the deviation is probably due to an incorrect
// line slope.
//根據(jù)注釋的意思是只有當樣本不被認為是極端異常值時才更新卡爾曼濾波器,言外之意就是網(wǎng)絡殘差值不能超過
// _numStdDevDelayOutlier * sqrt(_varNoise) = 30ms 默認值，隨_varNoise的大小變化而變化
if (fabs(deviation) < _numStdDevDelayOutlier * sqrt(_varNoise) ||
frameSizeBytes >
_avgFrameSize + _numStdDevFrameSizeOutlier * sqrt(_varFrameSize)) {
// Update the variance of the deviation from the line given by the Kalman
// filter.
EstimateRandomJitter(deviation, incompleteFrame);
// Prevent updating with frames which have been congested by a large frame,
// and therefore arrives almost at the same time as that frame.
// This can occur when we receive a large frame (key frame) which has been
// delayed. The next frame is of normal size (delta frame), and thus deltaFS
// will be << 0. This removes all frame samples which arrives after a key
// frame.
if ((!incompleteFrame || deviation >= 0.0) &&
static_cast<double>(deltaFS) > -0.25 * _maxFrameSize) {
// Update the Kalman filter with the new data
KalmanEstimateChannel(frameDelayMS, deltaFS);
}
} else { // 如果網(wǎng)絡殘差太大臼隔，說明噪聲偏移太大嘹裂，需要對測量噪聲進行校正，本次不進行卡爾曼預測和校正
int nStdDev =
(deviation >= 0) ? _numStdDevDelayOutlier : -_numStdDevDelayOutlier;
EstimateRandomJitter(nStdDev * sqrt(_varNoise), incompleteFrame);
}
// Post process the total estimated jitter
//6) 求得當前幀的jitterDelay最優(yōu)估計值
if (_startupCount >= kStartupDelaySamples) {
PostProcessEstimate();
} else {
_startupCount++;
}
}
-
_avgFrameSize表示平均幀大小摔握，使用滑動平均算法 -
_varFrameSize表示傳輸過程中幀大小的方差寄狼，象征了傳輸過程中每幀數(shù)據(jù)量的均勻性，_varFrameSize越大表示均勻性越差，所以在傳輸過程中泊愧，在發(fā)送端應當盡量使得每幀數(shù)據(jù)的數(shù)據(jù)量盡可能的接近（如盡量減少I幀的發(fā)送伊磺，按需所取） -
_maxFrameSize表示自會話以來最大幀的數(shù)據(jù)量大小 -
deltaFS表示當前幀和上一幀的數(shù)據(jù)量只差 -
deviation表示網(wǎng)絡殘差删咱，使用當前測量值減去估計值求得,代碼如下：
double VCMJitterEstimator::DeviationFromExpectedDelay(
int64_t frameDelayMS,
int32_t deltaFSBytes) const {
return frameDelayMS - (_theta[0] * deltaFSBytes + _theta[1]);
}
- 以上計算對應公式
3.6.3和3.6.4 - 上述第5步中根據(jù)當前網(wǎng)絡殘差是否在某范圍內來決定是只更新測量噪聲的方差還是即更新測量噪聲誤差也進行預測和校正工作
- 接下來著重分析
EstimateRandomJitter測量噪聲方差的計算和KalmanEstimateChannel函數(shù)的預測以及校正分析 - 同時以上
UpdateEstimate函數(shù)的過程中計算平均幀大小以及平均幀大小方差都使用了移動平均算法如下：
double avgFrameSize = _phi * _avgFrameSize + (1 - _phi) * frameSizeBytes;
_varFrameSize = VCM_MAX(
_phi * _varFrameSize + (1 - _phi) * (frameSizeBytes - avgFrameSize) *
(frameSizeBytes - avgFrameSize),1.0);
- 使用移動平均算法的原因是讓計算出來的值和前面幀的值有關聯(lián)屑埋，確保離散數(shù)據(jù)的平滑性，當_phi取1的時候就不平滑了
4.3）更新誤差方差
void VCMJitterEstimator::EstimateRandomJitter(double d_dT,
bool incompleteFrame) {
uint64_t now = clock_->TimeInMicroseconds();
//1)對幀率進行采樣統(tǒng)計
if (_lastUpdateT != -1) {
fps_counter_.AddSample(now - _lastUpdateT);
}
_lastUpdateT = now;
if (_alphaCount == 0) {
assert(false);
return;
}
//2) alpha = 399
double alpha =
static_cast<double>(_alphaCount - 1) / static_cast<double>(_alphaCount);
_alphaCount++;
if (_alphaCount > _alphaCountMax)
_alphaCount = _alphaCountMax;//_alphaCountMax = 400
// In order to avoid a low frame rate stream to react slower to changes,
// scale the alpha weight relative a 30 fps stream.
double fps = GetFrameRate();
if (fps > 0.0) {
double rate_scale = 30.0 / fps;
// At startup, there can be a lot of noise in the fps estimate.
// Interpolate rate_scale linearly, from 1.0 at sample #1, to 30.0 / fps
// at sample #kStartupDelaySamples.
if (_alphaCount < kStartupDelaySamples) {
rate_scale =
(_alphaCount * rate_scale + (kStartupDelaySamples - _alphaCount)) /
kStartupDelaySamples;//kStartupDelaySamples = 30
}
//alpha = pow(399/400,30.0 / fps)
alpha = pow(alpha, rate_scale);
}
double avgNoise = alpha * _avgNoise + (1 - alpha) * d_dT;
double varNoise =
alpha * _varNoise + (1 - alpha) * (d_dT - _avgNoise) * (d_dT - _avgNoise);
if (!incompleteFrame || varNoise > _varNoise) {
_avgNoise = avgNoise;
_varNoise = varNoise;
}
if (_varNoise < 1.0) {
// The variance should never be zero, since we might get stuck and consider
// all samples as outliers.
_varNoise = 1.0;
}
}
- 測量誤差的更新和幀率相關痰滋，
alpha為指數(shù)移動平均算法的指數(shù)加權系數(shù)收幀率大小的影響較大摘能，數(shù)學模型為遞減指數(shù)函數(shù)，當fps變低時敲街，rate_scale增大团搞，從而alpha會變小，最終avgNoise的值會增大聪富，意味著此時噪聲增大莺丑，實時幀率越接近30fps,效果越理想 -
_avgNoise為自會話以來的平均測量噪聲，_varNoise會噪聲方差墩蔓，變量d_dT表示網(wǎng)絡殘差 - 噪聲方差以及平均噪聲的更新都采用指數(shù)型移動平均算法求取梢莽，其特性是指數(shù)式遞減加權的移動平均，各數(shù)值的加權影響力隨時間呈指數(shù)式遞減奸披，時間越靠近當前時刻的數(shù)據(jù)加權影響力越大
double avgNoise = alpha * _avgNoise + (1 - alpha) * d_dT;
double varNoise =
alpha * _varNoise + (1 - alpha) * (d_dT - _avgNoise) * (d_dT - _avgNoise);
- 以上計算對應公式
3.5.3和3.5.4 - 當
_varNoise較大時昏名，說明此時網(wǎng)絡噪聲較大，此時后續(xù)的卡爾曼增益計算的結果會較小阵面，說明測量誤差較大轻局，使得當前幀的jitterDelay更趨近上一次的估計值，從而收當前網(wǎng)絡殘差的影響更小一些
4.4）Kalman Filter 預測及校正
- 結合第3章界的理論推導样刷，預測首先利用上一次的估計結果求得先驗估計誤差協(xié)方差
- 在
VCMJitterEstimator模塊中定義如下核心矩陣
double _thetaCov[2][2]; // Estimate covariance 先驗估計誤差協(xié)方差
double _Qcov[2][2]; // Process noise covariance 過程噪聲協(xié)方差對應Q向量
/**
@ frameDelayMS:為測量出來的當前幀和上一幀的幀間抖動延遲
@ deltaFSBytes：當前幀的數(shù)據(jù)量和上一幀的數(shù)據(jù)量之差
**/
void VCMJitterEstimator::KalmanEstimateChannel(int64_t frameDelayMS,
int32_t deltaFSBytes) {
double Mh[2];//P(i) = E[1/c(i) m(i)]
double hMh_sigma;
double kalmanGain[2];
double measureRes;
double t00, t01;
// Kalman filtering
/*1）計算先驗估計誤差協(xié)方差
結合公式3.4.1~3.4.4
e^-(i) = theta(i) - theta_hat^-(i)
P^-(i) = {e^-(i) * e^-(i)^T}
= E{(theta(i) - theta_hat^-(i)) * (theta(i) - theta_hat^-(i))^T}
= A * P(i-1) * A^T + Q
= P(i-1) + Q
= E(i-1) + Q
當前幀(i)的先驗估計誤差協(xié)防差 = 上一幀(i-1)的誤差協(xié)方差 + 過程噪聲協(xié)方差
*/
//Prediction
//M = M + Q = E(i-1) + 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];
/*
2) 校正 根據(jù)公式3.5.1~3.5.2計算卡爾曼增益
P^-(i) * H^T
k_bar(i) = ------------------------------------------------------
H * P^-(i) * H^T + R
P^-(i) * h_bar(i)
= ------------------------------------------------------ (3.5.1)
h_bar(i)^T * P^-(i) * h_bar(i) + R
( E(i-1) + Q(i) ) * h_bar(i)
= ------------------------------------------------------ (3.5.2)
var_v_hat(i) + h_bar(i)^T * (E(i-1) + Q(i)) * h_bar(i)
*/
// Kalman gain
// K = M*h'/(sigma2n + h*M*h') = M*h'/(1 + h*M*h') = M*h'/(var_v_hat(i) + h*M*h')
// h = [dFS 1] 其中dFS對應的入?yún)eltaFSBytes
// Mh = M*h' = _thetaCov[2][2] * [dFS 1]^
// hMh_sigma = h*M*h' + R = h_bar(i)^T * (E(i-1) + Q(i)) * h_bar(i) + R
Mh[0] = _thetaCov[0][0] * deltaFSBytes + _thetaCov[0][1];// 對應1/C(i) 信道傳輸速率的誤差協(xié)方差
Mh[1] = _thetaCov[1][0] * deltaFSBytes + _thetaCov[1][1];// 對應網(wǎng)絡排隊延遲m(i)的誤差協(xié)方差
// sigma weights measurements with a small deltaFS as noisy and
// measurements with large deltaFS as good
if (_maxFrameSize < 1.0) {
return;
}
//sigma為測量噪聲標準差的指數(shù)平均濾波仑扑，對應的是測量噪聲的協(xié)方差R
double sigma = (300.0 * exp(-fabs(static_cast<double>(deltaFSBytes)) /
(1e0 * _maxFrameSize)) +
1) *
sqrt(_varNoise);
if (sigma < 1.0) {
sigma = 1.0;
}
// hMh_sigma 對應H * P^-(i) * H^T = h_bar(i)^T * (E(i-1) + Q(i)) * h_bar(i) + R
// 對應公式(3.5.1)
//[dFS 1]^ * Mh = dFS * Mh[0] + Mh[1]
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;
}
//計算卡爾曼增益Mh / hMh_sigma
kalmanGain[0] = Mh[0] / hMh_sigma;
kalmanGain[1] = Mh[1] / hMh_sigma;
/*
3)根據(jù)公式3.6.1~3.6.4校正計算后驗估計值
theta_hat(i) = theta_hat^-(i) + k_bar(i) * (d(i) - H * theta_hat^-(i))
= theta_hat(i-1) + k_bar(i) * (d(i) - H * theta_hat(i-1)) (3.6.1)
= theta_hat(i-1) + k_bar(i) * d(i) - k_bar(i) * H * theta_hat(i-1)
= (1 - k_bar(i) * H) * theta_hat(i-1) + k_bar(i) * d(i) (3.6.2)
其中k_bar(i) ~ [0 ~ 1/H]
z(i) = d(i) - h_bar(i)^T * theta_hat(i-1) (3.6.3)
theta_hat(i) = theta_hat(i-1) + z(i) * k_bar(i) (3.6.4)
*/
// Correction
// theta = theta + K*(dT - h*theta)
// 計算網(wǎng)絡殘差，得到最優(yōu)估計值
measureRes = frameDelayMS - (deltaFSBytes * _theta[0] + _theta[1]);
_theta[0] += kalmanGain[0] * measureRes; //公式(3.6.4)
_theta[1] += kalmanGain[1] * measureRes; //公式(3.6.4)
if (_theta[0] < _thetaLow) {
_theta[0] = _thetaLow;
}
/**
4)根據(jù)公式3.7.1~3.7.4更新誤差協(xié)方差置鼻，為下一次預測提供最優(yōu)濾波器系數(shù)
e(i) = theta(i) - theta_hat(i) P(E(i)) ~ (0 , P) (3.7.1)
P(i) = E{e(i) * e(i)^T} (3.7.2)
= E{(theta(i) - theta_hat(i)) * (theta(i) - theta_hat(i))^T} (3.7.3)
= (I - k_bar(i) * H) * P^-(i) (3.7.4)
= (I - k_bar(i) * h_bar(i)^T) * (E(i-1) + Q(i)) (3.7.5)
= E(i)
*/
// 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);
}
以上分成4步完成預測以及校準過程镇饮，其中第4步涉及到舉證的減法以及乘法的運算
-
下面給出更新誤差協(xié)方差的矩陣運算公式
image-20210613181020401.png 到此為止卡爾曼濾波的一次更新以及迭代完成，總結其過程箕母，卡爾曼濾波主要是為了得到較為準確的信道傳輸速率的倒數(shù)
1/c(i)以及網(wǎng)絡排隊延遲m(i)最后通過如下計算
jitterDelay
double noiseThreshold = _noiseStdDevs * sqrt(_varNoise) - _noiseStdDevOffset;
jitterDelay = _theta[0] * (_maxFrameSize - _avgFrameSize) + NoiseThreshold();
其中
_varNoise為測量誤差储藐，通過計算當前幀的測量幀間延遲和上一次的卡爾曼最優(yōu)估計出來的傳輸速率的倒數(shù)1/c(i)以及網(wǎng)絡排隊延遲m(i)并計算出估計幀間延遲進行相減，得出網(wǎng)絡殘差依據(jù)網(wǎng)絡殘差更新測量誤差平均值以及平均方差
最終將最優(yōu)估計
jitterDelay作用到vcmTiming模塊用于計算最優(yōu)期望渲染時間至此代碼分析完畢
4.5）webrtc jitterdelay相關數(shù)據(jù)測試
-
通過實際調試嘶是，1080P@60fps 線上測試抓得如下數(shù)據(jù)
image-20210613173701391.png - 以上一共近2000幀數(shù)據(jù)钙勃，從上圖來看，在937幀數(shù)據(jù)左右聂喇，畫面從靜止編程動畫辖源，數(shù)據(jù)量瞬間增大，此時，平均幀大小的方差也隨之增大同木，從圖像右邊的測量噪聲的方差曲線來看浮梢，也發(fā)現(xiàn)測量噪聲的方差隨之增大
- 右下圖顯示的是實際測量的
frameDelayMs跛十、jitterMs(估計后)彤路、以及加上10ms操作延遲的變化曲線，由該曲線可以看出kalman filter最終確實能使實際的幀間抖動變得更加平滑芥映，但是從曲線來看洲尊，似乎準確性有待提升。 - 綜合分析奈偏，
VCMJitterEstimator模塊所得出的jitterDelay它的核心作用是給到vcmTiming模塊用于估計期望渲染時間坞嘀，為什么要這么做？ - 在實際的傳輸過程中若完全不做處理惊来，那么實際的幀間延遲就如上圖右下方藍色曲線所示丽涩，如果來一幀數(shù)據(jù)就送到解碼器讓解碼器做解碼并渲染，那么在數(shù)據(jù)量變化劇烈或突然出現(xiàn)網(wǎng)絡波動等一系列導致幀間延遲瞬時變化較大的時候可能會出現(xiàn)卡頓的問題裁蚁，并且視覺上也會有明顯的表象
總結：
- 本文通過分析
VCMJitterEstimator模塊的實現(xiàn)矢渊，學習卡爾曼濾波的應用場景，并深入分析其濾波步驟等 - 通過分析
VCMJitterEstimator模塊了解webrtc視頻接收流水線中jitterDelay的用途以及更新過程枉证，為進一步優(yōu)化卡頓做鋪墊
參考文獻
[1].A Google Congestion Control Algorithm for Real-Time Communication draft-alvestrand-rmcat-congestion-03
[2].Analysis and Design of the Google Congestion Control for Web Real-time Communication (WebRTC)
[3].DR_CAN 卡爾曼濾波視頻教程
[4].從放棄到精通矮男！卡爾曼濾波從理論到實踐
