| |
PA0RBC > DSP 15.12.96 21:59l 380 Lines 10205 Bytes #-10355 (0) @ EU
BID : 21471_PI8SAT
Read: DH9MAG DG3MJV DG4MLA DJ2TK DL1EEC DK5RAS GUEST
Subj: C-Calc. prgm DSK320C50
Path: DB0AAB<DB0KCP<DB0ULM<DB0LX<DB0RBS<DB0SEL<DB0ZDF<DB0AIS<DB0NDK<DB0ACH<
DB0RWI<DB0PKE<PI8DRS<PI8DAZ<PI8GCB<PI8WFL<PI8SAT
Sent: 961215/1647Z @:PI8SAT.#NH1.NLD.EU #:21471 [Den Helder] $:21471_PI8SAT
From: PA0RBC@PI8SAT.#NH1.NLD.EU
To : DSP@EU
Dr OM's,
Who can compile and send me back (in 7plus) this program typed in
by myself and written by John Wiseman KE3QG in the magazine
QEX july 1996 page 3 etc.
Excuses for type-errors.
Many thanks in advance fm Roeland, PA0RBC@PI8SAT.#NH1.NLD.EU
/*****************************************************************/
/* */
/* This C program calculates FIR filter coefficients for lowpass */
/* and bandpass configurations. Highpass filters may be designed*/
/* by modifying the function Coefficient() as shown. Floating */
/* point values my be used directly on the C30 DSP, and a */
/* conversion routine is given to allow hex output for the C50. */
/* This program compiles on a Sun workstation, and previous */
/* versions have been compiled on a PC with a Borland C compiler */
/* File names may have to be changed for your machine... */
/* */
/*****************************************************************/
#include <stdio.h>
#include <std1ib.h>
#include <math.h>
#define MAX 4096
void Conversion(void);
void Coefficeint(void);
void Window(void);
void Symmetric(void);
void FreqResponse(void);
void FreqRespons_16bit(void);
void ArrayMultiply(double WindoWt);
int m, 1, n2, i, j, TapNumber, FreqRespPoints, WindowType;
double q, q1, q2, am;
double pi;
double LPCoeffs[MAX], HPCoeffs[MAX], BPCoeffs[MAX], Coeffs_16bit[MAX];
/*******************************************************************/
/* */
/* The function Coefficient() calculates the tap values for */
/* specified filter parameters. The Fourier series method is */
/* used in this function. */
/* */
/*******************************************************************/
void Coefficient(void)
{
float LowerFC, UpperFC;
printf("\nls FIR Bandpass Filter Design With Windows\n";
printf("\nEnter: # of Taps, LowerFC, UpperFC, Window Type, # of Freq
Response Points\n\n");
scanf("%d%f%f%f%d%d", &TapNumber, &LowerFC, &UpperFC, &WindowType,
&FreqRespPoints);
n2=TapNumber/2;
m=(TapNumber+1)/2;
l=TabNumber-1;
pi=3.141592654;
am=(TapNumber+1.0)/2.0;
/* Calculate the FIR filter coefficients */
if (m==am)
{
LPCoeffs[m-1]=2.0*UpperFC;
HPCoeffs[m-1]=((2.0*LowerFC)+((0.5-(2.0*LowerFC))*2.0));
BPCoeffs[m-1]=(LPCoeffs[m-1]+HPCoeffs[m-1])-1.0;
}
for (j=0; j<n2; ++j)
{
q=pi*((j+1)-am);
LPCoefs[j]=(sin(UpperFC*2.0*q))/q;
HPCoefs[j]=(-(sin(LowerFC*2.0*q))/q);
/* HPCoeffs[j]=(pow(-1.0,j))*(-(sin(LowerFC*2.0*q))/q); */
/* Uncommenting this line and commenting out the previous one */
BPCoeffs[j]=LPCoeffs[j]+HPCoeffs[j];
}
return;
}
/*****************************************************************/
/* */
/* The function Window() calculates the window values for */
/* triangular, Hanning, Hamming, Blackman, or rectangular */
/* (no window - default) windows. */
/* */
/*****************************************************************/
void Window(void)
{
double WindowWt;
qp=pi/am;
q1=pi/(am-1.0);
switch (WindowType)
{
/* Multiply by a triangular window */
case 1:
for (j=0; j<am; ++j)
{
WindowWt=(j+1)/am;
ArrayMultiply(WindowWt);
}
break;
/* Multiply by a Hanning window */
case 2:
for (j=0; j<am; ++j)
{
WindowWt=0.5-(0.5*cos((j+1)*q));
ArrayMultiply(WindowWt);
}
break;
/* Multiply by a Hamming window */
case 3:
for (j=0; j<am; ++j)
{
WindowWt=0.54-(0.46*cos(j*q1));
ArrayMultiply(WindowWt);
break
*/ Multiply by a Blackman window */
case 4:
for (j=0; j<am; ++j)
{
WindowWt=0.42-(0.5*cos((j+1)*q))+(0.08*cos(2*(j+1)*q));
ArrayMultiply(WindowWt);
}
break;
}
return;
}
/**************************************************************/
/* */
/* The function Symmetric() calculates the symmetric portion */
/* of the filter coefficients, then writes these floating */
/* point values to the file "bpfilter.txt". These same values */
/* are written to the file "bpfilter.tfp", but each value is */
/* preceeded with the assembler directive .word for easy */
/* inclusion in the C30 assembly language code. */
/* */
/**************************************************************/
void Symmetric(void)
{
FILE *f;
FILE *ff;
f=fopen("bpfilter.ftp", "w");
ff=fopen("bpfilter.txt", "w");
/* Equate symmetric coefficients */
for (j=0; j<((Tapnumber-1)/2); 1-=2, ++j)
{
printf("%18f\n",BPCoeffs[j]);
fprintf(f, " .word %15.10f\n",BPCoeffs[j]);
fprintf(ff, "%15.10f\n", BPcoeffs[j]);
}
fclose(f);
fclose(ff);
return;
}
/**************************************************************/
/* */
/* The function FreqResponse() calculates the frequency */
/* response of the calculated FIR filter coefficients for */
/* response in dB. These values are printed on the screen, */
/* and c[j] is written to the file "freqresp.txt" for usage */
/* with plot programs. */
/* */
/**************************************************************/
void FreqResponse(void)
{
double bt;
double a[MAX], c[MAX];
FILE *f;
f=fopen("freqresp.txt", "w");
/* Calculate the frequency response */
printf("\nThe calculated frequency response */
printf("\nThe calculated frequency */
printf("\nThe calculated frequency response is:\n\n");
printf(" Magnitude Decibels\n\n");
q2=pi/FreqResPoints;
for (j=0; j<FreqRespPoints; ++j)
{
if (am==m)
(
bt=0.5*BPCoeffs(m-1);
)
for 9i=0; i<n2; ++i)
(
bt=bt+(BPCoeffs(i)*cos(q2*(am-(i+1))*j));
)
a(j)=2*bt;
/* Calculate frequency response in decibels */
c(j)=20.0*(log10(fabs(a(j))));
printf("%18f%25f\n", a[j], c[j]);
fprintf(f, "%10f\n", c[j];
}
fclose(f);
return;
{
/****************************************************************/
/* */
/* The function FreqResponse_16bit() calculates the FIR filter */
/* filter coefficients to 16-bit hex values (2's complement). */
/* These values are rounded to the nearest bit , then stored */
/* in the file "bpfilter.hex". Each coefficient is preceeded */
/* by the assembler directive .word for easy usage in C50 */
/* assembly code. */
/* */
/****************************************************************/
void Conversion(void)
{
long int hexnum;
FILE *ff;
ff=fopen("bpfilter.hex", "w");
for (j=0; j<TapNumber; ++j)
{
hexnum=(long int)(BPCoeffs[j]*(float)524288.0);
if ((hexnum & 15) >= 8)
{
hexnum += 16;
}
hexnum=((hexnum << 12) & 4294901760);
hexnum=(hexnum >> 16) &65535;
printf("%04x\n", hexnum);
fprintf(ff, " .word %04xh\n", hexnum);
Coeffs_16bit[j] = hexnum;
}
fclose(ff);
return;
}
/*****************************************************************/
/* */
/* The function FreqResponse_16bit() calculates the FIR filter */
/* frequency response for 16-bit resolution coefficients, for */
/* comparison with the "ideal" double precision floating point */
/* coefficients calculated elswhere. */
/* */
/*****************************************************************/
void FreqResponse_16bit(void)
{
double bt;
double a[MAX], c[MAX];
FILE *f;
f=fopen("freqresp_16bit.txt", "w");
/* Calculate the frequency response */
printf("\nThe calculated frequency response is:\n\n");
printf(" Magnitude Decibels\n\n");
q2=pi/FreqRespPoints;
for (j=0; j<FreqRespPoints; ++j)
{
if(Coeffs_16bit[j] >= 32768)
{
Coeffs_16bit[j] = (double)(Coeffs_16bit[j]/(float)32768.0);
Coeffs_16bit[j] -=2.0;
}
else
{
Coeffs_16bit[j] = (double)(Coeffs_16bit[j]/(float)32768.0);
}
}
for (j=0; j<FreqRespPoints; ++j)
{
if (am==m)
{
bt=0.5*Coeffs_16bit[m-1];
}
for (i=0; i<n2; ++i)
{
bt=0.5*Coeffs_16bit(m-1);
}
for (i=0; i<n2; ++i)
{
bt=bt+(Coeffs_16bit[i]*cos(q2*(am-(i+1))*j));
}
a[j]=2*bt;
/* Calculate frequency response in decibels */
C[j]=20.0*(log10(fabs(a[j])));
printf("%18f%25f\n", a[j],c[j]);
fprintf(f, "%10f\n", c[j]);
}
fclose(f);
return;
}
/**************************************************************/ */
/* */
/* The function ArrayMultiply() multiplies the FIR filter */
/* coefficients with the appropriate window values on a pair- */
/* wise basis. */
/* */
/**************************************************************/
void ArrayMultiply(double WindowWt)
{
BPCoeffs[j] *= WindowWt;
return;
}
/**************************************************************/
main()
{
Coefficient();
Window();
Symmetric();
FreqResponse();
Conversion();
FreqResponse_16bit();
return(0);
}
/***************************************************************/
Read previous mail | Read next mail
| |
