【计算机视觉】siamfc论文复现实现目标追踪
什么是目标跟踪
使用视频序列第一帧的图像(包括bounding box的位置)，来找出目标出现在后序帧位置的一种方法。
什么是孪生网络结构
孪生网络结构其思想是将一个训练样本(已知类别)和一个测试样本(未知类别)输入到两个CNN(这两个CNN往往是权值共享的)中，从而获得两个特征向量，然后通过计算这两个特征向量的的相似度，相似度越高表明其越可能是同一个类别。
给你一张我的正脸照（没有经过美颜处理的)，你该如何在人群中找到我呢？一种最直观的方案就是:“谁长得最像就是谁”。但是对于计算机来说，如何衡量“长得像”，并不是个简单的问题。这就涉及一种基本的运算——互相关(cross-correlation)。互相关运算可以用来度量两个信号之间的相似性。互相关得到的响应图中每个像素的响应高低代表着每个位置相似度的高低。
在目标领域中，最早利用这种思想的是SiamFC，其网络结构如上图。图中的φ就是CNN编码器，上下两个分支使用的CNN不仅结构相同，参数也是完全共享的（说白了就是同一个网络，并不存在孪生兄弟那样的设定)。z和x分别是要跟踪的目标模版图像(尺寸为127x127）和新的一帧中的搜索范围(尺寸为255x255)。二者经过同样的编码器后得到各自的特征图，对二者进行互相关运算后则会同样得到一个响应图(尺寸为17x17)，其每一个像素的值对应了x中与z等大的一个对应区域出现跟踪目标的概率。
互相关运算的步骤，像极了我们手里拿着一张目标的照片（模板图像），然后把这个照片按在需要寻找目标的图片上（搜索图像）进行移动，然后求重叠部分的相似度，从而找到这个目标，只不过为了计算机计算的方便，使用AlexNet对图像数据进行了编码/特征提取。
下面这个版本中有一些动图，还是会帮助理解的：https://github.com/rafellerc/Pytorch-SiamFC
SiamFC代码分析
我们对siamese的结构大致就讲完了，还有一些内容结合代码来讲，效果更好。
3.0 SiameseFC的结构
主要由以下几部分构成
- 特征提取网络AlexNet
- 互相关运算网络
3.0.1 图像特征提取网络
class AlexNet(nn.Module):def __init__(self):super(AlexNet, self).__init__()self.conv1 = nn.Sequential(nn.Conv2d(3, 96, 11, stride=2, padding=0),nn.BatchNorm2d(96),nn.ReLU(inplace=True),nn.MaxPool2d(kernel_size=3, stride=2, padding=0))self.conv2 = nn.Sequential(nn.Conv2d(96, 256, 5, stride=1, padding=0, groups=2),nn.BatchNorm2d(256),nn.ReLU(inplace=True),nn.MaxPool2d(kernel_size=3, stride=2, padding=0))self.conv3 = nn.Sequential(nn.Conv2d(256, 384, 3, stride=1, padding=0),nn.BatchNorm2d(384),nn.ReLU(inplace=True))self.conv4 = nn.Sequential(nn.Conv2d(384, 384, 3, stride=1, padding=0, groups=2),nn.BatchNorm2d(384),nn.ReLU(inplace=True))self.conv5 = nn.Sequential(nn.Conv2d(384, 256, 3, stride=1, padding=0, groups=2))def forward(self, x):conv1 = self.conv1(x)conv2 = self.conv2(conv1)conv3 = self.conv3(conv2)conv4 = self.conv4(conv3)conv5 = self.conv5(conv4)return conv5
3.0.2 互相关运算网络
class _corr(nn.Module):def __init__(self):super(_corr, self).__init__()#互相关运算，设batch_size=8def forward(self, z, x):kernel = z #[8,128,6,6]group = z.size(0) #8input = x.view(-1, group*x.size(1), x.size(2), x.size(3))#输出为[8,1,17,17], 那么反推input[8,128,22,22]，kernel[1,1024,6,6] group=128/1024？错误#所以先输出[1,8,17,17],再view变换维度成[8,1,17,17],那么input[1,1024,22,22],kernel[8,128,6,6],group=1024/128=8=batch_sizeresponse_maps = F.conv2d(input, kernel,groups=group)response_maps = response_maps.view(x.size(0),-1,response_maps.size(2), response_maps.size(3))return response_maps
至此，完成了对网络框架architecture的代码分析!!!
3.1 training
3.1.1图像预处理
小超up给出训练的框图如下。训练过程中，首先要获取训练数据集的所有视频序列(每个视频序列的所有帧)，我采用的是GOT-10k数据集训练；获取数据集之后进行图像预处理，对每一个视频序列抽取两帧图像并作数据增强处理(包括裁剪、resize等过程)，分别作为目标模板图像和搜索图像；把经过图像处理的所有图像对加载并以batch_size输入网络得到预测输出；建立标签和损失函数，损失函数的输入是预测输出，目标是标签；设置优化策略，梯度下降损失，最终得到网络模型。
先贴代码，再分析：
def train(data_dir, net_path=None,save_dir='pretrained'):#从文件中读取图像数据集seq_dataset = GOT10k(data_dir,subset='train',return_meta=False)#定义图像预处理方法transforms = SiamFCTransforms( exemplar_sz=cfg.exemplar_sz, #127instance_sz=cfg.instance_sz, #255context=cfg.context) #0.5#从读取的数据集每个视频序列配对训练图像并进行预处理，裁剪等train_dataset = GOT10kDataset(seq_dataset,transforms)
data_dir是存放GOT-10k数据集的文件路径，GOT-10k一共有9335个训练视频序列，seq_dataset返回的是所有视频序列的图片路径列表seq_dirs及对应groundtruth列表anno_files及一些其他信息，如下：
接下来是定义好图像预处理方法，在GOT10kDataset方法中对每个视频序列配对两帧图像，并使用定义好的图像处理方法，接下来直接进入该方法分析代码，GOT10kDataset的代码如下：
class GOT10kDataset(Dataset): #继承了torch.utils.data的Dataset类def __init__(self, seqs, transforms=None,pairs_per_seq=1):def __getitem__(self, index): #通过_sample_pair方法得到索引返回item=(z,x,box_z,box_x),然后经过transforms处理def __len__(self): #返回9335*pairs_per_seq对def _sample_pair(self, indices): #随机挑选两个索引，这里取的间隔不超过T=100def _filter(self, img0, anno, vis_ratios=None): #通过该函数筛选符合条件的有效索引val_indices
这里最重要的方法就是__getitem__，该方法最终返回处理后的图像，在内部首先调用了_sample_pair方法，用于提取两帧有效图片(有效的定义是图片目标的面积和高宽等有约束条件)的索引
因为这里我对网络在训练和测试时需要的数据集有一些疑问，通过源码以及网上博文做以下总结：
答：在训练过程中，一对图片（即模板图和搜索图）共同构成一个数据样本，而不是单独的每一帧。这意味着在Siamese网络的训练中，数据集是由多个这样的图像对组成的，每个图像对代表一个单独的训练样本。
在测试阶段，算法逐帧处理视频序列，每次只处理一个搜索图和对应的模板图。
def track(self, img_files, box, visualize=False):frame_num = len(img_files)boxes = np.zeros((frame_num, 4))boxes[0] = boxtimes = np.zeros(frame_num)for f, img_file in enumerate(img_files):img = ops.read_image(img_file)begin = time.time()if f == 0:self.init(img, box)else:boxes[f, :] = self.update(img)times[f] = time.time() - begin可见在测试阶段，在视频序列的第一帧选择目标区域作为模板图。这个模板图将用于整个视频序列的跟踪过程。
在得到这两帧图片和对应groundtruth之后通过定义好的transforms进行处理，transforms是SiamFCTransforms类的实例化对象，该类中主要继承了resize图片大小和各种裁剪方式等，如代码所示：
class SiamFCTransforms(object):def __init__(self, exemplar_sz=127, instance_sz=255, context=0.5):self.exemplar_sz = exemplar_szself.instance_sz = instance_szself.context = context#transforms_z/x是数据增强方法self.transforms_z = Compose([RandomStretch(), #随机resize图片大小,变化再[1 1.05]之内CenterCrop(instance_sz - 8), #中心裁剪 裁剪为255-8RandomCrop(instance_sz - 2 * 8), #随机裁剪 255-8->255-8-8CenterCrop(exemplar_sz), #中心裁剪 255-8-8->127ToTensor()]) #图片的数据格式从numpy转换成torch张量形式self.transforms_x = Compose([RandomStretch(), #s随机resize图片CenterCrop(instance_sz - 8), #中心裁剪 裁剪为255-8RandomCrop(instance_sz - 2 * 8), #随机裁剪 255-8->255-8-8ToTensor()]) #图片数据格式转化为torch张量def __call__(self, z, x, box_z, box_x): #z，x表示传进来的图像z = self._crop(z, box_z, self.instance_sz) #对z(x类似)图像 1、box转换(l,t,w,h)->(y,x,h,w)，并且数据格式转为float32,得到center[y,x],和target_sz[h,w]x = self._crop(x, box_x, self.instance_sz) #2、得到size=((h+(h+w)/2)*(w+(h+2)/2))^0.5*255(instance_sz)/127z = self.transforms_z(z) #3、进入crop_and_resize:传入z作为图片img，center，size，outsize=255(instance_sz),随机选方式填充，均值填充x = self.transforms_x(x) # 以center为中心裁剪一块边长为size大小的正方形框(注意裁剪时的padd边框填充问题)，再resize成out_size=255(instance_sz)return z, x
实例化对象后，直接从__call__开始运行代码，首先关注的应该是_crop函数，该函数将原始的两帧图片分别以目标为中心，裁剪一块包含上下文信息的patch，patch的边长定义如下：
式中，w、h分别表示目标的宽和高。下面具体讲里面的_crop函数：
def _crop(self, img, box, out_size):# convert box to 0-indexed and center based [y, x, h, w]box = np.array([box[1] - 1 + (box[3] - 1) / 2,box[0] - 1 + (box[2] - 1) / 2,box[3], box[2]], dtype=np.float32)center, target_sz = box[:2], box[2:]context = self.context * np.sum(target_sz)size = np.sqrt(np.prod(target_sz + context))size *= out_size / self.exemplar_szavg_color = np.mean(img, axis=(0, 1), dtype=float)interp = np.random.choice([cv2.INTER_LINEAR,cv2.INTER_CUBIC,cv2.INTER_AREA,cv2.INTER_NEAREST,cv2.INTER_LANCZOS4])patch = ops.crop_and_resize(img, center, size, out_size,border_value=avg_color, interp=interp)return patch
因为GOT-10k里面对于目标的bbox是以ltwh(即left, top, weight, height)形式给出的，上述代码一开始就先把输入的box变成center based，坐标形式变为[y, x, h, w]，结合下面这幅图就非常好理解
crop_and_resize：
def crop_and_resize(img, center, size, out_size,border_type=cv2.BORDER_CONSTANT,border_value=(0, 0, 0),interp=cv2.INTER_LINEAR):# convert box to corners (0-indexed)size = round(size) # the size of square cropcorners = np.concatenate((np.round(center - (size - 1) / 2),np.round(center - (size - 1) / 2) + size))corners = np.round(corners).astype(int)# pad image if necessarypads = np.concatenate((-corners[:2], corners[2:] - img.shape[:2]))npad = max(0, int(pads.max()))if npad > 0:img = cv2.copyMakeBorder(img, npad, npad, npad, npad,border_type, value=border_value)# crop image patchcorners = (corners + npad).astype(int)patch = img[corners[0]:corners[2], corners[1]:corners[3]]# resize to out_sizepatch = cv2.resize(patch, (out_size, out_size),interpolation=interp)return patch
在裁剪过程中会出现越界的情况，需要对原始图像边缘填充，填充值固定为图像的RGB均值，填充大小根据图像边缘越界最大值作为填充值，具体实现过程由以下代码完成。
# padding操作#corners表示目标的[ymin,xmin,ymax,xmax]pads = np.concatenate((-corners[:2], corners[2:] - img.shape[:2]))npad = max(0, int(pads.max())) #得到上下左右4个越界值中最大的与0对比，<0代表无越界if npad > 0:img = cv2.copyMakeBorder(img, npad, npad, npad, npad,cv2.BORDER_CONSTANT, value=img_average)
实验结果：
3.1.2加载训练数据、标签及损失函数
图像预处理完成后，得到了用与训练的9335对图像，将图像加载批量加载输入网络得到输出结果作为损失函数的input，损失函数的target是制定好的labels。
#加载训练数据集loader_dataset = DataLoader( dataset = train_dataset,batch_size=cfg.batch_size,shuffle=True,num_workers=cfg.num_workers,pin_memory=True,drop_last=True, )#初始化训练网络cuda = torch.cuda.is_available() #支持GPU为Truedevice = torch.device('cuda:0' if cuda else 'cpu') #cuda设备号为0model = AlexNet(init_weight=True)corr = _corr()model = model.to(device)corr = corr.to(device)# 设置损失函数和标签logist_loss = BalancedLoss()labels = _create_labels(size=[cfg.batch_size, 1, cfg.response_sz - 2, cfg.response_sz - 2])labels = torch.from_numpy(labels).to(device).float()
本小节主要讲网络输出的labels和损失函数，接下来只是小超up个人的一些理解，代码与论文理论部分形式不一致，但效果一样。先上图，论文中labels以及损失函数如下图：
然而代码中的labels值却是1和0，损失函数使用的是二值交叉熵损失函数F.binary_cross_entropy_with_logits，如下图推导所示，解释了为什么代码实现部分真正使用的labels值是1和0，而理论部分使用的是1和-1。
利用下面代码的这个_creat_labels方法可以得到标签。
def _create_labels(size):def logistic_labels(x, y, r_pos):# x^2+y^2<4 的位置设为为1，其他为0dist = np.sqrt(x ** 2 + y ** 2)labels = np.where(dist <= r_pos, #r_os=2np.ones_like(x), #np.ones_like(x),用1填充xnp.zeros_like(x)) #np.zeros_like(x),用0填充xreturn labels#获取标签的参数n, c, h, w = size # [8,1,15,15]x = np.arange(w) - (w - 1) / 2 #x=[-7 -6 ....0....6 7]y = np.arange(h) - (h - 1) / 2 #y=[-7 -6 ....0....6 7]x, y = np.meshgrid(x, y) #建立标签r_pos = cfg.r_pos / cfg.total_stride # 16/8labels = logistic_labels(x, y, r_pos)#重复batch_size个label，因为网络输出是batch_size张response maplabels = labels.reshape((1, 1, h, w)) #[1,1,15,15]labels = np.tile(labels, (n, c, 1, 1)) #将labels扩展[8,1,15,15]return labels
验证结果如下图，只截取了部分labels，得到的labels对应输入，大小都是[8,1,15,15]
if __name__ == '__main__':labels = _create_labels([8,1,15,15]) #返回的label.shape=(8,1,15,15)
其中关于np.tile、np.meshgrid、np.where函数的使用可以去看这篇博客，最后出来的一个batch下某一个通道下的label就是下面这样的
3.1.3 优化策略
这里主要说一下学习率lr，随着训练次数epoch增多而减小，具体值如下公式，式中，initial为初始学习率，gamma是定义的超参，epoch为训练次数。整个优化器及学习率调整实现代码如下：
#建立优化器，设置指数变化的学习率optimizer = optim.SGD(model.parameters(),lr=cfg.initial_lr, #初始化的学习率，后续会不断更新weight_decay=cfg.weight_decay, #λ=5e-4，正则化momentum=cfg.momentum) #v(now)=dx∗lr+v(last)∗momemtumgamma = np.power( #np.power(a,b) 返回a^bcfg.ultimate_lr / cfg.initial_lr,1.0 / cfg.epoch_num)lr_scheduler = ExponentialLR(optimizer, gamma) #指数形式衰减，lr=initial_lr*(gamma^epoch)=
3.1.4 模型的训练与保存
一切准备工作就绪后，就开始训练了。代码中设定epoch_num为50次，训练时密切加上model.train()，告诉网络处于训练状态，这样，网络运行时就会利用pytorch的自动求导机制求导；在测试时，改为model.eval()，关闭自动求导。模型训练的步骤如代码所示：
# loop over epochs
for epoch in range(self.cfg.epoch_num):# update lr at each epochself.lr_scheduler.step(epoch=epoch)# loop over dataloaderfor it, batch in enumerate(dataloader):loss = self.train_step(batch, backward=True)print('Epoch: {} [{}/{}] Loss: {:.5f}'.format(epoch + 1, it + 1, len(dataloader), loss))sys.stdout.flush()# save checkpointif not os.path.exists(save_dir):os.makedirs(save_dir)net_path = os.path.join(save_dir, 'siamfc_alexnet_e%d.pth' % (epoch + 1))torch.save(self.net.state_dict(), net_path)
至此此份repo的训练应该差不多结束了
3.2 testing
3.2.1 init(初始帧)
这一步目的是要得到6x6x128的目标模板feature map，在后续跟踪过程中，一直使用它作为卷积核不变，这也是一大SiamFC的缺点所在；init内还要初始化目标的中心所在位置、宽高、汉宁窗等后续跟踪使用。
#传入第一帧的gt和图片，初始化一些参数，计算一些之后搜索区域的中心等等def init(self, img, box):#设置成评估模式，测试模型时一开始要加这个，属于Pytorch，训练模型前，要加self.net.train()self.model.eval()#将原始的目标位置表示[l,t,w,h]->[center_y,center_x,h,w]yxhw = ltwh_to_yxhw(ltwh=box)self.center, self.target_sz = yxhw[:2], yxhw[2:]#创建汉宁窗口update使用self.response_upsz = cfg.response_up * cfg.response_sz # 16*17=272self.hann_window = creat_hanning_window(size=self.response_upsz)#三种尺度1.0375**(-1,0,1) 三种尺度self.scale_factors = three_scales()# patch边长,两种边长：目标模板图像z_sz和搜索图像x_zscontext = cfg.context * np.sum(self.target_sz) # 上下文信息(h+w)/2self.z_sz = np.sqrt(np.prod(self.target_sz + context)) # (h+(h+w)/2)*(w+(h+2)/2))^0.5self.x_sz = self.z_sz * cfg.instance_sz / cfg.exemplar_sz # (h+(h+w)/2)*(w+(h+2)/2))^0.5*255/127#图像的RGB均值，返回(aveR,aveG,aveB)self.avg_color = np.mean(img, axis=(0, 1))#裁剪一块以目标为中心，边长为z_sz大小的patch,然后将其resize成exemplar_sz的大小z = z_to127(img, self.center, self.z_sz, cfg.exemplar_sz, self.avg_color)z = torch.from_numpy(z).to( #torch.size=([1,3,127,127])self.device).permute(2, 0, 1).unsqueeze(0).float()self.kernel = self.model(z) #torch.size=([1,128,6,6])
follow代码，先设置model.eval()，原因在training过程中讲了，初始化目标的位置信息self.center和self.target_sz分别代表目标的中心位置和高宽；接着是汉宁窗口和尺度因子，汉宁窗口如下图展示的是17x17大小的，它与response map相乘会突出response map的中心并抑制边缘。窗口大小是17*16=272，而不是response map的大小17，这样做的好处是找到更加精确的位置，相比于在17x17大小的特征图上寻找最大值，272x272大小的特征图分辨率更高，位置更精确，所以在后续跟踪过程中，得到的17x17大小的response map会resize成272x272大小；代码使用的是SiamFC_3s，所以定义了三个不同的尺度；包含上下文信息的patch边长与training不同，训练时是将两张图片以相同的patch裁剪再resize成255x255x3大小，之后在一系列数据增强的过程中将其中一张图片裁剪成127x127x3作为目标模板图像。而test过程中，目标模板图像和搜索图像的patch边长不同，如代码所示，成比例255/127关系；最后就是对init传入的初始帧图片裁剪并输入网络得到跟踪过程使用的卷积核，大小为[1,128,6,6]。
3.2.2 update(后续帧)
Update就是tracking过程，传入的参数是后续帧图片，一共可以分为三个阶段：1、经过网络的正向推导得到response map；2、根据response map的最大值反推目标在原始图片中的位置；3、参数更新。
#传入后续帧，然后根据SiamFC跟踪流程返回目标的box坐标def update(self, img):self.model.eval()"""----------------正向推导得到response map最大值位置---------------------"""#三种patch边长 patch*3scalesx = x_to3s255(img,self.center,self.x_sz,self.scale_factors,cfg.instance_sz,self.avg_color)#numpy转为float的torch型张量x = torch.from_numpy(x).to(self.device).permute(0, 3, 1, 2).float()#[3,255,22,22]x = self.model(x)#得到三种尺度下的response mapresponses = self.corr(self.kernel, x) * cfg.out_reduce #[3，1，17，17]responses = responses.squeeze(1).cpu().numpy() #压缩为[3，17，17]并转为numpy作后续计算处理#将17x17大小的response map->[3,272,272]responses = map_to272(responses,out_size=self.response_upsz)#对尺度变化做出相应惩罚responses[:cfg.scale_num // 2] *= cfg.scale_penalty #response[0]*（0.9745惩罚项）responses[cfg.scale_num // 2 + 1:] *= cfg.scale_penalty #response[2]*（0.9745惩罚项）#找到最大值属于哪个response map，并把该response map赋给responsescale_id = np.argmax(np.amax(responses, axis=(1, 2))) #里面求得三个map的最大值 再对三个值求最大值 得到索引response = responses[scale_id] #[272，272]#一系列数据处理，重点在汉宁窗惩罚response = map_process(response,self.hann_window)loc = np.unravel_index(response.argmax(),response.shape) #unravel_index该函数可返回索引response.argmax()的元素的坐标，逐行拉伸，返回第几行第几个"""---------------由response map最大值位置反推目标在原图像的位置------------"""disp_in_response = np.array(loc) - (self.response_upsz - 1) / 2 #峰值点相对于response中心的位移disp = disp_in_response / 16disp = disp * 8disp = disp * self.x_sz * self.scale_factors[scale_id] / cfg.instance_szself.center += disp"""---------------参数更新------------"""
参考文档
siameseFC论文和代码解析
SiamFC 学习（论文、总结与分析）
siamfc-pytorch代码讲解（一）：backbone&head
siamfc-pytorch代码讲解（二）：train&siamfc
SiamFC代码分析(architecture、training、test)
siamfc前世今生
Siamese跟踪发展历程
视频推荐
目标跟踪零基础代码入门（一）：SiamFC_哔哩哔哩_bilibili