一、前言
最近想实现AVM拼接,看了不少博客和论文,不过比较愚钝,一直没能很好理解原理,尤其是怎么在实现时把下文式1与式2中Z1和Z2消除的,所以严谨的推导了一下对应的公式,如有不对,水平有限,烦请指出~
IPM变换(逆透视变换),顾名思义,即将正常透视效应消除的变换,变换结果为鸟瞰图BEV(俯视图);AVM环视拼接一般是将车身前后左右的四个鱼眼相机拼接成环视图(也是俯视),其使用的单应矩阵将四个相机转到同一个俯视坐标系。
投影变换的通俗理解就是:假设同一个相机分别在A、B两个不同位置,以不同的位姿拍摄同一个平面(重点是拍摄平面,例如桌面、墙面、地平面),生成了两张图象,这两张图象之间的关系就叫做投影变换。
二、公式推导
公式太难打了,全用word存了,在这直接截图
这也就解释了:为什么有的AVM算法是需要标定相机的内外参,而有的只提供单应矩阵H。
三、示例代码
#include <Eigen/Core>
#include <Eigen/Dense>
#include <iostream>
#include <opencv2/core/core.hpp>
#include <opencv2/core/eigen.hpp>
#include <opencv2/opencv.hpp>
#include <vector>void DistortEigenPoints(const Eigen::Vector3d &undis_pt, Eigen::Vector3d &dis_pt, double k1, double k2, double p1,double p2) {double x2 = undis_pt[0] * undis_pt[0];double y2 = undis_pt[1] * undis_pt[1];double r2 = x2 + y2;if (r2 == 0) return;double r4 = r2 * r2;double r6 = r2 * r4;double xy = undis_pt[0] * undis_pt[1];dis_pt[0] = undis_pt[0] * (1 + k1 * r2 + k2 * r4) + 2 * p1 * xy + p2 * (r2 + 2 * x2);dis_pt[1] = undis_pt[1] * (1 + k1 * r2 + k2 * r4) + p1 * (r2 + 2 * y2) + 2 * p2 * xy;return;
}float interp(const cv::Mat &img, float x, float y) {int ix = x;int iy = y;float dx = x - ix;float dy = y - iy;float ddx = 1.0f - dx;float ddy = 1.0f - dy;return ddx * ddy * img.data[iy * img.cols + ix] + ddx * dy * img.data[(iy + 1) * img.cols + ix] +dx * ddy * img.data[iy * img.cols + ix + 1] + dx * dy * img.data[(iy + 1) * img.cols + ix + 1];
}int main() {std::string img_path = "Mynt/00001772.jpg";Eigen::Matrix3d K_c;K_c << 1037.536376953125, 0, 600.182861328125, 0, 1038.532958984375, 358.40493774414062500, 0, 0, 1;double k1 = 0.03784750;double k2 = -0.051872;double p1 = -0.000938;double p2 = 0.000157;Eigen::Matrix3d R_gc;R_gc << 0.9962626012012678, -0.08634899803974432, -0.00216332734845117, -0.005878375212093168, -0.04279258802138065,-0.9990666840182884, 0.08617583276389137, 0.9953454902434454, -0.0431402468639808;Eigen::Vector3d t_gc;t_gc << 0.283289952815021, 1.136073800639606, -2.129837994618606;// 鸟瞰图设置为长宽20m,分辨率0.2mdouble W_m = 20;double H_m = 20;double dx = 0.02;double dy = 0.02;Eigen::Matrix3d K_g;K_g << 1.0 / dx, 0, W_m / (2 * dx), 0, -1.0 / dy, H_m / (2 * dy), 0, 0, 1;Eigen::Matrix3d J_gc;J_gc.block<3, 2>(0, 0) = R_gc.block<3, 2>(0, 0);J_gc.block<3, 1>(0, 2) = t_gc;//一定要带上图像类型,比如IMREAD_GRAYSCALE,默认的是IMREAD_COLORcv::Mat img = cv::imread(img_path, cv::IMREAD_GRAYSCALE);//一、未去畸变鸟瞰图Eigen::Matrix3d H = K_c * J_gc * K_g.inverse();cv::Mat trans_mat;cv::eigen2cv(H, trans_mat);cv::Mat bev_img(H_m / dy, W_m / dx, CV_8UC1);// warpPerspective是透视变换,所有trans_mat需要取逆cv::warpPerspective(img, bev_img, trans_mat.inv(), bev_img.size());//二、去畸变鸟瞰图Eigen::Matrix3d H_gc = J_gc * K_g.inverse();cv::Mat bev_img2(H_m / dy, W_m / dx, CV_8UC1);for (int row = 0; row < bev_img2.rows; row++) {for (int col = 0; col < bev_img2.cols; col++) {Eigen::Vector3d p_g(col, row, 1);Eigen::Vector3d P_c = H_gc * p_g;// 去掉在相机后面的点if (P_c[2] < 0) continue;P_c /= P_c[2];Eigen::Vector3d P_c_dis;DistortEigenPoints(P_c, P_c_dis, k1, k2, p1, p2);Eigen::Vector3d p_c = K_c * P_c_dis;// 在原图范围内的才能映射if (p_c[0] >= 0 && p_c[1] >= 0 && p_c[0] < img.cols - 1 && p_c[1] < img.rows - 1) {bev_img2.at<uchar>(row, col) = interp(img, p_c[0], p_c[1]);}}}cv::imshow("img", img);cv::imshow("bev_img", bev_img);cv::imshow("bev_img2", bev_img2);cv::waitKey(0);return 0;
}
效果展示
ps.
// 1. 读取图像cv::Mat img = cv::imread(img_path, cv::IMREAD_GRAYSCALE);
使用时一定要记得加对应图像类型 ,比如灰度的为cv::IMREAD_GRAYSCALE,否则读取的会为三通道的IMREAD_COLOR,导致整个转换出问题,查了个半死。。。
参考文献:
自动驾驶AVM环视算法--3D碗型投影模式的算法原理和代码实现_3d avm算法原理-CSDN博客
Fisheye Calibration Basics- MATLAB & Simulink- MathWorks 中国
AVM 环视拼接方法介绍
IPM 鸟瞰图公式转换与推导 - 古月居
https://blog.51cto.com/u_16099267/10196291
逆透视变换详解 及 代码实现(二)_uvlimitsp-CSDN博客
Online Camera Pose Optimization for the Surround-view System
IPM原理_ipm转换-CSDN博客
https://zhuanlan.zhihu.com/p/636990989