http://www.rosoo.net/a/201005/9521.html
#include <stdio.h>
#include <math.h>
#include "config.h"
#include "global.h"
/* private prototypes */
static void calc_actj _ANSI_ARGS_((unsigned char *frame));
static double var_sblk _ANSI_ARGS_((unsigned char *p, int lx));
/* rate control variables */
int Xi, Xp, Xb, r, d0i, d0p, d0b;
double avg_act;
static int R, T, d;
static double actsum;
static int Np, Nb, S, Q;
static int prev_mquant;
void rc_init_seq()
{
/* reaction parameter (constant)
* bit_rate , frame_rate从参数文档中在初始化的时候读取.
* r表示一个2帧图片的占用的带宽。也就是预测的buffer的大小。
*/
if (r==0) r = (int)floor(2.0*bit_rate/frame_rate + 0.5);
/* average activity
*avg_act 为average activity 对于一个宏块来说,它的activity值为4个8*8块中的最小值。它本
* 身是用在最近编码的归一化过程中的。400是设置的平均活动性,预设的,以后会修改。
* GOP的带宽R = G + R,其中G=bit_rate×N/picture_rate //N为(1+N(p)+N(b))GOP中帧的总数,
* 初始化GOP的时候会计算.自动纠错
*/
if (avg_act==0.0) avg_act = 400.0; /* remaining # of bits in GOP */
R = 0;
/* global complexity measure */ //
/*不理解这里。。这里设置全局复杂度,在预设I,P,B帧的带宽的时候用到。
* T(I)= max(R/[(1+(N(p)*X(p)/(X(i)*K(p))+(N(b)*X(b)/X(i)*K(b))] , bit_rate/(8*pic_rate));这个公式是求这个GOP中。
预设的I帧的带宽这里R只GOP的总带宽,N(p)是GOP中p帧的个数,X(p)值p帧的复杂度,X(i)指i帧的复杂度。K(p)=1.0 k(b)=1.4为经验值。(N(p)*X(p)/(X(i)*K(p) 表示GOP中所有的P帧相当于多少I帧(预测值)(N(b)*X(b)/X(i)*K(b))表示GOP中所有的B帧相当于多少I帧R/[(1+(N(p)*X(p)/(X(i)*K(p))+(N(b)*X(b)/X(i)*K(b))这个公式预测I帧占用的带宽。后面的 bit_rate/(8*pic_rate)表示I帧的带宽最小不能低于这个公式计算的带宽。。
*/
if (Xi==0) Xi = (int)floor(160.0*bit_rate/115.0 + 0.5); //I帧的复杂度
if (Xp==0) Xp = (int)floor( 60.0*bit_rate/115.0 + 0.5);
if (Xb==0) Xb = (int)floor( 42.0*bit_rate/115.0 + 0.5);
/* virtual buffer fullness */
if (d0i==0) d0i = (int)floor(10.0*r/31.0 + 0.5); //设置虚拟buffer的大小。每次计算宏块的量化参数的时候,根据这个值进行调整。这个值本身会调整。
if (d0p==0) d0p = (int)floor(10.0*r/31.0 + 0.5); //为P帧预留的虚拟buffer的大小。这个值在编码的时候,会减去前面已经编完宏块占用的buffer,根据剩下的buffer来计算Q(j)每个宏块的量化值。
if (d0b==0) d0b = (int)floor(1.4*10.0*r/31.0 + 0.5);
/*
if (d0i==0) d0i = (int)floor(10.0*r/(qscale_tab[0] ? 56.0 : 31.0) + 0.5);
if (d0p==0) d0p = (int)floor(10.0*r/(qscale_tab[1] ? 56.0 : 31.0) + 0.5);
if (d0b==0) d0b = (int)floor(1.4*10.0*r/(qscale_tab[2] ? 56.0 : 31.0) + 0.5);
*/
fprintf(statfile,"\nrate control: sequence initialization\n");
fprintf(statfile,
" initial global complexity measures (I,P,B): Xi=%d, Xp=%d, Xb=%d\n",
Xi, Xp, Xb);
fprintf(statfile," reaction parameter: r=%d\n", r);
fprintf(statfile,
" initial virtual buffer fullness (I,P,B): d0i=%d, d0p=%d, d0b=%d\n",
d0i, d0p, d0b);
fprintf(statfile," initial average activity: avg_act=%.1f\n", avg_act);
}
void rc_init_GOP(np,nb)
int np,nb;
{
R += (int) floor((1 + np + nb) * bit_rate / frame_rate + 0.5);
//给这个GOP的带宽,在一个图片编码结束后,这个R就会减去图片使用的带宽。
Np = fieldpic ? 2*np+1 : np;
Nb = fieldpic ? 2*nb : nb;
fprintf(statfile,"\nrate control: new group of pictures (GOP)\n");
fprintf(statfile," target number of bits for GOP: R=%d\n",R);
fprintf(statfile," number of P pictures in GOP: Np=%d\n",Np);
fprintf(statfile," number of B pictures in GOP: Nb=%d\n",Nb);
}
/* Note: we need to substitute K for the 1.4 and 1.0 constants -- this can
be modified to fit image content */
/* Step 1: compute target bits for current picture being coded */
void rc_init_pict(frame)
unsigned char *frame;
{
double Tmin;
switch (pict_type)
{
case I_TYPE:
T = (int) floor(R/(1.0+Np*Xp/(Xi*1.0)+Nb*Xb/(Xi*1.4)) + 0.5); //这里每编完一帧,R就减去编完帧使用的带宽。
d = d0i; //同时Np和Nb的个数也会减。由于一个GOP里只有一个I帧,I帧没有减。
break;
case P_TYPE:
T = (int) floor(R/(Np+Nb*1.0*Xb/(1.4*Xp)) + 0.5);
d = d0p;
break;
case B_TYPE:
T = (int) floor(R/(Nb+Np*1.4*Xp/(1.0*Xb)) + 0.5);
d = d0b;
break;
}
Tmin = (int) floor(bit_rate/(8.0*frame_rate) + 0.5); //最小的一帧图片带宽
if (T<Tmin)
T = Tmin;
S = bitcount(); //记录output的buffer里的bit个数,为了统计一帧编码出来后的bit个数
Q = 0;
calc_actj(frame); //计算每个宏块中的每个亮度块的最小方差存放到mbinfo中。
actsum = 0.0;
fprintf(statfile,"\nrate control: start of picture\n");
fprintf(statfile," target number of bits: T=%d\n",T);
}
static void calc_actj(frame)
unsigned char *frame;
{
int i,j,k;
unsigned char *p;
double actj,var;
k = 0;
for (j=0; j<height2; j+=16) //计算每个宏块中每个块的方差,取最小值放到宏块信息中。
for (i=0; i<width; i+=16)
{
p = frame + ((pict_struct==BOTTOM_FIELD)?width:0) + i + width2*j;
/* take minimum spatial activity measure of luminance blocks */
actj = var_sblk(p,width2); //计算块的方差
var = var_sblk(p+8,width2);
if (var<actj) actj = var;
var = var_sblk(p+8*width2,width2);
if (var<actj) actj = var;
var = var_sblk(p+8*width2+8,width2);
if (var<actj) actj = var; // 从actj = var_sblk(p,width2);开始,计算p开始,p+8开始 p+8*width2开始 p+8*width2+8开始的块的方差,得出最小值
if (!fieldpic && !prog_seq)
{
/* field */
var = var_sblk(p,width<<1);
if (var<actj) actj = var;
var = var_sblk(p+8,width<<1);
if (var<actj) actj = var;
var = var_sblk(p+width,width<<1);
if (var<actj) actj = var;
var = var_sblk(p+width+8,width<<1);
if (var<actj) actj = var;
}
actj+= 1.0; //
mbinfo[k++].act = actj;
}
}
void rc_update_pict()
{
double X;
S = bitcount() - S; /* total # of bits in picture */
R-= S; /* remaining # of bits in GOP */ //这里是计算GOP剩下的带宽
X = (int) floor(S*((0.5*(double)Q)/(mb_width*mb_height2)) + 0.5); //
d+= S - T; //T是这帧图片预测的带宽。S是实际编码输出的带宽。d的初值表示预测的(I/B/P)帧的buffer大小。
//这里表示当前帧编码后预测缓存的大小。因为编码上一帧时,可能占用当前缓存空间(S-T)。
avg_act = actsum/(mb_width*mb_height2); //平均活动因子,估算宏块活动的大小。这里将用来计算宏块的量化因子
switch (pict_type)
{
case I_TYPE:
Xi = X;
d0i = d;
break;
case P_TYPE:
Xp = X;
d0p = d;
Np--; //P帧的个数减一
break;
case B_TYPE:
Xb = X;
d0b = d;
Nb--;
break;
}
fprintf(statfile,"\nrate control: end of picture\n");
fprintf(statfile," actual number of bits: S=%d\n",S);
fprintf(statfile," average quantization parameter Q=%.1f\n",
(double)Q/(mb_width*mb_height2));
fprintf(statfile," remaining number of bits in GOP: R=%d\n",R);
fprintf(statfile,
" global complexity measures (I,P,B): Xi=%d, Xp=%d, Xb=%d\n",
Xi, Xp, Xb);
fprintf(statfile,
" virtual buffer fullness (I,P,B): d0i=%d, d0p=%d, d0b=%d\n",
d0i, d0p, d0b);
fprintf(statfile," remaining number of P pictures in GOP: Np=%d\n",Np);
fprintf(statfile," remaining number of B pictures in GOP: Nb=%d\n",Nb);
fprintf(statfile," average activity: avg_act=%.1f\n", avg_act);
}
/* compute initial quantization stepsize 步长 (at the beginning of picture) */
int rc_start_mb()
{
int mquant;//量化因子
if (q_scale_type) //这里是量化因子的两个表。由于量化参数的实际值是通过查表所得。两个表给出了两个不同的量化值。分别决定图片的质量。
{
mquant = (int) floor(2.0*d*31.0/r + 0.5); //计算量化因子,这是第一个宏块的量化因子,公式是d(j)*31/r
//其中d(j) = d(0) + B(j-1) - [T*(j-1) / MB_cnt]其中j表示第j个宏块。[]里的公式计算预测的前面所有编码宏块占用的带宽
/* clip mquant to legal (linear) range */ //B(j-1)表示前面编码宏块实际占有带宽
if (mquant<1)
mquant = 1;
if (mquant>112)
mquant = 112;
/* map to legal quantization level */
mquant = non_linear_mquant_table[map_non_linear_mquant[mquant]]; //编码存放的是量化步长,也就是量化表的下标。
}
else
{
mquant = (int) floor(d*31.0/r + 0.5);
mquant <<= 1;
/* clip mquant to legal (linear) range */
if (mquant<2)
mquant = 2;
if (mquant>62)
mquant = 62;
prev_mquant = mquant;
}
/*
fprintf(statfile,"rc_start_mb:\n");
fprintf(statfile,"mquant=%d\n",mquant);
*/
return mquant;
}
/* Step 2: measure virtual buffer - estimated buffer discrepancy */
int rc_calc_mquant(j)
int j;
{
int mquant;
double dj, Qj, actj, N_actj;
/* measure virtual buffer discrepancy from uniform distribution model */
dj = d + (bitcount()-S) - j*(T/(mb_width*mb_height2)); //d(j) = d(0) + B(j-1) - [T*(j-1) / MB_cnt]
/* scale against dynamic range of mquant and the bits/picture count */
Qj = dj*31.0/r;
/*Qj = dj*(q_scale_type ? 56.0 : 31.0)/r; */
actj = mbinfo[j].act;
actsum+= actj;
/* compute normalized activity */
N_actj = (2.0*actj+avg_act)/(actj+2.0*avg_act);
if (q_scale_type)
{
/* modulate mquant with combined buffer and local activity measures */
mquant = (int) floor(2.0*Qj*N_actj + 0.5);
/* clip mquant to legal (linear) range */
if (mquant<1)
mquant = 1;
if (mquant>112)
mquant = 112;
/* map to legal quantization level */
mquant = non_linear_mquant_table[map_non_linear_mquant[mquant]];
}
else
{
/* modulate mquant with combined buffer and local activity measures */
mquant = (int) floor(Qj*N_actj + 0.5);
mquant <<= 1;
/* clip mquant to legal (linear) range */
if (mquant<2)
mquant = 2;
if (mquant>62)
mquant = 62;
/* ignore small changes in mquant */
if (mquant>=8 && (mquant-prev_mquant)>=-4 && (mquant-prev_mquant)<=4)
mquant = prev_mquant;
prev_mquant = mquant;
}
Q+= mquant; /* for calculation of average mquant */
/*
fprintf(statfile,"rc_calc_mquant(%d): ",j);
fprintf(statfile,"bitcount=%d, dj=%f, Qj=%f, actj=%f, N_actj=%f, mquant=%d\n",
bitcount(),dj,Qj,actj,N_actj,mquant);
*/
return mquant;
}
/* compute variance of 8x8 block */
static double var_sblk(p,lx)
unsigned char *p;
int lx;
{
int i, j;
unsigned int v, s, s2;
s = s2 = 0;
for (j=0; j<8; j++)
{
for (i=0; i<8; i++)
{
v = *p++;
s+= v;
s2+= v*v;
}
p+= lx - 8;
}
return s2/64.0 - (s/64.0)*(s/64.0); //计算块的方差
}
/* VBV calculations
*
* generates warnings if underflow or overflow occurs
*/
/* vbv_end_of_picture
*
* - has to be called directly after writing picture_data()
* - needed for accurate VBV buffer overflow calculation
* - assumes there is no byte stuffing prior to the next start code
*/
static int bitcnt_EOP;
void vbv_end_of_picture()
{
bitcnt_EOP = bitcount();
bitcnt_EOP = (bitcnt_EOP + 7) & ~7; /* account for bit stuffing */
}
/* calc_vbv_delay
*
* has to be called directly after writing the picture start code, the
* reference point for vbv_delay
*/
void calc_vbv_delay()
{
double picture_delay;
static double next_ip_delay; /* due to frame reordering delay */
static double decoding_time;
/* number of 1/90000 s ticks until next picture is to be decoded */
if (pict_type == B_TYPE)
{
if (prog_seq)
{
if (!repeatfirst)
picture_delay = 90000.0/frame_rate; /* 1 frame */
else
{
if (!topfirst)
picture_delay = 90000.0*2.0/frame_rate; /* 2 frames */
else
picture_delay = 90000.0*3.0/frame_rate; /* 3 frames */
}
}
else
{
/* interlaced */
if (fieldpic)
picture_delay = 90000.0/(2.0*frame_rate); /* 1 field */
else
{
if (!repeatfirst)
picture_delay = 90000.0*2.0/(2.0*frame_rate); /* 2 flds */
else
picture_delay = 90000.0*3.0/(2.0*frame_rate); /* 3 flds */
}
}
}
else
{
/* I or P picture */
if (fieldpic)
{
if(topfirst==(pict_struct==TOP_FIELD))
{
/* first field */
picture_delay = 90000.0/(2.0*frame_rate);
}
else
{
/* second field */
/* take frame reordering delay into account */
picture_delay = next_ip_delay - 90000.0/(2.0*frame_rate);
}
}
else
{
/* frame picture */
/* take frame reordering delay into account*/
picture_delay = next_ip_delay;
}
if (!fieldpic || topfirst!=(pict_struct==TOP_FIELD))
{
/* frame picture or second field */
if (prog_seq)
{
if (!repeatfirst)
next_ip_delay = 90000.0/frame_rate;
else
{
if (!topfirst)
next_ip_delay = 90000.0*2.0/frame_rate;
else
next_ip_delay = 90000.0*3.0/frame_rate;
}
}
else
{
if (fieldpic)
next_ip_delay = 90000.0/(2.0*frame_rate);
else
{
if (!repeatfirst)
next_ip_delay = 90000.0*2.0/(2.0*frame_rate);
else
next_ip_delay = 90000.0*3.0/(2.0*frame_rate);
}
}
}
}
if (decoding_time==0.0)
{
/* first call of calc_vbv_delay */
/* we start with a 7/8 filled VBV buffer (12.5% back-off) */
picture_delay = ((vbv_buffer_size*16384*7)/8)*90000.0/bit_rate;
if (fieldpic)
next_ip_delay = (int)(90000.0/frame_rate+0.5);
}
/* VBV checks */
/* check for underflow (previous picture) */
if (!low_delay && (decoding_time < bitcnt_EOP*90000.0/bit_rate))
{
/* picture not completely in buffer at intended decoding time */
if (!quiet)
fprintf(stderr,"vbv_delay underflow! (decoding_time=%.1f, t_EOP=%.1f\n)",
decoding_time, bitcnt_EOP*90000.0/bit_rate);
}
/* when to decode current frame */
decoding_time += picture_delay;
/* warning: bitcount() may overflow (e.g. after 9 min. at 8 Mbit/s */
vbv_delay = (int)(decoding_time - bitcount()*90000.0/bit_rate);
/* check for overflow (current picture) */
if ((decoding_time - bitcnt_EOP*90000.0/bit_rate)
> (vbv_buffer_size*16384)*90000.0/bit_rate)
{
if (!quiet)
fprintf(stderr,"vbv_delay overflow!\n");
}
fprintf(statfile,
"\nvbv_delay=%d (bitcount=%d, decoding_time=%.2f, bitcnt_EOP=%d)\n",
vbv_delay,bitcount(),decoding_time,bitcnt_EOP);
if (vbv_delay<0)
{
if (!quiet)
fprintf(stderr,"vbv_delay underflow: %d\n",vbv_delay);
vbv_delay = 0;
}
{
if (!quiet)
fprintf(stderr,"vbv_delay overflow: %d\n",vbv_delay);
vbv_delay = 65535;
}
}