【实战】OpenCV+Python项目实战--光流估计

光流估计
光流是空间运动物体在观测成像平面上的像素运动的“瞬时速度”，根据各个像素点的速度矢量特征，可以对图像进行动态分析，例如目标跟踪。
要求如下：
（1）亮度恒定：同一点随着时间的变化，其亮度不会发生改变。

（2）小运动：随着时间的变化不会引起位置的剧烈变化，只有小运动情况下才能用前后帧之间单位位置变化引起的灰度变化去近似灰度对位置的偏导数。

（3）空间一致：一个场景上邻近的点投影到图像上也是邻近点，且邻近点速度一致。因为光流法基本方程约束只有一个，而要求x，y方向的速度，有两个未知变量。所以需要连立n多个方程求解。

Lucas-Kanade 算法
如何求解方程组呢？看起来一个像素点根本不够，在物体移动过程中还有哪些特性呢？

要保证u可逆，就ATA的特征值λ1和λ2都很大，角点刚好满足这个情况，所以一般都是利用角点来做。

常用函数：
cv2.goodFeaturesToTrack(old_gray, mask=None, **feature_params)

作用：用于获得光流估计所需要的角点
参数说明：old_gray表示输入图片，mask表示掩模，feature_params:maxCorners=100角点的最大个数（可以小于该数）,qualityLevel=0.3角点品质,minDistance=7即在这个范围内只存在一个品质最好的角点
返回值为(n, 1, 2)的矩阵

pl, st, err = cv2.calcOpticalFlowPyrLK(old_gray, frame_gray, p0,None, **lk_params)

用于获得光流检测后的角点位置

参数说明：pl表示光流检测后的角点位置，st表示是否是运动的角点，err表示是否出错，old_gray表示输入前一帧图片，frame_gray表示后一帧图片，p0表示需要检测的角点，lk_params：winSize表示选择多少个点进行u和v的求解，maxLevel表示空间金字塔的层数

2准备工作

2.1 numpy.random.randint

numpy.random.randint(low, high=None, size=None, dtype=‘l’)

函数的作用是，返回一个随机整型数，范围从低（包括）到高（不包括），即[low, high)。
如果没有写参数high的值，则返回[0,low)的值。

参数如下：

low: int 生成的数值最低要大于等于low。 （hign = None时，生成的数值要在[0, low)区间内）
high: int (可选) 如果使用这个值，则生成的数值在[low, high)区间。 size: int or tuple of ints(可选)
输出随机数的尺寸，比如size = (m * n* k)则输出同规模即m * n* k个随机数。默认是None的，仅仅返回满足要求的单一随机数。
dtype: dtype(可选)： 想要输出的格式。如int64、int等等

>>>np.random.randint(2, high=10, size=(2,3))
array([[6, 8, 7],
       [2, 5, 2]])

2.2 *args和 **kwargs

def aa(a, b, c):
    return a + b + c
# 1
print(aa(3, 4, 5))  # 返回值为：12

# 2
b = [1, 2, 3]
print(aa(*b))  # print(aa(b))  # 错误
# 返回值为：6
# “*”作用：拆开数列’b’的数值作为位置参数，并把这些位置参数传给函数’aa’来调用

# 3
b = [5, 3]
a = 1
print(aa(a, *b))  # 9


def xx(*args):
    print(*args)
    print(args)
def yy(*args):
    return args*2
xx(1, 2, 3)  # 1 2 3 和 (1, 2, 3)
print(yy(1, 2, 3))  # (1, 2, 3, 1, 2, 3)

def abc(x, *y):
    print(x)
    print(*y)
abc(1, 4, 2)  # 分别打印出1 和 4 2


# “**”作用：它unpack字典，并将字典中的数据项作为键值参数传给函数。
def xv(a, b, c):
    print(a)
    print(b)
    print(c)

xv(1, **{"b": 1, "c": 2})  # b和c都要对应函数


def xc(**c):
    print(c)
    # print(**c) 报错

xc(b=1, c=2)  # {'b': 1, 'c': 2}

d = dict(b=3, c=2)  # 相当于是{"b"=1, "c"=2}
xc(**d)  # {'b': 3, 'c': 2}
# print(**d)  # 单独的**d就会报错

2.3 cv2.calcOpticalFlowPyrLK()

参数：

• prevImage 前一帧图像

• nextImage 当前帧图像

• prevPts 待跟踪的特征点向量

• winSize 搜索窗口的大小

• maxLevel 最大的金字塔层数

返回：

• nextPts 输出跟踪特征点向量

• status 特征点是否找到，找到的状态为1，未找到的状态为0

• err表示是否出错

2.4 cv2.goodFeaturesToTrack()

cv2.goodFeaturesToTrack(old_gray, mask=None, maxCorners,qualityLevel=0.3,minDistance

这个函数可以帮我们使用 Shi-Tomasi 方法获取图像中 N 个最好的角点（如果你愿意的话也可以通过改变参数来使用 Harris 角点检测算法）。通常情况下，输入的应该是灰度图像。然后确定你想要检测到的角点数目。再设置角点的质量水平，0到 1之间。它代表了角点的最低质量，低于这个数的所有角点都会被忽略。最后在设置两个角点之间的最短欧式距离。根据这些信息，函数就能在图像上找到角点。所有低于质量水平的角点都会被忽略。然后再把合格角点按角点质量进行降序排列。函数会采用角点质量最高的那个角点（排序后的第一个），然后将它附近（最小距离之内）的角点都删掉。按着这样的方式最后返回N 个最佳角点。

3 代码实现

import numpy as np
import cv2

cap = cv2.VideoCapture('test.avi')

# 角点检测所需参数
feature_params = dict(maxCorners=100,
                      qualityLevel=0.3,
                      minDistance=7)

# lucas kanade参数
lk_params = dict(winSize=(15, 15),
                 maxLevel=2)  # 窗口大小为15*15，金字塔层数为2

# 随机颜色条
color = np.random.randint(0, 255, (100, 3))

# 拿到第一帧图像
ret, old_frame = cap.read()
old_gray = cv2.cvtColor(old_frame, cv2.COLOR_BGR2GRAY)
# 返回所有检测特征点，需要输入图像，角点最大数量（效率），品质因子（特征值越大的越好，来筛选）
# 距离相当于这区间有比这个角点强的，就不要这个弱的了
p0 = cv2.goodFeaturesToTrack(old_gray, mask=None, **feature_params)  # 寻找角点

# 创建一个mask
mask = np.zeros_like(old_frame)

while True:
    ret, frame = cap.read()  # 这个是取的第二帧图像，上面已经取出了第一帧图像
    frame_gray = cv2.cvtColor(frame, cv2.COLOR_BGR2GRAY)

    # 需要传入前一帧和当前图像以及前一帧检测到的角点
    p1, st, err = cv2.calcOpticalFlowPyrLK(old_gray, frame_gray, p0, None, **lk_params)

    # st=1表示
    good_new = p1[st == 1]
    print(p1.shape)  # (n,1,2)
    print(good_new.shape)  # (n, 2)
    good_old = p0[st == 1]

    # 绘制轨迹
    for i, (new, old) in enumerate(zip(good_new, good_old)):
        # new=[692.99805  83.00432]
        a, b = new.ravel()  # 或者[a, b] = new
        c, d = old.ravel()
        mask = cv2.line(mask, (a, b), (c, d), color[i].tolist(), 2)
        # python tolist()方法，将数组或者矩阵转换成列表
        frame = cv2.circle(frame, (a, b), 5, color[i].tolist(), -1)

    img = cv2.add(frame, mask)  # 这个相加不会超出边界

    cv2.imshow('frame', img)
    k = cv2.waitKey(150) & 0xff
    if k == 27:
        break

    # 更新
    old_gray = frame_gray.copy()
    p0 = good_new.reshape(-1, 1, 2)

cv2.destroyAllWindows()
cap.release()