2010-03-19 111 views
48

BOUNTY STATUS UPDATE:校正魚眼失真編程

I discovered how to map a linear lens,從destination座標source座標。

如何計算從魚眼到直線的中心徑向距離?

  • 1)。我真的很難扭轉它,並將源座標映射到目標座標。在我發佈的轉換函數風格的代碼中,反過來是什麼?

  • 2)。我也看到我的失真在一些鏡頭上不完美 - 大概是那些不是嚴格線性的。這些鏡頭的等效源和目標座標是多少?不僅僅是數學公式,請再次,更多的代碼...


Question as originally stated:

我有一些點,描述與魚眼鏡頭拍攝的照片位置。

我想將這些點轉換爲直線座標。我想要保持圖像的真實。

我發現this description如何產生魚眼效果,但沒有如何扭轉它。

還有一個blog post,它描述瞭如何使用工具來完成它;這些圖片是從該:

(1)SOURCEOriginal photo link

輸入:與魚眼失真原圖修正。

(2)DESTINATIONOriginal photo link

輸出:校正圖像(技術上也與透視校正,但是這是一個單獨的步驟)。

你如何計算從魚眼中心到直線的徑向距離?

我的函數存根看起來是這樣的:

Point correct_fisheye(const Point& p,const Size& img) { 
    // to polar 
    const Point centre = {img.width/2,img.height/2}; 
    const Point rel = {p.x-centre.x,p.y-centre.y}; 
    const double theta = atan2(rel.y,rel.x); 
    double R = sqrt((rel.x*rel.x)+(rel.y*rel.y)); 
    // fisheye undistortion in here please 
    //... change R ... 
    // back to rectangular 
    const Point ret = Point(centre.x+R*cos(theta),centre.y+R*sin(theta)); 
    fprintf(stderr,"(%d,%d) in (%d,%d) = %f,%f = (%d,%d)\n",p.x,p.y,img.width,img.height,theta,R,ret.x,ret.y); 
    return ret; 
} 

或者,我可以以某種方式將圖像從魚眼鏡頭轉換找到點之前直線,但我完全被OpenCV documentation迷惑。在OpenCV中有沒有一種簡單的方法來做到這一點,並且它的表現足夠好,可以用於實時視頻Feed嗎?

+0

我不太讓你在做什麼的C#代碼尋找。魚眼從一個球體映射到圖像平面。反向映射將從圖片返回到一個球體的權利?你在找什麼直線座標? – mtrw 2010-03-19 21:56:00

+0

@mtrw我的源圖像是魚眼失真,我想不破壞它 – Will 2010-03-20 06:18:37

+0

所以是http://photo.net/learn/fisheye/上的圖片你在找什麼? – mtrw 2010-03-20 20:46:32

回答

29

description you mention指出通過針孔照相機(一個不引入透鏡畸變)的突起由

R_u = f*tan(theta) 

並通過共同的魚眼鏡頭攝像機(即,失真)是投影模型通過

R_d = 2*f*sin(theta/2) 

你已經知道R_D和θ建模,如果你知道相機的焦距(用f表示)則修正圖像就等於在R_D和θ方面計算R_u。換句話說,

R_u = f*tan(2*asin(R_d/(2*f))) 

是你正在尋找的公式。估算焦距f可以通過校準相機或其他手段來解決,例如讓用戶提供關於圖像校正的好壞的反饋或使用來自原始場景的知識。

爲了解決使用OpenCV的相同問題,您必須獲取相機的固有參數和鏡頭失真係數。例如參見Learning OpenCV的第11章(不要忘記檢查correction)。然後你可以使用一個程序,這樣一個爲了扭轉鏡頭失真(與Python綁定OpenCV的書面):

#!/usr/bin/python 

# ./undistort 0_0000.jpg 1367.451167 1367.451167 0 0 -0.246065 0.193617 -0.002004 -0.002056 

import sys 
import cv 

def main(argv): 
    if len(argv) < 10: 
    print 'Usage: %s input-file fx fy cx cy k1 k2 p1 p2 output-file' % argv[0] 
    sys.exit(-1) 

    src = argv[1] 
    fx, fy, cx, cy, k1, k2, p1, p2, output = argv[2:] 

    intrinsics = cv.CreateMat(3, 3, cv.CV_64FC1) 
    cv.Zero(intrinsics) 
    intrinsics[0, 0] = float(fx) 
    intrinsics[1, 1] = float(fy) 
    intrinsics[2, 2] = 1.0 
    intrinsics[0, 2] = float(cx) 
    intrinsics[1, 2] = float(cy) 

    dist_coeffs = cv.CreateMat(1, 4, cv.CV_64FC1) 
    cv.Zero(dist_coeffs) 
    dist_coeffs[0, 0] = float(k1) 
    dist_coeffs[0, 1] = float(k2) 
    dist_coeffs[0, 2] = float(p1) 
    dist_coeffs[0, 3] = float(p2) 

    src = cv.LoadImage(src) 
    dst = cv.CreateImage(cv.GetSize(src), src.depth, src.nChannels) 
    mapx = cv.CreateImage(cv.GetSize(src), cv.IPL_DEPTH_32F, 1) 
    mapy = cv.CreateImage(cv.GetSize(src), cv.IPL_DEPTH_32F, 1) 
    cv.InitUndistortMap(intrinsics, dist_coeffs, mapx, mapy) 
    cv.Remap(src, dst, mapx, mapy, cv.CV_INTER_LINEAR + cv.CV_WARP_FILL_OUTLIERS, cv.ScalarAll(0)) 
    # cv.Undistort2(src, dst, intrinsics, dist_coeffs) 

    cv.SaveImage(output, dst) 


if __name__ == '__main__': 
    main(sys.argv) 

還要注意的是OpenCV的使用非常不同的鏡頭失真模型的一個在網絡您鏈接到的頁面。

+0

這依賴於有權訪問相機。我不知道,我只是看一段視頻。 此外,它發生在我身上,攝像機是大規模生產的,並且不能單獨確定地有多少變化?原文中鏈接的工具不需要有人用棋盤站在鏡頭前面! – Will 2010-03-21 20:28:37

+2

只需調整變焦,攝像機參數就會在同一臺攝像機中有所不同。此外,您可以依靠自動校準技術而不是使用棋盤。 無論如何,我編輯了我的答案,以解決問題的第一部分,爲您提供您正在尋找的公式。 – jmbr 2010-03-21 21:58:16

+0

謝謝jmbr!世界正在開始變得有意義。但是我實際上得不到R_u = f * tan(2 * asin(R_d /(2 * f)))給我任何東西,但是NaN。我沒有三角學的知識,它現在正在阻止我:( – Will 2010-03-22 06:22:40

3

如果你認爲你的公式是準確的,你可以COMPUT一個精確的公式與三角函數,像這樣:

Rin = 2 f sin(w/2) -> sin(w/2)= Rin/2f 
Rout= f tan(w)  -> tan(w)= Rout/f 

(Rin/2f)^2 = [sin(w/2)]^2 = (1 - cos(w))/2 -> cos(w) = 1 - 2(Rin/2f)^2 
(Rout/f)^2 = [tan(w)]^2 = 1/[cos(w)]^2 - 1 

-> (Rout/f)^2 = 1/(1-2[Rin/2f]^2)^2 - 1 

然而,隨着@jmbr說,實際的攝像機失真將取決於鏡頭和變焦上。而不是依賴於一個固定的公式,你可能想嘗試一個多項式展開:

Rout = Rin*(1 + A*Rin^2 + B*Rin^4 + ...) 

通過調整第一個A,然後高階係數,可以計算出任何合理的本地功能(擴展的形式利用問題的對稱性)。特別是,應該有可能計算初始係數以接近上面的理論函數。

此外,爲了獲得好的結果,您需要使用插值濾鏡生成校正圖像。只要失真不是太大,您可以使用您用來線性重新縮放圖像的那種過濾器,而不會出現太多問題。

編輯:按照您的要求,爲上述公式等同比例因子:

(Rout/f)^2 = 1/(1-2[Rin/2f]^2)^2 - 1 
-> Rout/f = [Rin/f] * sqrt(1-[Rin/f]^2/4)/(1-[Rin/f]^2/2) 

如果打印上面的公式旁邊譚(琳/ F),你可以看到,他們都非常相似形狀。基本上,在sin(w)與w大不相同之前,切線的變形變得嚴重。

逆式應該是這樣的:

Rin/f = [Rout/f]/sqrt(sqrt(([Rout/f]^2+1) * (sqrt([Rout/f]^2+1) + 1)/2) 
7

(原始海報,提供一種替代)

下面的函數映射目的地(直線的)座標來源(魚眼失真)座標。 通過試錯(我會很感激扭轉它的幫助)

我得到了這一點:我沒有從根本上把握爲什麼這個代碼工作,解釋並提高了精度讚賞

def dist(x,y): 
    return sqrt(x*x+y*y) 

def correct_fisheye(src_size,dest_size,dx,dy,factor): 
    """ returns a tuple of source coordinates (sx,sy) 
     (note: values can be out of range)""" 
    # convert dx,dy to relative coordinates 
    rx, ry = dx-(dest_size[0]/2), dy-(dest_size[1]/2) 
    # calc theta 
    r = dist(rx,ry)/(dist(src_size[0],src_size[1])/factor) 
    if 0==r: 
     theta = 1.0 
    else: 
     theta = atan(r)/r 
    # back to absolute coordinates 
    sx, sy = (src_size[0]/2)+theta*rx, (src_size[1]/2)+theta*ry 
    # done 
    return (int(round(sx)),int(round(sy))) 

當與3.0的因子所使用的,它成功地undistorts用作示例的圖像(I製成在質量插補沒有嘗試):

死鏈接

(這是從博客文章,比較:)

Using Panotools

+0

你的代碼工作的原因是你正在通過因子theta(現在是一個比例,而不是一個角度)縮放(rx,ry)。如果原始鏡頭有維基文章說)從視角到圖像偏移量的「線性映射」,我相信atan(r)/ r是正確的 – comingstorm 2010-03-24 16:07:21

+1

映射的反向應該是以tan(r')/ r' ,其中r'是從未失真圖像中心開始的半徑 – comingstorm 2010-03-24 16:27:21

+1

這兩個工作的原因是,如果你有矢量v'= v * k(| v |),並且你想| v'| = f(| v |),則取第一個方程的絕對值:| v'| = | v | * k(| v |)= f(| v |) - |)= F(| V |)/ | V | – comingstorm 2010-03-24 16:32:32

2

我盲目IMPL從here發表的公式,所以我不能保證它會做你需要的。

使用auto_zoom可獲取zoom參數的值。


def dist(x,y): 
    return sqrt(x*x+y*y) 

def fisheye_to_rectilinear(src_size,dest_size,sx,sy,crop_factor,zoom): 
    """ returns a tuple of dest coordinates (dx,dy) 
     (note: values can be out of range) 
crop_factor is ratio of sphere diameter to diagonal of the source image""" 
    # convert sx,sy to relative coordinates 
    rx, ry = sx-(src_size[0]/2), sy-(src_size[1]/2) 
    r = dist(rx,ry) 

    # focal distance = radius of the sphere 
    pi = 3.1415926535 
    f = dist(src_size[0],src_size[1])*factor/pi 

    # calc theta 1) linear mapping (older Nikon) 
    theta = r/f 

    # calc theta 2) nonlinear mapping 
    # theta = asin (r/(2 * f)) * 2 

    # calc new radius 
    nr = tan(theta) * zoom 

    # back to absolute coordinates 
    dx, dy = (dest_size[0]/2)+rx/r*nr, (dest_size[1]/2)+ry/r*nr 
    # done 
    return (int(round(dx)),int(round(dy))) 


def fisheye_auto_zoom(src_size,dest_size,crop_factor): 
    """ calculate zoom such that left edge of source image matches left edge of dest image """ 
    # Try to see what happens with zoom=1 
    dx, dy = fisheye_to_rectilinear(src_size, dest_size, 0, src_size[1]/2, crop_factor, 1) 

    # Calculate zoom so the result is what we wanted 
    obtained_r = dest_size[0]/2 - dx 
    required_r = dest_size[0]/2 
    zoom = required_r/obtained_r 
    return zoom 
3

我發現這個PDF文件,我已經證明了數學是正確的(除了線vd = *xd**fv+v0 which should say vd = **yd**+fv+v0)。

http://perception.inrialpes.fr/CAVA_Dataset/Site/files/Calibration_OpenCV.pdf

它不使用所有最新的合作efficients是OpenCV的可用的,但我相信它可以很容易適應。

double k1 = cameraIntrinsic.distortion[0]; 
double k2 = cameraIntrinsic.distortion[1]; 
double p1 = cameraIntrinsic.distortion[2]; 
double p2 = cameraIntrinsic.distortion[3]; 
double k3 = cameraIntrinsic.distortion[4]; 
double fu = cameraIntrinsic.focalLength[0]; 
double fv = cameraIntrinsic.focalLength[1]; 
double u0 = cameraIntrinsic.principalPoint[0]; 
double v0 = cameraIntrinsic.principalPoint[1]; 
double u, v; 


u = thisPoint->x; // the undistorted point 
v = thisPoint->y; 
double x = (u - u0)/fu; 
double y = (v - v0)/fv; 

double r2 = (x*x) + (y*y); 
double r4 = r2*r2; 

double cDist = 1 + (k1*r2) + (k2*r4); 
double xr = x*cDist; 
double yr = y*cDist; 

double a1 = 2*x*y; 
double a2 = r2 + (2*(x*x)); 
double a3 = r2 + (2*(y*y)); 

double dx = (a1*p1) + (a2*p2); 
double dy = (a3*p1) + (a1*p2); 

double xd = xr + dx; 
double yd = yr + dy; 

double ud = (xd*fu) + u0; 
double vd = (yd*fv) + v0; 

thisPoint->x = ud; // the distorted point 
thisPoint->y = vd; 
3

我拿了JMBR所做的,並基本上扭轉了它。他將畸變圖像的半徑(Rd,即與圖像中心的距離(以像素爲單位))找出Ru,即未失真圖像的半徑的公式。

你想走另一條路。對於未失真(處理後的圖像)中的每個像素,您想知道失真圖像中相應的像素是什麼。 換句話說,給(xu,yu) - >(xd,yd)。然後,將失真圖像中的每個像素與失真圖像中的相應像素進行替換。

從JMBR做的地方開始,我做了相反的事情,找到Rd作爲Ru的函數。我得到:

Rd = f * sqrt(2) * sqrt(1 - 1/sqrt(r^2 +1)) 

其中f是以像素爲單位的焦距(我會在後面解釋),和r = Ru/f

我的相機的焦距是2.5毫米。我的CCD上每個像素的大小是6平方米。因此f是2500/6 = 417個像素。這可以通過試驗和錯誤找到。

Finding Rd允許您使用極座標查找變形圖像中的對應像素。

從中心點每個像素的角度是相同的:

theta = arctan((yu-yc)/(xu-xc))其中XC,YC是中心點。

然後,

xd = Rd * cos(theta) + xc 
yd = Rd * sin(theta) + yc 

確保你知道哪個象限你英寸

下面是我用

public class Analyzer 
{ 
     private ArrayList mFisheyeCorrect; 
     private int mFELimit = 1500; 
     private double mScaleFESize = 0.9; 

     public Analyzer() 
     { 
      //A lookup table so we don't have to calculate Rdistorted over and over 
      //The values will be multiplied by focal length in pixels to 
      //get the Rdistorted 
      mFisheyeCorrect = new ArrayList(mFELimit); 
      //i corresponds to Rundist/focalLengthInPixels * 1000 (to get integers) 
      for (int i = 0; i < mFELimit; i++) 
      { 
       double result = Math.Sqrt(1 - 1/Math.Sqrt(1.0 + (double)i * i/1000000.0)) * 1.4142136; 
       mFisheyeCorrect.Add(result); 
      } 
     } 

     public Bitmap RemoveFisheye(ref Bitmap aImage, double aFocalLinPixels) 
     { 
      Bitmap correctedImage = new Bitmap(aImage.Width, aImage.Height); 
      //The center points of the image 
      double xc = aImage.Width/2.0; 
      double yc = aImage.Height/2.0; 
      Boolean xpos, ypos; 
      //Move through the pixels in the corrected image; 
      //set to corresponding pixels in distorted image 
      for (int i = 0; i < correctedImage.Width; i++) 
      { 
       for (int j = 0; j < correctedImage.Height; j++) 
       { 
        //which quadrant are we in? 
        xpos = i > xc; 
        ypos = j > yc; 
        //Find the distance from the center 
        double xdif = i-xc; 
        double ydif = j-yc; 
        //The distance squared 
        double Rusquare = xdif * xdif + ydif * ydif; 
        //the angle from the center 
        double theta = Math.Atan2(ydif, xdif); 
        //find index for lookup table 
        int index = (int)(Math.Sqrt(Rusquare)/aFocalLinPixels * 1000); 
        if (index >= mFELimit) index = mFELimit - 1; 
        //calculated Rdistorted 
        double Rd = aFocalLinPixels * (double)mFisheyeCorrect[index] 
             /mScaleFESize; 
        //calculate x and y distances 
        double xdelta = Math.Abs(Rd*Math.Cos(theta)); 
        double ydelta = Math.Abs(Rd * Math.Sin(theta)); 
        //convert to pixel coordinates 
        int xd = (int)(xc + (xpos ? xdelta : -xdelta)); 
        int yd = (int)(yc + (ypos ? ydelta : -ydelta)); 
        xd = Math.Max(0, Math.Min(xd, aImage.Width-1)); 
        yd = Math.Max(0, Math.Min(yd, aImage.Height-1)); 
        //set the corrected pixel value from the distorted image 
        correctedImage.SetPixel(i, j, aImage.GetPixel(xd, yd)); 
       } 
      } 
      return correctedImage; 
     } 
}