ページ 11

ガウス消去法の射影変換(座標変換の計算について)

Posted: 2014年2月26日(水) 14:32
by もう六時G
現在、ガウス消去法による射影変換を勉強しています。
パラメータの計算までは求められているのですが、変換前の座標と変換パラメータを用いて変換後の座標を計算するところで詰まっています。
以下、現在書いているコードです。効率の良し悪しはとりあえず放っておいてます。
最後の変換の計算部分についてよくないところ、間違っているところを指摘頂ければと思います。

よろしくお願い致します。

コード:

/* Gaussの単純消去法(部分ピボット選択) */
#include <stdio.h>
#include <stdlib.h>
#include <math.h>

#define X_SIZE 5000		/*画像の横サイズを定義*/
#define Y_SIZE 5000		/*画像の縦サイズを定義*/
#define MAXN 32

unsigned char screen[Y_SIZE][X_SIZE][3],img_work[Y_SIZE][X_SIZE][3],projection_work[Y_SIZE][X_SIZE][3];

int main(void)
{

/*--------------------【画像取得プログラム】用変数型--------------------*/

	FILE *fp;	/*ファイルポインタの宣言*/
	unsigned char header[54];			/*unsigned:符号無しの整数型、char:文字データ、header[54]:0~53の1次配列*/
									
	char input_file[128],gray_file[128],projection_file[128];	/*ファイル名の配列を定義*/


/*--------------------【射影変換プログラム】用変数型--------------------*/

	int n=8, m, i, j, k, l, p, imax, x1, x2, x3, x0, y1, y2, y3, y0, u1, u2, u3, u0, v1, v2, v3, v0, U, V, z;
	double b, t, t1, a[MAXN][MAXN + 1];
	double x[4],y[4],u[4],v[4];




/*--------------------ここから【画像取得プログラム】--------------------*/



/*----------画像ファイルを読み込む----------*/



	printf("処理を行う画像ファイル名:(拡張子:bmp)\n---");
	scanf("%s",input_file);					/*入力ファイル名を指定*/
	printf("処理を行う画像ファイル:「 %s 」.\n\n",input_file);
	fp=fopen(input_file,"rb"); 	/* 処理する画像をオープンする	
					Windowsビットマップ形式 X_SIZE*Y_SIZEピクセル,24ビットカラー、rb:バイナリデータの読み込みに使用 */
	if((fp=fopen(input_file,"rb"))==NULL){		/*ファイルが開けない場合のエラー表示*/
		printf("画像ファイルが正常に開けませんでした.\n");
		exit(1);
	}
	fread(header,1,54,fp); /* ヘッダ(54バイト)を飛ばす */
	fread(screen,1,X_SIZE*Y_SIZE*3,fp); /* 残りはデータ(最下行から順に入る) */
	fclose(fp);

/*----------カラー画像をグレースケールに変換する----------*/


	int X,Y,col;
	for(X=0;X<X_SIZE;X++){
		for(Y=0;Y<Y_SIZE;Y++){
			for(col=0;col<=2;col++){
				img_work[Y][X][col]=0.1145*screen[Y][X][0]
					+0.5866*screen[Y][X][1]+0.2989*screen[Y][X][2];	
			}
		}
	}

/*----------グレースケール画像を保存する----------*/


	printf("グレースケールにする画像ファイル名を入力してください.(拡張子:bmp)\n---");
	scanf("%s",gray_file);		/*グレースケールとするファイル名を指定*/
	fp=fopen(gray_file,"wb");	/*処理した画像を書き込むファイルをオープンする*/
	if((fp=fopen(gray_file,"wb"))==NULL){		/*ファイルが開けない場合のエラー表示*/
		printf("画像ファイルが正常に開けませんでした.\n");
		exit(1);
	}
	fwrite(header,1,54,fp); /* ヘッダ */
	fwrite(img_work,1,X_SIZE*Y_SIZE*3,fp); /* データ */
	fclose(fp);
	printf("グレースケールの画像ファイル「 %s 」を保存しました.\n\n",gray_file);




/*--------------------ここまで【画像取得プログラム】--------------------*/














/*--------------------ここから【射影変換プログラム】--------------------*/


   /* 行列の大きさをセット */
  printf("ガウスの消去法により射影変換のパラメータを求めます。\n");
  printf("↓ 8元連立方程式の結果 ↓\n\n");




//一般画像なら・・・ →〈1〉
//テンプレ画像なら・・・ →〈2〉


/*-----------↓ 〈1〉 ↓-----------*/
/*
  printf("変換前後の座標をそれぞれ入力してください。(左上から時計回りの順)\n");

	for(j = 0; j<4 ; j++){
		printf("x%d=", j+1);
			scanf("%lf", &x[j]);
		printf("y%d=", j+1);
			scanf("%lf", &y[j]);
		}
	for(m = 0; m<4 ; m++){
		printf("u%d=", m+1);
			scanf("%lf", &u[m]);
		printf("v%d=", m+1);
			scanf("%lf", &v[m]);
		}

	for(j = 0; j<4 ; j++){
		printf("x%d=%lf\n",j,x[j]);
		printf("y%d=%lf\n",j,y[j]);
		}
	for(j = 0; j<4 ; j++){
		printf("u%d=%lf\n",j,u[j]);
		printf("v%d=%lf\n",j,v[j]);
		}
*/

/*-----------↑ 〈1〉 ↑-----------*/


/*-----------↓ 〈2〉 ↓-----------*/

x[0]=229;
y[0]=84;

x[1]=333;
y[1]=84;

x[2]=438;
y[2]=0;

x[3]=0;
y[3]=0;


u[0]=0;
v[0]=700;

u[1]=203;
v[1]=700;

u[2]=203;
v[2]=0;

u[3]=0;
v[3]=0;



a[0][0]=x[0];
a[0][1]=y[0];
a[0][2]=1;
a[0][3]=0;
a[0][4]=0;
a[0][5]=0;
a[0][6]=-x[0]*u[0];
a[0][7]=-y[0]*u[0];
a[0][8]=u[0];

a[1][0]=x[1];
a[1][1]=y[1];
a[1][2]=1;
a[1][3]=0;
a[1][4]=0;
a[1][5]=0;
a[1][6]=-x[1]*u[1];
a[1][7]=-y[1]*u[1];
a[1][8]=u[1];

a[2][0]=x[2];
a[2][1]=y[2];
a[2][2]=1;
a[2][3]=0;
a[2][4]=0;
a[2][5]=0;
a[2][6]=-x[2]*u[2];
a[2][7]=-y[2]*u[2];
a[2][8]=u[2];

a[3][0]=x[3];
a[3][1]=y[3];
a[3][2]=1;
a[3][3]=0;
a[3][4]=0;
a[3][5]=0;
a[3][6]=-x[3]*u[3];
a[3][7]=-y[3]*u[3];
a[3][8]=u[3];

a[4][0]=0;
a[4][1]=0;
a[4][2]=0;
a[4][3]=x[0];
a[4][4]=y[0];
a[4][5]=1;
a[4][6]=-x[0]*v[0];
a[4][7]=-y[0]*v[0];
a[4][8]=v[0];

a[5][0]=0;
a[5][1]=0;
a[5][2]=0;
a[5][3]=x[1];
a[5][4]=y[1];
a[5][5]=1;
a[5][6]=-x[1]*v[1];
a[5][7]=-y[1]*v[1];
a[5][8]=v[1];

a[6][0]=0;
a[6][1]=0;
a[6][2]=0;
a[6][3]=x[2];
a[6][4]=y[2];
a[6][5]=1;
a[6][6]=-x[2]*v[2];
a[6][7]=-y[2]*v[2];
a[6][8]=v[2];

a[7][0]=0;
a[7][1]=0;
a[7][2]=0;
a[7][3]=x[3];
a[7][4]=y[3];
a[7][5]=1;
a[7][6]=-x[3]*v[3];
a[7][7]=-y[3]*v[3];
a[7][8]=v[3];

/*-----------↑ 〈2〉 ↑-----------*/



  for (k = 0; k < n; k++) {
    /* ピポットを走査する */
    t = 0; imax = k;
    for (l = k; l < n; l++) {
      if (fabs(a[l][k]) > t) {
        imax = l;
        t = a[l][k];
      }
    }
    /* 見つかったピボット行をk行と入れ替える */
    for (l = 0; l < n + 1; l++) {
      t = a[k][l]; 
		a[k][l] = a[imax][l]; 
		a[imax][l] = t;	// t ---> t1
    }
    /* 後は同じ */
    b = a[k][k];
    for (j = 0; j < n + 1; j++) a[k][j] /= b;	// a /= b  [a = a/b]
    
    /* ピポット行よりも下の行について */
    for (i = k + 1; i < n; i++) {
      /* ピポットと同じ列を0に */
      b = -a[i][k];
      for (j = 0; j < n + 1; j++) a[i][j] += b * a[k][j];
    }
  }

  /* 次に後置代入 */
  for (k = n - 2; k >= 0; k--)
    for (j = k + 1; j < n; j++)
      a[k][n] -= a[k][j] * a[j][n];
  
  /* 答を表示 */
  printf("Gauss単純消去法(部分ピボット選択)による解\n");
	for (i = 0; i < n; i++){
   	printf("--- x[%d]=%lf\n", i, a[i][n]);
	}

/*--------------------画像の変換処理--------------------*/   
//  (ここまでは一応正しい出力がされていることを確認済みです。)


//  (ここからの計算がうまくいきません。)

	for(X=0;X<X_SIZE;X++){
		for(Y=0;Y<Y_SIZE;Y++){
			for(col=0;col<=2;col++){
				U=((int)a[0][8]*X+(int)a[1][8]*Y+(int)a[2][8])/((int)a[6][8]*X+(int)a[7][8]*Y+1);
				V=((int)a[3][8]*X+(int)a[4][8]*Y+(int)a[5][8])/((int)a[6][8]*X+(int)a[7][8]*Y+1);
				projection_work[Y][X][col]=img_work[V][U][col];
			}
		}
	}


	printf("射影変換した画像ファイル名を入力してください.(拡張子:bmp)\n---");
	scanf("%s",projection_file);		/*グレースケールとするファイル名を指定*/
	fp=fopen(projection_file,"wb");	/*処理した画像を書き込むファイルをオープンする*/
	
	if((fp=fopen(projection_file,"wb"))==NULL){		/*ファイルが開けない場合のエラー表示*/
		printf("画像ファイルが正常に開けませんでした.\n");
		exit(1);
	}
	fwrite(header,1,54,fp); /* ヘッダ */
	fwrite(projection_work,1,X_SIZE*Y_SIZE*3,fp); /* データ */
	fclose(fp);
	printf("射影変換した画像ファイル「 %s 」を保存しました.\n\n",projection_file);


/*--------------------ここまで【射影変換プログラム】--------------------*/


  return 0;
}

Re: ガウス消去法の射影変換(座標変換の計算について)

Posted: 2014年2月26日(水) 15:04
by usao
303,304行を

コード:

//before
U=((int)a[0][8]*X+(int)a[1][8]*Y+(int)a[2][8])/((int)a[6][8]*X+(int)a[7][8]*Y+1);
V=((int)a[3][8]*X+(int)a[4][8]*Y+(int)a[5][8])/((int)a[6][8]*X+(int)a[7][8]*Y+1);

コード:

//after
U=(int)(   (a[0][8]*X+a[1][8]*Y+a[2][8]) / (a[6][8]*X+a[7][8]*Y+1)   );
V=(int)(   (a[3][8]*X+a[4][8]*Y+a[5][8]) / (a[6][8]*X+a[7][8]*Y+1)   );
としたら,なにかしら改善しませんか?

あと,この計算は colに関するforの中に入れなくてもいいような.

Re: ガウス消去法の射影変換(座標変換の計算について)

Posted: 2014年2月26日(水) 15:12
by usao
それと,対応点の座標次第では(U,V)が画像領域外になる可能性もあると思うので
ちゃんと判定した方がよいのではないでしょうか.

Re: ガウス消去法の射影変換(座標変換の計算について)

Posted: 2014年2月26日(水) 15:14
by もう六時G
回答ありがとうございます。

指摘頂いた通り直して実行してみたところ、Segmentation faultが出てしまいました。

書き忘れていたのですが、
質問に書いてあるプログラムは一応実行できますが、射影変換されるはずの画像がおかしくなってしまうといった状態です。

Re: ガウス消去法の射影変換(座標変換の計算について)

Posted: 2014年2月26日(水) 16:40
by usao
>対応点の座標次第では(U,V)が画像領域外になる可能性もある
Segmentation faultは, これのせいではなく,ということでしょうか?


あと,これって,画像変換の式が
 変換後座標(U,V) = f( 変換前座標(X,Y) )
なのではないでしょうか.
(あなたのコードをまともに見たわけではないので,ここらへんがどうなっているのかわかりませんが)
もし仮にそうだとしたら,
>projection_work[Y][X][col]=img_work[V][col];
これが違うかも.
オフトピック
関連して,この手の画像変換を行う際は,普通は(?)逆方向に
 変換前座標(X,Y)  = f_inv( 変換後座標(U,V) )
という計算をおこなって,(U,V)に対応する(X,Y)を求めます.
そうでないと,結果画像が穴あきになってしまうので.

Re: ガウス消去法の射影変換(座標変換の計算について)

Posted: 2014年2月27日(木) 18:37
by もう六時G
回答ありがとうございます。

指摘して頂いたことも参考にいろいろ試してみました。


質問文プログラムの303,304行目について

コード:

U=(a[0][8]*X+a[1][8]*Y+a[2][8]) / (a[6][8]*X+a[7][8]*Y+1);
V=(a[3][8]*X+a[4][8]*Y+a[5][8]) / (a[6][8]*X+a[7][8]*Y+1);
としてよいか試すために以下のプログラムを実行しました。

コード:

/* Gaussの単純消去法(部分ピボット選択) */
#include <stdio.h>
#include <stdlib.h>
#include <math.h>

#define X_SIZE 3000		/*画像の横サイズを定義*/
#define Y_SIZE 3000		/*画像の縦サイズを定義*/
#define MAXN 32

unsigned char screen[Y_SIZE][X_SIZE][3],img_work[Y_SIZE][X_SIZE][3],projection_work[Y_SIZE][X_SIZE][3];

int main(void)
{

/*--------------------【画像取得プログラム】用変数型--------------------*/

	FILE *fp;	/*ファイルポインタの宣言*/
	unsigned char header[54];			/*unsigned:符号無しの整数型、char:文字データ、header[54]:0~53の1次配列*/
									
	char input_file[128],gray_file[128],projection_file[128];	/*ファイル名の配列を定義*/


/*--------------------【射影変換プログラム】用変数型--------------------*/

	int n=8, m, i, j, k, l, p, imax, x1, x2, x3, x0, y1, y2, y3, y0, u1, u2, u3, u0, v1, v2, v3, v0, U1, V1, U2, V2, U3, V3, U4, V4, z;
	double b, t, t1, a[MAXN][MAXN + 1];
	double x[4],y[4],u[4],v[4];




/*--------------------ここから【画像取得プログラム】--------------------*/



/*----------画像ファイルを読み込む----------*/



	printf("処理を行う画像ファイル名:(拡張子:bmp)\n---");
	scanf("%s",input_file);					/*入力ファイル名を指定*/
	printf("処理を行う画像ファイル:「 %s 」.\n\n",input_file);
	fp=fopen(input_file,"rb"); 	/* 処理する画像をオープンする	
					Windowsビットマップ形式 X_SIZE*Y_SIZEピクセル,24ビットカラー、rb:バイナリデータの読み込みに使用 */
	if((fp=fopen(input_file,"rb"))==NULL){		/*ファイルが開けない場合のエラー表示*/
		printf("画像ファイルが正常に開けませんでした.\n");
		exit(1);
	}
	fread(header,1,54,fp); /* ヘッダ(54バイト)を飛ばす */
	fread(screen,1,X_SIZE*Y_SIZE*3,fp); /* 残りはデータ(最下行から順に入る) */
	fclose(fp);

/*----------カラー画像をグレースケールに変換する----------*/


	int X,Y,col;
	for(X=0;X<X_SIZE;X++){
		for(Y=0;Y<Y_SIZE;Y++){
			for(col=0;col<=2;col++){
				img_work[Y][X][col]=0.1145*screen[Y][X][0]
					+0.5866*screen[Y][X][1]+0.2989*screen[Y][X][2];	
			}
		}
	}

/*----------グレースケール画像を保存する----------*/


//	printf("グレースケールにする画像ファイル名を入力してください.(拡張子:bmp)\n---");
//	scanf("%s",gray_file);		/*グレースケールとするファイル名を指定*/
//	fp=fopen(gray_file,"wb");	/*処理した画像を書き込むファイルをオープンする*/
//	if((fp=fopen(gray_file,"wb"))==NULL){		/*ファイルが開けない場合のエラー表示*/
//		printf("画像ファイルが正常に開けませんでした.\n");
//		exit(1);
//	}
//	fwrite(header,1,54,fp); /* ヘッダ */
//	fwrite(img_work,1,X_SIZE*Y_SIZE*3,fp); /* データ */
//	fclose(fp);
//	printf("グレースケールの画像ファイル「 %s 」を保存しました.\n\n",gray_file);




/*--------------------ここまで【画像取得プログラム】--------------------*/














/*--------------------ここから【射影変換プログラム】--------------------*/


   /* 行列の大きさをセット */
  printf("ガウスの消去法により射影変換のパラメータを求めます。\n");
  printf("↓ 8元連立方程式の結果 ↓\n\n");




//一般画像なら		・・・ →〈1〉
//テンプレ画像なら・・・ →〈2〉


/*-----------↓ 〈1〉 ↓-----------*/
/*
  printf("変換前後の座標をそれぞれ入力してください。(左上から時計回りの順)\n");

	for(j = 0; j<4 ; j++){
		printf("x%d=", j+1);
			scanf("%lf", &x[j]);
		printf("y%d=", j+1);
			scanf("%lf", &y[j]);
		}
	for(m = 0; m<4 ; m++){
		printf("u%d=", m+1);
			scanf("%lf", &u[m]);
		printf("v%d=", m+1);
			scanf("%lf", &v[m]);
		}

	for(j = 0; j<4 ; j++){
		printf("x%d=%lf\n",j,x[j]);
		printf("y%d=%lf\n",j,y[j]);
		}
	for(j = 0; j<4 ; j++){
		printf("u%d=%lf\n",j,u[j]);
		printf("v%d=%lf\n",j,v[j]);
		}
*/

/*-----------↑ 〈1〉 ↑-----------*/


/*-----------↓ 〈2〉 ↓-----------*/

x[0]=229;
y[0]=84;

x[1]=333;
y[1]=84;

x[2]=438;
y[2]=0;

x[3]=0;
y[3]=0;


u[0]=0;
v[0]=700;

u[1]=203;
v[1]=700;

u[2]=203;
v[2]=0;

u[3]=0;
v[3]=0;



a[0][0]=x[0];
a[0][1]=y[0];
a[0][2]=1;
a[0][3]=0;
a[0][4]=0;
a[0][5]=0;
a[0][6]=-x[0]*u[0];
a[0][7]=-y[0]*u[0];
a[0][8]=u[0];

a[1][0]=x[1];
a[1][1]=y[1];
a[1][2]=1;
a[1][3]=0;
a[1][4]=0;
a[1][5]=0;
a[1][6]=-x[1]*u[1];
a[1][7]=-y[1]*u[1];
a[1][8]=u[1];

a[2][0]=x[2];
a[2][1]=y[2];
a[2][2]=1;
a[2][3]=0;
a[2][4]=0;
a[2][5]=0;
a[2][6]=-x[2]*u[2];
a[2][7]=-y[2]*u[2];
a[2][8]=u[2];

a[3][0]=x[3];
a[3][1]=y[3];
a[3][2]=1;
a[3][3]=0;
a[3][4]=0;
a[3][5]=0;
a[3][6]=-x[3]*u[3];
a[3][7]=-y[3]*u[3];
a[3][8]=u[3];

a[4][0]=0;
a[4][1]=0;
a[4][2]=0;
a[4][3]=x[0];
a[4][4]=y[0];
a[4][5]=1;
a[4][6]=-x[0]*v[0];
a[4][7]=-y[0]*v[0];
a[4][8]=v[0];

a[5][0]=0;
a[5][1]=0;
a[5][2]=0;
a[5][3]=x[1];
a[5][4]=y[1];
a[5][5]=1;
a[5][6]=-x[1]*v[1];
a[5][7]=-y[1]*v[1];
a[5][8]=v[1];

a[6][0]=0;
a[6][1]=0;
a[6][2]=0;
a[6][3]=x[2];
a[6][4]=y[2];
a[6][5]=1;
a[6][6]=-x[2]*v[2];
a[6][7]=-y[2]*v[2];
a[6][8]=v[2];

a[7][0]=0;
a[7][1]=0;
a[7][2]=0;
a[7][3]=x[3];
a[7][4]=y[3];
a[7][5]=1;
a[7][6]=-x[3]*v[3];
a[7][7]=-y[3]*v[3];
a[7][8]=v[3];

/*-----------↑ 〈2〉 ↑-----------*/



  for (k = 0; k < n; k++) {
    /* ピポットを走査する */
    t = 0; imax = k;
    for (l = k; l < n; l++) {
      if (fabs(a[l][k]) > t) {
        imax = l;
        t = a[l][k];
      }
    }


    /* 見つかったピボット行をk行と入れ替える */
    for (l = 0; l < n + 1; l++) {
      t = a[k][l]; 
		a[k][l] = a[imax][l]; 
		a[imax][l] = t;	// t ---> t1
    }

    /* 後は同じ */
    b = a[k][k];
    for (j = 0; j < n + 1; j++) {
		a[k][j] /= b;	// a /= b  [a = a/b]
    }

    /* ピポット行よりも下の行について */
    for (i = k + 1; i < n; i++) {
      /* ピポットと同じ列を0に */
      b = -a[i][k];
      for (j = 0; j < n + 1; j++) {
			a[i][j] += b * a[k][j];
		}
    }
  }

  /* 次に後置代入 */
  for (k = n - 2; k >= 0; k--)
    for (j = k + 1; j < n; j++)
      a[k][n] -= a[k][j] * a[j][n];
  
  /* 答を表示 */
  printf("Gauss単純消去法(部分ピボット選択)による解\n");
	for (i = 0; i < n; i++){
   	printf("--- x[%d]=%lf\n", i, a[i][n]);
	}

/*--------------------画像の変換処理--------------------*/

//	for(Y=0;Y<Y_SIZE;Y++){
//		for(X=0;X<X_SIZE;X++){
			for(col=0;col<=2;col++){


				U1=(a[0][8]*229+a[1][8]*84+a[2][8]) / (a[6][8]*229+a[7][8]*84+1);
				V1=(a[3][8]*229+a[4][8]*84+a[5][8]) / (a[6][8]*229+a[7][8]*84+1);
				U2=(a[0][8]*333+a[1][8]*84+a[2][8]) / (a[6][8]*333+a[7][8]*84+1);
				V2=(a[3][8]*333+a[4][8]*84+a[5][8]) / (a[6][8]*333+a[7][8]*84+1);
				U3=(a[0][8]*438+a[1][8]*0+a[2][8]) / (a[6][8]*438+a[7][8]*0+1);
				V3=(a[3][8]*438+a[4][8]*0+a[5][8]) / (a[6][8]*438+a[7][8]*0+1);
				U4=(a[0][8]*438+a[1][8]*0+a[2][8]) / (a[6][8]*438+a[7][8]*0+1);
				V4=(a[3][8]*438+a[4][8]*0+a[5][8]) / (a[6][8]*438+a[7][8]*0+1);
			//	U=(a[0][8]*(double)X+a[1][8]*(double)Y+a[2][8]) / (a[6][8]*(double)X+a[7][8]*(double)Y+1);
			//	V=(a[3][8]*(double)X+a[4][8]*(double)Y+a[5][8]) / (a[6][8]*(double)X+a[7][8]*(double)Y+1);
			//	U=((int)a[0][8]*X+(int)a[1][8]*Y+a[2][8]) / ((int)a[6][8]*X+(int)a[7][8]*Y+1);
			//	V=((int)a[3][8]*X+(int)a[4][8]*Y+a[5][8]) / ((int)a[6][8]*X+(int)a[7][8]*Y+1);
				printf("U1=%d\n",U1);
				printf("V1=%d\n",V1);
				printf("U2=%d\n",U2);
				printf("V2=%d\n",V2);
				printf("U3=%d\n",U3);
				printf("V3=%d\n",V3);
				printf("U4=%d\n",U4);
				printf("V4=%d\n",V4);
			//	projection_work[Y][X][col]=img_work[V][U][col];			//Is this line wrong ?
				projection_work[V1][U1][col]=img_work[Y][X][col];		//alternative
				projection_work[V2][U2][col]=img_work[Y][X][col];		//alternative
				projection_work[V3][U3][col]=img_work[Y][X][col];		//alternative
				projection_work[V4][U4][col]=img_work[Y][X][col];		//alternative
			}
//		}
//	}


//	printf("射影変換した画像ファイル名を入力してください.(拡張子:bmp)\n---");
//	scanf("%s",projection_file);		/*グレースケールとするファイル名を指定*/
//	fp=fopen(projection_file,"wb");	/*処理した画像を書き込むファイルをオープンする*/
	
//	if((fp=fopen(projection_file,"wb"))==NULL){		/*ファイルが開けない場合のエラー表示*/
//		printf("画像ファイルが正常に開けませんでした.\n");
//		exit(1);
//	}
//	fwrite(header,1,54,fp); /* ヘッダ */
//	fwrite(projection_work,1,X_SIZE*Y_SIZE*3,fp); /* データ */
//	fclose(fp);
//	printf("射影変換した画像ファイル「 %s 」を保存しました.\n\n",projection_file);


/*--------------------ここまで【射影変換プログラム】--------------------*/


  return 0;
}


これであれば変換前座標に対応した変換後座標が出るようです。
しかし、質問文のプログラムの303,304行目を実際に変えると segmentation fault が出てしまいます。

自分なりに配列の範囲?を変えるなどして調べましたが、原因がわかりません。

また、305行目

コード:

projection_work[Y][X][col]=img_work[V][U][col];
の直し方もよくわからないので教えていただければと思います。

よろしくお願い致します。

Re: ガウス消去法の射影変換(座標変換の計算について)

Posted: 2014年2月27日(木) 20:26
by usao
>これであれば変換前座標に対応した変換後座標が出るようです。

ということは,計算内容はこう
>変換後座標(U,V) = f( 変換前座標(X,Y) )
だということですね.
(X,Y)が入力で,(U,V)が出力という形.

で,少なくとも,最初に指定した対応点については,改善した計算式コードによって
(この入出力関係の上での)正しい出力が得られることが確認できた,と.


とりあえずsegmentation faultと言われる理由を解明するために,

コード:

for( Y=0; Y<Y_SIZE; Y++ ){
  for( X=0; X<X_SIZE; X++ ){
    U = ...;
    V = ...;
    ( X,Y )→( U,V ) という出力をしてみる
  }
}
みたいなことをされてみてはいかがでしょうか.
UやVが負の値だったり,U>=X_SIZE とか V>=Y_SIZE とかになっている箇所があるのではないでしょうか.


>また、305行目
>projection_work[Y][X][col]=img_work[V][col];
>の直し方もよくわからないので教えていただければと思います。

入力画像が img_work で,出力画像が projection_work なのであれば,当然ながら,
img_workにはX,Yを,projection_workにはU,Vを使うのだと思いますが,どうなんでしょう.