#define pi 3.141593
#define f 3 //原正弦波频率
#define cx pi/3 //原正弦波初相,每次采样的时候是随机的
double out[10000]; //存插值后各点值用的,大点好
double in[10000]; //采样值缓存,大点好
//使用Kaiser窗截断,并使用参数9
//num:采样点数
//insert_num:每2个采样点插入的点数
//T:周期数
//start_phase:生成的正弦波的起始相位
//a:生成正弦波的a系数
//b:生成正弦波的b系数
void test_sinc(int sample_num, int insert_num, double T, double start_phase, double a, double b)
{
insert_num ++;
for(int m = 0;m < sample_num;m++){//采样
in[m] = a*sin( (2*pi*m*T)/(sample_num) + start_phase) + b;
}
for(int u = 0, i =0;u< sample_num;u++){
out[i++] = in[u];
for(int k = 1;k < insert_num;k++){
double sum = 0;
double t = (double)u + (double)k/((double)insert_num);
for(int m = 0;m < sample_num;m++){
double tmp = pi*(t -(double)m);
sum += in[m]*sin(tmp)/tmp;
}
out[i++] = sum;
}
}
}
//生成hn系数数组的算法,len是数组总长度,insert_num是插入点数
void SincHn(double *data, int len, int insert_num)
{
int end = len /2;
const int start = -end;
if(len%2 != 0)
end++;
insert_num ++;
for(int i = start; i < end; i++){
if(i == 0){
*data++ = 1;
continue;
}
double tmp = pi * (double)i/ (double)insert_num;
*data++ = sin(tmp)/tmp;
}
}
unsigned long Fact(unsigned long k)
{
unsigned long sum = 1;
for(unsigned long i = 1; i <= k; i++)
sum *= i;
return sum;
}
double besseli_0(double x, int len)
{
double sum = 0;
for(int i = 0; i < len; i++){
sum += std::pow(-1, i) * std::pow(x/2, 2*i)/(Fact((unsigned long)i) * Fact((unsigned long)i+1));
}
return sum;
}
void Kaiser(double *data, int N, double beta)
{
double theta;
for( int n = 0; n < N; n++){
theta = beta * std::sqrt( 1 - std::pow( ( (2 * (double)n/(N -1)) - 1),2 ) );
double res = besseli_0(theta, 10) / besseli_0(beta, 10);
if(res < 0)
res = 0;
*data++ = res;
}
}
一键复制
编辑
Web IDE
原始数据
按行查看
历史