/*
July 15, 2015
Iowa Hills Software LLC
http://www.iowahills.com
If you find a problem with this code, please leave us a note on:
http://www.iowahills.com/feedbackcomments.html
Source: ~Projects\Common\BasicFIRFilterCode.cpp
This generic FIR filter code is described in most textbooks.
e.g. Discrete Time Signal Processing, Oppenheim and Shafer
A nice paper on this topic is:
http://dea.brunel.ac.uk/cmsp/Home_Saeed_Vaseghi/Chapter05-DigitalFilters.pdf
This code first generates either a low pass, high pass, band pass, or notch
impulse response for a rectangular window. It then applies a window to this
impulse response.
There are several windows available, including the Kaiser, Sinc, Hanning,
Blackman, and Hamming. Of these, the Kaiser and Sinc are probably the most useful
for FIR filters because their sidelobe levels can be controlled with the Beta parameter.
This is a typical function call:
BasicFIR(FirCoeff, NumTaps, PassType, OmegaC, BW, wtKAISER, Beta);
BasicFIR(FirCoeff, 33, LPF, 0.2, 0.0, wtKAISER, 3.2);
33 tap, low pass, corner frequency at 0.2, BW=0 (ignored in the low pass code),
Kaiser window, Kaiser Beta = 3.2
These variables should be defined similar to this:
double FirCoeff[MAXNUMTAPS];
int NumTaps; NumTaps can be even or odd, but must be less than the FirCoeff array size.
TPassTypeName PassType; PassType is an enum defined in the header file. LPF, HPF, BPF, or NOTCH
double OmegaC 0.0 < OmegaC < 1.0 The filters corner freq, or center freq if BPF or NOTCH
double BW 0.0 < BW < 1.0 The filters band width if BPF or NOTCH
TWindowType WindowType; WindowType is an enum defined in the header to be one of these.
wtNONE, wtKAISER, wtSINC, wtHANNING, .... and others.
double Beta; 0 <= Beta <= 10.0 Beta is used with the Kaiser, Sinc, and Sine windows only.
It controls the transition BW and sidelobe level of the filters.
If you want to use it, Kaiser originally defined Beta as follows.
He derived its value based on the desired sidelobe level, dBAtten.
double dBAtten, Beta, Beta1=0.0, Beta2=0.0;
if(dBAtten < 21.0)dBAtten = 21.0;
if(dBAtten > 50.0)Beta1 = 0.1102 * (dBAtten - 8.7);
if(dBAtten >= 21.0 && dBAtten <= 50.0) Beta2 = 0.5842 * pow(dBAtten - 21.0, 0.4) + 0.07886 * (dBAtten - 21.0);
Beta = Beta1 + Beta2;
*/
//---------------------------------------------------------------------------
#pragma hdrstop // for Embarcadero's C++ Builder
#include "WindowedFIRFilterWebCode.h"
#include
#include // defines new(std::nothrow)
#include // for Embarcadero's ShowMessage function.
#pragma package(smart_init) // for C++ Builder
#define M_2PI 6.28318530717958647692 // M_PI should be in the math.h file
//---------------------------------------------------------------------------
// This first calculates the impulse response for a rectangular window.
// It then applies the windowing function of choice to the impulse response.
void BasicFIR(double *FirCoeff, int NumTaps, TPassTypeName PassType, double OmegaC, double BW, TWindowType WindowType, double WinBeta)
{
int j;
double Arg, OmegaLow, OmegaHigh;
switch(PassType)
{
case LPF:
for(j=0; j -1.0E-5 && x < 1.0E-5)return(1.0);
return(sin(x)/x);
}
//---------------------------------------------------------------------------
// These are the various windows definitions. These windows can be used for either
// FIR filter design or with an FFT for spectral analysis.
// Sourced verbatim from: ~MyDocs\Code\Common\FFTFunctions.cpp
// For definitions, see this article: http://en.wikipedia.org/wiki/Window_function
// This function has 6 inputs
// Data is the array, of length N, containing the data to to be windowed.
// This data is either a FIR filter sinc pulse, or the data to be analyzed by an fft.
// WindowType is an enum defined in the header file, which is at the bottom of this file.
// e.g. wtKAISER, wtSINC, wtHANNING, wtHAMMING, wtBLACKMAN, ...
// Alpha sets the width of the flat top.
// Windows such as the Tukey and Trapezoid are defined to have a variably wide flat top.
// As can be seen by its definition, the Tukey is just a Hanning window with a flat top.
// Alpha can be used to give any of these windows a partial flat top, except the Flattop and Kaiser.
// Alpha = 0 gives the original window. (i.e. no flat top)
// To generate a Tukey window, use a Hanning with 0 < Alpha < 1
// To generate a Bartlett window (triangular), use a Trapezoid window with Alpha = 0.
// Alpha = 1 generates a rectangular window in all cases. (except the Flattop and Kaiser)
// Beta is used with the Kaiser, Sinc, and Sine windows only.
// These three windows are primarily used for FIR filter design, not spectral analysis.
// In FIR filter design, Beta controls the filter's transition bandwidth and the sidelobe levels.
// The code ignores Beta except in the Kaiser, Sinc, and Sine window cases.
// UnityGain controls whether the gain of these windows is set to unity.
// Only the Flattop window has unity gain by design. The Hanning window, for example, has a gain
// of 1/2. UnityGain = true will set the gain of all these windows to 1.
// Then, when the window is applied to a signal, the signal's energy content is preserved.
// Don't use this with FIR filter design however. Since most of the enegy in an FIR sinc pulse
// is in the middle of the window, the window needs a peak amplitude of one, not unity gain.
// Setting UnityGain = true will simply cause the resulting FIR filter to have excess gain.
// If using these windows for FIR filters, start with the Kaiser, Sinc, or Sine windows and
// adjust Beta for the desired transition BW and sidelobe levels (set Alpha = 0).
// While the FlatTop is an excellent window for spectral analysis, don't use it for FIR filter design.
// It has a peak amplitude of ~ 4.7 which causes the resulting FIR filter to have about this much gain.
// It works poorly for FIR filters even if you adjust its peak amplitude.
// The Trapezoid also works poorly for FIR filter design.
// If using these windows with an fft for spectral analysis, start with the Hanning, Gauss, or Flattop.
// When choosing a window for spectral analysis, you must trade off between resolution and amplitude accuracy.
// The Hanning has the best resolution while the Flatop has the best amplitude accuracy.
// The Gauss is midway between these two for both accuracy and resolution.
// These three were the only windows available in the HP 89410A Vector Signal Analyzer. Which is to say,
// unless you have specific windowing requirements, use one of these 3 for general purpose signal analysis.
// Set UnityGain = true when using any of these windows for spectral analysis to preserve the signal's enegy level.
void WindowData(double *Data, int N, TWindowType WindowType, double Alpha, double Beta, bool UnityGain)
{
if(WindowType == wtNONE) return;
int j, M, TopWidth;
double dM, *WinCoeff;
if(WindowType == wtKAISER || WindowType == wtFLATTOP )Alpha = 0.0;
if(Alpha < 0.0)Alpha = 0.0;
if(Alpha > 1.0)Alpha = 1.0;
if(Beta < 0.0)Beta = 0.0;
if(Beta > 10.0)Beta = 10.0;
WinCoeff = new(std::nothrow) double[N+2];
if(WinCoeff == NULL)
{
ShowMessage("Failed to allocate memory in FFTFunctions::WindowFFTData() ");
return;
}
TopWidth = (int)( Alpha * (double)N );
if(TopWidth%2 != 0)TopWidth++;
if(TopWidth > N)TopWidth = N;
M = N - TopWidth;
dM = M + 1;
// Calculate the window for N/2 points, then fold the window over (at the bottom).
// TopWidth points will be set to 1.
if(WindowType == wtKAISER)
{
double Arg;
for(j=0; j 0. Cannot be applied to a Kaiser or Flat Top.
if(WindowType != wtKAISER && WindowType != wtFLATTOP)
{
for(j=M/2; j