论文名称：Visual Object Tracking using Adaptive Correlation Filters
文献地址：https://www.cs.colostate.edu/~vision/publications/bolme_cvpr10.pdf
源码地址：https://github.com/xingqing45678/Mosse_CF
由于基于CNNs深度学习在单目标跟踪方法的参数量和计算量都较大，难以与目标检测算法一起移植到嵌入式中。同时，受CVPR2020 AutoTrack的影响，开始从基于传统CF，DCF思想角度入手，对相关滤波CF的鼻祖MOSSE进行攻读。
本文主要介绍相关滤波系列算法的开篇——MOSSE基本原理及其python代码实现流程。
由于论文中涉及到大量傅里叶变换的公式，可能显得MOSSE晦涩难懂，本篇尽量淡化傅里叶变换的公式内容，突出MOSSE在跟踪过程中所做的工作。
相关滤波
相关滤波（CF）源于信号处理领域，有这么一句话"两个信号越相似，其相关值越高。在跟踪，就是找到与跟踪目标响应最大的项" 。通过后续代码的解读，会体会到相关值越高在MOSSE中的作用。
相关操作
假设有两个信号f和g，则两个信号的相关性（correlation）为：
其中，f∗表示f的复共轭...（说好的淡化公式的... 看不明白不要紧，只要明白卷积就可以往下看）
对于图像而言，相关性可以理解为相关核在图像上进行点积操作，如下图所示。图像的相关公式可以表达为g = f ⓧ h，具体而言，对应元素相乘在求和：
从上图的结果可以看出，输出正好是相关核旋转了180°。
具体步骤：
a. 相关核，中心分别位于输入图像的f(i, j)像素上进行滑动；（会补零）
b. 利用上式进行点积操作，并求和，得到输出图像对应的像素值g(i, j);
c. 在输出图像上进行滑动，重复b的操作，求出输出图像上的所有像素值。
你可能会觉得相关操作和卷积操作一样，但二者实际上是有区别的。
卷积操作
图像卷积操作的公式：g = f ★ h，具体而言，对应到每个像素的表达式如下所示，可以看出，卷积操作满足交换律。
In Convolution operation, the kernel is first flipped by an angle of 180 degrees and is then applied to the image.
也就是说， 卷积操作与相关操作的差异在于卷积将核kernel旋转了180度（顺逆都一样）；相关可以反应两个信号的相似度，卷积不可以；卷积满足交换律，相关不可以；卷积可以直接通过卷积定理（时域上的卷积等于频域上的乘积）来加速运算，相关不可以。
参考文章：https://towardsdatascience.com/convolution-vs-correlation-af868b6b4fb5
MOSSE
基本思想
MOSSE的基本思想：以上一帧目标位置为bbox，在当前帧图像上截取目标，使用相关核找到当前截取图像上的最大响应，这里的最大响应就是当前目标的中心，然后更新当前帧目标的中心，并通过第一帧给出的bbox的宽高为宽高，生成此时物体的框体。（这里有一个假设，就是相邻帧物体不会发生过大的位移，即使用上一帧的bbox截取当前帧的图像，仍可以获得到目标的中心） --- 可以看出，MOSSE对于物体大小的变化，以及物体的消失没有丝毫抵抗能力。
具体操作流程
通过基本思想可以看出，需要首先获得一个可以在目标中心获得最大响应的相关核h。即下述公式操作后，g中最大值是物体的中心。
由于上式在计算机中计算量较大，作者对上式采用快速傅里叶变换(FFT):
便得到论文中给出的
因此，跟踪的任务便是找到相应的H*
在实际跟踪中，需要考虑目标外观的变化等因素的影响，需要进行一个优化的求解，即输出的响应与期望的响应越接近越好：
该优化的求解可以根据偏导数进行求取，求取结果可以直接表达为：
以上的操作均是在跟踪第一帧完成的。也就是说，相关核H的确定是通过第一帧标出的目标在线得出的。而上式中的i，对应到代码中即为若干次图像的仿射变换。而对应的输出响应的期望G*是高斯核，后续会结合代码具体解读。
为了提升算法的鲁棒性，应对目标在跟踪过程中的外观变换，需要对核进行在线的更新，更新策略如下所示：
值得一提的是，F、G、H的尺寸都是一致的。并且，在工程实践中确实将H分解为A，B两部分，并进行迭代更新。
代码解读
初始化
初始化应用于第一帧的功能实现，包括初始bbox信息的获取（中心点、宽、高、特征），以及核的生成。其中，包含图像的预处理和准备工作：
- 灰度化（输入特征为灰度，特征表达能力有限）
- 归一化
- 余弦窗的生成
- 期望输出响应G的生成（高斯核，同时让图片的左边缘与右边缘连接，上边缘与下边缘连接）
之后，便是最终要的核的生成。代码如下：
1 def init(self,first_frame,bbox): 2 if len(first_frame.shape)!=2: 3 assert first_frame.shape[2] == 3 4 first_frame=cv2.cvtColor(first_frame, cv2.COLOR_BGR2GRAY) 5 first_frame=first_frame.astype(np.float32)/255 # 归一化 6 x,y,w,h=tuple(bbox) # x,y为bbox的左下角点，和宽高 7 self._center=(x+w/2, y+h/2) # 中心点坐标 8 self.w,self.h=w,h 9 w,h=int(round(w)), int(round(h)) # 取整 10 self.cos_window=cos_window((w,h)) # 生成宽高的一维hanning窗 后再进行矩阵乘法 形成h*w的hanning矩阵 11 self._fi=cv2.getRectSubPix(first_frame,(w,h), self._center) # 根据宽高和中心点的位置截取图像/特征 12 self._G=np.fft.fft2(gaussian2d_labels((w,h), self.sigma)) # 高斯相应图 13 self.crop_size=(w,h) # 裁剪尺寸 14 self._Ai=np.zeros_like(self._G) 15 self._Bi=np.zeros_like(self._G) 16 for _ in range(8): 17 fi=self._rand_warp(self._fi) # 对裁剪的图像进行8次仿射变化(旋转)，从而求得H* 18 Fi=np.fft.fft2(self._preprocessing(fi,self.cos_window)) # 进行余弦窗处理 19 self._Ai += self._G*np.conj(Fi) # conj共轭负数 20 self._Bi += Fi*np.conj(Fi) # 求解H*所需
其中，代码第4-5行：图像转化为灰度图；
代码第6-9行：获取第一帧图像的信息，将左下角点，宽高信息转化为中心点坐标。（bbox是 通过cv2.selectROI()获得到目标的左下角点坐标和宽高信息。）
代码第10行：获取余弦窗；
代码第11行：通过中心点坐标和宽高在第一帧图像上截取目标图像；
代码第12行：获得高斯响应作为期望输出G。
代码第16-20行：通过8次仿射变换，获取相关核H。
余弦窗的生成
余弦窗通过汉宁窗生成。
1 def cos_window(sz): 2 cos_window = np.hanning(int(sz[1]))[:, np.newaxis].dot(np.hanning(int(sz[0]))[np.newaxis, :]) 3 return cos_window
汉宁窗的模样
通过下图余弦矩阵的颜色分布可以看出，余弦窗具有两边小中间大的特性，并且四周的值趋近于零，有利于突出靠近中心的目标。
期望响应输出G的生成
G相当于是真值，但此处不需要人为的标注。而是通过高斯函数生成即可。高斯函数产生的最大值在图像的中心。生成后，需要进行傅里叶变换转化到频域。
1 def gaussian2d_labels(sz,sigma): 2 w, h=sz 3 xs, ys = np.meshgrid(np.arange(w), np.arange(h)) # 生成两个矩阵 一个矩阵各行是0-w-1 一个矩阵各列是0-h-1 4 center_x, center_y = w / 2, h / 2 5 dist = ((xs - center_x) ** 2 + (ys - center_y) ** 2) / (sigma**2) 6 labels = np.exp(-0.5*dist) 7 return labels
为什么可以用高斯函数产生的值来作为期望的响应输出呢？这是由于高斯函数自带中间大两边小的特性，十分符合物体中心是最大值的要求。并且在第一帧框定目标时，也会标出十分贴合的bbox，bbox的中心可以近似为物体的中心。
随机仿射变换
通过随机数的生成控制仿射变换的力度。
1 def _rand_warp(self,img): # 输入是裁剪下的图形 2 h, w = img.shape[:2] 3 C = .1 4 ang = np.random.uniform(-C, C) 5 c, s = np.cos(ang), np.sin(ang) 6 W = np.array([[c + np.random.uniform(-C, C), -s + np.random.uniform(-C, C), 0], 7 [s + np.random.uniform(-C, C), c + np.random.uniform(-C, C), 0]]) 8 center_warp = np.array([[w / 2], [h / 2]]) 9 tmp = np.sum(W[:, :2], axis=1).reshape((2, 1)) 10 W[:, 2:] = center_warp - center_warp * tmp # 仿射矩阵 11 warped = cv2.warpAffine(img, W, (w, h), cv2.BORDER_REFLECT) # 进行仿射变换 W是仿射矩阵 12 return warped
在线更新
在线更新核心操作在于相关核的更新（下述代码第6行）以及影响相关核的两个因素的更新（下述代码22-23行）。
从代码14-18行以及返回值可以看出，在线更新的过程只有中心点位置的更新，宽高不更新，如果物体的大小发生明显的变化，很难自适应。
1 def update(self,current_frame,vis=False): 2 if len(current_frame.shape)!=2: 3 assert current_frame.shape[2]==3 4 current_frame=cv2.cvtColor(current_frame,cv2.COLOR_BGR2GRAY) 5 current_frame=current_frame.astype(np.float32)/255 6 Hi=self._Ai/self._Bi # Ht kernel 7 fi=cv2.getRectSubPix(current_frame,(int(round(self.w)),int(round(self.h))),self._center) # 裁剪 8 fi=self._preprocessing(fi,self.cos_window) 9 Gi=Hi*np.fft.fft2(fi) # fft 10 gi=np.real(np.fft.ifft2(Gi)) # 返回负数类型参数的实部， 以及二维傅里叶反变换 11 if vis is True: 12 self.score=gi 13 curr=np.unravel_index(np.argmax(gi, axis=None),gi.shape) # 找到最大那个相应的位置 14 dy,dx=curr[0]-(self.h/2),curr[1]-(self.w/2) # 中心点位置移动量 15 x_c,y_c=self._center 16 x_c+=dx 17 y_c+=dy # 中心点更新 18 self._center=(x_c,y_c) 19 fi=cv2.getRectSubPix(current_frame,(int(round(self.w)),int(round(self.h))),self._center) 20 fi=self._preprocessing(fi,self.cos_window) 21 Fi=np.fft.fft2(fi) 22 self._Ai=self.interp_factor*(self._G*np.conj(Fi))+(1-self.interp_factor)*self._Ai 23 self._Bi=self.interp_factor*(Fi*np.conj(Fi))+(1-self.interp_factor)*self._Bi 24 return [self._center[0]-self.w/2,self._center[1]-self.h/2,self.w,self.h]
缺点
- 输入的特征为单通道灰度图像，特征表达能力有限
- 没有尺度更新，对于尺度变化的跟踪目标不敏感