commit dd0f54d3bc783eda8af59155f995f333bf2c5601 Author: Jens Ahrensfeld Date: Sat Jul 19 07:44:42 2014 +0000 Initial import git-svn-id: http://moon:8086/svn/software/trunk/libsrc/iir@1 b431acfa-c32f-4a4a-93f1-934dc6c82436 diff --git a/Audio-EQ-Cookbook.txt b/Audio-EQ-Cookbook.txt new file mode 100755 index 0000000..9ecb0c2 --- /dev/null +++ b/Audio-EQ-Cookbook.txt @@ -0,0 +1,280 @@ + + Cookbook formulae for audio EQ biquad filter coefficients +---------------------------------------------------------------------------- + by Robert Bristow-Johnson + + +All filter transfer functions were derived from analog prototypes (that +are shown below for each EQ filter type) and had been digitized using the +Bilinear Transform. BLT frequency warping has been taken into account for +both significant frequency relocation (this is the normal "prewarping" that +is necessary when using the BLT) and for bandwidth readjustment (since the +bandwidth is compressed when mapped from analog to digital using the BLT). + +First, given a biquad transfer function defined as: + + b0 + b1*z^-1 + b2*z^-2 + H(z) = ------------------------ (Eq 1) + a0 + a1*z^-1 + a2*z^-2 + +This shows 6 coefficients instead of 5 so, depending on your architechture, +you will likely normalize a0 to be 1 and perhaps also b0 to 1 (and collect +that into an overall gain coefficient). Then your transfer function would +look like: + + (b0/a0) + (b1/a0)*z^-1 + (b2/a0)*z^-2 + H(z) = --------------------------------------- (Eq 2) + 1 + (a1/a0)*z^-1 + (a2/a0)*z^-2 + +or + + 1 + (b1/b0)*z^-1 + (b2/b0)*z^-2 + H(z) = (b0/a0) * --------------------------------- (Eq 3) + 1 + (a1/a0)*z^-1 + (a2/a0)*z^-2 + + +The most straight forward implementation would be the "Direct Form 1" +(Eq 2): + + y[n] = (b0/a0)*x[n] + (b1/a0)*x[n-1] + (b2/a0)*x[n-2] + - (a1/a0)*y[n-1] - (a2/a0)*y[n-2] (Eq 4) + +This is probably both the best and the easiest method to implement in the +56K and other fixed-point or floating-point architechtures with a double +wide accumulator. + + + +Begin with these user defined parameters: + + Fs (the sampling frequency) + + f0 ("wherever it's happenin', man." Center Frequency or + Corner Frequency, or shelf midpoint frequency, depending + on which filter type. The "significant frequency".) + + dBgain (used only for peaking and shelving filters) + + Q (the EE kind of definition, except for peakingEQ in which A*Q is + the classic EE Q. That adjustment in definition was made so that + a boost of N dB followed by a cut of N dB for identical Q and + f0/Fs results in a precisely flat unity gain filter or "wire".) + + _or_ BW, the bandwidth in octaves (between -3 dB frequencies for BPF + and notch or between midpoint (dBgain/2) gain frequencies for + peaking EQ) + + _or_ S, a "shelf slope" parameter (for shelving EQ only). When S = 1, + the shelf slope is as steep as it can be and remain monotonically + increasing or decreasing gain with frequency. The shelf slope, in + dB/octave, remains proportional to S for all other values for a + fixed f0/Fs and dBgain. + + + +Then compute a few intermediate variables: + + A = sqrt( 10^(dBgain/20) ) + = 10^(dBgain/40) (for peaking and shelving EQ filters only) + + w0 = 2*pi*f0/Fs + + cos(w0) + sin(w0) + + alpha = sin(w0)/(2*Q) (case: Q) + = sin(w0)*sinh( ln(2)/2 * BW * w0/sin(w0) ) (case: BW) + = sin(w0)/2 * sqrt( (A + 1/A)*(1/S - 1) + 2 ) (case: S) + + FYI: The relationship between bandwidth and Q is + 1/Q = 2*sinh(ln(2)/2*BW*w0/sin(w0)) (digital filter w BLT) + or 1/Q = 2*sinh(ln(2)/2*BW) (analog filter prototype) + + The relationship between shelf slope and Q is + 1/Q = sqrt((A + 1/A)*(1/S - 1) + 2) + + 2*sqrt(A)*alpha = sin(w0) * sqrt( (A^2 + 1)*(1/S - 1) + 2*A ) + is a handy intermediate variable for shelving EQ filters. + + +Finally, compute the coefficients for whichever filter type you want: + (The analog prototypes, H(s), are shown for each filter + type for normalized frequency.) + + +LPF: H(s) = 1 / (s^2 + s/Q + 1) + + b0 = (1 - cos(w0))/2 + b1 = 1 - cos(w0) + b2 = (1 - cos(w0))/2 + a0 = 1 + alpha + a1 = -2*cos(w0) + a2 = 1 - alpha + + + +HPF: H(s) = s^2 / (s^2 + s/Q + 1) + + b0 = (1 + cos(w0))/2 + b1 = -(1 + cos(w0)) + b2 = (1 + cos(w0))/2 + a0 = 1 + alpha + a1 = -2*cos(w0) + a2 = 1 - alpha + + + +BPF: H(s) = s / (s^2 + s/Q + 1) (constant skirt gain, peak gain = Q) + + b0 = sin(w0)/2 = Q*alpha + b1 = 0 + b2 = -sin(w0)/2 = -Q*alpha + a0 = 1 + alpha + a1 = -2*cos(w0) + a2 = 1 - alpha + + +BPF: H(s) = (s/Q) / (s^2 + s/Q + 1) (constant 0 dB peak gain) + + b0 = alpha + b1 = 0 + b2 = -alpha + a0 = 1 + alpha + a1 = -2*cos(w0) + a2 = 1 - alpha + + + +notch: H(s) = (s^2 + 1) / (s^2 + s/Q + 1) + + b0 = 1 + b1 = -2*cos(w0) + b2 = 1 + a0 = 1 + alpha + a1 = -2*cos(w0) + a2 = 1 - alpha + + + +APF: H(s) = (s^2 - s/Q + 1) / (s^2 + s/Q + 1) + + b0 = 1 - alpha + b1 = -2*cos(w0) + b2 = 1 + alpha + a0 = 1 + alpha + a1 = -2*cos(w0) + a2 = 1 - alpha + + + +peakingEQ: H(s) = (s^2 + s*(A/Q) + 1) / (s^2 + s/(A*Q) + 1) + + b0 = 1 + alpha*A + b1 = -2*cos(w0) + b2 = 1 - alpha*A + a0 = 1 + alpha/A + a1 = -2*cos(w0) + a2 = 1 - alpha/A + + + +lowShelf: H(s) = A * (s^2 + (sqrt(A)/Q)*s + A)/(A*s^2 + (sqrt(A)/Q)*s + 1) + + b0 = A*( (A+1) - (A-1)*cos(w0) + 2*sqrt(A)*alpha ) + b1 = 2*A*( (A-1) - (A+1)*cos(w0) ) + b2 = A*( (A+1) - (A-1)*cos(w0) - 2*sqrt(A)*alpha ) + a0 = (A+1) + (A-1)*cos(w0) + 2*sqrt(A)*alpha + a1 = -2*( (A-1) + (A+1)*cos(w0) ) + a2 = (A+1) + (A-1)*cos(w0) - 2*sqrt(A)*alpha + + + +highShelf: H(s) = A * (A*s^2 + (sqrt(A)/Q)*s + 1)/(s^2 + (sqrt(A)/Q)*s + A) + + b0 = A*( (A+1) + (A-1)*cos(w0) + 2*sqrt(A)*alpha ) + b1 = -2*A*( (A-1) + (A+1)*cos(w0) ) + b2 = A*( (A+1) + (A-1)*cos(w0) - 2*sqrt(A)*alpha ) + a0 = (A+1) - (A-1)*cos(w0) + 2*sqrt(A)*alpha + a1 = 2*( (A-1) - (A+1)*cos(w0) ) + a2 = (A+1) - (A-1)*cos(w0) - 2*sqrt(A)*alpha + + + + + +FYI: The bilinear transform (with compensation for frequency warping) +substitutes: + + 1 1 - z^-1 + (normalized) s <-- ----------- * ---------- + tan(w0/2) 1 + z^-1 + + and makes use of these trig identities: + + sin(w0) 1 - cos(w0) + tan(w0/2) = ------------- (tan(w0/2))^2 = ------------- + 1 + cos(w0) 1 + cos(w0) + + + resulting in these substitutions: + + + 1 + cos(w0) 1 + 2*z^-1 + z^-2 + 1 <-- ------------- * ------------------- + 1 + cos(w0) 1 + 2*z^-1 + z^-2 + + + 1 + cos(w0) 1 - z^-1 + s <-- ------------- * ---------- + sin(w0) 1 + z^-1 + + 1 + cos(w0) 1 - z^-2 + = ------------- * ------------------- + sin(w0) 1 + 2*z^-1 + z^-2 + + + 1 + cos(w0) 1 - 2*z^-1 + z^-2 + s^2 <-- ------------- * ------------------- + 1 - cos(w0) 1 + 2*z^-1 + z^-2 + + + The factor: + + 1 + cos(w0) + ------------------- + 1 + 2*z^-1 + z^-2 + + is common to all terms in both numerator and denominator, can be factored + out, and thus be left out in the substitutions above resulting in: + + + 1 + 2*z^-1 + z^-2 + 1 <-- ------------------- + 1 + cos(w0) + + + 1 - z^-2 + s <-- ------------------- + sin(w0) + + + 1 - 2*z^-1 + z^-2 + s^2 <-- ------------------- + 1 - cos(w0) + + + In addition, all terms, numerator and denominator, can be multiplied by a + common (sin(w0))^2 factor, finally resulting in these substitutions: + + + 1 <-- (1 + 2*z^-1 + z^-2) * (1 - cos(w0)) + + s <-- (1 - z^-2) * sin(w0) + + s^2 <-- (1 - 2*z^-1 + z^-2) * (1 + cos(w0)) + + 1 + s^2 <-- 2 * (1 - 2*cos(w0)*z^-1 + z^-2) + + + The biquad coefficient formulae above come out after a little + simplification. diff --git a/Debug/iir.obj b/Debug/iir.obj new file mode 100755 index 0000000..5bd1650 Binary files /dev/null and b/Debug/iir.obj differ diff --git a/Debug/iir.pch b/Debug/iir.pch new file mode 100755 index 0000000..c99853a Binary files /dev/null and b/Debug/iir.pch differ diff --git a/Debug/vc60.idb b/Debug/vc60.idb new file mode 100755 index 0000000..5cdb016 Binary files /dev/null and b/Debug/vc60.idb differ diff --git a/Debug/vc60.pdb b/Debug/vc60.pdb new file mode 100755 index 0000000..2841187 Binary files /dev/null and b/Debug/vc60.pdb differ diff --git a/Iir.c b/Iir.c new file mode 100755 index 0000000..b07d9a9 --- /dev/null +++ b/Iir.c @@ -0,0 +1,301 @@ +/*************************************************************************/ +/* iir.c */ +/*************************************************************************/ +#include "stdio.h" +#include "stdlib.h" +#include "math.h" +#include "iir.h" + +/*************************************************************************/ +/* Global Variables */ +/*************************************************************************/ +const char *filterTypeString[] = +{ + "Unknown filter type", + "Butterworth-Lowpass", + "Butterworth-Highpass", + "Butterworth-Bandpass", + "Butterworth-Bandstop", + "Peaking-EQ", + "Low Shelving-EQ", + "High Shelving-EQ" +}; + +/******************************************************************************/ +void IIRCalcFilterCoeff(struct _sIIRCoeff *pCoeff, iir_float_t fa, iir_float_t fg, iir_float_t q, unsigned order, unsigned filterType) +{ + unsigned p; + iir_float_t qp; + +// IIRInit(pCoeff, order); + + for(p=0; p < order/2;p++) + { + qp = q * IIRCalcQp(p+1, order); + IIRCalcPartFilterCoeff2(&pCoeff[p], 1.0, fa, fg, qp, filterType); + } +} + +int IIRCalcPartFilterCoeff1(struct _sIIRCoeff *pCoeff, iir_float_t fa, iir_float_t fg, iir_float_t Qi, unsigned filterType) +{ + iir_float_t K, a0; + iir_float_t alpha, omega, ks, kc; + unsigned error; + + omega = (iir_float_t)(2*iir_pi*fg/fa); + ks = (iir_float_t)sin(omega); + kc = (iir_float_t)cos(omega); + alpha = 0.5f*ks /Qi; + + K = IIRBilTrans(fg, fa); + a0 = K/Qi + 1; + + switch(filterType) + { + case IIR_FILTERTYPE_LOWPASS: + + pCoeff->ak0 = 1.0f; + pCoeff->ak1 = (1 - K/Qi)/a0; + pCoeff->ak2 = 0.0; + + pCoeff->bk0 = 1.0f/a0; + pCoeff->bk1 = 1.0f/a0; + pCoeff->bk2 = 0.0; + + break; + + case IIR_FILTERTYPE_HIGHPASS: + + pCoeff->ak0 = 1.0f; + pCoeff->ak1 = (1 - K/Qi) /a0; + pCoeff->ak2 = 0.0; + + pCoeff->bk0 = 1.0f*K /a0; + pCoeff->bk1 = -1.0f*K /a0; + pCoeff->bk2 = 0.0; + + break; + + default: + error = -1; + break; + } + return error; +} + +int IIRCalcPartFilterCoeff2(struct _sIIRCoeff *pCoeff, iir_float_t A, iir_float_t fa, iir_float_t fg, iir_float_t qp, unsigned filterType) +{ + iir_float_t a0; + iir_float_t alpha, omega, ks, kc; + unsigned error; + + omega = (iir_float_t)(2*iir_pi*fg/fa); + ks = (iir_float_t)sin(omega); + kc = (iir_float_t)cos(omega); + alpha = 0.5f*ks /qp; + + error = 0; + switch(filterType) + { + case IIR_FILTERTYPE_LOWPASS: + + a0 = 1 + alpha; + pCoeff->ak0 = 1.0f; + pCoeff->ak1 = -2.0f*kc /a0; + pCoeff->ak2 = (1 - alpha) /a0; + + pCoeff->bk0 = 0.5f*(1 - kc) /a0; + pCoeff->bk1 = (1 - kc) /a0; + pCoeff->bk2 = 0.5f*(1 - kc) /a0; + + break; + + case IIR_FILTERTYPE_HIGHPASS: + + a0 = 1 + alpha; + pCoeff->ak0 = 1.0f; + pCoeff->ak1 = -2.0f*kc /a0; + pCoeff->ak2 = (1 - alpha) /a0; + + pCoeff->bk0 = 0.5f*(1 + kc) /a0; + pCoeff->bk1 = -(1 + kc) /a0; + pCoeff->bk2 = 0.5f*(1 + kc) /a0; + + break; + + case IIR_FILTERTYPE_BANDPASS: + + a0 = 1 + alpha; + pCoeff->ak0 = 1.0f; + pCoeff->ak1 = -2.0f*kc /a0; + pCoeff->ak2 = (1 - alpha) /a0; + + pCoeff->bk0 = alpha /a0; + pCoeff->bk1 = 0; + pCoeff->bk2 = -alpha /a0; + + break; + + case IIR_FILTERTYPE_BANDSTOP: + + a0 = 1 + alpha; + pCoeff->ak0 = 1.0f; + pCoeff->ak1 = -2.0f*kc /a0; + pCoeff->ak2 = (1 - alpha) /a0; + + pCoeff->bk0 = 1.0f /a0; + pCoeff->bk1 = -2.0f*kc /a0; + pCoeff->bk2 = 1.0f /a0; + + break; + + case IIR_FILTERTYPE_PEAKING: + + a0 = 1 + (alpha/A); + pCoeff->ak0 = 1.0f; + pCoeff->ak1 = -2.0f*kc /a0; + pCoeff->ak2 = (1 - (alpha/A)) /a0; + + pCoeff->bk0 = (1 + (alpha*A)) /a0; + pCoeff->bk1 = -2.0f*kc /a0; + pCoeff->bk2 = (1 - (alpha*A)) /a0; + + break; + + default: + error = -1; + break; + } + return error; +} + +void IIR(struct _sIIRCoeff *pCoeff, iir_float_t *xn, iir_float_t *yn, unsigned order, unsigned numPoints) +{ + iir_float_t xp, yp; + unsigned i, p; + unsigned numSec = order/2; + + for (i=0; iorder = order; + pObj->pX = (iir_float_t*)malloc((order+1)*sizeof(iir_float_t)); + pObj->pY = (iir_float_t*)malloc((order+1)*sizeof(iir_float_t)); + memset(pObj->pX, 0, (order+1)*sizeof(iir_float_t)); + memset(pObj->pY, 0, (order+1)*sizeof(iir_float_t)); +} + +void IIR_lin_free(iir_lin_t *pObj) +{ + if (pObj->pX) + free(pObj->pX); + + if (pObj->pY) + free(pObj->pY); + + pObj->order = 0; + +} + +iir_float_t IIR_lin_process(iir_lin_t *pObj, iir_float_t *pB, iir_float_t *pA, iir_float_t x) +{ + unsigned i; + iir_float_t y; + + if (!pObj->order) + return 0; + + for (i=pObj->order; i >= 1; i--) + pObj->pX[i] = pObj->pX[i-1]; + + for (i=pObj->order; i >= 1; i--) + pObj->pY[i] = pObj->pY[i-1]; + + pObj->pX[0] = x; + + y = 0; + for (i=0; i <= pObj->order; i++) + y += pObj->pX[i]*pB[i]; + + for (i=1; i <= pObj->order; i++) + y -= pObj->pY[i]*pA[i]; + + pObj->pY[0] = y; + + return y; + +} + diff --git a/Iir.c.double b/Iir.c.double new file mode 100755 index 0000000..129e24d --- /dev/null +++ b/Iir.c.double @@ -0,0 +1,225 @@ +/*************************************************************************/ +/* iir.c +/*************************************************************************/ +#include "stdio.h" +#include "math.h" +#include "iir.h" + +/*************************************************************************/ +/* Global Variables +/*************************************************************************/ +const char *filterTypeString[] = +{ + "Unknown filter type", + "Butterworth-Lowpass", + "Butterworth-Highpass", + "Butterworth-Bandpass", + "Butterworth-Bandstop", + "Peaking-EQ", + "Low Shelving-EQ", + "High Shelving-EQ" +}; + +/******************************************************************************/ +void IIRCalcFilterCoeff(struct _sIIRCoeff *pCoeff, double fa, double fg, double q, unsigned order, unsigned filterType) +{ + unsigned p; + double qp; + + IIRInit(pCoeff, order); + + for(p=0; p < order/2;p++) + { + qp = q * IIRCalcQp(p+1, order); + IIRCalcPartFilterCoeff2(&pCoeff[p], 1.0, fa, fg, qp, filterType); + } +} + +int IIRCalcPartFilterCoeff1(struct _sIIRCoeff *pCoeff, double fa, double fg, double Qi, unsigned filterType) +{ + double K, a0; + double alpha, omega, ks, kc; + unsigned error; + + omega = 2*pi*fg/fa; + ks = sin(omega); + kc = cos(omega); + alpha = 0.5*ks /Qi; + + K = IIRBilTrans(fg, fa); + a0 = K/Qi + 1; + + switch(filterType) + { + case IIR_FILTERTYPE_LOWPASS: + + pCoeff->ak0 = 1.0; + pCoeff->ak1 = (1 - K/Qi)/a0; + pCoeff->ak2 = 0.0; + + pCoeff->bk0 = 1.0/a0; + pCoeff->bk1 = 1.0/a0; + pCoeff->bk2 = 0.0; + + break; + + case IIR_FILTERTYPE_HIGHPASS: + + pCoeff->ak0 = 1.0; + pCoeff->ak1 = (1 - K/Qi) /a0; + pCoeff->ak2 = 0.0; + + pCoeff->bk0 = 1.0*K /a0; + pCoeff->bk1 = -1.0*K /a0; + pCoeff->bk2 = 0.0; + + break; + + default: + error = -1; + break; + } + return error; +} + +int IIRCalcPartFilterCoeff2(struct _sIIRCoeff *pCoeff, double A, double fa, double fg, double qp, unsigned filterType) +{ + double a0; + double alpha, omega, ks, kc; + unsigned error; + + omega = 2*pi*fg/fa; + ks = sin(omega); + kc = cos(omega); + alpha = 0.5*ks /qp; + + error = 0; + switch(filterType) + { + case IIR_FILTERTYPE_LOWPASS: + + a0 = 1 + alpha; + pCoeff->ak0 = 1.0; + pCoeff->ak1 = -2.0*kc /a0; + pCoeff->ak2 = (1 - alpha) /a0; + + pCoeff->bk0 = 0.5*(1 - kc) /a0; + pCoeff->bk1 = (1 - kc) /a0; + pCoeff->bk2 = 0.5*(1 - kc) /a0; + + break; + + case IIR_FILTERTYPE_HIGHPASS: + + a0 = 1 + alpha; + pCoeff->ak0 = 1.0; + pCoeff->ak1 = -2.0*kc /a0; + pCoeff->ak2 = (1 - alpha) /a0; + + pCoeff->bk0 = 0.5*(1 + kc) /a0; + pCoeff->bk1 = -(1 + kc) /a0; + pCoeff->bk2 = 0.5*(1 + kc) /a0; + + break; + + case IIR_FILTERTYPE_BANDPASS: + + a0 = 1 + alpha; + pCoeff->ak0 = 1.0; + pCoeff->ak1 = -2.0*kc /a0; + pCoeff->ak2 = (1 - alpha) /a0; + + pCoeff->bk0 = alpha /a0; + pCoeff->bk1 = 0; + pCoeff->bk2 = -alpha /a0; + + break; + + case IIR_FILTERTYPE_BANDSTOP: + + a0 = 1 + alpha; + pCoeff->ak0 = 1.0; + pCoeff->ak1 = -2.0*kc /a0; + pCoeff->ak2 = (1 - alpha) /a0; + + pCoeff->bk0 = 1.0 /a0; + pCoeff->bk1 = -2.0*kc /a0; + pCoeff->bk2 = 1.0 /a0; + + break; + + case IIR_FILTERTYPE_PEAKING: + + a0 = 1 + (alpha/A); + pCoeff->ak0 = 1.0; + pCoeff->ak1 = -2.0*kc /a0; + pCoeff->ak2 = (1 - (alpha/A)) /a0; + + pCoeff->bk0 = (1 + (alpha*A)) /a0; + pCoeff->bk1 = -2.0*kc /a0; + pCoeff->bk2 = (1 - (alpha*A)) /a0; + + break; + + default: + error = -1; + break; + } + return error; +} + +void IIR(struct _sIIRCoeff *pCoeff, double *xn, double *yn, unsigned order, unsigned numPoints) +{ + double xp, yp; + unsigned i, p; + unsigned numSec = order/2; + + for (i=0; i +#include +#include +#include +#include +#include "iir.h" + +/******************************************************************************/ +const char *filterTypeString[] = +{ + "Unknown filter type", + "Butterworth-Lowpass", + "Butterworth-Highpass", + "Butterworth-Bandpass", + "Butterworth-Bandstop" +}; + +/******************************************************************************/ +int IIRCalcPartFilterCoeff1(CIIRCoeff *pCoeff, double fg, double fa, double Qi, unsigned filterType) +{ + double K, A0; + unsigned error; + + error = 0; + switch(filterType) + { + case IIR_FILTERTYPE_LOWPASS: + + K = IIRBilTrans(fg, fa); + A0 = 1.0 /(K/Qi + 1); + + pCoeff->m_pak[0] = 1.0; + pCoeff->m_pak[1] = (1 - K/Qi) * A0; + + pCoeff->m_pbk[0] = 1.0;// * A0; + pCoeff->m_pbk[1] = 1.0;// * A0; + pCoeff->m_aScale = A0; + break; + + case IIR_FILTERTYPE_HIGHPASS: + + K = IIRBilTrans(fg, fa); + A0 = 1.0 /(K/Qi + 1); + + pCoeff->m_pak[0] = 1.0; + pCoeff->m_pak[1] = (1 - K/Qi) * A0; + + pCoeff->m_pbk[0] = 1.0*K;// * A0; + pCoeff->m_pbk[1] = -1.0*K;// * A0; + pCoeff->m_aScale = A0; + break; + + default: + error = -1; + break; + } + return error; +} + +int IIRCalcPartFilterCoeff2(CIIRCoeff *pCoeff, double fg, double fa, double Qi, unsigned filterType) +{ + double K, KK, A0,B0; + double alpha, omega, sn, cs; + unsigned error; + + K = IIRBilTrans(fg, fa); + KK = K*K; + A0 = 1.0 /(1 + KK + K/Qi); + + omega = 2*pi*fg/fa; + sn = sin(omega); + cs = cos(omega); + alpha = 0.5*sn /Qi; + + error = 0; + switch(filterType) + { + case IIR_FILTERTYPE_LOWPASS: + + pCoeff->m_pak[0] = 1.0 ; + pCoeff->m_pak[1] = 2 *(1 - KK) * A0; + pCoeff->m_pak[2] = (1 + KK - K/Qi) * A0; + + pCoeff->m_pbk[0] = 1.0*A0; + pCoeff->m_pbk[1] = 2.0*A0; + pCoeff->m_pbk[2] = 1.0*A0; + + pCoeff->m_aScale = A0; + + break; + + case IIR_FILTERTYPE_HIGHPASS: + + pCoeff->m_pak[0] = 1.0; + pCoeff->m_pak[1] = 2 *(1 - KK) * A0; + pCoeff->m_pak[2] = (KK - K/Qi + 1) * A0; + + pCoeff->m_pbk[0] = 1.0*KK * A0; + pCoeff->m_pbk[1] = -2.0*KK * A0; + pCoeff->m_pbk[2] = 1.0*KK * A0; + pCoeff->m_aScale = A0; + break; + + default: + error = -1; + break; + } + return error; +} + +double IIRBilTrans(double fg, double fa) +{ + return 1.0/(tan(pi*fg/fa)); +} + +void IIR(double *xn, double *yn, CIIRCoeff *pCoeff, unsigned numPoints) +{ + unsigned n, k; + double y1, y2; + + for (n=0; n < numPoints; n++) + { + y1 = 0; + y2 = 0; + for (k=0; k <= pCoeff->m_Nb; k++) + { + y1 = y1 + (pCoeff->m_pbk[k] * xn[pCoeff->m_Nb-k+n]); + if (!_finite(y1)) + printf("\nException: MATH ERROR!\n"); + } +// y1 = y1 *pCoeff->m_aScale; + + for (k=1; k <= pCoeff->m_Na; k++) + { + y2 = y2 - (pCoeff->m_pak[k] * yn[pCoeff->m_Na-k+n]); + if (!_finite(y2)) + printf("\nException: MATH ERROR!\n"); + } +// y2 = y2 *pCoeff->m_bScale; + + yn[pCoeff->m_Nb+n] = (y1 + y2); + if (!_finite(yn[pCoeff->m_Nb+n])) + printf("\nException: MATH ERROR!\n"); + + } +} + +int IIRCalcFilterCoeff(double fg, double fa, double Qi, unsigned N, CIIRCoeff *pCoeff, unsigned filterType) +{ + unsigned p, order_ap, order_bp; + div_t result; + unsigned numEvenFilterParts, filterCnt; + + CIIRCoeff Temp(N, N); + CIIRCoeff Coeff1(1,1); + CIIRCoeff Coeff2(2,2); + + double Qp; + FILE *pFile; + + pFile = fopen("filter.out","w"); + fprintf(pFile,"IIR-Filter Version 1.0\n"); + fprintf(pFile,"Filter Coefficients for %u-Order-%s, Qi = %4.2f\n",N, filterTypeString[filterType],Qi); + fprintf(pFile,"fg = %9.2f Hz\nfa = %9.2f Hz\n",fg, fa); + + p = 1; + filterCnt = 1; + + switch (N) + { + case 0: + break; + + default: + order_ap = 2; + order_bp = 2; + result = div(N,2); + numEvenFilterParts = result.quot; + + if (result.rem != 0) + { + IIRCalcPartFilterCoeff1(&Coeff1, fg, fa, 1.0*Qi, filterType); + fprintf(pFile,"\n1.Partfilter Np = 1, Qp = 1.00\n"); + IIRPrintCoeff(pFile,&Coeff1, 1); + filterCnt++; + } + + while (p <= numEvenFilterParts) + { + Qp = IIRCalcQp(p, N); + if (p == 1) + { + IIRCalcPartFilterCoeff2(&Temp, fg, fa, Qp*Qi, filterType); + fprintf(pFile,"\n%u.Partfilter Np = 2, Qp = %4.2f\n", filterCnt, Qp); + IIRPrintCoeff(pFile, &Temp, 2); + filterCnt++; + p++; + if (numEvenFilterParts > 1) + continue; + memcpy(pCoeff->m_pak, Temp.m_pak, (order_ap+1)*sizeof(double)); + memcpy(pCoeff->m_pbk, Temp.m_pbk, (order_bp+1)*sizeof(double)); + continue; + } + + IIRCalcPartFilterCoeff2(&Coeff2, fg, fa, Qp*Qi, filterType); + fprintf(pFile,"\n%u.Partfilter Np = 2, Qp = %4.2f\n", filterCnt, Qp); + IIRPrintCoeff(pFile,&Coeff2, 2); + order_ap = IIRMulPolynom(Coeff2.m_pak, 2, Temp.m_pak, order_ap, pCoeff->m_pak); + order_bp = IIRMulPolynom(Coeff2.m_pbk, 2, Temp.m_pbk, order_bp, pCoeff->m_pbk); + memcpy(Temp.m_pak, pCoeff->m_pak, (order_ap+1)*sizeof(double)); + memcpy(Temp.m_pbk, pCoeff->m_pbk, (order_bp+1)*sizeof(double)); + filterCnt++; + p++; + } + + if (result.rem != 0) + { + if (result.quot == 0) + { + memcpy(pCoeff->m_pak, Coeff1.m_pak, 2*sizeof(double)); + memcpy(pCoeff->m_pbk, Coeff1.m_pbk, 2*sizeof(double)); + } + else + { + order_ap = IIRMulPolynom(Coeff1.m_pak, 1, Temp.m_pak, order_ap, pCoeff->m_pak); + order_bp = IIRMulPolynom(Coeff1.m_pbk, 1, Temp.m_pbk, order_bp, pCoeff->m_pbk); + } + } + + fprintf(pFile,"\n\nResulted Filter N = %u\n",N); + IIRPrintCoeff(pFile,pCoeff, N); +// ScaleCoeff(pCoeff); + fprintf(pFile,"\n\nNormalized Filterkernel N = %u\n",N); + IIRPrintCoeff(pFile,pCoeff, N); + break; + } + + fclose(pFile); +return 0; +} + + +unsigned IIRMulPolynom(double *pA, unsigned orderA, double *pB, unsigned orderB, double *pProduct) +{ + unsigned cntA, cntB, newOrder; + + newOrder = orderA+orderB; + memset(pProduct, 0, (newOrder+1)*sizeof(double)); + + for (cntA=0; cntA <= orderA; cntA++) + { + for (cntB=0; cntB <= orderB; cntB++) + pProduct[cntA+cntB] += pA[cntA] * pB[cntB]; + } + + return newOrder; +} + +double IIRCalcQp(unsigned p, unsigned N) +{ + return 1.0/(2*sin(pi*(2*p-1)/(2*N))); +} + +void IIRPrintCoeff(FILE *pFile, CIIRCoeff *pCoeff, unsigned N) +{ + unsigned k; + + + for (k=0; k <= N; k++) + { + fprintf(pFile,"a[%2u] = %9.6g, b[%2u] = %9.6g\n",k,pCoeff->m_pak[k],k,pCoeff->m_pbk[k]); + } + fprintf(pFile,"aScale = %9.6g, bScale = %9.6g\n\n",pCoeff->m_aScale, pCoeff->m_bScale); + + fprintf(pFile,";DSP Coefficients\ncoef\n"); + for (k=N; k > 0; k--) + { + fprintf(pFile,"\tdc\t%9.7g\t; a%u\n",pCoeff->m_pak[k]/2.0,k); + } + for (k=N; k > 0; k--) + { + fprintf(pFile,"\tdc\t%9.7g\t; b%u\n",pCoeff->m_pbk[k]/2.0,k); + } + +} + + +void ScaleCoeff(CIIRCoeff *pCoeff) +{ + unsigned i; + double val; + + val=0; + for(i=0; i <= pCoeff->m_Na; i++) + val=MaxMag(val, pCoeff->m_pak[i]); + + for(i=0; i <= pCoeff->m_Na; i++) + pCoeff->m_pak[i] /= val; + + pCoeff->m_aScale = val; + + val=0; + for(i=0; i <= pCoeff->m_Nb; i++) + val=MaxMag(val, pCoeff->m_pbk[i]); + + for(i=0; i <= pCoeff->m_Nb; i++) + pCoeff->m_pbk[i] /= val; + + pCoeff->m_bScale = val; + +} + + +double MinMag(double val1, double val2) +{ + if(fabs(val1) < fabs(val2)) + return fabs(val1); + + return fabs(val2); +} + +double MaxMag(double val1, double val2) +{ + if(fabs(val1) > fabs(val2)) + return fabs(val1); + + return fabs(val2); +} diff --git a/Iir.h b/Iir.h new file mode 100755 index 0000000..b997844 --- /dev/null +++ b/Iir.h @@ -0,0 +1,91 @@ +/******************************************************************************/ +/* iir.h */ +/******************************************************************************/ +#ifndef IIR_H +#define IIR_H + +#include + +#ifndef iir_pi +#define iir_pi 3.1415926535897932384626433832795 +#endif + +#define IIR_FILTERTYPE_UNKNOWN 0x00000000 +#define IIR_FILTERTYPE_LOWPASS 0x00000001 +#define IIR_FILTERTYPE_HIGHPASS 0x00000002 +#define IIR_FILTERTYPE_BANDPASS 0x00000003 +#define IIR_FILTERTYPE_BANDSTOP 0x00000004 +#define IIR_FILTERTYPE_PEAKING 0x00000005 +#define IIR_FILTERTYPE_LOWSHELF 0x00000006 +#define IIR_FILTERTYPE_HIGHSHELF 0x00000007 + + +/******************************************************************************/ +#ifndef iir_float_t +#define iir_float_t float +#endif + +typedef struct _sComplex +{ + iir_float_t pRealData, pImagData; +} Complex; + +typedef struct _sIIRCoeff +{ + + iir_float_t ak0, ak1, ak2; + iir_float_t bk0, bk1, bk2; + iir_float_t xn1, xn2; + iir_float_t yn1, yn2; + +}IIRCOEFF; + +typedef struct _sIIRParam +{ + /* General Params */ + iir_float_t fg, Qf; + + /* for shelving EQs */ + iir_float_t beta; + + /* for peaking and shelving EQs */ + iir_float_t A; +}IIRPARAM; + +typedef struct _siir_lin_t +{ + unsigned order; + iir_float_t *pX, *pY; +} iir_lin_t; + + +/******************************************************************************/ +#ifdef __cplusplus +extern "C" { +#endif + +void IIRInit(struct _sIIRCoeff *pCoeff, unsigned order); +int IIRCalcPartFilterCoeff1(struct _sIIRCoeff *pCoeff, iir_float_t fa, iir_float_t fg, iir_float_t qp, unsigned filterType); +int IIRCalcPartFilterCoeff2(struct _sIIRCoeff *pCoeff, iir_float_t A, iir_float_t fa, iir_float_t fg, iir_float_t qp, unsigned filterType); +void IIRCalcFilterCoeff(struct _sIIRCoeff *pCoeff, iir_float_t fa, iir_float_t fg, iir_float_t q, unsigned order, unsigned filterType); + +iir_float_t IIRBilTrans(iir_float_t fg, iir_float_t fa); +iir_float_t IIRCalcQp(unsigned p, unsigned N); +void IIR(struct _sIIRCoeff *pCoeff, iir_float_t *xn, iir_float_t *yn, unsigned order, unsigned numPoints); +iir_float_t IIRS(struct _sIIRCoeff *pCoeff, iir_float_t xn, unsigned order); + +void IIRSSE(struct _sIIRCoeff *pCoeff, iir_float_t *xn, iir_float_t *yn, unsigned order, unsigned numPoints); +void IIRPrintCoeff(FILE *pFile, struct _sIIRCoeff *pCoeff, unsigned N); +void ScaleCoeff(struct _sIIRCoeff *pCoeff); +iir_float_t MinMag(iir_float_t val1, iir_float_t val2); +iir_float_t MaxMag(iir_float_t val1, iir_float_t val2); + +void IIR_lin_init(iir_lin_t *pObj, unsigned order); +void IIR_lin_free(iir_lin_t *pObj); +iir_float_t IIR_lin_process(iir_lin_t *pObj, iir_float_t *pB, iir_float_t *pA, iir_float_t x); + +#ifdef __cplusplus +} +#endif +#endif // IIR_H +/******************************************************************************/ diff --git a/Iir.h.double b/Iir.h.double new file mode 100755 index 0000000..e90ed61 --- /dev/null +++ b/Iir.h.double @@ -0,0 +1,69 @@ +/******************************************************************************/ +/* iir.h +/******************************************************************************/ +#ifndef IIR_H +#define IIR_H + +#define pi 3.1415926535897932384626433832795 + +#define IIR_FILTERTYPE_UNKNOWN 0x00000000 +#define IIR_FILTERTYPE_LOWPASS 0x00000001 +#define IIR_FILTERTYPE_HIGHPASS 0x00000002 +#define IIR_FILTERTYPE_BANDPASS 0x00000003 +#define IIR_FILTERTYPE_BANDSTOP 0x00000004 +#define IIR_FILTERTYPE_PEAKING 0x00000005 +#define IIR_FILTERTYPE_LOWSHELF 0x00000006 +#define IIR_FILTERTYPE_HIGHSHELF 0x00000007 + + +/******************************************************************************/ +typedef struct _sComplex +{ + double pRealData, pImagData; +} Complex; + +typedef struct _sIIRCoeff +{ + + double ak0, ak1, ak2; + double bk0, bk1, bk2; + double xn1, xn2; + double yn1, yn2; + +}IIRCOEFF; + +typedef struct _sIIRParam +{ + /* General Params */ + double fg, Qf; + + /* for shelving EQs */ + double beta; + + /* for peaking and shelving EQs */ + double A; +}IIRPARAM; + +/******************************************************************************/ +#ifdef __cplusplus +extern "C" { +#endif + +void IIRInit(struct _sIIRCoeff *pCoeff, unsigned order); +int IIRCalcPartFilterCoeff1(struct _sIIRCoeff *pCoeff, double fa, double fg, double qp, unsigned filterType); +int IIRCalcPartFilterCoeff2(struct _sIIRCoeff *pCoeff, double A, double fa, double fg, double qp, unsigned filterType); +void IIRCalcFilterCoeff(struct _sIIRCoeff *pCoeff, double fa, double fg, double q, unsigned order, unsigned filterType); + +double IIRBilTrans(double fg, double fa); +double IIRCalcQp(unsigned p, unsigned N); +void IIR(struct _sIIRCoeff *pCoeff, double *xn, double *yn, unsigned order, unsigned numPoints); +void IIRPrintCoeff(FILE *pFile, struct _sIIRCoeff *pCoeff, unsigned N); +void ScaleCoeff(struct _sIIRCoeff *pCoeff); +double MinMag(double val1, double val2); +double MaxMag(double val1, double val2); + +#ifdef __cplusplus +} +#endif +#endif // IIR_H +/******************************************************************************/ diff --git a/Iir.opt b/Iir.opt new file mode 100755 index 0000000..d73589b Binary files /dev/null and b/Iir.opt differ diff --git a/Iir2.c b/Iir2.c new file mode 100755 index 0000000..a887245 --- /dev/null +++ b/Iir2.c @@ -0,0 +1,295 @@ +/*************************************************************************/ +/* iir.c */ +/*************************************************************************/ +#include "stdio.h" +#include "stdlib.h" +#include "math.h" +#include "iir2.h" + +/*************************************************************************/ +/* Global Variables */ +/*************************************************************************/ +const char *filterTypeString[] = +{ + "Unknown filter type", + "Butterworth-Lowpass", + "Butterworth-Highpass", + "Butterworth-Bandpass", + "Butterworth-Bandstop", + "Peaking-EQ", + "Low Shelving-EQ", + "High Shelving-EQ" +}; + +/******************************************************************************/ +void IIRCalcFilterCoeff(iir_coef_t *pCoeff, iir_float_t fa, iir_float_t fg, iir_float_t q, unsigned order, unsigned filterType) +{ + unsigned p; + iir_float_t qp; + +// IIRInit(pCoeff, order); + + for(p=0; p < order/2;p++) + { + qp = q * IIRCalcQp(p+1, order); + IIRCalcPartFilterCoeff2(&pCoeff[p], 1.0, fa, fg, qp, filterType); + } +} + +int IIRCalcPartFilterCoeff1(iir_coef_t *pCoeff, iir_float_t fa, iir_float_t fg, iir_float_t Qi, unsigned filterType) +{ + iir_float_t K, a0; + iir_float_t alpha, omega, ks, kc; + unsigned error; + + omega = (iir_float_t)(2*pi*fg/fa); + ks = (iir_float_t)sin(omega); + kc = (iir_float_t)cos(omega); + alpha = 0.5f*ks /Qi; + + K = IIRBilTrans(fg, fa); + a0 = K/Qi + 1; + + switch(filterType) + { + case IIR_FILTERTYPE_LOWPASS: + + pCoeff->ak0 = 1.0f; + pCoeff->ak1 = (1 - K/Qi)/a0; + pCoeff->ak2 = 0.0; + + pCoeff->bk0 = 1.0f/a0; + pCoeff->bk1 = 1.0f/a0; + pCoeff->bk2 = 0.0; + + break; + + case IIR_FILTERTYPE_HIGHPASS: + + pCoeff->ak0 = 1.0f; + pCoeff->ak1 = (1 - K/Qi) /a0; + pCoeff->ak2 = 0.0; + + pCoeff->bk0 = 1.0f*K /a0; + pCoeff->bk1 = -1.0f*K /a0; + pCoeff->bk2 = 0.0; + + break; + + default: + error = (unsigned)-1; + break; + } + return error; +} + +int IIRCalcPartFilterCoeff2(iir_coef_t *pCoeff, iir_float_t A, iir_float_t fa, iir_float_t fg, iir_float_t qp, unsigned filterType) +{ + iir_float_t a0; + iir_float_t alpha, omega, ks, kc; + unsigned error; + + omega = (iir_float_t)(2*pi*fg/fa); + ks = (iir_float_t)sin(omega); + kc = (iir_float_t)cos(omega); + alpha = 0.5f*ks /qp; + + error = 0; + switch(filterType) + { + case IIR_FILTERTYPE_LOWPASS: + + a0 = 1 + alpha; + pCoeff->ak0 = 1.0f; + pCoeff->ak1 = -2.0f*kc /a0; + pCoeff->ak2 = (1 - alpha) /a0; + + pCoeff->bk0 = 0.5f*(1 - kc) /a0; + pCoeff->bk1 = (1 - kc) /a0; + pCoeff->bk2 = 0.5f*(1 - kc) /a0; + + break; + + case IIR_FILTERTYPE_HIGHPASS: + + a0 = 1 + alpha; + pCoeff->ak0 = 1.0f; + pCoeff->ak1 = -2.0f*kc /a0; + pCoeff->ak2 = (1 - alpha) /a0; + + pCoeff->bk0 = 0.5f*(1 + kc) /a0; + pCoeff->bk1 = -(1 + kc) /a0; + pCoeff->bk2 = 0.5f*(1 + kc) /a0; + + break; + + case IIR_FILTERTYPE_BANDPASS: + + a0 = 1 + alpha; + pCoeff->ak0 = 1.0f; + pCoeff->ak1 = -2.0f*kc /a0; + pCoeff->ak2 = (1 - alpha) /a0; + + pCoeff->bk0 = alpha /a0; + pCoeff->bk1 = 0; + pCoeff->bk2 = -alpha /a0; + + break; + + case IIR_FILTERTYPE_BANDSTOP: + + a0 = 1 + alpha; + pCoeff->ak0 = 1.0f; + pCoeff->ak1 = -2.0f*kc /a0; + pCoeff->ak2 = (1 - alpha) /a0; + + pCoeff->bk0 = 1.0f /a0; + pCoeff->bk1 = -2.0f*kc /a0; + pCoeff->bk2 = 1.0f /a0; + + break; + + case IIR_FILTERTYPE_PEAKING: + + a0 = 1 + (alpha/A); + pCoeff->ak0 = 1.0f; + pCoeff->ak1 = -2.0f*kc /a0; + pCoeff->ak2 = (1 - (alpha/A)) /a0; + + pCoeff->bk0 = (1 + (alpha*A)) /a0; + pCoeff->bk1 = -2.0f*kc /a0; + pCoeff->bk2 = (1 - (alpha*A)) /a0; + + break; + + default: + error = -1; + break; + } + return error; +} + +void IIR(iir_state_t *pState, iir_coef_t *pCoeff, iir_float_t *xn, iir_float_t *yn, unsigned order, unsigned numPoints) +{ + iir_float_t xp, yp; + unsigned i, p; + unsigned numSec = order/2; + + for (i=0; iorder = order; + pObj->pX = (iir_float_t*)malloc((order+1)*sizeof(iir_float_t)); + pObj->pY = (iir_float_t*)malloc((order+1)*sizeof(iir_float_t)); + memset(pObj->pX, 0, (order+1)*sizeof(iir_float_t)); + memset(pObj->pY, 0, (order+1)*sizeof(iir_float_t)); +} + +void IIR_lin_free(iir_lin_t *pObj) +{ + if (pObj->pX) + free(pObj->pX); + + if (pObj->pY) + free(pObj->pY); + + pObj->order = 0; + +} + +iir_float_t IIR_lin_process(iir_lin_t *pObj, iir_float_t *pB, iir_float_t *pA, iir_float_t x) +{ + unsigned i; + iir_float_t y; + + if (!pObj->order) + return 0; + + for (i=pObj->order; i >= 1; i--) + pObj->pX[i] = pObj->pX[i-1]; + + for (i=pObj->order; i >= 1; i--) + pObj->pY[i] = pObj->pY[i-1]; + + pObj->pX[0] = x; + + y = 0; + for (i=0; i <= pObj->order; i++) + y += pObj->pX[i]*pB[i]; + + for (i=1; i <= pObj->order; i++) + y -= pObj->pY[i]*pA[i]; + + pObj->pY[0] = y; + + return y; + +} + diff --git a/Iir2.h b/Iir2.h new file mode 100755 index 0000000..9718dcd --- /dev/null +++ b/Iir2.h @@ -0,0 +1,97 @@ +/******************************************************************************/ +/* iir.h */ +/******************************************************************************/ +#ifndef IIR_H +#define IIR_H + +#include + +#ifndef pi +#define pi 3.1415926535897932384626433832795 +#endif + +#define IIR_FILTERTYPE_UNKNOWN 0x00000000 +#define IIR_FILTERTYPE_LOWPASS 0x00000001 +#define IIR_FILTERTYPE_HIGHPASS 0x00000002 +#define IIR_FILTERTYPE_BANDPASS 0x00000003 +#define IIR_FILTERTYPE_BANDSTOP 0x00000004 +#define IIR_FILTERTYPE_PEAKING 0x00000005 +#define IIR_FILTERTYPE_LOWSHELF 0x00000006 +#define IIR_FILTERTYPE_HIGHSHELF 0x00000007 + + +/******************************************************************************/ +#ifndef iir_float_t +#define iir_float_t float +#endif + +typedef struct _sComplex +{ + iir_float_t pRealData, pImagData; +} Complex; + +typedef struct _siir_coef_t +{ + + iir_float_t ak0, ak1, ak2; + iir_float_t bk0, bk1, bk2; + +} iir_coef_t; + +typedef struct _siir_state_t +{ + + iir_float_t xn1, xn2; + iir_float_t yn1, yn2; + +} iir_state_t; + +typedef struct _sIIRParam +{ + /* General Params */ + iir_float_t fg, Qf; + + /* for shelving EQs */ + iir_float_t beta; + + /* for peaking and shelving EQs */ + iir_float_t A; +}IIRPARAM; + +typedef struct _siir_lin_t +{ + unsigned order; + iir_float_t *pX, *pY; +} iir_lin_t; + + +/******************************************************************************/ +#ifdef __cplusplus +extern "C" { +#endif + +void IIRInit(iir_state_t *pState, unsigned order); +int IIRCalcPartFilterCoeff1(iir_coef_t *pCoeff, iir_float_t fa, iir_float_t fg, iir_float_t qp, unsigned filterType); +int IIRCalcPartFilterCoeff2(iir_coef_t *pCoeff, iir_float_t A, iir_float_t fa, iir_float_t fg, iir_float_t qp, unsigned filterType); +void IIRCalcFilterCoeff(iir_coef_t *pCoeff, iir_float_t fa, iir_float_t fg, iir_float_t q, unsigned order, unsigned filterType); + +iir_float_t IIRBilTrans(iir_float_t fg, iir_float_t fa); +iir_float_t IIRCalcQp(unsigned p, unsigned N); +void IIR(iir_state_t *pState, iir_coef_t *pCoeff, iir_float_t *xn, iir_float_t *yn, unsigned order, unsigned numPoints); +iir_float_t IIRS(iir_state_t *pState, iir_coef_t *pCoeff, iir_float_t xn, unsigned order); + +void IIRSSE(iir_coef_t *pCoeff, iir_float_t *xn, iir_float_t *yn, unsigned order, unsigned numPoints); +void IIRPrintCoeff(FILE *pFile, iir_coef_t *pCoeff, unsigned N); +void ScaleCoeff(iir_coef_t *pCoeff); +iir_float_t MinMag(iir_float_t val1, iir_float_t val2); +iir_float_t MaxMag(iir_float_t val1, iir_float_t val2); + +void IIR_lin_init(iir_lin_t *pObj, unsigned order); +void IIR_lin_free(iir_lin_t *pObj); +iir_float_t IIR_lin_process(iir_lin_t *pObj, iir_float_t *pB, iir_float_t *pA, iir_float_t x); + +#ifdef __cplusplus +} +#endif +#endif // IIR_H +/******************************************************************************/ diff --git a/Iir_h.bak b/Iir_h.bak new file mode 100755 index 0000000..b67dea7 --- /dev/null +++ b/Iir_h.bak @@ -0,0 +1,68 @@ +/******************************************************************************/ + +/******************************************************************************/ +#define pi 3.1415926535897932384626433832795 + +#define IIR_FILTERTYPE_UNKNOWN 0x00000000 +#define IIR_FILTERTYPE_LOWPASS 0x00000001 +#define IIR_FILTERTYPE_HIGHPASS 0x00000002 +#define IIR_FILTERTYPE_BANDPASS 0x00000003 +#define IIR_FILTERTYPE_BANDSTOP 0x00000004 + +/******************************************************************************/ +typedef struct _sComplex +{ + double pRealData, pImagData; +} Complex; + +class CIIRCoeff +{ +public: + CIIRCoeff(unsigned Na=0, unsigned Nb=0) + { + Init(Na, Nb); + } + ~CIIRCoeff() + { + if(m_pak!=0) + delete [] m_pak; + + if(m_pbk!=0) + delete [] m_pbk; + } + void Init(unsigned Na=0, unsigned Nb=0) + { + m_pak= 0; + m_pbk= 0; + m_Na = Na; + m_Nb = Nb; + m_aScale = 1.0; + m_bScale = 1.0; + + if(Na!=0) + m_pak = new double[Na+1]; + + if(Nb!=0) + m_pbk = new double[Nb+1]; + } + + double *m_pak, m_aScale; + unsigned m_Na; + double *m_pbk, m_bScale; + unsigned m_Nb; +}; + +/******************************************************************************/ +int IIRCalcPartFilterCoeff1(CIIRCoeff *pCoeff, double fg, double fa, double Qi, unsigned filterType); +int IIRCalcPartFilterCoeff2(CIIRCoeff *pCoeff, double fg, double fa, double Qi, unsigned filterType); +int IIRCalcFilterCoeff(double fg, double fa, double Qi, unsigned N, CIIRCoeff *pCoeff, unsigned filterType); +double IIRBilTrans(double fg, double fa); +double IIRCalcQp(unsigned p, unsigned N); +unsigned IIRMulPolynom(double *pA, unsigned orderA, double *pB, unsigned orderB, double *pProduct); +void IIR(double *xn, double *yn, CIIRCoeff *pCoeff, unsigned numPoints); +void IIRPrintCoeff(FILE *pFile, CIIRCoeff *pCoeff, unsigned N); +void ScaleCoeff(CIIRCoeff *pCoeff); +double MinMag(double val1, double val2); +double MaxMag(double val1, double val2); + +/******************************************************************************/ diff --git a/Release/iir.obj b/Release/iir.obj new file mode 100755 index 0000000..911fe34 Binary files /dev/null and b/Release/iir.obj differ diff --git a/Release/iir.pch b/Release/iir.pch new file mode 100755 index 0000000..174c66e Binary files /dev/null and b/Release/iir.pch differ diff --git a/Release/vc60.idb b/Release/vc60.idb new file mode 100755 index 0000000..8a0ca79 Binary files /dev/null and b/Release/vc60.idb differ diff --git a/iir.dsp b/iir.dsp new file mode 100755 index 0000000..f7d98fc --- /dev/null +++ b/iir.dsp @@ -0,0 +1,100 @@ +# Microsoft Developer Studio Project File - Name="iir" - Package Owner=<4> +# Microsoft Developer Studio Generated Build File, Format Version 6.00 +# ** NICHT BEARBEITEN ** + +# TARGTYPE "Win32 (x86) Static Library" 0x0104 + +CFG=iir - Win32 Debug +!MESSAGE Dies ist kein gültiges Makefile. Zum Erstellen dieses Projekts mit NMAKE +!MESSAGE verwenden Sie den Befehl "Makefile exportieren" und führen Sie den Befehl +!MESSAGE +!MESSAGE NMAKE /f "iir.mak". +!MESSAGE +!MESSAGE Sie können beim Ausführen von NMAKE eine Konfiguration angeben +!MESSAGE durch Definieren des Makros CFG in der Befehlszeile. Zum Beispiel: +!MESSAGE +!MESSAGE NMAKE /f "iir.mak" CFG="iir - Win32 Debug" +!MESSAGE +!MESSAGE Für die Konfiguration stehen zur Auswahl: +!MESSAGE +!MESSAGE "iir - Win32 Release" (basierend auf "Win32 (x86) Static Library") +!MESSAGE "iir - Win32 Debug" (basierend auf "Win32 (x86) Static Library") +!MESSAGE + +# Begin Project +# PROP AllowPerConfigDependencies 0 +# PROP Scc_ProjName "" +# PROP Scc_LocalPath "" +CPP=cl.exe +RSC=rc.exe + +!IF "$(CFG)" == "iir - Win32 Release" + +# PROP BASE Use_MFC 0 +# PROP BASE Use_Debug_Libraries 0 +# PROP BASE Output_Dir "Release" +# PROP BASE Intermediate_Dir "Release" +# PROP BASE Target_Dir "" +# PROP Use_MFC 0 +# PROP Use_Debug_Libraries 0 +# PROP Output_Dir "Release" +# PROP Intermediate_Dir "Release" +# PROP Target_Dir "" +# ADD BASE CPP /nologo /W3 /GX /O2 /D "WIN32" /D "NDEBUG" /D "_MBCS" /D "_LIB" /YX /FD /c +# ADD CPP /nologo /W3 /GX /O2 /I "../../include" /D "WIN32" /D "NDEBUG" /D "_MBCS" /D "_LIB" /YX /FD /c +# ADD BASE RSC /l 0x407 /d "NDEBUG" +# ADD RSC /l 0x407 /d "NDEBUG" +BSC32=bscmake.exe +# ADD BASE BSC32 /nologo +# ADD BSC32 /nologo +LIB32=link.exe -lib +# ADD BASE LIB32 /nologo +# ADD LIB32 /nologo /out:"..\..\lib\release\iir.lib" + +!ELSEIF "$(CFG)" == "iir - Win32 Debug" + +# PROP BASE Use_MFC 0 +# PROP BASE Use_Debug_Libraries 1 +# PROP BASE Output_Dir "Debug" +# PROP BASE Intermediate_Dir "Debug" +# PROP BASE Target_Dir "" +# PROP Use_MFC 0 +# PROP Use_Debug_Libraries 1 +# PROP Output_Dir "Debug" +# PROP Intermediate_Dir "Debug" +# PROP Target_Dir "" +# ADD BASE CPP /nologo /W3 /Gm /GX /ZI /Od /D "WIN32" /D "_DEBUG" /D "_MBCS" /D "_LIB" /YX /FD /GZ /c +# ADD CPP /nologo /W3 /Gm /GX /ZI /Od /I "../../include" /D "WIN32" /D "_DEBUG" /D "_MBCS" /D "_LIB" /YX /FD /GZ /c +# ADD BASE RSC /l 0x407 /d "_DEBUG" +# ADD RSC /l 0x407 /d "_DEBUG" +BSC32=bscmake.exe +# ADD BASE BSC32 /nologo +# ADD BSC32 /nologo +LIB32=link.exe -lib +# ADD BASE LIB32 /nologo +# ADD LIB32 /nologo /out:"..\..\lib\debug\iir.lib" + +!ENDIF + +# Begin Target + +# Name "iir - Win32 Release" +# Name "iir - Win32 Debug" +# Begin Group "Quellcodedateien" + +# PROP Default_Filter "cpp;c;cxx;rc;def;r;odl;idl;hpj;bat" +# Begin Source File + +SOURCE=.\iir.cpp +# End Source File +# End Group +# Begin Group "Header-Dateien" + +# PROP Default_Filter "h;hpp;hxx;hm;inl" +# Begin Source File + +SOURCE=..\..\Include\Iir.h +# End Source File +# End Group +# End Target +# End Project diff --git a/iir.dsw b/iir.dsw new file mode 100755 index 0000000..34721a7 --- /dev/null +++ b/iir.dsw @@ -0,0 +1,29 @@ +Microsoft Developer Studio Workspace File, Format Version 6.00 +# WARNUNG: DIESE ARBEITSBEREICHSDATEI DARF NICHT BEARBEITET ODER GELÖSCHT WERDEN! + +############################################################################### + +Project: "iir"=.\iir.dsp - Package Owner=<4> + +Package=<5> +{{{ +}}} + +Package=<4> +{{{ +}}} + +############################################################################### + +Global: + +Package=<5> +{{{ +}}} + +Package=<3> +{{{ +}}} + +############################################################################### + diff --git a/iir.ncb b/iir.ncb new file mode 100755 index 0000000..c176375 Binary files /dev/null and b/iir.ncb differ diff --git a/iir.plg b/iir.plg new file mode 100755 index 0000000..dc251bf --- /dev/null +++ b/iir.plg @@ -0,0 +1,28 @@ + + +
+
Erstellungsprotokoll

+

+--------------------Konfiguration: iir - Win32 Debug--------------------
+

+
Befehlszeilen

+Erstellen der temporären Datei  "E:\WIN95\TEMP\RSP4230.TMP" mit Inhalten
+[
+/nologo /MLd /W3 /Gm /GX /ZI /Od /I "../../include" /D "WIN32" /D "_DEBUG" /D "_MBCS" /D "_LIB" /Fp"Debug/iir.pch" /YX /Fo"Debug/" /Fd"Debug/" /FD /GZ /c 
+"G:\work\Develop\MSVC\LIBSRC\IIR\iir.cpp"
+]
+Creating command line "cl.exe @E:\WIN95\TEMP\RSP4230.TMP" 
+Erstellen der Befehlzeile "link.exe -lib /nologo /out:"..\..\lib\debug\iir.lib"  .\Debug\iir.obj "
+
Ausgabefenster

+Kompilierung läuft...
+iir.cpp
+g:\work\develop\msvc\libsrc\iir\iir.cpp(63) : warning C4101: 'B0' : Unreferenzierte lokale Variable
+Bibliothek wird erstellt...
+
+
+
+
Ergebnisse

+iir.lib - 0 Fehler, 1 Warnung(en)
+
+ +