iIR数字滤波器设计在工程实践中具有十分重要的作用，使用C语言实现对深化理解滤波器设计的理论和功能具有深远的意义，凡是选择本项目并能达到项目基本要求可以获得90分以上的成绩。
要求：报告中必须有能够实现以下设计指标的数字低通滤波器的证明：通带和阻带截止频率分别为0.2π和0.35π，通带最大衰减为1dB，阻带最小衰减为10dB。
提示：理论知识参考《数字信号处理》教材。
这是c语言课程设计,写一个文章来记录一下,方便以后查漏补缺,
基础知识
低通滤波器的主要性能指标有以下几个：通带截止频率fp、阻带截止频率fs、通带衰减 ( Ap)、阻带衰减 ( As) 以及归一化频率时需要用到的-3dB的转折频率fc。
根据给定的参数设计模拟滤波器，然后进行变数变换，求取数字滤波器的方法，称为滤波器的间接设计。做为数字滤波器的设计基础的模拟滤波器，称之为原型滤波器。
模拟滤波器的设计
巴特沃斯滤波器因其在通带内的幅值特性具有最大平坦的特性而闻名，是四种经典滤波器中最简单的，巴特沃斯滤波器只需要两个参数表征，滤波器的阶数N和-3dB处的截止频率
。其幅度平方函数为：
N是滤波器的阶数，从幅度平方函数可以看出，N阶滤波器有2N个极点，而且这2N个极点均布在一个圆上，圆的半径为，称之为巴特沃斯圆，巴特沃斯滤波器系统是一个线性系统，要使其稳定，其极点必须位于S平面的左半平面，所以取左半平面内的N个极点作为滤波器的极点，滤波器就是稳定的了，求出极点之后，计算模拟滤波器的系数as、bs，然后通过双线性变换（不懂得自行查书）由模拟域到数字域，求出系数az和bz 。最后通过差分方程就可以计算滤波结果了。
次数N的计算为
巴特沃斯低通滤波器的传递函数，可由其振幅特性的分母多项式求得。其分母多项式
根据S解开，可以得到极点。这里，为了方便处理，我们分为两种情况去解这个方程。当N为偶数的时候，
这里，使用了欧拉公式
。同样的，当N为奇数的时候，
同样的，这里也使用了欧拉公式。归纳以上，极点的解为
上式所求得的极点，是在s平面内，在半径为Ωc的圆上等间距的点，其数量为2N个。为了使得其IIR滤波器稳定，那么，只能选取极点在S平面左半平面的点。选定了稳定的极点之后，其模拟滤波器的传递函数就可由下式求得。
变换为Z域
然后通过差分方程就可以计算滤波结果了
C语言的实现
1.求阶数
代码:
N = ceil(0.5*( log10 (( pow (10, Stopband_attenuation/10) - 1)/
( pow (10, Passband_attenuation/10) - 1)) / log10 (Stopband/Passband) ));
2.求极点
for(k = 0;k <= ((2*N)-1) ; k++)
{
if(Cutoff*cos((2*k+N-1)*(pi/(2*N))) < 0.0)
{
poles[count].x = -Cutoff*cos((2*k+N-1)*(pi/(2*N)));
poles[count].y = -Cutoff*sin((2*k+N-1)*(pi/(2*N)));
count++;
if (count == N) break;
}
}
这里面 poles 是复数类型_complex
_Complex是C99新增的关键字，表示一种基本数据类型——复数
该类型的出现主要是为了解决工程和数学计算上很多涉及到复数计算的情况
3.求模拟滤波器系数
计算出稳定的极点之后，就可以进行传递函数的计算了。传递的函数的计算，就像下式一样
这里，为了得到模拟滤波器的系数，需要将分母乘开。很显然，这里的极点不一定是整数，或者来说，这里的乘开需要做复数运算。其复数的乘法代码如下，
//复数乘
_complex ComplexMul(_complex a, _complex b)
{
_complex c;
c.x = a.x*b.x - a.y*b.y;
c.y = a.y*b.x + a.x*b.y;
return c;
}
有了乘法代码之后，我们现在简单的情况下，看看其如何计算其滤波器系数。我们做如下假设`
这个时候，其传递函数为
将其乘开，其大致的关系就像下图所示一样。
计算的关系一目了然，这样的话，实现就简单多了。高阶的情况下也一样，重复这种计算就可以了。其代码为
Res[0].x = poles[0].x;
Res[0].y = poles[0].y;
Res[1].x = 1;
Res[1].y= 0;
for(count_1 = 0;count_1 < N-1;count_1++)//N个极点相乘次数
{
for(count = 0;count <= count_1 + 2;count++)
{
if(0 == count)
{
Res_Save[count] = ComplexMul( Res[count], poles[count_1+1] );
}
else if((count_1 + 2) == count)
{
Res_Save[count].x += Res[count - 1].x;
Res_Save[count].y += Res[count - 1].y;
}
else
{
Res_Save[count] = ComplexMul( Res[count], poles[count_1+1] );
Res_Save[count].x += Res[count - 1].x;
Res_Save[count].y += Res[count - 1].y;
}
}
for(count = 0;count <= N;count++)//Res[i]=a(i),i越大次数越高
{
Res[count].x = Res_Save[count].x;
Res[count].y = Res_Save[count].y;
*(a + N - count) = Res[count].x;
}
}
*(b+N) = *(a+N);
到此，我们就可以得到一个模拟滤波器巴特沃斯低通滤波器了。
4.S变Z
公式
代码
int Count = 0,Count_1 = 0,Count_2 = 0,Count_Z = 0;
double *Res, *Res_Save;
Res = new double[N+1]();
Res_Save = new double[N+1]();
memset(Res, 0, sizeof(double)*(N+1));
memset(Res_Save, 0, sizeof(double)*(N+1));
for(Count_Z = 0;Count_Z <= N;Count_Z++)
{
*(az+Count_Z) = 0;
*(bz+Count_Z) = 0;
}
for(Count = 0;Count<=N;Count++)
{
for(Count_Z = 0;Count_Z <= N;Count_Z++)
{
Res[Count_Z] = 0;
Res_Save[Count_Z] = 0;
}
Res_Save [0] = 1;
for(Count_1 = 0; Count_1 < N-Count;Count_1++)//计算（1-Z^-1）^N-Count的系数,
{ //Res_Save[]=Z^-1多项式的系数，从常数项开始
for(Count_2 = 0; Count_2 <= Count_1+1;Count_2++)
{
if(Count_2 == 0)
{
Res[Count_2] += Res_Save[Count_2];
}
else if((Count_2 == (Count_1+1))&&(Count_1 != 0))
{
Res[Count_2] += -Res_Save[Count_2 - 1];
}
else
{
Res[Count_2] += Res_Save[Count_2] - Res_Save[Count_2 - 1];
}
}
for(Count_Z = 0;Count_Z<= N;Count_Z++)
{
Res_Save[Count_Z] = Res[Count_Z] ;
Res[Count_Z] = 0;
}
}
for(Count_1 = (N-Count); Count_1 < N;Count_1++)//计算(1-Z^-1)^N-Count*（1+Z^-1）^Count的系数,
{ //Res_Save[]=Z^-1多项式的系数，从常数项开始
for(Count_2 = 0; Count_2 <= Count_1+1;Count_2++)
{
if(Count_2 == 0)
{
Res[Count_2] += Res_Save[Count_2];
}
else if((Count_2 == (Count_1+1))&&(Count_1 != 0))
{
Res[Count_2] += Res_Save[Count_2 - 1];
}
else
{
Res[Count_2] += Res_Save[Count_2] + Res_Save[Count_2 - 1];
}
}
for(Count_Z = 0;Count_Z<= N;Count_Z++)
{
Res_Save[Count_Z] = Res[Count_Z] ;
Res[Count_Z] = 0;
}
}
for(Count_Z = 0;Count_Z<= N;Count_Z++)
{
*(az+Count_Z) += pow(2,N-Count) * (*(as+Count)) * Res_Save[Count_Z];
*(bz+Count_Z) += (*(bs+Count)) * Res_Save[Count_Z];
}
}//最外层for循环
for(Count_Z = N;Count_Z >= 0;Count_Z--)
{
*(bz+Count_Z) = (*(bz+Count_Z))/(*(az+0));
*(az+Count_Z) = (*(az+Count_Z))/(*(az+0));
}
5.差分方程计算滤波结果
采用滤波器直接II型结构，可以减少一半的中间缓存内存，具体代码如下：
double FiltButter(double *pdAz, //滤波器参数表1
double *pdBz, //滤波器参数表2
int nABLen, //参数序列的长度
double dDataIn,//输入数据
double *pdBuf) //数据缓冲区
{
int i;
int nALen;
int nBLen;
int nBufLen;
double dOut;
if( nABLen<1 )return 0.0;
//根据参数,自动求取序列有效长度
nALen = nABLen;
for( i=nABLen-1; i; --i )
{
if( *(pdAz+i) != 0.0 )//从最后一个系数判断是否为0
{
nALen = i+1;
break;
}
}
//printf("%lf ", nALen);
if( i==0 ) nALen = 0;
nBLen = nABLen;
for( i=nABLen-1; i; --i )
{
if( *(pdBz+i) != 0.0 )
{
nBLen = i+1;
break;
}
}
//printf("%lf ", nBLen);
if( i==0 ) nBLen = 0;
//计算缓冲区有效长度
nBufLen = nALen;
if( nALen < nBLen)
nBufLen = nBLen;
//滤波: 与系数a卷乘
dOut = ( *pdAz ) * dDataIn; // a(0) * x(i)
for( i=1; i<nALen; i++) // a(i) * w(n-i),i=1toN
{
dOut -= *(pdAz+i) * *(pdBuf + (nBufLen - 1) - i);
}
//卷乘结果保存为缓冲序列的最后一个
*(pdBuf + nBufLen - 1) = dOut;
//滤波: 与系数b卷乘
dOut = 0.0;
for( i=0; i<nBLen; i++) // b(i) * w(n-i)
{
dOut += *(pdBz+i) * *(pdBuf + (nBufLen - 1) - i);
}
//丢弃缓冲序列中最早的一个数, 最后一个数清零
for( i=0; i<nBufLen-1; i++)
{
*(pdBuf + i) = *(pdBuf + i + 1);
}
*(pdBuf + nBufLen - 1) = 0;
//返回输出值
return dOut;
}
完整程序
#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <iostream>
#include <malloc.h>
#include <graphics.h>
#include<conio.h>//getch()函数头文件
using namespace std;
#define pi ((double)3.1415926)
struct DESIGN_SPECIFICATION
{
double Passband;
double Stopband;
double Passband_attenuation;
double Stopband_attenuation;
};
//复数乘
_complex ComplexMul(_complex a, _complex b)
{
_complex c;
c.x = a.x*b.x - a.y*b.y;
c.y = a.y*b.x + a.x*b.y;
return c;
}
//复数除
_complex ComplexDiv(_complex a, _complex b)
{
_complex c;
double dTmp = b.x * b.x + b.y * b.y;
//除数太小,防溢出
if (dTmp < 1.0E-300) {
c.x = 1.0E300;
c.y = 1.0E300;
return c;
}
c.x = (a.x * b.x + a.y * b.y) / dTmp;
c.y = (a.x * b.y - a.y * b.x) / dTmp;
return c;
}
//计算滤波器阶数N
int Buttord(double Passband,
double Stopband,
double Passband_attenuation,
double Stopband_attenuation)
{
int N;
printf("Wp = %lf [rad/sec] \n", Passband);
printf("Ws = %lf [rad/sec] \n", Stopband);
printf("Ap = %lf [dB] \n", Passband_attenuation);
printf("As = %lf [dB] \n", Stopband_attenuation);
printf("--------------------------------------------------------\n");
N = ceil(0.5*(log10((pow(10, Stopband_attenuation / 10) - 1) /
(pow(10, Passband_attenuation / 10) - 1)) / log10(Stopband / Passband)));
return (int)N;
}
//计算S平面的滤波器系数
int Butter(int N, double Cutoff, double *a, double *b)
{
double dk = 0;
int k = 0;
int count = 0, count_1 = 0;
_complex *poles, *Res, *Res_Save;
poles = new _complex[N]();
Res = new _complex[N + 1]();
Res_Save = new _complex[N + 1]();
memset(poles, 0, sizeof(_complex)*(N));
memset(Res, 0, sizeof(_complex)*(N + 1));
memset(Res_Save, 0, sizeof(_complex)*(N + 1));
if ((N % 2) == 0) dk = 0.5;
else dk = 0;
for (k = 0; k <= ((2 * N) - 1); k++)
{
if (Cutoff*cos((2 * k + N - 1)*(pi / (2 * N))) < 0.0)
{
poles[count].x = -Cutoff * cos((2 * k + N - 1)*(pi / (2 * N)));
poles[count].y = -Cutoff * sin((2 * k + N - 1)*(pi / (2 * N)));
count++;
if (count == N) break;
}
}
printf("Pk = \n");
for (count = 0; count < N; count++)
{
printf("(%lf) + (%lf i) \n", -poles[count].x, -poles[count].y);
}
printf("--------------------------------------------------------\n");
Res[0].x = poles[0].x;
Res[0].y = poles[0].y;
Res[1].x = 1;
Res[1].y = 0;
for (count_1 = 0; count_1 < N - 1; count_1++)//N个极点相乘次数
{
for (count = 0; count <= count_1 + 2; count++)
{
if (0 == count)
{
Res_Save[count] = ComplexMul(Res[count], poles[count_1 + 1]);
}
else if ((count_1 + 2) == count)
{
Res_Save[count].x += Res[count - 1].x;
Res_Save[count].y += Res[count - 1].y;
}
else
{
Res_Save[count] = ComplexMul(Res[count], poles[count_1 + 1]);
Res_Save[count].x += Res[count - 1].x;
Res_Save[count].y += Res[count - 1].y;
}
}
for (count = 0; count <= N; count++)//Res[i]=a(i),i越大次数越高
{
Res[count].x = Res_Save[count].x;
Res[count].y = Res_Save[count].y;
*(a + N - count) = Res[count].x;
}
}
*(b + N) = *(a + N);
//*(b+N) = pow(Cutoff,N);
//------------------------打印 bs as ---------------------------------//
printf("bs = [");
for (count = 0; count <= N; count++)
{
printf("%lf ", *(b + count));
}
printf(" ] \n");
printf("as = [");
for (count = 0; count <= N; count++)
{
printf("%lf ", *(a + count));
}
printf(" ] \n");
printf("--------------------------------------------------------\n");
delete[]poles;
delete[]Res;
delete[]Res_Save;
return (int)1;
}
//双线性Z变换，S->Z
int Bilinear(int N,
double *as, double *bs,
double *az, double *bz)
{
int Count = 0, Count_1 = 0, Count_2 = 0, Count_Z = 0;
double *Res, *Res_Save;
Res = new double[N + 1]();
Res_Save = new double[N + 1]();
memset(Res, 0, sizeof(double)*(N + 1));
memset(Res_Save, 0, sizeof(double)*(N + 1));
for (Count_Z = 0; Count_Z <= N; Count_Z++)
{
*(az + Count_Z) = 0;
*(bz + Count_Z) = 0;
}
for (Count = 0; Count <= N; Count++)
{
for (Count_Z = 0; Count_Z <= N; Count_Z++)
{
Res[Count_Z] = 0;
Res_Save[Count_Z] = 0;
}
Res_Save[0] = 1;
for (Count_1 = 0; Count_1 < N - Count; Count_1++)//计算（1-Z^-1）^N-Count的系数,
{ //Res_Save[]=Z^-1多项式的系数，从常数项开始
for (Count_2 = 0; Count_2 <= Count_1 + 1; Count_2++)
{
if (Count_2 == 0)
{
Res[Count_2] += Res_Save[Count_2];
}
else if ((Count_2 == (Count_1 + 1)) && (Count_1 != 0))
{
Res[Count_2] += -Res_Save[Count_2 - 1];
}
else
{
Res[Count_2] += Res_Save[Count_2] - Res_Save[Count_2 - 1];
}
}
for (Count_Z = 0; Count_Z <= N; Count_Z++)
{
Res_Save[Count_Z] = Res[Count_Z];
Res[Count_Z] = 0;
}
}
for (Count_1 = (N - Count); Count_1 < N; Count_1++)//计算(1-Z^-1)^N-Count*（1+Z^-1）^Count的系数,
{ //Res_Save[]=Z^-1多项式的系数，从常数项开始
for (Count_2 = 0; Count_2 <= Count_1 + 1; Count_2++)
{
if (Count_2 == 0)
{
Res[Count_2] += Res_Save[Count_2];
}
else if ((Count_2 == (Count_1 + 1)) && (Count_1 != 0))
{
Res[Count_2] += Res_Save[Count_2 - 1];
}
else
{
Res[Count_2] += Res_Save[Count_2] + Res_Save[Count_2 - 1];
}
}
for (Count_Z = 0; Count_Z <= N; Count_Z++)
{
Res_Save[Count_Z] = Res[Count_Z];
Res[Count_Z] = 0;
}
}
for (Count_Z = 0; Count_Z <= N; Count_Z++)
{
*(az + Count_Z) += pow(2, N - Count) * (*(as + Count)) * Res_Save[Count_Z];
*(bz + Count_Z) += (*(bs + Count)) * Res_Save[Count_Z];
}
}//最外层for循环
for (Count_Z = N; Count_Z >= 0; Count_Z--)
{
*(bz + Count_Z) = (*(bz + Count_Z)) / (*(az + 0));
*(az + Count_Z) = (*(az + Count_Z)) / (*(az + 0));
}
//------------------------display---------------------------------//
printf("bz = [");
for (Count_Z = 0; Count_Z <= N; Count_Z++)
{
printf("%lf ", *(bz + Count_Z));
}
printf(" ] \n");
printf("az = [");
for (Count_Z = 0; Count_Z <= N; Count_Z++)
{
printf("%lf ", *(az + Count_Z));
}
printf(" ] \n");
printf("--------------------------------------------------------\n");
delete[]Res;
delete[]Res_Save;
return (int)1;
}
// 执行滤波
// 返回值: 滤波处理后的数
double FiltButter(double *pdAz, //滤波器参数表1
double *pdBz, //滤波器参数表2
int nABLen, //参数序列的长度
double dDataIn,//输入数据
double *pdBuf) //数据缓冲区
{
int i;
int nALen;
int nBLen;
int nBufLen;
double dOut;
if (nABLen < 1)return 0.0;
//根据参数,自动求取序列有效长度
nALen = nABLen;
for (i = nABLen - 1; i; --i)
{
if (*(pdAz + i) != 0.0)//从最后一个系数判断是否为0
{
nALen = i + 1;
break;
}
}
//printf("%lf ", nALen);
if (i == 0) nALen = 0;
nBLen = nABLen;
for (i = nABLen - 1; i; --i)
{
if (*(pdBz + i) != 0.0)
{
nBLen = i + 1;
break;
}
}
//printf("%lf ", nBLen);
if (i == 0) nBLen = 0;
//计算缓冲区有效长度
nBufLen = nALen;
if (nALen < nBLen)
nBufLen = nBLen;
//滤波: 与系数a卷乘
dOut = (*pdAz) * dDataIn; // a(0) * x(i)
for (i = 1; i < nALen; i++) // a(i) * w(n-i),i=1toN
{
dOut -= *(pdAz + i) * *(pdBuf + (nBufLen - 1) - i);
}
//卷乘结果保存为缓冲序列的最后一个
*(pdBuf + nBufLen - 1) = dOut;
//滤波: 与系数b卷乘
dOut = 0.0;
for (i = 0; i < nBLen; i++) // b(i) * w(n-i)
{
dOut += *(pdBz + i) * *(pdBuf + (nBufLen - 1) - i);
}
//丢弃缓冲序列中最早的一个数, 最后一个数清零
for (i = 0; i < nBufLen - 1; i++)
{
*(pdBuf + i) = *(pdBuf + i + 1);
}
*(pdBuf + nBufLen - 1) = 0;
//返回输出值
return dOut;
}
int main(void)
{
int count;
double Fcutoff;
struct DESIGN_SPECIFICATION IIR_Filter;
IIR_Filter.Passband = (double)((pi*2)/10); //通带截止频率[rad]
IIR_Filter.Stopband = (double)(( pi*35)/100); //阻带截止频率[rad]
IIR_Filter.Passband_attenuation = 1; //通带最大衰减[dB]
IIR_Filter.Stopband_attenuation = 10; //阻带最小衰减[dB]
int N;
IIR_Filter.Passband = 2 * tan((IIR_Filter.Passband) / 2); //[rad/sec]
IIR_Filter.Stopband = 2 * tan((IIR_Filter.Stopband) / 2); //[rad/sec]
N = Buttord(IIR_Filter.Passband,
IIR_Filter.Stopband,
IIR_Filter.Passband_attenuation,
IIR_Filter.Stopband_attenuation);
printf("N: %d \n", N);
Fcutoff = IIR_Filter.Passband / (pow((pow(10, 0.1*IIR_Filter.Passband_attenuation) - 1), 1.0 / (2 * N)));
printf("Wc: %lf \n", Fcutoff);
printf("--------------------------------------------------------\n");
double *as = new double[N + 1]();
double *bs = new double[N + 1]();
memset(as, 0, sizeof(double)*(N + 1));
memset(bs, 0, sizeof(double)*(N + 1));
Butter(N, Fcutoff, as, bs);//计算模拟滤波器系数
double *az;
double *bz;
az = new double[N + 1]();
bz = new double[N + 1]();
memset(az, 0, sizeof(double)*(N + 1));
memset(bz, 0, sizeof(double)*(N + 1));
Bilinear(N, as, bs, az, bz); //模拟变为数字滤波器
//生成叠加信号
printf(" 请输入信号1频率 单位 Hz \n");
double w1, w2;
scanf("%lf", &w1);
printf(" 请输入信号2频率 单位 Hz \n");
scanf("%lf", &w2);
printf("采样频率为 10kHz \n");
FILE* Input_Data1;
Input_Data1 = fopen("E:\\yssj.txt", "w");
printf("生成原始叠加信号并且保存到E:\\yssj.txt\n");
double t, s1;
for (int i = 0; i < 10000; i++)//4秒，产生更多数据
{
t = i / 10000.0; //采样频率
s1 = sin(2 * pi *w1 * t) + sin(2 * pi *w2* t);//设定频率为w1Hz+w2hz
fprintf(Input_Data1, "%lf,\n", s1);
}
fclose(Input_Data1);
char s[100];
int length;
double *data_Buffer;
data_Buffer = (double *)malloc(sizeof(double)*(N + 1));
memset(data_Buffer, 0, sizeof(double)*(N + 1));
double Output = 0;
FILE* Input_Data;
FILE* Output_Data;
Input_Data = fopen("E:\\yssj.txt", "r");
Output_Data = fopen("E:\\lbh.txt", "w");
length = N + 1;
long DataCnt = 0;
for (long j = 0; !feof(Input_Data); j++)
{
fscanf(Input_Data, "%s\n", &s);
DataCnt++;
}
// cout<<DataCnt<<"\n";
rewind(Input_Data);
double *D = new double[DataCnt]();
for (int i = 0; i < DataCnt; i++)
{
fscanf(Input_Data, "%s\n", &s);
D[i] = atof(s);
Output = FiltButter(az,bz,length,D[i],data_Buffer);
fprintf(Output_Data, "%lf\n", Output);
//cout<<Output<<"\n";
}
free(data_Buffer);
fclose(Input_Data);
fclose(Output_Data);
cout << "按任意键绘制图形" << "\n";
_getch();
Input_Data = fopen("E:\\yssj.txt", "r");
Input_Data1 = fopen("E:\\lbh.txt", "r");
char s11[100];
long DataCnt1 = 0;
for (long j = 0; !feof(Input_Data); j++)
{
fscanf(Input_Data, "%s\n", &s);
DataCnt++;
}
for (long j = 0; !feof(Input_Data1); j++)
{
fscanf(Input_Data1, "%s\n", &s11);
DataCnt1++;
}
rewind(Input_Data);
rewind(Input_Data1);
double y1=0, x1=0;
double *D1 = new double[DataCnt]();
initgraph(800, 1080); // 初始化640x480的绘图屏幕
for (int i = 0; i < DataCnt; i++)
{
fscanf(Input_Data, "%s\n", &s);
D1[i] = atof(s);
s1 = D1[i];
t = i / 10000.0; //采样频率
line(x1 * 50000, 200 + y1 * 100, t * 50000, 200 + s1 * 100);
x1 = t;
y1 = s1;
}
y1 = 0;
x1 = 0;
for (int i = 0; i < DataCnt1; i++)
{
fscanf(Input_Data1, "%s\n", &s);
D1[i] = atof(s);
s1 = D1[i];
t = i / 10000.0; //采样频率
line(x1 * 50000, 600 + y1 * 100, t * 50000, 600 + s1 * 100);
x1 = t;
y1 = s1;
}
// printf("Finish \n");
delete[]D;
delete[]as;
delete[]bs;
delete[]az;
delete[]bz;
_getch();
closegraph();
fclose(Input_Data);
fclose(Input_Data1);
return (int)0;
}
运行结果
参考资料：http://blog.csdn.net/zhoufan900428/article/details/9069475
具体源码 : https://download.csdn.net/download/qq_37131451/13767299
3271
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
