最近做一个项目,用到广角镜头。畸变较大,所以就研究了一下畸变修正算法。主要是读了读 Carsten Steger 等所著 Machine Vision Algorithms and Applications 一书 3.9 节。把里面的算法实现了一下。这本书里写的修正方法应该算是最简单的了。只有一个参数 k。k > 0 修正桶形畸变,k < 0 修正枕形畸变。
下面把代码贴上来,里面插值算法只实现了最临近插值和双线性插值。一般来收双线性插值出来的效果就足够好了,所以就没实现更复杂的插值算法。
#ifndef RDCORRECTION_H
#define RDCORRECTION_H
#include <QImage>
/**
* @brief The RDCorrection class 镜头的畸变修正
* 只考虑了镜头的径向畸变,也就是桶形畸变和枕形畸变。
* 使用的模型比较简单,只有一个参数 k > 0 修正桶形畸变,k < 0 修正枕形畸变。
* 修正后的图像会有一定的放缩。四周也有可能出现黑边或者切掉一些。
* 这里假设径向畸变可以表示为:
* \bar u = \frac{2}{1+\sqrt{1- 4k (u^2+v^2)}} u
* \bar v = \frac{2}{1+\sqrt{1- 4k (u^2+v^2)}} v
* 那么修正公式为:
* u = \frac{1}{1 + k (\bar u^2 + \bar v^2)} \bar u
* v = \frac{1}{1 + k (\bar u^2 + \bar v^2)} \bar v
* 这个公式参考了 Carsten Steger 等所著 Machine Vision Algorithms and Applications 一书 3.9 节。
*/
class RDCorrection
{
public:
RDCorrection();
/**
* @brief map_nearest 镜头畸变修正,输出图像与输入图像同样大小。采用最临近插值
* @param in 输入的图像,灰度图或者 RGB32 类型
* @param k 修正系数,k > 0 修正桶形畸变,k < 0 修正枕形畸变。
* @return 修正后的图像
*/
QImage map_nearest(const QImage &in, double k);
/**
* @brief map_nearest 镜头畸变修正,输出图像与输入图像同样大小。采用双线性插值
* @param in 输入的图像,灰度图或者 RGB32 类型
* @param k 修正系数,k > 0 修正桶形畸变,k < 0 修正枕形畸变。
* @return 修正后的图像
*/
QImage map_bilinear(const QImage &in, double k);
/**
* @brief map_nearest 镜头畸变修正,采用最临近插值。输出的图像可以与输入图像不相同的尺寸。
* @param in 输入图像。
* @param k 修正系数,k > 0 修正桶形畸变,k < 0 修正枕形畸变。
* @param out 输出图像。传入这个参数时要求 format 与输入图像相同。
* @return true 表示成功。
*/
bool map_nearest(const QImage &in, double k, QImage &out);
/**
* @brief map_nearest 镜头畸变修正,采用双线性插值。输出的图像可以与输入图像不相同的尺寸。
* @param in 输入图像。
* @param k 修正系数,k > 0 修正桶形畸变,k < 0 修正枕形畸变。
* @param out 输出图像。传入这个参数时要求 format 与输入图像相同。
* @return true 表示成功。
*/
bool map_bilinear(const QImage &in, double k, QImage &out);
void setBackground(QRgb bg);
private:
bool map_rgb32_nearest(const QImage &in, double k, QImage &out);
bool map_rgb32_bilinear(const QImage &in, double k, QImage &out);
bool map_gray_nearest(const QImage &in, double k, QImage &out);
bool map_gray_bilinear(const QImage &in, double k, QImage &out);
QImage map_gray_nearest(const QImage &in, double k);
QImage map_gray_bilinear(const QImage &in, double k);
QImage map_rgb32_nearest(const QImage &in, double k);
QImage map_rgb32_bilinear(const QImage &in, double k);
inline uchar bilinear(uchar in00, uchar in01, uchar in10, uchar in11, double x, double y)
{
return in00 + (in01 - in00) * x + (in10 - in00) * y + (in11 - in01 - in10 + in00) * x * y;
}
inline QRgb bilinear(QRgb in00, QRgb in01, QRgb in10, QRgb in11, double x, double y)
{
// f= b1 + b2 * x + b3 * y + b4 * x * y
// b1 = in(0,0)
// b2 = in(1,0) - in(0,0)
// b3 = in(0,1) - in(0,0)
// b4 = in(1,1) - in(1,0) - in(0,1) + in(0,0)
int c00 = qRed(in00);
int c01 = qRed(in01);
int c10 = qRed(in10);
int c11 = qRed(in11);
int red = c00 + (c01 - c00) * x + (c10 - c00) * y + (c11 - c01 - c10 + c00) * x * y;
// int red = (1 - x) * (1 - y) * c00 + (1 - x) * y * c01 + (1 - y) * x * c10 + x * y * c11;
c00 = qGreen(in00);
c01 = qGreen(in01);
c10 = qGreen(in10);
c11 = qGreen(in11);
int green = c00 + (c01 - c00) * x + (c10 - c00) * y + (c11 - c01 - c10 + c00) * x * y;//
//int green = (1 - x) * (1 - y) * c00 + (1 - x) * y * c01 + (1 - y) * x * c10 + x * y * c11;
c00 = qBlue(in00);
c01 = qBlue(in01);
c10 = qBlue(in10);
c11 = qBlue(in11);
int blue = c00 + (c01 - c00) * x + (c10 - c00) * y + (c11 - c01 - c10 + c00) * x * y;
//int blue = (1 - x) * (1 - y) * c00 + (1 - x) * y * c01 + (1 - y) * x * c10 + x * y * c11;
return qRgb(red, green, blue);
}
QRgb m_bg_rgb32;
uchar m_bg_gray;
};
#endif // RDCORRECTION_H
#include "rdcorrection.h"
RDCorrection::RDCorrection()
{
m_bg_rgb32 = qRgb(0, 0, 0);
m_bg_gray = 0;
}
QImage RDCorrection::map_bilinear(const QImage &in, double k)
{
if(in.isNull())
{
return QImage();
}
switch (in.format())
{
case QImage::Format_RGB32:
case QImage::Format_ARGB32:
case QImage::Format_ARGB32_Premultiplied:
return map_rgb32_bilinear(in, k);
break;
case QImage::Format_Grayscale8:
case QImage::Format_Indexed8:
return map_gray_bilinear(in, k);
break;
default:
break;
}
return QImage();
}
QImage RDCorrection::map_nearest(const QImage &in, double k)
{
if(in.isNull())
{
return QImage();
}
switch (in.format())
{
case QImage::Format_RGB32:
case QImage::Format_ARGB32:
case QImage::Format_ARGB32_Premultiplied:
return map_rgb32_nearest(in, k);
break;
case QImage::Format_Grayscale8:
case QImage::Format_Indexed8:
return map_gray_nearest(in, k);
break;
default:
break;
}
return QImage();
}
bool RDCorrection::map_nearest(const QImage &in, double k, QImage &out)
{
if(in.isNull() || out.isNull())
{
return false;
}
switch (in.format())
{
case QImage::Format_RGB32:
case QImage::Format_ARGB32:
case QImage::Format_ARGB32_Premultiplied:
return map_rgb32_nearest(in, k, out);
break;
case QImage::Format_Grayscale8:
case QImage::Format_Indexed8:
return map_gray_nearest(in, k, out);
break;
default:
break;
}
return false;
}
bool RDCorrection::map_bilinear(const QImage &in, double k, QImage &out)
{
if(in.isNull() || out.isNull())
{
return false;
}
switch (in.format())
{
case QImage::Format_RGB32:
case QImage::Format_ARGB32:
case QImage::Format_ARGB32_Premultiplied:
return map_rgb32_bilinear(in, k, out);
break;
case QImage::Format_Grayscale8:
case QImage::Format_Indexed8:
return map_gray_bilinear(in, k, out);
break;
default:
break;
}
return false;
}
bool RDCorrection::map_rgb32_bilinear(const QImage &in, double k, QImage &out)
{
if(out.isNull() || in.isNull() || out.format() != in.format())
{
return false;
}
int cols_out = out.width();
int rows_out = out.height();
int cols_in = in.width();
int rows_in = in.height();
double xc_out = cols_out / 2.0;
double yc_out = rows_out / 2.0;
double xc_in = cols_in / 2.0;
double yc_in = rows_in / 2.0;
const uchar * source = in.bits();
int lineByte = in.bytesPerLine();
for(int y = 0; y < rows_out; y ++)
{
QRgb * dest = (QRgb * )out.scanLine(y);
for(int x = 0; x < cols_out; x ++)
{
double xx = (x - xc_out) / 1000.0;
double yy = (y - yc_out) / 1000.0;
double scale = 1 / (1 + k * (xx * xx + yy * yy));
xx = xx * scale * 1000.0 + xc_in;
yy = yy * scale * 1000.0 + yc_in;
int ix = xx;
int iy = yy;
if(iy < 0 || iy >= rows_in || ix < 0 || ix >= cols_in)
{
dest[x] = m_bg_rgb32;
}
else
{
QRgb * line0 = (QRgb *) (source + iy * lineByte);
QRgb * line1 = (QRgb *) (source + iy * lineByte + lineByte);
dest[x] = bilinear(line0[ix], line0[ix + 1], line1[ix], line1[ix + 1], xx - ix, yy - iy);
}
}
}
return true;
}
bool RDCorrection::map_rgb32_nearest(const QImage &in, double k, QImage &out)
{
if(out.isNull() || in.isNull() || out.format() != in.format())
{
return false;
}
int cols_out = out.width();
int rows_out = out.height();
int cols_in = in.width();
int rows_in = in.height();
double xc_out = cols_out / 2.0;
double yc_out = rows_out / 2.0;
double xc_in = cols_in / 2.0;
double yc_in = rows_in / 2.0;
const uchar * source = in.bits();
int lineByte = in.bytesPerLine();
for(int y = 0; y < rows_out; y ++)
{
QRgb * dest = (QRgb * )out.scanLine(y);
for(int x = 0; x < cols_out; x ++)
{
double xx = (x - xc_out) / 1000.0;
double yy = (y - yc_out) / 1000.0;
double scale = 1 / (1 + k * (xx * xx + yy * yy));
xx = xx * scale * 1000.0 + xc_in;
yy = yy * scale * 1000.0 + yc_in;
int ix = xx;
int iy = yy;
if(iy < 0 || iy >= rows_in || ix < 0 || ix >= cols_in)
{
dest[x] = m_bg_rgb32;
}
else
{
QRgb * line = (QRgb *) (source + iy * lineByte);
dest[x] = line[ix];
}
}
}
return true;
}
bool RDCorrection::map_gray_nearest(const QImage &in, double k, QImage &out)
{
if(out.isNull() || in.isNull() || out.format() != in.format())
{
return false;
}
int cols_out = out.width();
int rows_out = out.height();
int cols_in = in.width();
int rows_in = in.height();
double xc_out = cols_out / 2.0;
double yc_out = rows_out / 2.0;
double xc_in = cols_in / 2.0;
double yc_in = rows_in / 2.0;
const uchar * source = in.bits();
int lineByte = in.bytesPerLine();
for(int y = 0; y < rows_out; y ++)
{
uchar * dest = (uchar * )out.scanLine(y);
for(int x = 0; x < cols_out; x ++)
{
double xx = (x - xc_out) / 1000.0;
double yy = (y - yc_out) / 1000.0;
double scale = 1 / (1 + k * (xx * xx + yy * yy));
xx = xx * scale * 1000.0 + xc_in;
yy = yy * scale * 1000.0 + yc_in;
int ix = xx;
int iy = yy;
if(iy < 0 || iy >= rows_in || ix < 0 || ix >= cols_in)
{
dest[x] = m_bg_gray;
}
else
{
uchar * line0 = (uchar *) (source + iy * lineByte);
uchar * line1 = (uchar *) (source + iy * lineByte + lineByte);
dest[x] = bilinear(line0[ix], line0[ix + 1], line1[ix], line1[ix + 1], xx - ix, yy - iy);
}
}
}
return true;
}
bool RDCorrection::map_gray_bilinear(const QImage &in, double k, QImage &out)
{
if(out.isNull() || in.isNull() || out.format() != in.format())
{
return false;
}
int cols_out = out.width();
int rows_out = out.height();
int cols_in = in.width();
int rows_in = in.height();
double xc_out = cols_out / 2.0;
double yc_out = rows_out / 2.0;
double xc_in = cols_in / 2.0;
double yc_in = rows_in / 2.0;
const uchar * source = in.bits();
int lineByte = in.bytesPerLine();
for(int y = 0; y < rows_out; y ++)
{
uchar * dest = (uchar * )out.scanLine(y);
for(int x = 0; x < cols_out; x ++)
{
double xx = (x - xc_out) / 1000.0;
double yy = (y - yc_out) / 1000.0;
double scale = 1 / (1 + k * (xx * xx + yy * yy));
xx = xx * scale * 1000.0 + xc_in;
yy = yy * scale * 1000.0 + yc_in;
int ix = xx;
int iy = yy;
if(iy < 0 || iy >= rows_in || ix < 0 || ix >= cols_in)
{
dest[x] = m_bg_gray;
}
else
{
uchar * line = (uchar *) (source + iy * lineByte);
dest[x] = line[ix];
}
}
}
return true;
}
QImage RDCorrection::map_gray_nearest(const QImage &in, double k)
{
QImage out(in.width(), in.height(), in.format());
map_rgb32_nearest(in, k, out);
return out;
}
QImage RDCorrection::map_gray_bilinear(const QImage &in, double k)
{
QImage out(in.width(), in.height(), in.format());
map_gray_bilinear(in, k, out);
return out;
}
QImage RDCorrection::map_rgb32_nearest(const QImage &in, double k)
{
QImage out(in.width(), in.height(), in.format());
map_gray_nearest(in, k, out);
return out;
}
QImage RDCorrection::map_rgb32_bilinear(const QImage &in, double k)
{
QImage out(in.width(), in.height(), in.format());
map_rgb32_bilinear(in, k, out);
return out;
}
下面是测试用例。