正弦信号采样

[jsjygm]

给定正弦信号的周期T，采样频率P，如何获取采样的值。用c。
如果要32个点，这样行不：

for(int i=0;i<32;i++)
a=sin(2*3.14*T*i*P);

[就想叫yoko]

http://wenwen.soso.com/z/q125166245.htm

[ryfdizuo]

	const int len = 50;
	const int amplitude = 30;
	const int PathSegment = 20;
	for (int i=0; i<=PathSegment; i++)
	{
		float cy = amplitude * sin( i * 1.0f/PathSegment * PI );
		float cz = i * 1.0f/PathSegment * len;

		path.push_back( RVec3(0.0f, cy, cz) );
	}

上面代码把正弦曲线0-π区间映射到 0-len长度，振幅0-1映射到0-amplitude

[ryfdizuo]

PathSegment是采样点的个数。

[jsjygm]

[code=C/C++][path.push_back]

[pathuang68]

[Quote=引用楼主 jsjygm 的回复:]
给定正弦信号的周期T，采样频率P，如何获取采样的值。用c。
如果要32个点，这样行不：

for(int i=0;i<32;i++)
a=sin(2*3.14*T*i*P);
[/Quote]
不行。

[bayuan0414]

const double pi = 3.14;
int i; 
int Ta;//信号频率
int fs;//采样频率
for (i=0; i<32; i++) {
    printf("x(%d) = %f\n", i, sin(2*pi/Ta/fs*i));
}

[bayuan0414]

[Quote=引用 8 楼 jsjygm 的回复:]

还有，那个正弦信号应该是无限长的，也就是x是无穷的。
[/Quote]
信号是无限长的，但是你只需要取32个点。
x(t)=A*sin(w*t)=A*sin(2*pi/Ta*t) 其中Ta是信号周期
x(n)=A*sin(2*pi/Ta*T*n) 其中T是采样周期（#6的注释有误）
在信号周期和采样周期固定的情况下，一个信号周期内能采多少个点基本是固定的（新序列x(n)是否是周期的以及周期是多少取决于Ta和T）。
我理解你这里的32个采样点就是简单的从无限长的采样点中连续取32个点。

[jsjygm]

[Quote=引用 2 楼 dizuo 的回复:]
C/C++ code
const int len = 50;
const int amplitude = 30;
const int PathSegment = 20;
for (int i=0; i<=PathSegment; i++)
{
float cy = amplitude * sin( i * 1.0f/PathSeg……
[/Quote]
path.push_back( RVec3(0.0f, cy, cz) );
是什么意思，我在msdn上也没找到相应的解释呀。

[jsjygm]

还有，那个正弦信号应该是无限长的，也就是x是无穷的。

[jsjygm]

我们老让做一个fft,我的代码是这样的：

#include <stdio.h>
#include <math.h>
#include <stdlib.h>

#define N 32
typedef struct
{
	double real;/*实部*/
	double img;/*虚部*/
}complex;
void fft(); /*快速傅立叶变换*/
void ifft(); /*快速傅立叶逆变换*/
void initW(); /*计算旋转因子*/
void change(); /*变址运算*/
void add(complex ,complex ,complex *); /*复数加*/
void mul(complex ,complex ,complex *); /*复数乘*/
void sub(complex ,complex ,complex *); /*复数减*/
void divi(complex ,complex ,complex *);/*复数除*/
void output(); /*输出结果*/
complex x[N],*W;/*输出序列的值*/
double xk[N];
//int 32=0;/*序列大小*/
double PI;



int main()
{
	int ch;
	int method;
	int i,fc;
	double fre;
	const int len = 50;
	const int amplitude = 32;
	const int PathSegment = 32;
	while(1)
	{
		printf("\n	1 进行fft运算     2 退出：");
		scanf("%d",&ch);
		if(ch==2)
			exit(1);
		
		
		PI=atan(1)*4;    /*pi等于4乘以1.0的正切值*/

		printf("采样周期选择：1 0.000625  2 0.005  3 0.0046875  4 0.004:");
		scanf("%d",&fc);
		switch(fc)
		{
		case 1:fre=0.000625; break;
		case 2:fre=0.005;    break;
		case 3:fre=0.0046875;break;
		case 4:fre=0.004;    break;
		default: exit(1);
		}
		
		
		/*采样*/
		for(i=0;i<32;i++)
		{
			x[i].img=0;
			x[i].real=sin(2*PI*50*fre*i);
			xk[i]=x[i].real;
//			printf("%f	",xk[i]);
		}     



/*for (i=0; i<=PathSegment; i++)
	{
		float cy = amplitude * sin( i * 1.0f/PathSegment * PI );
		float cz = i * 1.0f/PathSegment * len;
			
		::path.push_back( RVec3(0.0f, cy, cz) );
		xk[i]=x[i]=cy;
		}*/




		initW();
		
		printf("Use FFT(1) or IFFT(2):"); //选择变换方式
		scanf("%d",&method);
		if(method==1)
			fft();
		else
			ifft();
		output();     //输出结果
		
	}
	return 0;
}




/*进行基2 fft运算*/
void fft()
{
	int i=0,j=0,k=0,l=0;
	complex up,down,product;
	change();
	for(i=0;i< log(32)/log(2) ;i++) /*一级蝶形运算*/
	{
//		printf("%f %f    ",W[i].img,W[i].real);
		l=1<<i;
		for(j=0;j<32;j+= 2*l ) /*一组蝶形运算*/
		{
			for(k=0;k<l;k++) /*一个蝶形运算*/
			{
				mul(x[j+k+l],W[32*k/2/l],&product);
				add(x[j+k],product,&up);
				sub(x[j+k],product,&down);
				x[j+k]=up;
				x[j+k+l]=down;
			}
//			printf("		%.4f\n",xk[i]);
		}
		
	}
}

void ifft()
{
	int i=0,j=0,k=0,l=32;
	complex up,down;
	for(i=0;i< (int)( log(32)/log(2) );i++) /*一级蝶形运算*/
	{
		l/=2;
		for(j=0;j<32;j+= 2*l ) /*一组蝶形运算*/
		{
			for(k=0;k<l;k++) /*一个蝶形运算*/
			{
				add(x[j+k],x[j+k+l],&up);
				up.real/=2;up.img/=2;
				sub(x[j+k],x[j+k+l],&down);
				down.real/=2;down.img/=2;
				divi(down,W[32*k/2/l],&down);
				x[j+k]=up;
				x[j+k+l]=down;
			}
		}
	}
	change();
}




/*计算旋转因子*/
void initW()
{
	int i;
	W=(complex *)malloc(sizeof(complex) * 32);
	for(i=0;i<32;i++)
	{
		W[i].real=cos(2*PI/32*i);
		W[i].img=-1*sin(2*PI/32*i);
//		printf("%f   %f",W[i].img,W[i].real);
	}
}






/*变址运算，将x（n）码位倒置*/
void change()
{
	complex temp;
	unsigned short i=0,j=0,k=0;
	double t;

	for(i=0;i<32;i++)
	{
		k=i;j=0;
		t=(log(32)/log(2));  //计算log以2为底32的对数
		while( (t--)>0 )
		{
			j=j<<1;   //左移
			j|=(k & 1);    //j=j|(k&1)
			k=k>>1;   //右移
		}
		if(j>i)
		{
			temp=x[i]; 
			x[i]=x[j];
			x[j]=temp;
		}
	}
	for(i=0;i<32;i++)
	{
		xk[i]=x[i].real;
	}
}



/*输出结果*/
void output() 
{
	int i;
	printf("fft计算结果如下：\n");
	for(i=0;i<32;i++)
	{
		printf("%4f	",xk[i]);

		printf("	%f",x[i].real);
		if(x[i].img>=0.0001)
			printf("+%4fj\n",x[i].img);
		else if(fabs(x[i].img)<0.0001)
			printf("\n");
		else
			printf("%4fj\n",x[i].img);
		
	}
}
void add(complex a,complex b,complex *c)
{
	c->real=a.real+b.real;
	c->img=a.img+b.img;
}
void mul(complex a,complex b,complex *c)
{
	c->real=a.real*b.real - a.img*b.img;
	c->img=a.real*b.img + a.img*b.real;
}
void sub(complex a,complex b,complex *c)
{
	c->real=a.real-b.real;
	c->img=a.img-b.img;
}
void divi(complex a,complex b,complex *c)
{
	c->real=( a.real*b.real+a.img*b.img )/(b.real*b.real+b.img*b.img);
	c->img=( a.img*b.real-a.real*b.img)/(b.real*b.real+b.img*b.img);
}

前面输入的是采样经过变换之后的值，后面是fft值，结果为什么没有虚部？