現在地 >> メニュー >> サンプルコード::OpenCV >> OpenCV::パワースペクトルの描画

関連 OpenCV::離散フーリエ変換::線形性


問題


入力データを離散フーリエ変換して、そのパワースペクトルを計算し棒グラフで描画せよ。
その他条件:
 わかりやすいように、棒グラフ両端の色を変える。

答え

#include <iostream>
#include <cstdio>

#include <cv.h>
#include <highgui.h>

/********* 各種外部変数 ***********/
const CvSize WindowSize = {300,300};

double data[][2]=
{
{0,0}, //0+0i
{0,0},
{1,0},  //1+0i
{1,0},
{1,0},
{1,0},
{0,0},
{0,0},
};

const int ELEMENT = sizeof(data)/sizeof(data[0]);//入力データ数の計算

/*********** プロトタイプ宣言 ***********/
void Show_Result(const CvMat &dataM);
double GET_MAX(const double *power_spectrum);


/************** ここからメイン関数 **********************/
int main()
{


	CvMat dataM = cvMat (1,ELEMENT, CV_64FC2, &data[0][0]); //縦1,横ELEMENT

	puts("------------- 離散フーリエの計算結果 -------------");
	cvDFT(&dataM,&dataM,CV_DXT_FORWARD,dataM.rows);
	for(int loop = 0;loop < ELEMENT; ++loop)
	{
		printf("%f,%f \n",dataM.data.db[loop*2],dataM.data.db[loop*2+1]);
	}


	Show_Result(dataM); //グラフにして表示


	return EXIT_SUCCESS;
}


/**********変換結果を出力*************/
void Show_Result(const CvMat &dataM)
{

	using std::cout;


	IplImage *imgA = cvCreateImage(WindowSize,IPL_DEPTH_8U,3);//ウィンドウの初期化
	cvSet (imgA, cvScalarAll (255), 0);


	double Power_Spectrum[ELEMENT]; //パワースペクトル用配列

	puts("----------- パワースペクトル ------------------");
	for(int loop = 0; loop < ELEMENT; ++loop) //パワースペクトルの計算
	{
		Power_Spectrum[loop] = pow(dataM.data.db[loop*2],2) + pow(dataM.data.db[loop*2+1],2); 
		cout << Power_Spectrum[loop] << "\n";
	}



	double Max = GET_MAX(&Power_Spectrum[0]);  //パワースペクトルの最大値を取得



	/*** 出力用にスケール変換 ***/
	for(int loop = 0; loop < ELEMENT; ++loop)
	{
		Power_Spectrum[loop] = (Power_Spectrum[loop]/Max) * imgA->height;
	}


	//横幅決定
	int unit_width = cvRound( imgA->width / static_cast<double>(ELEMENT) );
	

	/*右端の余白からシフト分を計算*/
	int space_shift = cvRound*1*0.5);


	/**** 描画 ***/
	CvScalar color={0,255,0};
	bool flag = false;
	for (int j = 0; j < ELEMENT; j++)
	{
		if (j == 0 || j == (ELEMENT - 1) )
		{
			flag = true; //両端を赤く塗るフラグ
		}

		cvLine(imgA,cvPoint(j * unit_width+space_shift , imgA->height),cvPoint (j * unit_width +space_shift ,
			imgA->height - cvRound (Power_Spectrum[j])), flag ? CV_RGB(255,0,0) : color ,3);

		flag = false;
	}



	cvNamedWindow("after DFT",CV_WINDOW_AUTOSIZE);
	cvShowImage("after DFT",imgA);

	cvWaitKey(0);

	cvReleaseImage(& imgA);
	cvDestroyWindow("after DFT");

}


/*** 配列の最大値を取得 ***/
double GET_MAX(const double *power_spectrum)
{
	double temp = power_spectrum[0];
	for(int loop = 1; loop < ELEMENT; ++loop)
	{
		if(temp < power_spectrum[loop])
		{
			temp = power_spectrum[loop];
		}

	}

	return temp;

}

別ver.:実数のみの入力データを用意し、そこからデータを作成する場合。


(上の場合、すでに実部、虚部のデータを用意していた)


#include <iostream>
#include <cstdio>

#include <cv.h>
#include <highgui.h>


/********* 各種外部変数 ***********/
const CvSize WindowSize = {300,300};
double input[]={0,0,1,1,1,1,0,0};
const int ELEMENT = sizeof(input)/sizeof(input[0]);


/*********** プロトタイプ宣言 ***********/
CvMat *GetMat(double *data);
void Show_Result(const CvMat *dataM);
double GET_MAX(const double *power_spectrum);


/************** ここからメイン関数 **********************/
int main()
{

	CvMat *dataM = GetMat(input);


	puts("------------- 離散フーリエの計算結果 -------------");
	cvDFT(dataM,dataM,CV_DXT_FORWARD,dataM->rows);
	for(int loop = 0;loop < ELEMENT; ++loop)
	{
		printf("%f,%f \n",dataM->data.db[loop*2],dataM->data.db[loop*2+1]);
	}


	Show_Result(dataM); //グラフにして表示

	cvReleaseMat(&dataM);

	return EXIT_SUCCESS;
}


/********** データから計算用に行列を作成 ***********/
CvMat *GetMat(double *data)
{
	CvMat *complex = cvCreateMat(1,ELEMENT, CV_64FC2);
	CvMat Re = cvMat(1,ELEMENT, CV_64FC1,&input[0]);
	CvMat *Im = cvCreateMat(1,ELEMENT, CV_64FC1);
	cvSetZero(Im);

	cvMerge(&Re,Im,NULL,NULL,complex);
	cvReleaseMat(&Im);
	return complex;
}


/**********変換結果を出力*************/
void Show_Result(const CvMat *dataM)
{

	using std::cout;


	IplImage *imgA = cvCreateImage(WindowSize,IPL_DEPTH_8U,3);//ウィンドウの初期化
	cvSet (imgA, cvScalarAll (255), 0);

	double Power_Spectrum[ELEMENT];//パワースペクトル用配列


	puts("----------- パワースペクトル ------------------");
	for(int loop = 0; loop < ELEMENT; ++loop) //パワースペクトルの計算
	{
		Power_Spectrum[loop] = pow(dataM->data.db[loop*2],2) + pow(dataM->data.db[loop*2+1],2); 
		cout << Power_Spectrum[loop] << "\n";
	}



	double Max = GET_MAX(&Power_Spectrum[0]);  //パワースペクトルの最大値を取得



	/*** 出力用にスケール変換 ***/
	for(int loop = 0; loop < ELEMENT; ++loop)
	{
		Power_Spectrum[loop] = (Power_Spectrum[loop]/Max) * imgA->height;
	}


	//横幅決定
	int unit_width = cvRound( imgA->width / static_cast<double>(ELEMENT) );
	

	/*右端の余白からシフト分を計算*/
	int space_shift = cvRound*2*0.5);


	/**** 描画 ***/
	CvScalar color={0,255,0};
	bool flag = false;
	for (int j = 0; j < ELEMENT; j++)
	{
		if (j == 0 || j == (ELEMENT - 1) )
		{
			flag = true; //両端を赤く塗るフラグ
		}

		cvLine(imgA,cvPoint(j * unit_width+space_shift , imgA->height),cvPoint (j * unit_width +space_shift ,
			imgA->height - cvRound (Power_Spectrum[j])), flag ? CV_RGB(255,0,0) : color ,3);

		flag = false;
	}


	cvNamedWindow("after DFT",CV_WINDOW_AUTOSIZE);
	cvShowImage("after DFT",imgA);

	cvWaitKey(0);

	cvReleaseImage(& imgA);
	cvDestroyWindow("after DFT");

}


/*** 配列の最大値を取得 ***/
double GET_MAX(const double *power_spectrum)
{
	double temp = power_spectrum[0];
	for(int loop = 1; loop < ELEMENT; ++loop)
	{
		if(temp < power_spectrum[loop])
		{
			temp = power_spectrum[loop];
		}

	}

	return temp;

}


関連 OpenCV::離散フーリエ変換::線形性

目次

― その他 ―

Wiki内検索

計測中...(07.10.8〜)

Save The World






▲よろしければ広告のクリックもお願いします


▲ランキングに参加しました

管理人/副管理人のみ編集できます