Matlabでレンズ歪み補正・レクティフィケーション（Lens distortion correction, Rectification）を効率的に処理する方法

先日，Matlabで2次元データ配列である画像のアフィン変換，射影変換にコーディングテクニックをご紹介しました．

 http://miraimage-lab.hatenablog.com/entry/2014/07/04/165947

アフィン変換，射影変換であれば，Image Processing Toolboxのimwarpを利用可能なようです．私は使ったことがありませんが．

http://www.mathworks.co.jp/jp/help/images/performing-general-2-d-spatial-transformations.html

meshgridを使えば，アフィン変換や射影変換だけでなく，レンズ歪みの補正などの非線形な画像変換も簡単に実装できます！！ステレオカメラの画像の平行化処理などにも使えます．

for文を使ってコーディングすると物凄く重たくなるので，Matlab初心者の方は参考にしてみてください．

-------------------------------------------------------------------------------------------------------------------------

  % 画像を読み込み
  im=imread('test.jpg');
  if size(im, 3) > 1
      im=rgb2gray(im);
  end
 
  % 画像サイズ
  width = size(im,2); height = size(im,1);  
 
  % 内部パラメータ: Intrinsic matrix (補正前)
  P_bf = [ 880 0   width/2-10;
           0   880 height/2+10;
           0   0   1];
  % 内部パラメータ: Intrinsic matrix (補正後)
  P_af = [ 880 0   width/2;
           0   880 height/2;
           0   0   1];
  % レンズ歪みパラメータ: Distortion parameter
  k1=-0.3; k2=0.075; k3=0; p1=-0.0004; p2=0.0002;
             
  % 処理領域の各画素の座標をmeshgridで一括生成
  [xim,yim] = meshgrid(1:width, 1:height);
  points = [xim(:)';yim(:)';ones(1, size(xim(:),1))];

  % 変換後の座標を一気に計算 ※※ forループを使わない ※※
  points1 = inv(P_af)*points;
  r2 = points1(1,:).^2 + points1(2,:).^2;
  k_radial =  1 + k1 * r2 + k2 * r2.^2 + k3 * r2.^3;
  delta_x = [2*p1*points1(1,:).*points1(2,:) + p2*(r2 + 2*points1(1,:).^2);
  p1 * (r2 + 2*points1(2,:).^2)+2*p2*points1(1,:).*points1(2,:)];
  points2(1,:) = points1(1,:).*k_radial + delta_x(1,:);
  points2(2,:) = points1(2,:).*k_radial + delta_x(2,:);
  points2(3,:) = ones(1, size(xim(:),1));
  points3 = P_bf*points2;
 
  % 画像と同じように2次元配列に変換
  x_warped = reshape(points3(1,:),[height width]);
  y_warped = reshape(points3(2,:),[height width]);

  % 元の画像の各画素の座標をmeshgridで一括生成
  [xim,yim] = meshgrid(1:size(im,2),1:size(im,1));

  % 各画素の輝度値を取得．下記の例ではバイリニア('linear')
  im_warped = interp2(xim,yim,double(im),x_warped,y_warped,'linear');

  % 結果を表示
  imshow(uint8(im_warped));

-------------------------------------------------------------------------------------------------------------------------

※内部パラメータ、レンズ歪みパラメータは、カメラキャリブレーション（Camera Calibration）で事前に求めてください．

※ステレオカメラ（Stereo Camera）、ステレオビジョン（Stereo Vision）

※平行化処理: レクティフィケーション（Rectification）