通过对电场做DFT方法计算S参数的C语言代码

////////////////////////////////////////////////////////////////////////// // 计算S参数的代码 // obs-port1.txt : 记录输入端口入射波的数据 // obs-total.txt : 记录输入端口总场的数据 // obs-port2.txt : 记录输出端口的电场数据 // 2010年12月3日 ////////////////////////////////////////////////////////////////////////// #include <complex> using namespace std; typedef complex<double> Complex ; #define PI 3.14159265358979323846 #define TN 4000 //时间采样点数 #define FN 300 //频率采样点数 ////////////////////////////////////////////////////////////////////////// // 计算s参数 ////////////////////////////////////////////////////////////////////////// void main() { printf("/n Calculate S Parameter Start../n"); double dx = 0.4064E-3 ; // //0.402*14 //0.402*6 = double dy = 0.4233E-3 ; //*5 double dz = 0.265E-3 ; //0.265e-3; double cc = 2.99792458e8; double dt = 0.99 / (cc * sqrt(1.0 / (dx * dx) + 1.0 / (dy * dy) + 1.0 / (dz * dz))); int i; double dB; double inc_t[TN]; //入射场时域数值 double tot_t[TN]; //总场时域数值 double tra_t[TN]; //传输场时域数值 Complex inc_f[FN]; Complex ref_f[FN]; Complex tot_f[FN]; Complex tra_f[FN]; Complex s11,s21; FILE* inc = fopen("obs-port1.txt","r");//打开1端口入射文件 FILE* tot = fopen("obs-total.txt","r");//打开1端口总场文件 FILE* tra = fopen("obs-port2.txt","r");//打开2端口电场文件 if (!inc) { printf("obs-port1.txt lost!"); return; } if (!tot) { printf("obs-total.txt lost!"); return; } if (!tra) { printf("obs-port2.txt"); return; } FILE* sf11 = fopen("result_s11.txt","w"); fprintf(sf11,"Freq(GHz) /t S11(real) /t S11(imag) /t S11(abs) /t S11(dB)/n"); FILE* sf21 = fopen("result_s21.txt","w"); fprintf(sf21,"Freq(GHz) /t S21(real) /t S21(imag) /t S21(abs) /t S21(dB)/n"); double df = 0.1e9;//频率间隔100MHz for(i=0;i<FN;i++){ inc_f[i] = Complex(0.,0.); ref_f[i] = Complex(0.,0.); tot_f[i] = Complex(0.,0.); tra_f[i] = Complex(0.,0.); } for(i=0;i<TN;i++){ fscanf(inc,"%lf",&inc_t[i]); fscanf(tot,"%lf",&tot_t[i]); fscanf(tra,"%lf",&tra_t[i]); } //对电场做离散傅里叶变换 for(i=0;i<FN;i++){ for(int n=0;n<TN;n++){ inc_f[i] += inc_t[n]*exp(Complex(0,-2*PI*df*i*dt*n)); tot_f[i] += tot_t[n]*exp(Complex(0,-2*PI*df*i*dt*n)); tra_f[i] += tra_t[n]*exp(Complex(0,-2*PI*df*i*dt*n)); } //S11计算 s11 = (tot_f[i]-inc_f[i]) / inc_f[i]; dB = 20*log10(abs(s11)); fprintf(sf11,"%lf /t %lf /t %lf /t %lf /t %lf /n" , (df*i/1e9) , real(s11), imag(s11), abs(s11), dB); //S21计算 s21 = tra_f[i] / inc_f[i] ; dB = 20*log10(abs(s21)); fprintf(sf21,"%lf /t %lf /t %lf /t %lf /t %lf /n" , (df*i/1e9) , real(s21), imag(s21), abs(s21), dB); } fclose(inc); fclose(tot); fclose(sf11); fclose(sf21); printf("/n Calculate S Parameter Finished ../n/n"); return; }