【C#】画像処理 ~ 幾何補正 / 2次元アフィン変換編 ~

■ 幾何補正(geometric correction)

 * 画像に対する歪みやズレを取り除く処理

■ サンプル

 * 以下で実装してみる
  + 補正方法 : 2次元アフィン変換
  + 補間方法 : 最近傍補間(ニアレストネイバー)

# 「2次元アフィン変換」「最近傍補間(ニアレストネイバー)」については、以下の関連記事を参照のこと
https://blogs.yahoo.co.jp/dk521123/38073336.html
https://blogs.yahoo.co.jp/dk521123/38080694.html
 * 以下のサイトを参考に実装してみる
http://mf-atelier.sakura.ne.jp/mf-atelier/modules/tips/program/algorithm/a10.html
using System;
using System.Collections.Generic;
using System.Drawing;
using System.Windows.Forms;

namespace SampleForm
{
  public partial class Form1 : Form
  {
    private const int DataNumberForAffineCoefficients = 6;

    private Bitmap originalBitmap = null;

    public Form1()
    {
      InitializeComponent();
    }
    private void Form1_Load(object sender, EventArgs e)
    {
      // 画像ファイルのImageオブジェクトを作成する
      this.originalBitmap = new Bitmap(@"20161215052204.gif");
      this.pictureBox1.Image = this.originalBitmap;
    }


    private void button1_Click(object sender, EventArgs e)
    {
      var inputs = new List<PointF>()
      {
        new PointF(0, 0),
        new PointF(this.originalBitmap.Width, this.originalBitmap.Height),
        new PointF(this.originalBitmap.Width, 0),
      };
      var outputs = new List<PointF>()
      {
        new PointF(0, 0),
        new PointF((float)(this.originalBitmap.Width * 1.5), (float)(this.originalBitmap.Height * 1.5)),
        new PointF((float)(this.originalBitmap.Width * 1.5), 0),
      };

      var affineCoefficients = CalculateAffineCoefficients(inputs, outputs);
      if (affineCoefficients == null)
      {
        MessageBox.Show("Error...", "Error", MessageBoxButtons.OK, MessageBoxIcon.Error);
        return;
      }

      // 逆変換
      var transferedBitmap = new Bitmap(this.originalBitmap.Width, this.originalBitmap.Height);
      for (int x = 0; x < this.originalBitmap.Width; ++x)
      {
        for (int y = 0; y < this.originalBitmap.Height; ++y)
        {
          Point correctedPoint = GetCorrectPointByAffineTransfer(affineCoefficients, x, y);

          Color color;
          if (correctedPoint.X >= 0 && correctedPoint.X < this.originalBitmap.Width &&
            correctedPoint.Y >= 0 && correctedPoint.Y < this.originalBitmap.Height)
          {
            color = this.originalBitmap.GetPixel(correctedPoint.X, correctedPoint.Y);
          }
          else
          {
            color = Color.Black;
          }

          transferedBitmap.SetPixel(x, y, color);
        }
      }

      this.pictureBox2.Image = transferedBitmap;
    }

    /// <summary>
    /// アフィン変換係数を求める
    /// </summary>
    /// <param name="observedPoints"></param>
    /// <param name="transferedPoints"></param>
    /// <returns></returns>
    private double[] CalculateAffineCoefficients(List<PointF> observedPoints, List<PointF> transferedPoints)
    {
      if (observedPoints.Count < 3 || observedPoints.Count != transferedPoints.Count)
      {
        // 3点未満の際は、計算不可なので、nullを返す
        return null;
      }

      var matrixValues = new double[DataNumberForAffineCoefficients / 2][]
      {
        new double[] { 0.0, 0.0, 0.0 },
        new double[] { 0.0, 0.0, 0.0 },
        new double[] { 0.0, 0.0, 0.0 }
      };
      var vectors = new double[DataNumberForAffineCoefficients] { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 };

      int index = 0;
      foreach (var observedPoint in observedPoints)
      {
        var observedX = (double)observedPoint.X;
        var observedY = (double)observedPoint.Y;
        var transferedX = (double)transferedPoints[index].X;
        var transferedY = (double)transferedPoints[index].Y;

        matrixValues[0][0] += observedX * observedX;
        var xy = observedX * observedY;
        matrixValues[0][1] += xy;
        matrixValues[1][0] += xy;
        matrixValues[0][2] += observedX;
        matrixValues[2][0] += observedX;
        matrixValues[1][1] += observedY * observedY;
        matrixValues[1][2] += observedY;
        matrixValues[2][1] += observedY;
        matrixValues[2][2] += 1;

        vectors[0] += observedX * transferedX;
        vectors[1] += observedY * transferedX;
        vectors[2] += transferedX;
        vectors[3] += observedX * transferedY;
        vectors[4] += observedY * transferedY;
        vectors[5] += transferedY;

        index++;
      }

      // 逆行列を計算する
      if (!CalculateInverseMatrix(matrixValues, out double[][] inverseMatrixValues))
      {
        // 逆行列算出に失敗したため、null
        return null;
      }

      var affineCoefficients = new double[DataNumberForAffineCoefficients];
      affineCoefficients[0] = inverseMatrixValues[0][0] * vectors[0] + inverseMatrixValues[0][1] * vectors[1] + inverseMatrixValues[0][2] * vectors[2];
      affineCoefficients[1] = inverseMatrixValues[1][0] * vectors[0] + inverseMatrixValues[1][1] * vectors[1] + inverseMatrixValues[1][2] * vectors[2];
      affineCoefficients[2] = inverseMatrixValues[2][0] * vectors[0] + inverseMatrixValues[2][1] * vectors[1] + inverseMatrixValues[2][2] * vectors[2];
      affineCoefficients[3] = inverseMatrixValues[0][0] * vectors[3] + inverseMatrixValues[0][1] * vectors[4] + inverseMatrixValues[0][2] * vectors[5];
      affineCoefficients[4] = inverseMatrixValues[1][0] * vectors[3] + inverseMatrixValues[1][1] * vectors[4] + inverseMatrixValues[1][2] * vectors[5];
      affineCoefficients[5] = inverseMatrixValues[2][0] * vectors[3] + inverseMatrixValues[2][1] * vectors[4] + inverseMatrixValues[2][2] * vectors[5];

      return affineCoefficients;
    }

    // 逆行列算出。詳細は、以下の関連記事を参照
    // https://blogs.yahoo.co.jp/dk521123/38068366.html
    private static bool CalculateInverseMatrix(
      double[][] targetMatrixes, out double[][] outInverseMatrixes)
    {
      var dataNumber = targetMatrixes.GetLength(0);
      outInverseMatrixes = new double[dataNumber][];
      //単位行列を作る
      for (int i = 0; i < dataNumber; i++)
      {
        outInverseMatrixes[i] = new double[dataNumber];
        for (int j = 0; j < dataNumber; j++)
        {
          outInverseMatrixes[i][j] = (i == j) ? 1.0f : 0.0f;
        }
      }

      // 掃き出し法
      for (int i = 0; i < dataNumber; i++)
      {
        if (targetMatrixes[i][i] == 0)
        {
          return false;
        }
        var tempValue = 1.0 / targetMatrixes[i][i];
        for (int j = 0; j < dataNumber; j++)
        {
          targetMatrixes[i][j] *= tempValue;
          outInverseMatrixes[i][j] *= tempValue;
        }
        for (int j = 0; j < dataNumber; j++)
        {
          if (i != j)
          {
            tempValue = targetMatrixes[j][i];
            for (int k = 0; k < dataNumber; k++)
            {
              targetMatrixes[j][k] -= targetMatrixes[i][k] * tempValue;
              outInverseMatrixes[j][k] -= outInverseMatrixes[i][k] * tempValue;
            }
          }
        }
      }
      return true;
    }

    /// <summary>
    /// アフィン変換による補正座標を返す
    /// </summary>
    /// <param name="affineCoefficients">アフィン係数</param>
    /// <param name="targetPoint">対象座標</param>
    /// <returns>アフィン変換による補正座標</returns>
    private Point GetCorrectPointByAffineTransfer(double[] affineCoefficients, int x, int y)
    {
      var correctedX = affineCoefficients[0] * x + affineCoefficients[1] * y + affineCoefficients[2];
      var correctedY = affineCoefficients[3] * x + affineCoefficients[4] * y + affineCoefficients[5];

      // 座標にするので四捨五入
      return new Point(
        (int)Math.Round(correctedX, MidpointRounding.AwayFromZero),
        (int)Math.Round(correctedY, MidpointRounding.AwayFromZero));
    }
  }
}


関連記事

C# / 画像処理

画像処理 ~ 輝度変更 (明るさ) ~
https://blogs.yahoo.co.jp/dk521123/37844934.html
画像処理 ~ コントラスト ~
https://blogs.yahoo.co.jp/dk521123/37852480.html
画像処理 ~ シャープネス ~
https://blogs.yahoo.co.jp/dk521123/37837353.html
画像処理 ~ アフィン変換 ~
https://blogs.yahoo.co.jp/dk521123/38061211.html
画像処理 ~ アフィン変換で任意角度の回転を自作する ~
https://blogs.yahoo.co.jp/dk521123/38093149.html
画像処理 ~ 幾何補正 / 擬似2次元アフィン変換編 ~
https://blogs.yahoo.co.jp/dk521123/38106095.html

その他

画像処理 ~ 幾何補正について ~
https://blogs.yahoo.co.jp/dk521123/38073336.html
画像処理 ~ 補間方法について ~
https://blogs.yahoo.co.jp/dk521123/38080694.html
C# で 行列計算を自作する ~ 逆行列編 ~
https://blogs.yahoo.co.jp/dk521123/38068366.html