#include<iostream>
#include<cmath>
#include<fstream>
using namespace std;
const double pi=4*atan(1);
void butterworthlowpass(double,double,double,double,double);
int main()
{
	//數字低通濾波器性能要求
    double f,fc,fst;//f:抽樣頻率；fc:通帶截止頻率；fst:阻帶截止頻率
    double q1,q2;//q1:通帶允許的最大衰減；q2:阻帶允許的最小衰減
	cout<<"please input the parameter:"<<endl;
	cin>>f>>fc>>fst>>q1>>q2;
	butterworthlowpass( f,fc,fst, q1, q2);
    return 0;}

void butterworthlowpass(double f,double fc,double fst,double q1,double q2)
{	double Wc,Wst;//數字頻率
     Wc=2*pi*fc/f;
	 Wst=2*pi*fst/f;
	//轉換為模擬低通濾波器性能要求
    double Uc,Ust,Usc,Ksc;
     Uc=2*f*tan(Wc/2);
     Ust=2*f*tan(Wst/2);
     Usc=Ust/Uc;
     Ksc=sqrt((pow(10,0.1*q1)-1)/(pow(10,0.1*q2)-1));
    double log(double num);
    double v;
    v=-log(Ksc)/log(Usc);//計算階數
    double ceil(double x);
   int N=ceil(v);//階數圓整
    double Up;
    Up=Uc*pow((pow(10,0.1*q1)-1),-(1/(2*N)));
    double a,b[20],d;//模擬低通濾波器系統函數參數
	a=1;
	d=1;
	int k;
	for(k=1;k<=(N/2);k++)
	{b[k]=-2*Up*cos(pi/2+(2*k-1)*pi/(2*N));
	}

	double c;
     c=2*f;
	double A0[20],A1[50],A2[50],B0[50],B1[50],B2[50],R[50];//數字低通濾波器系統函數參數
	ofstream outfile("f1.dat",ios::out);//定義文件流對象，打開磁盤文件“f1.dat"
	if(!outfile)
	{ 
	   cerr<<"open error!"<<endl;
	exit(1);
	}
	for (k=1;k<=N/2;k++)
	{
		R[k]=c*c+c*b[k]+Up*Up;
		A0[k]=Up*Up/R[k];
		A1[k]=2*A0[k];
		A2[k]=A0[k];
		B0[k]=1;
		B1[k]=(2*Up*Up-2*c*c)/R[k];
		B2[k]=(Up*Up-b[k]*c+c*c)/R[k];
	}
cout<<v<<endl;
		if (N%2==0)
		{
			for (k=1;k<=N/2;k++)
				outfile<<endl<<A0[k]<<" "<<A1[k]<<" "<<A2[k]<<" "<<B0[k]<<" "<<B1[k]<<" "<<B2[k]<<endl;
				}
		else 
		{
			outfile<<endl<<1/(1+c)<<" "<<1/(1+c)<<" "<<1<<" "<<(1-c)/(1+c)<<endl;
for (k=1;k<=N/2;k++)
				outfile<<A0[k]<<" "<<A1[k]<<" "<<A2[k]<<" "<<B0[k]<<" "<<B1[k]<<" "<<B2[k]<<endl;
		}
		
	outfile.close();	}