这段代码实现了基于运动补偿的前景检测算法。

主要功能包括：

• 运动补偿模块：使用基于网格的 KLT 特征跟踪算法计算两帧之间的运动，然后通过单应性变换实现帧间运动补偿。
• 前景检测模块：结合两帧运动补偿结果，通过帧间差分计算前景掩码。
• 异常处理：添加了图像加载检查和异常捕获，提高了代码的健壮性。
• 路径处理：自动创建保存目录，避免因目录不存在导致的错误。

使用时需要提供三帧连续图像：两个参考帧和当前帧。代码会计算出前景掩码并保存为图像文件。

import cv2
import numpy as np
import os
import sysdef motion_compensate(frame1, frame2):"""使用基于网格的KLT特征跟踪实现两帧之间的运动补偿参数:frame1: 前一帧图像(BGR格式)frame2: 当前帧图像(BGR格式)返回:compensated: 运动补偿后的图像mask: 补偿区域的掩码avg_dst: 平均运动距离motion_x: x方向平均运动量motion_y: y方向平均运动量homography_matrix: 单应性变换矩阵"""# 设置LK光流参数lk_params = dict(winSize=(15, 15), maxLevel=3,criteria=(cv2.TERM_CRITERIA_EPS | cv2.TERM_CRITERIA_COUNT, 30, 0.003))# 图像预处理和网格点生成width = frame2.shape[1]height = frame2.shape[0]scale = 2  # 放大图像以获得更精确的跟踪# 调整图像大小以提高特征点检测精度frame1_grid = cv2.resize(frame1, (960 * scale, 540 * scale), dst=None, interpolation=cv2.INTER_CUBIC)frame2_grid = cv2.resize(frame2, (960 * scale, 540 * scale), dst=None, interpolation=cv2.INTER_CUBIC)width_grid = frame2_grid.shape[1]height_grid = frame2_grid.shape[0]gridSizeW = 32 * 2  # 网格宽度gridSizeH = 24 * 2  # 网格高度# 生成网格点作为特征点p1 = []grid_numW = int(width_grid / gridSizeW - 1)grid_numH = int(height_grid / gridSizeH - 1)for i in range(grid_numW):for j in range(grid_numH):# 将点放置在每个网格中心point = (np.float32(i * gridSizeW + gridSizeW / 2.0), np.float32(j * gridSizeH + gridSizeH / 2.0))p1.append(point)p1 = np.array(p1)pts_num = grid_numW * grid_numHpts_prev = p1.reshape(pts_num, 1, 2)# 计算光流pts_cur, st, err = cv2.calcOpticalFlowPyrLK(frame1_grid, frame2_grid, pts_prev, None, **lk_params)# 选择跟踪成功的点good_new = pts_cur[st == 1]  # 当前帧中的跟踪点good_old = pts_prev[st == 1]  # 前一帧中的跟踪点# 计算运动距离和位移motion_distance = []translate_x = []translate_y = []for i, (new, old) in enumerate(zip(good_new, good_old)):a, b = new.ravel()c, d = old.ravel()motion_distance0 = np.sqrt((a - c) * (a - c) + (b - d) * (b - d))# 过滤异常大的运动值if motion_distance0 > 50:continuetranslate_x0 = a - ctranslate_y0 = b - dmotion_distance.append(motion_distance0)translate_x.append(translate_x0)translate_y.append(translate_y0)motion_dist = np.array(motion_distance)motion_x = np.mean(np.array(translate_x)) if translate_x else 0motion_y = np.mean(np.array(translate_y)) if translate_y else 0avg_dst = np.mean(motion_dist) if motion_dist.size > 0 else 0# 计算单应性变换矩阵if len(good_old) < 15:# 点太少时使用近似恒等变换homography_matrix = np.array([[0.999, 0, 0], [0, 0.999, 0], [0, 0, 1]])else:# 使用RANSAC算法估计单应性矩阵homography_matrix, status = cv2.findHomography(good_new, good_old, cv2.RANSAC, 3.0)# 应用单应性变换进行运动补偿compensated = cv2.warpPerspective(frame1, homography_matrix, (width, height), flags=cv2.INTER_LINEAR + cv2.WARP_INVERSE_MAP)# 生成掩码以指示变换区域vertex = np.array([[0, 0], [width, 0], [width, height], [0, height]], dtype=np.float32).reshape(-1, 1, 2)homo_inv = np.linalg.inv(homography_matrix)vertex_trans = cv2.perspectiveTransform(vertex, homo_inv)vertex_transformed = np.array(vertex_trans, dtype=np.int32).reshape(1, 4, 2)im = np.zeros(frame1.shape[:2], dtype='uint8')cv2.polylines(im, vertex_transformed, 1, 255)cv2.fillPoly(im, vertex_transformed, 255)mask = 255 - imreturn compensated, mask, avg_dst, motion_x, motion_y, homography_matrixdef FD_mask(lastFrame1, lastFrame2, currentFrame, save_path='data/mask.jpg'):"""使用两帧运动补偿计算前景掩码参数:lastFrame1: 第一参考帧(BGR格式)lastFrame2: 第二参考帧(BGR格式)currentFrame: 当前帧(BGR格式)save_path: 结果掩码保存路径"""# 图像预处理:高斯模糊和灰度转换lastFrame1 = cv2.GaussianBlur(lastFrame1, (11, 11), 0)lastFrame1 = cv2.cvtColor(lastFrame1, cv2.COLOR_BGR2GRAY)lastFrame2 = cv2.GaussianBlur(lastFrame2, (11, 11), 0)lastFrame2 = cv2.cvtColor(lastFrame2, cv2.COLOR_BGR2GRAY)currentFrame = cv2.GaussianBlur(currentFrame, (11, 11), 0)currentFrame = cv2.cvtColor(currentFrame, cv2.COLOR_BGR2GRAY)# 计算第一参考帧到第二参考帧的运动补偿img_compensate1, mask1, avg_dist1, motion_x1, motion_y1, homo_matrix = motion_compensate(lastFrame1, lastFrame2)frameDiff1 = cv2.absdiff(lastFrame2, img_compensate1)# 计算当前帧到第二参考帧的运动补偿img_compensate2, mask2, avg_dist2, motion_x2, motion_y2, homo_matrix2 = motion_compensate(currentFrame, lastFrame2)frameDiff2 = cv2.absdiff(lastFrame2, img_compensate2)# 融合两个差分结果frameDiff = (frameDiff1 + frameDiff2) / 2# 保存结果os.makedirs(os.path.dirname(save_path), exist_ok=True)cv2.imwrite(save_path, frameDiff)print(f'前景掩码已保存至: {save_path}')return frameDiffif __name__ == "__main__":# 示例:加载三帧图像并计算前景掩码# 请确保这些图像存在，或者修改为您自己的图像路径try:lastFrame1 = cv2.imread('data/Test_images/images/phantom05_0600.jpg')lastFrame3 = cv2.imread('data/Test_images/images/phantom05_0602.jpg')currentFrame = cv2.imread('data/Test_images/images/phantom05_0604.jpg')if lastFrame1 is None or lastFrame3 is None or currentFrame is None:print("错误: 无法加载图像，请检查文件路径!")else:FD_mask(lastFrame1, lastFrame3, currentFrame)except Exception as e:print(f"程序执行出错: {e}")    

#include <opencv2/opencv.hpp>
#include <vector>
#include <iostream>
#include <numeric>
#include <cmath>using namespace cv;
using namespace std;struct MotionCompensationResult {Mat compensated;Mat mask;float avg_dst;float motion_x;float motion_y;Mat homography_matrix;
};MotionCompensationResult motion_compensate(const Mat& frame1, const Mat& frame2) {// KLT 跟踪参数TermCriteria term_criteria(TermCriteria::EPS | TermCriteria::COUNT, 30, 0.003);Size win_size(15, 15);int max_level = 3;// 图像缩放int scale = 2;Mat frame1_grid, frame2_grid;resize(frame1, frame1_grid, Size(960 * scale, 540 * scale), 0, 0, INTER_CUBIC);resize(frame2, frame2_grid, Size(960 * scale, 540 * scale), 0, 0, INTER_CUBIC);// 创建网格点int gridSizeW = 32 * 2;int gridSizeH = 24 * 2;int grid_numW = static_cast<int>(frame2_grid.cols / gridSizeW - 1);int grid_numH = static_cast<int>(frame2_grid.rows / gridSizeH - 1);vector<Point2f> p1;for (int i = 0; i < grid_numW; i++) {for (int j = 0; j < grid_numH; j++) {p1.push_back(Point2f(i * gridSizeW + gridSizeW / 2.0f, j * gridSizeH + gridSizeH / 2.0f));}}int pts_num = grid_numW * grid_numH;Mat pts_prev = Mat(p1).reshape(2, pts_num);// 计算光流vector<Point2f> pts_cur;vector<uchar> status;vector<float> err;calcOpticalFlowPyrLK(frame1_grid, frame2_grid, pts_prev, pts_cur, status, err, win_size, max_level, term_criteria);// 筛选好点vector<Point2f> good_new, good_old;for (size_t i = 0; i < status.size(); i++) {if (status[i]) {good_new.push_back(pts_cur[i]);good_old.push_back(p1[i]);}}// 计算运动距离和位移vector<float> motion_distance;vector<float> translate_x, translate_y;for (size_t i = 0; i < good_new.size(); i++) {float dx = good_new[i].x - good_old[i].x;float dy = good_new[i].y - good_old[i].y;float dist = sqrt(dx * dx + dy * dy);if (dist > 50) continue;motion_distance.push_back(dist);translate_x.push_back(dx);translate_y.push_back(dy);}// 计算平均值float avg_dst = 0, motion_x = 0, motion_y = 0;if (!motion_distance.empty()) {avg_dst = accumulate(motion_distance.begin(), motion_distance.end(), 0.0f) / motion_distance.size();}if (!translate_x.empty()) {motion_x = accumulate(translate_x.begin(), translate_x.end(), 0.0f) / translate_x.size();motion_y = accumulate(translate_y.begin(), translate_y.end(), 0.0f) / translate_y.size();}// 计算单应性矩阵Mat homography_matrix;if (good_old.size() < 15) {homography_matrix = (Mat_<double>(3, 3) << 0.999, 0, 0, 0, 0.999, 0, 0, 0, 1);} else {homography_matrix = findHomography(good_new, good_old, RANSAC, 3.0);}// 运动补偿Mat compensated;warpPerspective(frame1, compensated, homography_matrix, Size(frame1.cols, frame1.rows), INTER_LINEAR + WARP_INVERSE_MAP);// 计算掩膜vector<Point2f> vertex;vertex.push_back(Point2f(0, 0));vertex.push_back(Point2f(frame1.cols, 0));vertex.push_back(Point2f(frame1.cols, frame1.rows));vertex.push_back(Point2f(0, frame1.rows));Mat vertex_mat = Mat(vertex).reshape(2);Mat homo_inv = homography_matrix.inv();vector<Point2f> vertex_trans;perspectiveTransform(vertex_mat, vertex_trans, homo_inv);vector<Point> vertex_transformed;for (const auto& pt : vertex_trans) {vertex_transformed.push_back(Point(static_cast<int>(pt.x), static_cast<int>(pt.y)));}Mat mask = Mat::zeros(frame1.size(), CV_8UC1);vector<vector<Point>> contours;contours.push_back(vertex_transformed);polylines(mask, contours, true, Scalar(255), 1);fillPoly(mask, contours, Scalar(255));mask = 255 - mask;return {compensated, mask, avg_dst, motion_x, motion_y, homography_matrix};
}void FD_mask(const Mat& lastFrame1, const Mat& lastFrame2, const Mat& currentFrame, const string& save_path = "mask.jpg") {// 图像预处理Mat lastFrame1_gray, lastFrame2_gray, currentFrame_gray;GaussianBlur(lastFrame1, lastFrame1_gray, Size(11, 11), 0);GaussianBlur(lastFrame2, lastFrame2_gray, Size(11, 11), 0);GaussianBlur(currentFrame, currentFrame_gray, Size(11, 11), 0);cvtColor(lastFrame1_gray, lastFrame1_gray, COLOR_BGR2GRAY);cvtColor(lastFrame2_gray, lastFrame2_gray, COLOR_BGR2GRAY);cvtColor(currentFrame_gray, currentFrame_gray, COLOR_BGR2GRAY);// 第一组运动补偿auto result1 = motion_compensate(lastFrame1_gray, lastFrame2_gray);Mat frameDiff1;absdiff(lastFrame2_gray, result1.compensated, frameDiff1);// 第二组运动补偿auto result2 = motion_compensate(currentFrame_gray, lastFrame2_gray);Mat frameDiff2;absdiff(lastFrame2_gray, result2.compensated, frameDiff2);// 计算最终差分Mat frameDiff;frameDiff1.convertTo(frameDiff1, CV_32F);frameDiff2.convertTo(frameDiff2, CV_32F);frameDiff = (frameDiff1 + frameDiff2) / 2;frameDiff.convertTo(frameDiff, CV_8U);// 保存结果imwrite(save_path, frameDiff);cout << "done!" << endl;
}int main() {// 读取图像Mat lastFrame1 = imread("data/Test_images/images/phantom05_0600.jpg");Mat lastFrame3 = imread("data/Test_images/images/phantom05_0602.jpg");Mat currentFrame = imread("data/Test_images/images/phantom05_0604.jpg");/*Mat lastFrame1 = imread("250514_430.bmp");Mat lastFrame3 = imread("250514_463.bmp");Mat currentFrame = imread("250514_490.bmp");*/// 检查图像是否成功加载if (lastFrame1.empty() || lastFrame3.empty() || currentFrame.empty()) {cout << "无法加载图像!" << endl;return -1;}// 执行帧间差分FD_mask(lastFrame1, lastFrame3, currentFrame);return 0;
}