DSP

基于opencv的一维Fourier变换

2019-07-13 19:17发布

此程序主要实现了fft、ifft计算。原本是用于cpu雷达信号脉冲压缩的C++实现。程序调试后无错误,可以运行,程序需要读取磁盘中的.dat文件实现数据导入。 调试环境:VS2015+opencv //2016.11.2 @Xidian //DSP雷达信号处理2 //调用opencv中dft()函数实现fft和ifft #include #include #include using namespace std; using namespace cv; const int N = 7680; const int M = 8192; void main() { double count1 = getTickCount(); //-------------------对回波信号echo进行8192点FFT----------------- FILE* fp; fp=fopen ("echo.dat", "r"); float temp[7680 * 2] = {0}; vector<float> echo; // for (int i = 0; i < 7680*2; i++) { fscanf(fp, "%f", &temp[i]); echo.push_back(temp[i]); } fclose(fp); //echo信号的实部、虚部 Mat_<float> echo_real = Mat::zeros(1, 8192, CV_32F); Mat_<float> echo_img = Mat::zeros(1, 8192, CV_32F); vector<float>hh[2]; for (int i = 0; i < 7680 * 2; i += 2) { hh[0].push_back(echo[i]); hh[1].push_back(echo[i + 1]); } for (int i = 0; i < 7680; i++) { echo_real.at<float>(i) = hh[0][i]; echo_img.at<float>(i) = hh[1][i]; } //将echo的实、虚部组合为复数complex vectorfloat>> data; data.push_back(echo_real); data.push_back(echo_img); Mat complex; merge(data, complex); //计算FFT dft(complex, complex); //将echo信号的8192点fft写入.dat文件 vectorfloat>> echo_fft; split(complex, echo_fft); FILE* fp1; FILE* fp2; fp1 = fopen("echo_fft_real.dat", "w"); fp2 = fopen("echo_fft_img.dat", "w"); for (int i = 0; i < 8192; i++) { //将echo进行FFT后的实部写入文本 fprintf(fp1, "%f ", echo_fft[0].at<float>(i)); //将echo进行FFT后的虚部写入文本 fprintf(fp2, "%f ", echo_fft[1].at<float>(i)); } //-----------------对脉压系数coeff进行8192点FFT----------------- //Mat_ coeff_real = Mat::zeros(1, 8192,CV_32F); //coeff的实部 //Mat_ coeff_img = Mat::zeros(1, 8192, CV_32F); //coeff的虚部 //Mat_ coeff_complex; //coeff的复数 //Mat_ coeff_fft = Mat::zeros(1, 8192, CV_32F); //coeff的8192点FFT // //FILE* fp3; //fp3 = fopen("coeff.dat", "r"); //float temp3[7680 * 2]; //for (int i = 0; i <7680*2; i++) //{ //读取coeff.dat文本数据 // fscanf(fp3, "%f", &temp3[i]); //} //vector GG[2]; //for (int i = 0; i < 7680*2; i+=2) //{ // GG[0].push_back(temp3[i]); // GG[1].push_back(temp3[i + 1]); //} //for (int i = 0; i < 7680; i++) //{ // coeff_real.at(i) = GG[0][i]; // coeff_img.at(i) = GG[1][i]; //} //vector>coeff; //coeff.push_back(coeff_real); //coeff.push_back(coeff_img); ////组合复数 //merge(coeff, coeff_complex); ////FFT计算 //dft(coeff_complex, coeff_complex); // //----------------------读取coeff_fft并保存为复数--------------------------- Mat_<float> coeff_fftReal2 = Mat::zeros(1, 8192,CV_32F); //coeff_fft的实部 Mat_<float> coeff_fftImg2 = Mat::zeros(1, 8192, CV_32F); //coeff_fft的虚部 Mat coeff_fftComplex2; //coeff_fft的复数 FILE* fp4; fp4 = fopen("coeff_fft.dat", "r"); float temp4[7680 * 2] = {0}; for (int i = 0; i < 7680*2; i++) { fscanf(fp4, "%f", &temp4[i]); } vector<float> JJ[2]; for (int i = 0; i < 7680*2; i+=2) { JJ[0].push_back(temp4[i]); JJ[1].push_back(temp4[i+1]); } for (int i = 0; i < 7680; i++) { coeff_fftReal2.at<float>(i) = JJ[0][i]; coeff_fftImg2.at<float>(i) = JJ[1][i]; } vectorfloat>>coeff_fft; coeff_fft.push_back(coeff_fftReal2); coeff_fft.push_back(coeff_fftImg2); //组合为复数 merge(coeff_fft, coeff_fftComplex2); //------------------频域相乘----------------------------- Mat_<float> pc_freq0Img; Mat_<float> pc_freq0Real; Mat pc_freq0Complex; vectorfloat>> pc_freq0; for (int i = 0; i < 8192; i++) { Vec2f a = complex.at(i); Vec2f b = coeff_fftComplex2.at(i); float real = a[0] * b[0] - a[1] * b[1]; float img = a[0] * b[1] + a[1] * b[0]; img = -1.0*img; /* cout << "real=" << real << endl; cout << "imag=" << img << endl;*/ pc_freq0Real.push_back(real); pc_freq0Img.push_back(img); //共轭 } //组合为复数 pc_freq0.push_back(pc_freq0Real); pc_freq0.push_back(pc_freq0Img); merge(pc_freq0, pc_freq0Complex); /*cout << pc_freq0Complex << endl;*/ //---------------------IFFT计算----------------------- //直接计算ifft dft(pc_freq0Complex, pc_freq0Complex, DFT_SCALE); //???与matlab计算结果相差太大 /*cout << pc_freq0Complex << endl;*/ //用fft计算ifft //step1 计算fft //dft(pc_freq0Complex, pc_freq0Complex); ////step2 共轭、乘以1/N //for (int i = 0; i < 8192; i++) //{ // Vec2f bb = pc_freq0Complex.at(i); // float real2 = bb[0]; // float imag2 = bb[1]; // imag2 = imag2*(-1.0); // pc_freq0Complex.at(i) = (1.0/8192.0)*Vec2f(real2, imag2); //} // //cout << pc_freq0Complex << endl; //------------------计算结果写入.dat文件(包含暂态点)------------------------ FILE* fp5; FILE* fp6; fp5 = fopen("pc_freq1_real.dat","w"); fp6 = fopen("pc_freq1_imag.dat","w"); for (int i = 0; i < 8192; i++) { Vec2f c = pc_freq0Complex.at(i); /*cout << "pc_freq0.real=" << c[0]<fprintf(fp5, "%f ",c[0]); fprintf(fp6, "%f ", c[1]); } fclose(fp5); fclose(fp6); double count2 = getTickCount(); double freq = getTickFrequency(); double time = (count2 - count1) / freq; cout << "运行时间=" << time << endl; cout << "Successful!" << endl; system("PAUSE"); }