  /*
  
    COPYRIGHT 1992  Steven D. Swift
    Rights to copy are granted provided that the
    program is neither sold nor used for commercial
    purposes and that this notice stay with all
    copies.
    
history:

    rev 0:    December 19, 1991    Steven D. Swift, P.E.
                                Seattle, Washington 
    rev 1:    Feb 22, 1993
              Fixed bug in unequal term section (output) that
              caused bad L values for k=m

    rev 2:    March 15, 1993
              Found more bugs:  a:  center capacitor for N=5 missing
              Adjusted so normalized values are correct for
              unequal loads.
              Adjust for correct unequal load values
              (note: this uses Bartlett's Bisection theorem and
 	      appears to be bogus)
              Added test for odd Order.		

    rev 3:    March 23, 1993
              added minimum order calculation routine

****************************************************
        cheb.c
        
this program calculates the component values for
nth odd order chebychev low pass filters.
No restrictions are placed on the values and
the program does not check for standard values.
Although the output is in a form most appropriate for
filters that end up with L's in uH and C's in pF.

The formulas used in this program and the examples
used to test it came from:

	Electronic Filter Design Handbook, 2nd ed
        Williams and Taylor, p 241ff.

Required inputs are:

        Fc:            3db frequency, in MHz
        Rs:            source impedance, in ohms
        Rl:            load impedance, in ohms
        n:             order of desired filter, must be odd

If cutoff frequency and stopband edge are known then the program
will find the minimum odd order required.

The outputs are:

        component value, inductor or capacitor indicated by
        units printed. (and order if program calculates it)

--------------------------------------------------------------
***********************************************************/

#include <stdio.h>
#include <math.h>

    double A[25];                /* dummy coef register            */
    double B[25];                /* dummy coef register            */
    double G[25];                /* dummy coef register            */
    int M; /* minimum order required, always odd */
    
main() 
{
    double pi = 3.141592654; 
    double sin();                /* need sine to calculate coef     */
    double pow();  
    int k = 0;                    /* counter for coeficients        */
    double Fc ;            /* 3dB cutoff frequency    hertz        */
    double Wc ;            /* 3dB cutoff frequency    radians        */
    double Fs ;            /* 3dB stopband frequency    hertz        */
    double Ws ;            /* 3dB stopband frequency    radians        */
    double Rs = 50.0;            /* source and load impedance    */
    double Rl = 50.0;            /* source and load impedance    */
    double Ripple ;  /* ripple allowed in passband, in dB */
    int n = 3;                    /* order of filter                */
    double meg = 1.0e6;
    double cosh();	/* hyperbolic functions needed */
    double tanh();
    double sinh();
    double sqrt();
    double ceil();  /*ceiling function */
    double beta;
    double Y,At;
    double log();
    double Acosh();
    double Asinh();
    double A1;
    int isodd();
    int order();
    double epsilon;

/* input user chosen variables first */

    printf("\nCHEBYCHEV ODD ORDER LPF PROGRAM\n");
    printf("This version assumes capacitive input and output\n");

    printf("Input allowed passband ripple in dB:  ");
    scanf("%lf",&Ripple);

    printf("Input cutoff frequency (MHz):  ");
    scanf("%lf",&Fc);
    Wc = 2.0*pi*Fc*meg;
    Fc = Fc*meg;

    printf("\nInput the order of the filter (25 max)\n");
    printf("Input a 1 if you want program to find order:  ");
    scanf("%d", &n);

    if( n != 1) {
    if( isodd(n) == 0 ){
	printf("\nOrder not odd!\n");
	exit();
	}
    if(n>25){
	printf("\nOrder too large!\n");
	exit();
	}
               }

    if( n == 1 ) {
    printf("\nInput stopband edge frequency (MHz):  ");
    scanf("%lf",&Fs);
    Ws = 2.0*pi*Fs*meg;

    printf("\nInput stopband edge attenuation (dB):  ");
    scanf("%lf",&At);

    n = order(Wc,Ws,At,Ripple);

}
printf("Odd order of %d chosen\n", n);

    init(n);                /* initialize arrays to all zeroes */

    printf("Input source resistance (ohms):  ");
    scanf("%lf",&Rs);

    printf("Input load resistance (ohms):  ");
    scanf("%lf",&Rl);
/* 
     calculate various factors needed in later equations 
*/
    epsilon = sqrt( ( pow(10.0, Ripple/10.0) -1.0 )) ;
    beta = log( 1.0/(tanh(Ripple/17.37)));
    Y = sinh(beta/(2.0*n));
    A1 = Asinh(1.0/epsilon)/n;

/*
   calculate normalized filter coeficients and store them in
   the array g[]
*/
    A[1] = sin ( (pi/(2.0*n)));  /* find A1 */
    B[1] = Y*Y + sin(pi/n)*(sin(pi/n));

    for(k = 2; k <= n; k++)
        {	
	A[k] = sin((2.0*k-1.0)*pi/(2.0*n));
	B[k] = Y*Y + sin(k*pi/n)*(sin(k*pi/n));
	}

    G[1] = (2.0*A[1]*cosh(A1))/Y;

    for(k = 2; k <= n; k++)
        {
	G[k] = (4.0*A[k-1]*A[k]*cosh(A1)*cosh(A1))/(B[k-1]*G[k-1]);
	}

   output(n,Fc,Rs,Rl);

}

/*  initialize array variables  */

init(n)
int n;
{
    int i;
    for (i=0; i<=n; i++)
        {
        G[i] = 0.0;
        A[i] = 0.0;
        B[i] = 0.0;
        }
 }

/* output():  prints out the results of the program */

output(n,f,Ri,Ro)

double f;
double Ri;
double Ro;
int n;
{
    int k;
    int m;
    double w;
    w = 2.0*3.141592654*f;
   printf("\nNormalized Values for %d order CHEBYCHEV filter @ %.2e Hz\n",n,f);
   printf("\nNormalized Values assume Equal Termination\n");

    m = (int) ceil( ( (double) n )/2.0 ) ;
    
    for (k=1; k<=n; k++)
        {
            printf("  #: %d        %.4lf\n",k,G[k]);
        }

    k = 2;
   while (k <= n)
    {
 if (k < m) 
 printf("Inductor L%d = %.4lf uH \n",k,1.0e6*G[k]*Ri/w );
 if (k == m)
printf("Inductor L%d = %.4lf uH\n",k,1.0e6*((G[k]*Ri/(2.0*w)+G[k]*Ro/(2.0*w))));
 if (k > m) 
 printf("Inductor L%d = %.4lf uH \n",k,1.0e6*G[k]*Ro/w );
 k = k+2;
    }

    k = 1;
   while (k <= n)
    {
 if (k < m) 
     printf("Capacitor C%d = %.2lf pF \n",k,1.0e12*G[k]/(Ri*w) );
 if (k == m)
printf("Capacitor C%d = %.2lf pF \n",k,1.0e12*((G[k]/(2.0*Ri*w))+(G[k]/(2.0*Ri*w))));
 if (k > m) 
     printf("Capacitor C%d = %.2lf pF \n",k,1.0e12*G[k]/(Ro*w) );
     k = k+2;
    }
}    

int isodd(p)
int p;
{
return(p % 2);
/* out = 0 if not odd */
}

int order(wc,ws,a,r) 
    double wc,ws; 
    double a,r;
{
    double ep2;
    double N;
    double Acosh();
    double pow();
    double sqrt();
    double ceil();
    int M;
	ep2 =  ( pow(10.0, r/10.0 ) ) - 1.0 ;
	N  = ( Acosh(sqrt(((pow(10.0,a/10.0)) -1.0)/ep2)))/Acosh(ws/wc);
	M = (int) ceil( N );
	if ( M % 2 == 0 ) M = M+1;
        return(M);
}

double Asinh(x)
double x;
{
	double log();
	double sqrt();
	return(log(x + sqrt(x*x+1.0)));
}

double Acosh(x)
double x;
{
 	double log();
	double sqrt();
	return(log(x + sqrt(x*x-1.0)));
}
