Cv运动分析与对象跟踪（续）

全文再续，书接上一回

目录：

4 光流

4.1 CalcOpticalFlowHS
4.2 CalcOpticalFlowLK
4.3 CalcOpticalFlowBM
4.4 CalcOpticalFlowPyrLK
5 预估器
5.1 CvKalman
5.2 CreateKalman
5.3 ReleaseKalman
5.4 KalmanPredict
5.5 KalmanCorrect
5.6 CvConDensation
5.7 CreateConDensation
5.8 ReleaseConDensation
5.9 ConDensInitSampleSet

5.10 ConDensUpdateByTime

4.光流
4.1）CalcOpticalFlowHS：计算两幅图像的光流

void cvCalcOpticalFlowHS( const CvArr* prev, const CvArr* curr, int use_previous,CvArr* velx, CvArr* vely, double lambda,CvTermCriteria criteria );
prev
第一幅图像, 8-比特, 单通道.
curr
第二幅图像, 8-比特, 单通道.
use_previous
使用以前的 (输入) 速度域
velx
光流的水平部分，与输入图像大小一样, 32-比特，浮点数, 单通道.
vely
光流的垂直部分，与输入图像大小一样, 32-比特， 浮点数, 单通道.
lambda
Lagrangian 乘子
criteria
速度计算的终止条件
函数 cvCalcOpticalFlowHS 为输入图像的每一个象素计算光流，使用 Horn & Schunck 算法 [Horn81].

4.2）CalcOpticalFlowLK：计算两幅图像的光流

void cvCalcOpticalFlowLK( const CvArr* prev, const CvArr* curr, CvSize win_size,CvArr* velx, CvArr* vely );
prev
第一幅图像, 8-比特, 单通道.
curr
第二幅图像, 8-比特, 单通道.
win_size
用来归类象素的平均窗口尺寸 （Size of the averaging window used for grouping pixels）
velx
光流的水平部分，与输入图像大小一样, 32-比特， 浮点数, 单通道.
vely
光流的垂直部分，与 输入图像大小一样, 32-比特， 浮点数, 单通道.
函数 cvCalcOpticalFlowLK 为输入图像的每一个象素计算光流，使用 Lucas & Kanade 算法 [Lucas81].

4.3）CalcOpticalFlowBM：用块匹配方法计算两幅图像的光流

void cvCalcOpticalFlowBM( const CvArr* prev, const CvArr* curr, CvSize block_size,CvSize shift_size, CvSize max_range, int use_previous,CvArr* velx, CvArr* vely );
prev
第一幅图像, 8-比特, 单通道.
curr
第二幅图像, 8-比特, 单通道.
block_size
比较的基本块尺寸
shift_size
块坐标的增量
max_range
块周围象素的扫描邻域的尺寸
use_previous
使用以前的 (输入) 速度域
velx
光流的水平部分，尺寸为 floor((prev->width - block_size.width)/shiftSize.width) × floor((prev->height - block_size.height)/shiftSize.height) , 32-比特，浮点数, 单通道.
vely
光流的垂直部分，与 velx 大小一样，32-比特，浮点数, 单通道.
函数 cvCalcOpticalFlowBM 为重叠块 block_size.width×block_size.height 中的每一个象素计算光流，因此其速度域小于整个图像的速度域。对每一个在图像 prev 中的块，函数试图在 curr 中某些原始块或其偏移 (velx(x0,y0),vely(x0,y0)) 块的邻域里寻找类似的块，如同在前一个函数调用中所计算的类似(如果 use_previous=1)

4.4）CalcOpticalFlowPyrLK：计算一个稀疏特征集的光流，使用金字塔中的迭代 Lucas-Kanade 方法

void cvCalcOpticalFlowPyrLK( const CvArr* prev, const CvArr* curr, CvArr* prev_pyr, CvArr* curr_pyr,
                             const CvPoint2D32f* prev_features, CvPoint2D32f* curr_features,
                             int count, CvSize win_size, int level, char* status,
                             float* track_error, CvTermCriteria criteria, int flags );
prev
在时间 t 的第一帧
curr
在时间 t + dt 的第二帧
prev_pyr
第一帧的金字塔缓存. 如果指针非 NULL , 则缓存必须有足够的空间来存储金字塔从层 1 到层 #level 的内容。尺寸 (image_width+8)*image_height/3 比特足够了
curr_pyr
与 prev_pyr 类似， 用于第二帧
prev_features
需要发现光流的点集
curr_features
包含新计算出来的位置的 点集
count
特征点的数目
win_size
每个金字塔层的搜索窗口尺寸
level
最大的金字塔层数。如果为 0 , 不使用金字塔 (即金字塔为单层), 如果为 1 , 使用两层，下面依次类推。
status
数组。如果对应特征的光流被发现，数组中的每一个元素都被设置为 1， 否则设置为 0。
error
双精度数组，包含原始图像碎片与移动点之间的差。为可选参数，可以是 NULL .
criteria
准则，指定在每个金字塔层，为某点寻找光流的迭代过程的终止条件。
flags
其它选项：
CV_LKFLOW_PYR_A_READY , 在调用之前，第一帧的金字塔已经准备好
CV_LKFLOW_PYR_B_READY , 在调用之前，第二帧的金字塔已经准备好
CV_LKFLOW_INITIAL_GUESSES , 在调用之前，数组 B 包含特征的初始坐标 （Hunnish: 在本节中没有出现数组 B，不知是指的哪一个）
函数 cvCalcOpticalFlowPyrLK 实现了金字塔中 Lucas-Kanade 光流计算的稀疏迭代版本 ([Bouguet00])。 它根据给出的前一帧特征点坐标计算当前视频帧上的特征点坐标。 函数寻找具有子象素精度的坐标值。

两个参数 prev_pyr 和 curr_pyr 都遵循下列规则： 如果图像指针为 0, 函数在内部为其分配缓存空间，计算金字塔，然后再处理过后释放缓存。 否则，函数计算金字塔且存储它到缓存中，除非设置标识 CV_LKFLOW_PYR_A[B]_READY 。 图像应该足够大以便能够容纳 Gaussian 金字塔数据。 调用函数以后，金字塔被计算而且相应图像的标识可以被设置，为下一次调用准备就绪 (比如：对除了第一个图像的所有图像序列，标识 CV_LKFLOW_PYR_A_READY 被设置).

5.预估器
5.1）CvKalman：Kalman 滤波器状态
typedef struct CvKalman
{
    int MP;                     /* 测量向量维数 */
    int DP;                     /* 状态向量维数 */
    int CP;                     /* 控制向量维数 */
    /* 向后兼容字段 */
#if 1
    float* PosterState;         /* =state_pre->data.fl */
    float* PriorState;          /* =state_post->data.fl */
    float* DynamMatr;           /* =transition_matrix->data.fl */
    float* MeasurementMatr;     /* =measurement_matrix->data.fl */
    float* MNCovariance;        /* =measurement_noise_cov->data.fl */
    float* PNCovariance;        /* =process_noise_cov->data.fl */
    float* KalmGainMatr;        /* =gain->data.fl */
    float* PriorErrorCovariance;/* =error_cov_pre->data.fl */
    float* PosterErrorCovariance;/* =error_cov_post->data.fl */
    float* Temp1;               /* temp1->data.fl */
    float* Temp2;               /* temp2->data.fl */
#endif
    CvMat* state_pre;           /* 预测状态 (x'(k)): 
                                    x(k)=A*x(k-1)+B*u(k) */
    CvMat* state_post;          /* 矫正状态 (x(k)):
                                    x(k)=x'(k)+K(k)*(z(k)-H*x'(k)) */
    CvMat* transition_matrix;   /* 状态传递矩阵 state transition matrix (A) */
    CvMat* control_matrix;      /* 控制矩阵 control matrix (B)
                                   (如果没有控制，则不使用它)*/
    CvMat* measurement_matrix;  /* 测量矩阵 measurement matrix (H) */
    CvMat* process_noise_cov;   /* 过程噪声协方差矩阵
                                        process noise covariance matrix (Q) */
    CvMat* measurement_noise_cov; /* 测量噪声协方差矩阵
                                          measurement noise covariance matrix (R) */
    CvMat* error_cov_pre;       /* 先验误差计协方差矩阵
                                        priori error estimate covariance matrix (P'(k)):
                                     P'(k)=A*P(k-1)*At + Q)*/
    CvMat* gain;                /* Kalman 增益矩阵 gain matrix (K(k)):
                                    K(k)=P'(k)*Ht*inv(H*P'(k)*Ht+R)*/
    CvMat* error_cov_post;      /* 后验错误估计协方差矩阵
                                        posteriori error estimate covariance matrix (P(k)):
                                     P(k)=(I-K(k)*H)*P'(k) */
    CvMat* temp1;               /* 临时矩阵 temporary matrices */
    CvMat* temp2;
    CvMat* temp3;
    CvMat* temp4;
    CvMat* temp5;
}
CvKalman;
结构 CvKalman 用来保存 Kalman 滤波器状态。它由函数 cvCreateKalman 创建，由函数f cvKalmanPredict 和 cvKalmanCorrect 更新，由 cvReleaseKalman 释放. 通常该结构是为标准 Kalman 所使用的 (符号和公式都借自非常优秀的 Kalman 教程 [Welch95]):
系统运动方程：x_k=A\cdot x_{k-1}+ B\cdot u_k+w_k
系统观测方程：z_k=H\cdot x_k + v_k
其中:
xk(xk − 1) - 系统在时刻 k (k-1) 的状态向量 （state of the system at the moment k (k-1)）
zk - 在时刻 k 的系统状态测量向量 （measurement of the system state at the moment k）
uk - 应用于时刻 k 的外部控制 (external control applied at the moment k)
wk 和 vk 分别为正态分布的运动和测量噪声
p(w) ~ N(0,Q)
p(v) ~ N(0,R),
即,
Q - 运动噪声的相关矩阵，常量或变量
R - 测量噪声的相关矩阵，常量或变量
对标准 Kalman 滤波器，所有矩阵: A, B, H, Q 和 R 都是通过 cvCreateKalman 在分配结构 CvKalman 时初始化一次。但是，同样的结构和函数，通过在当前系统状态邻域中线性化扩展 Kalman 滤波器方程，可以用来模拟扩展 Kalman 滤波器，在这种情况下， A, B, H (也许还有 Q 和 R) 在每一步中都被更新。

5.2）CreateKalman：分配 Kalman 滤波器结构
CvKalman* cvCreateKalman( int dynam_params, int measure_params, int control_params=0 );
dynam_params
状态向量维数
measure_params
测量向量维数
control_params
控制向量维数
函数 cvCreateKalman 分配 CvKalman 以及它的所有矩阵和初始参数

5.3）ReleaseKalman：释放 Kalman 滤波器结构
void cvReleaseKalman( CvKalman** kalman );
kalman
指向 Kalman 滤波器结构的双指针
函数 cvReleaseKalman 释放结构 CvKalman 和里面所有矩阵

5.4）KalmanPredict：估计后来的模型状态
const CvMat* cvKalmanPredict( CvKalman* kalman, const CvMat* control=NULL );
#define cvKalmanUpdateByTime cvKalmanPredict
kalman
Kalman 滤波器状态
control
控制向量 (uk), 如果没有外部控制 (control_params=0) 应该为 NULL
函数 cvKalmanPredict 根据当前状态估计后来的随机模型状态，并存储于 kalman->state_pre:
x'_k=A \cdot x_{k-1}+ B \cdot u_k
P'_k=A \cdot P_{k-1} \cdot A^T + Q,
其中
x'k 是预测状态 (kalman->state_pre),
xk − 1 是前一步的矫正状态 (kalman->state_post)，应该在开始的某个地方初始化，即缺省为零向量，
uk 是外部控制(control 参数),
P'k 是先验误差相关矩阵 (kalman->error_cov_pre)
Pk − 1 是前一步的后验误差相关矩阵(kalman->error_cov_post),应该在开始的某个地方初始化，即缺省为单位矩阵.
函数返回估计得到的状态值

5.5）KalmanCorrect：调节模型状态
const CvMat* cvKalmanCorrect( CvKalman* kalman, const CvMat* measurement );
#define cvKalmanUpdateByMeasurement cvKalmanCorrect
kalman
被更新的 Kalman 结构的指针
measurement
指向测量向量的指针，向量形式为 CvMat
函数 cvKalmanCorrect 在给定的模型状态的测量基础上，调节随机模型状态：
K_k=P'_k\cdot H^T \cdot (H \cdot P'_k \cdot H^T+R)^{-1}
x_k=x'_k+K_k \cdot (z_k-H \cdot x'_k)
P_k=(I-K_k \cdot H)\cdot P'_k
其中
zk - 给定测量(mesurement parameter)
Kk - Kalman "增益" 矩阵
函数存储调节状态到 kalman->state_post 中并且输出时返回它。

5.6）CvConDensation：ConDensaation 状态
typedef struct CvConDensation
{
    int MP;     // 测量向量的维数： Dimension of measurement vector
    int DP;     // 状态向量的维数： Dimension of state vector
    float* DynamMatr;       // 线性动态系统矩阵：Matrix of the linear Dynamics system
    float* State;           // 状态向量： Vector of State
    int SamplesNum;         // 粒子数： Number of the Samples
    float** flSamples;      // 粒子向量数组： array of the Sample Vectors
    float** flNewSamples;   // 粒子向量临时数组： temporary array of the Sample Vectors
    float* flConfidence;    // 每个粒子的置信度(译者注：也就是粒子的权值)：Confidence for each Sample
    float* flCumulative;    // 权值的累计： Cumulative confidence
    float* Temp;            // 临时向量：Temporary vector
    float* RandomSample;    // 用来更新粒子集的随机向量： RandomVector to update sample set
    CvRandState* RandS;     // 产生随机向量的结构数组： Array of structures to generate random vectors
} CvConDensation;
结构 CvConDensation中条件概率密度传播(译者注：粒子滤波的一种特例)（Con-Dens-Ation: 单词 CONditional DENSity propagATION 的缩写）跟踪器的状态。该算法描述可参考 http://www.dai.ed.ac.uk/CVonline/LOCAL_COPIES/ISARD1/condensation.html

5.7）CreateConDensation：分配 ConDensation 滤波器结构
CvConDensation* cvCreateConDensation( int dynam_params, int measure_params, int sample_count );
dynam_params
状态向量的维数
measure_params
测量向量的维数
sample_count
粒子数
函数 cvCreateConDensation 创建结构 CvConDensation 并且返回结构指针。

5.8）ReleaseConDensation：释放 ConDensation 滤波器结构
void cvReleaseConDensation( CvConDensation** condens );
condens
要释放结构的双指针
函数 cvReleaseConDensation 释放结构 CvConDensation (见cvConDensation) 并且清空所有事先被开辟的内存空间。

5.9）ConDensInitSampleSet：初始化 ConDensation 算法中的粒子集
void cvConDensInitSampleSet( CvConDensation* condens, CvMat* lower_bound, CvMat* upper_bound );
condens
需要初始化的结构指针
lower_bound
每一维的下界向量
upper_bound
每一维的上界向量
函数 cvConDensInitSampleSet 在指定区间内填充结构 CvConDensation 中的样例数组。

5.10）ConDensUpdateByTime：估计下个模型状态
void cvConDensUpdateByTime( CvConDensation* condens );
condens
要更新的结构指针
函数 cvConDensUpdateByTime 从当前状态估计下一个随机模型状态。

（完）