Radial Distortion Correction - 삭제예정

[한창호]

1. Radial distortion 이란

렌즈 왜곡(radial distortion)은 다음과 같은 형태의 종류가 있다. (pincushion, barrel)

Radial distortion(방사왜곡)

Tangential distortion(접선왜곡)

합쳐서

         K 값은 barrel 에서는 음수(-)로 되며, pincusion 에서는 양수(+)가 된다.

         r값은 normalize된 값으로 1이하의 값이어야 한다.

(출처: https://en.wikipedia.org/wiki/Distortion_%28optics%29#Software_correction)

2. Radial distortion을 만들어 보자

1) 먼저 원 영상 Bmp파일을 읽어 화면에 띄우는 프로그램을 만든다.

샘플 bmp 읽고쓰기 소스 : SampleBmp

축소/확대 (양선형보간법) 영상처리 소스 ch08.zip  (Visual C++ 6.0, OpenCV 사용없이 순수 영상처리 코드)

2) 왜곡 함수를 만들어 영상을 왜곡시켜보자.

정상적인 입력 영상에 왜곡 기법을 적용하고 양선형 보간법을 이용하여 보정된 영상을 출력한다. k1, k2, k3에 따라 어떻게 변하는지 관찰한다. 여기서 r값은 정규화된 값이다.

 void UndistortImg(float k1,float k2,float k3,  BYTE* inImg, BYTE* outImg, int width, int height)

{

  int xc=width/2, yc=height/2; // 중심좌표

  for (int y = 0; y<height; y++)

    for (int x = 0; x<width; x++)

    {

      // ------------------------------------------------------------

      // 왜곡 함수를 이용해서 영상을 왜곡시킨다.

      // k1,k2,k3 적용하여 radial distortion을 구현한다.

      float r = sqrt( (float)(y - yc)*(y - yc) / (2 * yc*yc) +

                     (float)(x - xc)*(x - xc) / (2 * xc*xc) );

      float dist = (1 + k1*r*r + k2*r*r*r*r + k3*r*r*r*r*r*r);

      float f_yu = (y-yc) * dist + yc; // 보정된 좌표 yu

      float f_xu = (x-xc) * dist + xc; // 보정된 좌표 xu

      // ------------------------------------------------------------

      // 이동한 ,  값을 구한다.

      int i_yu = floor(f_yu); //예: floor(2.8)=2.0

      int i_xu = floor(f_xu);

      float dy = f_yu – i_yu; // beta value

      float dx = f_xu – i_xu; // alpha value

      //원이미지의 범위를 벗어나는 경우 0값 할당 

      if (i_yu<0 || i_yu>height-1 || i_xu<0 || i_xu>width-1) {

        outImg[ y*width + x ] = 0;

      }

      //원 이미지의 범위 내에 존재 양선형 보간 수행 

      else {

        I1 = (float)inImg[i_yu*width + i_xu]; // (x, y)

        I2 = (float)inImg[i_yu*width + (i_xu+1)]; // (x, y+1)

        I3 = (float)inImg[(i_yu+1)*width + i_xu]; // (x+1, y)

        I4 = (float)inImg[(i_yu+1)*width + (i_xu+1)]; // (x+1, y+1)

        //양선형 보간을 통한 새로운 밝기값 계산

        newValue = (BYTE)(I1*(1-dx)*(1-dy) + I2*dx*(1-dy) + I3*(1-dx)*dy + I4*dx*dy);

        outImg[ y*width + x ] = newValue;

      }

    }

}

다음 그림은 pincushion 모양의 왜곡 영상을 구현한 결과이다.

다음 그림은 barrel 모양의 왜곡을 구현한 결과이다.

3. Radial distortion 구하기

[문제] 다음과 같은 왜곡영상이 주어졌을 때 보정을 해보자. (사진을 다운받아 실험해보자)

[풀이]

1. 왜곡 영상에서 점들의 중심좌표를 얻는다. (왼쪽그림의 왜곡 영상으로부터 구한다. (x,y))

2. 보정된 영상의 점들의 중심좌표를 얻는다. (임의로 계산할 수 있다. (x^,y^))

3. 방사왜곡 algorithm으로 구한다.

Radial distortion을 구하는 방법에 대해 요약해 보면 다음과 같다.

식(1)에서 (x^,y^)은 왜곡이 제거된 좌표이고, (x,y)는 왜곡 영상 좌표이다. (xc,yc)는 중심 좌표이다. 단 여기서 사용하는 좌표는 normalized된 좌표를 사용한다.

식(3)의 L(r)은 왜곡을 제거하는 undistortion model을 정의하는 함수이다.

모든 점을 만족하는 L(r) 함수를 구하기 위해 k 계수들을 구하면 된다. 또한 k 계수는 최소자승법과 역행렬을 이용하여 구할 수 있다.

위 수식은 아래 curve fitting 소스를 이용하여 구할 수 있다.

구해진 L(r)을 이용하여 왜곡이 제거된 영상을 구한다.

4. 참고

distortion 소스코드: Algebraic Lens Distortion Model Estimation

pdf: Algebraic Lens Distortion Model Estimation.pdf

src: algebraic_lens_distortion_model_estimation_basic.tar.gz

curve fitting 소스: curve fitting(least square) source code

영상 왜곡보정 예제 Matlab 코드   http://darkpgmr.tistory.com/31

function undistort_image(Id, Iu, fx, fy, cx ,cy, skew_c, k1, k2, k3, p1, p2)

    w: image width of Iu

    h: image height of Iu

    for y = 0:h-1,

        for x = 0:w-1,

            y_nu = (y-cy)/fy;

            x_nu = (x-cx)/fx – skew_c*y_nu;

            ru2 = x_nu*x_nu + y_nu*y_nu; // ru2 = ru*ru

            radial_d = 1 + k1*ru2 + k2*ru2*ru2 + k3*ru2*ru2*ru2;

            x_nd = radial_d*x_nu + 2*p1*x_nu*y_nu + p2*(ru2 + 2*x_nu*x_nu);

            y_nd = radial_d*y_nu + p1*(ru2 + 2*y_nu*y_nu) + 2*p2*x_nu*y_nu;

            x_pd = fx*(x_nd + skew_c*y_nd) + cx;

            y_pd = fy*y_nd + cy;

            Iu(x, y) = Id(x_pd, y_pd);

        end;

    end;    

return Iu;

영상좌표 왜곡보정 (예제 코드)

function undistort(point p_d)

      p_d’ = normalize(p_d);

      p_u’ = p_d’;

      while (1)

            err = distort_normal(p_u’) – p_d’;

            p_u’ = p_u’ – err;

            if(err<err_threshold) break;

      end;

      p_u = denormalize(p_u’);

return p_u;

function normalize(x, y)

      y_n = (y-cy)/fy;

      x_n = (x-cx)/fx – skew_c*y_n;

return (x_n, y_n);

function denormalize(x, y)

      x_p = fx*(x + skew_c*y) + cx;

      y_p = fy*y + cy;

return (x_p, y_p);

function distort_normal(x, y)

      r2 = x*x + y*y;

      radial_d = 1 + k1*r2 + k2*r2*r2 + k3*r2*r2*r2;

      x_d = radial_d*x + 2*p1*x*y + p2*(r2 + 2*x*x);

      y_d = radial_d*y + p1*(r2 + 2*y*y) + 2*p2*x*y;

return (x_d, y_d);

참고문서:

카메라 보정 소개  camcalib.ppt

lecture9.pdf

Calibration_Ahmed Elgammal.pdf

rob2-08-camera-calibration.pdf