mirror of
https://github.com/f4exb/sdrangel.git
synced 2024-12-23 01:55:48 -05:00
Fix and refactor fir filter
This commit is contained in:
parent
b9adbfb1d8
commit
d7f8208814
@ -22,8 +22,7 @@
|
||||
#include "dsp/nco.h"
|
||||
#include "dsp/interpolator.h"
|
||||
#include "dsp/agc.h"
|
||||
#include "dsp/bandpass.h"
|
||||
#include "dsp/lowpass.h"
|
||||
#include "dsp/firfilter.h"
|
||||
#include "dsp/phaselockcomplex.h"
|
||||
#include "audio/audiofifo.h"
|
||||
#include "util/movingaverage.h"
|
||||
|
@ -23,7 +23,7 @@
|
||||
#include "dsp/channelsamplesink.h"
|
||||
#include "dsp/nco.h"
|
||||
#include "dsp/interpolator.h"
|
||||
#include "dsp/lowpass.h"
|
||||
#include "dsp/firfilter.h"
|
||||
#include "dsp/movingaverage.h"
|
||||
#include "dsp/fftfilt.h"
|
||||
#include "dsp/phaselock.h"
|
||||
|
@ -22,8 +22,7 @@
|
||||
#include "dsp/phasediscri.h"
|
||||
#include "dsp/nco.h"
|
||||
#include "dsp/interpolator.h"
|
||||
#include "dsp/lowpass.h"
|
||||
#include "dsp/bandpass.h"
|
||||
#include "dsp/firfilter.h"
|
||||
#include "dsp/afsquelch.h"
|
||||
#include "dsp/afsquelch.h"
|
||||
#include "audio/audiofifo.h"
|
||||
|
@ -24,8 +24,7 @@
|
||||
#include "dsp/phasediscri.h"
|
||||
#include "dsp/nco.h"
|
||||
#include "dsp/interpolator.h"
|
||||
#include "dsp/lowpass.h"
|
||||
#include "dsp/bandpass.h"
|
||||
#include "dsp/firfilter.h"
|
||||
#include "dsp/afsquelch.h"
|
||||
#include "dsp/agc.h"
|
||||
#include "dsp/ctcssdetector.h"
|
||||
|
@ -23,7 +23,7 @@
|
||||
#include "dsp/channelsamplesink.h"
|
||||
#include "dsp/nco.h"
|
||||
#include "dsp/interpolator.h"
|
||||
#include "dsp/lowpass.h"
|
||||
#include "dsp/firfilter.h"
|
||||
#include "util/movingaverage.h"
|
||||
#include "dsp/fftfilt.h"
|
||||
#include "dsp/phasediscri.h"
|
||||
|
@ -27,8 +27,7 @@
|
||||
#include "dsp/interpolator.h"
|
||||
#include "util/movingaverage.h"
|
||||
#include "dsp/agc.h"
|
||||
#include "dsp/bandpass.h"
|
||||
#include "dsp/lowpass.h"
|
||||
#include "dsp/firfilter.h"
|
||||
#include "dsp/phaselockcomplex.h"
|
||||
#include "dsp/freqlockcomplex.h"
|
||||
#include "util/message.h"
|
||||
|
@ -27,7 +27,7 @@
|
||||
#include "dsp/phasediscri.h"
|
||||
#include "dsp/movingaverage.h"
|
||||
#include "dsp/agc.h"
|
||||
#include "dsp/bandpass.h"
|
||||
#include "dsp/firfilter.h"
|
||||
#include "util/udpsinkutil.h"
|
||||
#include "audio/audiofifo.h"
|
||||
|
||||
|
@ -29,9 +29,7 @@
|
||||
#include "dsp/nco.h"
|
||||
#include "dsp/ncof.h"
|
||||
#include "dsp/interpolator.h"
|
||||
#include "dsp/lowpass.h"
|
||||
#include "dsp/bandpass.h"
|
||||
#include "dsp/highpass.h"
|
||||
#include "dsp/firfilter.h"
|
||||
#include "dsp/raisedcosine.h"
|
||||
#include "dsp/fmpreemphasis.h"
|
||||
#include "util/lfsr.h"
|
||||
|
@ -27,8 +27,7 @@
|
||||
#include "dsp/nco.h"
|
||||
#include "dsp/ncof.h"
|
||||
#include "dsp/interpolator.h"
|
||||
#include "dsp/lowpass.h"
|
||||
#include "dsp/bandpass.h"
|
||||
#include "dsp/firfilter.h"
|
||||
#include "dsp/filterrc.h"
|
||||
#include "util/movingaverage.h"
|
||||
#include "dsp/cwkeyer.h"
|
||||
|
@ -29,9 +29,7 @@
|
||||
#include "dsp/nco.h"
|
||||
#include "dsp/ncof.h"
|
||||
#include "dsp/interpolator.h"
|
||||
#include "dsp/lowpass.h"
|
||||
#include "dsp/bandpass.h"
|
||||
#include "dsp/highpass.h"
|
||||
#include "dsp/firfilter.h"
|
||||
#include "dsp/raisedcosine.h"
|
||||
#include "dsp/fmpreemphasis.h"
|
||||
#include "util/lfsr.h"
|
||||
|
@ -101,6 +101,7 @@ set(sdrbase_SOURCES
|
||||
dsp/filtermbe.cpp
|
||||
dsp/filerecord.cpp
|
||||
dsp/filerecordinterface.cpp
|
||||
dsp/firfilter.cpp
|
||||
dsp/fmpreemphasis.cpp
|
||||
dsp/freqlockcomplex.cpp
|
||||
dsp/interpolator.cpp
|
||||
@ -108,7 +109,6 @@ set(sdrbase_SOURCES
|
||||
dsp/glspectrumsettings.cpp
|
||||
dsp/hbfilterchainconverter.cpp
|
||||
dsp/hbfiltertraits.cpp
|
||||
dsp/lowpass.cpp
|
||||
dsp/mimochannel.cpp
|
||||
dsp/nco.cpp
|
||||
dsp/ncof.cpp
|
||||
@ -245,6 +245,7 @@ set(sdrbase_HEADERS
|
||||
dsp/filtermbe.h
|
||||
dsp/filerecord.h
|
||||
dsp/filerecordinterface.h
|
||||
dsp/firfilter.h
|
||||
dsp/fmpreemphasis.h
|
||||
dsp/freqlockcomplex.h
|
||||
dsp/gfft.h
|
||||
@ -266,7 +267,6 @@ set(sdrbase_HEADERS
|
||||
dsp/inthalfbandfiltersti.h
|
||||
dsp/kissfft.h
|
||||
dsp/kissengine.h
|
||||
dsp/lowpass.h
|
||||
dsp/mimochannel.h
|
||||
dsp/misc.h
|
||||
dsp/movingaverage.h
|
||||
|
@ -1,134 +0,0 @@
|
||||
#ifndef INCLUDE_BANDPASS_H
|
||||
#define INCLUDE_BANDPASS_H
|
||||
|
||||
#define _USE_MATH_DEFINES
|
||||
#include <math.h>
|
||||
#include "dsp/dsptypes.h"
|
||||
|
||||
// #undef M_PI
|
||||
// #define M_PI 3.14159265358979323846
|
||||
|
||||
template <class Type> class Bandpass {
|
||||
public:
|
||||
Bandpass() : m_ptr(0) { }
|
||||
|
||||
void create(int nTaps, double sampleRate, double lowCutoff, double highCutoff)
|
||||
{
|
||||
std::vector<Real> taps_lp;
|
||||
std::vector<Real> taps_hp;
|
||||
double wcl = 2.0 * M_PI * lowCutoff;
|
||||
double Wcl = wcl / sampleRate;
|
||||
double wch = 2.0 * M_PI * highCutoff;
|
||||
double Wch = wch / sampleRate;
|
||||
int i;
|
||||
|
||||
// check constraints
|
||||
if (!(nTaps & 1))
|
||||
{
|
||||
qDebug("Bandpass filter has to have an odd number of taps");
|
||||
nTaps++;
|
||||
}
|
||||
|
||||
// make room
|
||||
m_samples.resize(nTaps);
|
||||
|
||||
for (int i = 0; i < nTaps; i++) {
|
||||
m_samples[i] = 0;
|
||||
}
|
||||
|
||||
m_ptr = 0;
|
||||
m_taps.resize(nTaps / 2 + 1);
|
||||
taps_lp.resize(nTaps / 2 + 1);
|
||||
taps_hp.resize(nTaps / 2 + 1);
|
||||
|
||||
// generate Sinc filter core
|
||||
for (i = 0; i < nTaps / 2 + 1; i++)
|
||||
{
|
||||
if (i == (nTaps - 1) / 2)
|
||||
{
|
||||
taps_lp[i] = Wch / M_PI;
|
||||
taps_hp[i] = -(Wcl / M_PI);
|
||||
}
|
||||
else
|
||||
{
|
||||
taps_lp[i] = sin(((double)i - ((double)nTaps - 1.0) / 2.0) * Wch) / (((double)i - ((double)nTaps - 1.0) / 2.0) * M_PI);
|
||||
taps_hp[i] = -sin(((double)i - ((double)nTaps - 1.0) / 2.0) * Wcl) / (((double)i - ((double)nTaps - 1.0) / 2.0) * M_PI);
|
||||
}
|
||||
}
|
||||
|
||||
taps_hp[(nTaps - 1) / 2] += 1;
|
||||
|
||||
// apply Hamming window and combine lowpass and highpass
|
||||
for (i = 0; i < nTaps / 2 + 1; i++)
|
||||
{
|
||||
taps_lp[i] *= 0.54 + 0.46 * cos((2.0 * M_PI * ((double)i - ((double)nTaps - 1.0) / 2.0)) / (double)nTaps);
|
||||
taps_hp[i] *= 0.54 + 0.46 * cos((2.0 * M_PI * ((double)i - ((double)nTaps - 1.0) / 2.0)) / (double)nTaps);
|
||||
m_taps[i] = -(taps_lp[i]+taps_hp[i]);
|
||||
}
|
||||
|
||||
m_taps[(nTaps - 1) / 2] += 1;
|
||||
|
||||
// normalize
|
||||
Real sum = 0;
|
||||
|
||||
for (i = 0; i < (int)m_taps.size() - 1; i++) {
|
||||
sum += m_taps[i] * 2;
|
||||
}
|
||||
|
||||
sum += m_taps[i] - 1;
|
||||
|
||||
for (i = 0; i < (int)m_taps.size(); i++) {
|
||||
m_taps[i] /= sum;
|
||||
}
|
||||
}
|
||||
|
||||
Type filter(Type sample)
|
||||
{
|
||||
Type acc = 0;
|
||||
int a = m_ptr;
|
||||
int b = a - 1;
|
||||
int i, n_taps, size;
|
||||
|
||||
m_samples[m_ptr] = sample;
|
||||
size = m_samples.size(); // Valgrind optim (2)
|
||||
|
||||
while (b < 0) {
|
||||
b += size;
|
||||
}
|
||||
|
||||
n_taps = m_taps.size() - 1; // Valgrind optim
|
||||
|
||||
for(i = 0; i < n_taps; i++)
|
||||
{
|
||||
acc += (m_samples[a] + m_samples[b]) * m_taps[i];
|
||||
a++;
|
||||
|
||||
while (a >= size) {
|
||||
a -= size;
|
||||
}
|
||||
|
||||
b--;
|
||||
|
||||
while(b < 0) {
|
||||
b += size;
|
||||
}
|
||||
}
|
||||
|
||||
acc += m_samples[a] * m_taps[i];
|
||||
|
||||
m_ptr++;
|
||||
|
||||
while (m_ptr >= size) {
|
||||
m_ptr -= size;
|
||||
}
|
||||
|
||||
return acc;
|
||||
}
|
||||
|
||||
private:
|
||||
std::vector<Real> m_taps;
|
||||
std::vector<Type> m_samples;
|
||||
int m_ptr;
|
||||
};
|
||||
|
||||
#endif // INCLUDE_BANDPASS_H
|
52
sdrbase/dsp/firfilter.cpp
Normal file
52
sdrbase/dsp/firfilter.cpp
Normal file
@ -0,0 +1,52 @@
|
||||
///////////////////////////////////////////////////////////////////////////////////
|
||||
// Copyright (C) 2020 kasper93 //
|
||||
// written by Kacper Michajłow and Edouard Griffiths //
|
||||
// //
|
||||
// This program is free software; you can redistribute it and/or modify //
|
||||
// it under the terms of the GNU General Public License as published by //
|
||||
// the Free Software Foundation as version 3 of the License, or //
|
||||
// (at your option) any later version. //
|
||||
// //
|
||||
// This program is distributed in the hope that it will be useful, //
|
||||
// but WITHOUT ANY WARRANTY; without even the implied warranty of //
|
||||
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the //
|
||||
// GNU General Public License V3 for more details. //
|
||||
// //
|
||||
// You should have received a copy of the GNU General Public License //
|
||||
// along with this program. If not, see <http://www.gnu.org/licenses/>. //
|
||||
///////////////////////////////////////////////////////////////////////////////////
|
||||
|
||||
#include "firfilter.h"
|
||||
|
||||
void FirFilterGenerators::generateLowPassFilter(int nTaps, double sampleRate, double cutoff, std::vector<Real> &taps)
|
||||
{
|
||||
if (!(nTaps & 1))
|
||||
{
|
||||
qDebug("Filter has to have an odd number of taps");
|
||||
nTaps++;
|
||||
}
|
||||
|
||||
double Wc = (2.0 * M_PI * cutoff) / sampleRate;
|
||||
int halfTaps = nTaps / 2 + 1;
|
||||
taps.resize(halfTaps);
|
||||
|
||||
for (int i = 0; i < halfTaps; ++i)
|
||||
{
|
||||
if (i == halfTaps - 1)
|
||||
{
|
||||
taps[i] = Wc / M_PI;
|
||||
}
|
||||
else
|
||||
{
|
||||
int n = i - (nTaps - 1) / 2;
|
||||
taps[i] = sin(n * Wc) / (n * M_PI);
|
||||
}
|
||||
}
|
||||
|
||||
// Blackman window
|
||||
for (int i = 0; i < halfTaps; i++)
|
||||
{
|
||||
int n = i - (nTaps - 1) / 2;
|
||||
taps[i] *= 0.42 + 0.5 * cos((2.0 * M_PI * n) / nTaps) + 0.08 * cos((4.0 * M_PI * n) / nTaps);
|
||||
}
|
||||
}
|
153
sdrbase/dsp/firfilter.h
Normal file
153
sdrbase/dsp/firfilter.h
Normal file
@ -0,0 +1,153 @@
|
||||
///////////////////////////////////////////////////////////////////////////////////
|
||||
// Copyright (C) 2020 kasper93 //
|
||||
// written by Kacper Michajłow and Edouard Griffiths //
|
||||
// //
|
||||
// This program is free software; you can redistribute it and/or modify //
|
||||
// it under the terms of the GNU General Public License as published by //
|
||||
// the Free Software Foundation as version 3 of the License, or //
|
||||
// (at your option) any later version. //
|
||||
// //
|
||||
// This program is distributed in the hope that it will be useful, //
|
||||
// but WITHOUT ANY WARRANTY; without even the implied warranty of //
|
||||
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the //
|
||||
// GNU General Public License V3 for more details. //
|
||||
// //
|
||||
// You should have received a copy of the GNU General Public License //
|
||||
// along with this program. If not, see <http://www.gnu.org/licenses/>. //
|
||||
///////////////////////////////////////////////////////////////////////////////////
|
||||
|
||||
#pragma once
|
||||
|
||||
#define _USE_MATH_DEFINES
|
||||
#include <math.h>
|
||||
#include "dsp/dsptypes.h"
|
||||
#include "export.h"
|
||||
|
||||
class FirFilterGenerators
|
||||
{
|
||||
public:
|
||||
static void generateLowPassFilter(int nTaps, double sampleRate, double cutoff, std::vector<Real> &taps);
|
||||
};
|
||||
|
||||
template <class Type>
|
||||
class FirFilter
|
||||
{
|
||||
public:
|
||||
Type filter(Type sample)
|
||||
{
|
||||
Type acc = 0;
|
||||
int n_samples = m_samples.size();
|
||||
int n_taps = m_taps.size() - 1;
|
||||
int a = m_ptr;
|
||||
int b = a == 0 ? n_samples - 1 : a - 1;
|
||||
|
||||
m_samples[m_ptr] = sample;
|
||||
|
||||
for (size_t i = 0; i < n_taps; ++i)
|
||||
{
|
||||
acc += (m_samples[a++] + m_samples[b--]) * m_taps[i];
|
||||
|
||||
if (a == n_samples) {
|
||||
a = 0;
|
||||
}
|
||||
|
||||
if (b == -1) {
|
||||
b = n_samples - 1;
|
||||
}
|
||||
}
|
||||
|
||||
acc += m_samples[a] * m_taps[n_taps];
|
||||
|
||||
if (++m_ptr == n_samples) {
|
||||
m_ptr = 0;
|
||||
}
|
||||
|
||||
return acc;
|
||||
}
|
||||
|
||||
protected:
|
||||
void init(int nTaps)
|
||||
{
|
||||
m_ptr = 0;
|
||||
m_samples.resize(nTaps);
|
||||
|
||||
for (int i = 0; i < nTaps; i++) {
|
||||
m_samples[i] = 0;
|
||||
}
|
||||
}
|
||||
|
||||
void normalize(Real sum_fix = 0.0)
|
||||
{
|
||||
Real sum = 0;
|
||||
size_t i;
|
||||
|
||||
for (i = 0; i < m_taps.size() - 1; ++i) {
|
||||
sum += m_taps[i] * 2.0;
|
||||
}
|
||||
|
||||
sum += m_taps[i] + sum_fix;
|
||||
|
||||
for (i = 0; i < m_taps.size(); ++i) {
|
||||
m_taps[i] /= sum;
|
||||
}
|
||||
}
|
||||
|
||||
protected:
|
||||
std::vector<Real> m_taps;
|
||||
std::vector<Type> m_samples;
|
||||
size_t m_ptr;
|
||||
};
|
||||
|
||||
template <class T>
|
||||
struct Lowpass : public FirFilter<T>
|
||||
{
|
||||
public:
|
||||
void create(int nTaps, double sampleRate, double cutoff)
|
||||
{
|
||||
this->init(nTaps);
|
||||
FirFilterGenerators::generateLowPassFilter(nTaps, sampleRate, cutoff, this->m_taps);
|
||||
this->normalize();
|
||||
}
|
||||
};
|
||||
|
||||
template <class T>
|
||||
struct Bandpass : public FirFilter<T>
|
||||
{
|
||||
void create(int nTaps, double sampleRate, double lowCutoff, double highCutoff)
|
||||
{
|
||||
this->init(nTaps);
|
||||
FirFilterGenerators::generateLowPassFilter(nTaps, sampleRate, highCutoff, this->m_taps);
|
||||
std::vector<Real> highPass;
|
||||
FirFilterGenerators::generateLowPassFilter(nTaps, sampleRate, lowCutoff, highPass);
|
||||
|
||||
for (size_t i = 0; i < highPass.size(); ++i) {
|
||||
highPass[i] = -highPass[i];
|
||||
}
|
||||
|
||||
highPass[highPass.size() - 1] += 1;
|
||||
|
||||
for (size_t i = 0; i < this->m_taps.size(); ++i) {
|
||||
this->m_taps[i] = -(this->m_taps[i] + highPass[i]);
|
||||
}
|
||||
|
||||
this->m_taps[this->m_taps.size() - 1] += 1;
|
||||
this->normalize(-1.0);
|
||||
}
|
||||
};
|
||||
|
||||
template <class T>
|
||||
struct Highpass : public FirFilter<T>
|
||||
{
|
||||
void create(int nTaps, double sampleRate, double cutoff)
|
||||
{
|
||||
this->init(nTaps);
|
||||
FirFilterGenerators::generateLowPassFilter(nTaps, sampleRate, cutoff, this->m_taps);
|
||||
|
||||
for (size_t i = 0; i < this->m_taps.size(); ++i) {
|
||||
this->m_taps[i] = -this->m_taps[i];
|
||||
}
|
||||
|
||||
this->m_taps[this->m_taps.size() - 1] += 1;
|
||||
this->normalize(-1.0);
|
||||
}
|
||||
};
|
@ -1,114 +0,0 @@
|
||||
#ifndef INCLUDE_HIGHPASS_H
|
||||
#define INCLUDE_HIGHPASS_H
|
||||
|
||||
#define _USE_MATH_DEFINES
|
||||
#include <math.h>
|
||||
#include "dsp/dsptypes.h"
|
||||
|
||||
template <class Type> class Highpass {
|
||||
public:
|
||||
Highpass() { }
|
||||
|
||||
void create(int nTaps, double sampleRate, double cutoff)
|
||||
{
|
||||
double wc = 2.0 * M_PI * cutoff;
|
||||
double Wc = wc / sampleRate;
|
||||
int i;
|
||||
|
||||
// check constraints
|
||||
if (!(nTaps & 1))
|
||||
{
|
||||
qDebug("Highpass filter has to have an odd number of taps");
|
||||
nTaps++;
|
||||
}
|
||||
|
||||
// make room
|
||||
m_samples.resize(nTaps);
|
||||
|
||||
for (int i = 0; i < nTaps; i++) {
|
||||
m_samples[i] = 0;
|
||||
}
|
||||
|
||||
m_ptr = 0;
|
||||
m_taps.resize(nTaps / 2 + 1);
|
||||
|
||||
// generate Sinc filter core for lowpass but inverting every other tap for highpass keeping center tap
|
||||
for (i = 0; i < nTaps / 2 + 1; i++)
|
||||
{
|
||||
if (i == (nTaps - 1) / 2) {
|
||||
m_taps[i] = -(Wc / M_PI);
|
||||
} else {
|
||||
m_taps[i] = -sin(((double)i - ((double)nTaps - 1.0) / 2.0) * Wc) / (((double)i - ((double)nTaps - 1.0) / 2.0) * M_PI);
|
||||
}
|
||||
}
|
||||
|
||||
m_taps[(nTaps - 1) / 2] += 1;
|
||||
|
||||
// apply Hamming window
|
||||
for (i = 0; i < nTaps / 2 + 1; i++) {
|
||||
m_taps[i] *= 0.54 + 0.46 * cos((2.0 * M_PI * ((double)i - ((double)nTaps - 1.0) / 2.0)) / (double)nTaps);
|
||||
}
|
||||
|
||||
// normalize
|
||||
Real sum = 0;
|
||||
|
||||
for (i = 0; i < (int)m_taps.size() - 1; i++) {
|
||||
sum += m_taps[i] * 2;
|
||||
}
|
||||
|
||||
sum += m_taps[i] - 1;
|
||||
|
||||
for (i = 0; i < (int)m_taps.size(); i++) {
|
||||
m_taps[i] /= sum;
|
||||
}
|
||||
}
|
||||
|
||||
Type filter(Type sample)
|
||||
{
|
||||
Type acc = 0;
|
||||
int a = m_ptr;
|
||||
int b = a - 1;
|
||||
int i, n_taps, size;
|
||||
|
||||
m_samples[m_ptr] = sample;
|
||||
size = m_samples.size(); // Valgrind optim (2)
|
||||
|
||||
while (b < 0) {
|
||||
b += size;
|
||||
}
|
||||
|
||||
n_taps = m_taps.size() - 1; // Valgrind optim
|
||||
|
||||
for (i = 0; i < n_taps; i++)
|
||||
{
|
||||
acc += (m_samples[a] + m_samples[b]) * m_taps[i];
|
||||
a++;
|
||||
|
||||
while (a >= size) {
|
||||
a -= size;
|
||||
}
|
||||
|
||||
b--;
|
||||
|
||||
while (b < 0) {
|
||||
b += size;
|
||||
}
|
||||
}
|
||||
|
||||
acc += m_samples[a] * m_taps[i];
|
||||
m_ptr++;
|
||||
|
||||
while (m_ptr >= size) {
|
||||
m_ptr -= size;
|
||||
}
|
||||
|
||||
return acc;
|
||||
}
|
||||
|
||||
private:
|
||||
std::vector<Real> m_taps;
|
||||
std::vector<Type> m_samples;
|
||||
int m_ptr;
|
||||
};
|
||||
|
||||
#endif // INCLUDE_HIGHPASS_H
|
@ -1,6 +0,0 @@
|
||||
#include <stdio.h>
|
||||
#include <QtGlobal>
|
||||
#define _USE_MATH_DEFINES
|
||||
#include <math.h>
|
||||
#include <vector>
|
||||
#include "dsp/lowpass.h"
|
@ -1,112 +0,0 @@
|
||||
#ifndef INCLUDE_LOWPASS_H
|
||||
#define INCLUDE_LOWPASS_H
|
||||
|
||||
#define _USE_MATH_DEFINES
|
||||
#include <math.h>
|
||||
#include "dsp/dsptypes.h"
|
||||
|
||||
template <class Type> class Lowpass {
|
||||
public:
|
||||
Lowpass() : m_ptr(0) { }
|
||||
|
||||
void create(int nTaps, double sampleRate, double cutoff)
|
||||
{
|
||||
double wc = 2.0 * M_PI * cutoff;
|
||||
double Wc = wc / sampleRate;
|
||||
int i;
|
||||
|
||||
// check constraints
|
||||
if (!(nTaps & 1))
|
||||
{
|
||||
qDebug("Lowpass filter has to have an odd number of taps");
|
||||
nTaps++;
|
||||
}
|
||||
|
||||
// make room
|
||||
m_samples.resize(nTaps);
|
||||
|
||||
for (int i = 0; i < nTaps; i++) {
|
||||
m_samples[i] = 0;
|
||||
}
|
||||
|
||||
m_ptr = 0;
|
||||
m_taps.resize(nTaps / 2 + 1);
|
||||
|
||||
// generate Sinc filter core
|
||||
for (i = 0; i < nTaps / 2 + 1; i++)
|
||||
{
|
||||
if(i == (nTaps - 1) / 2) {
|
||||
m_taps[i] = Wc / M_PI;
|
||||
} else {
|
||||
m_taps[i] = sin(((double)i - ((double)nTaps - 1.0) / 2.0) * Wc) / (((double)i - ((double)nTaps - 1.0) / 2.0) * M_PI);
|
||||
}
|
||||
}
|
||||
|
||||
// apply Hamming window
|
||||
for (i = 0; i < nTaps / 2 + 1; i++) {
|
||||
m_taps[i] *= 0.54 + 0.46 * cos((2.0 * M_PI * ((double)i - ((double)nTaps - 1.0) / 2.0)) / (double)nTaps);
|
||||
}
|
||||
|
||||
// normalize
|
||||
Real sum = 0;
|
||||
|
||||
for (i = 0; i < (int)m_taps.size() - 1; i++) {
|
||||
sum += m_taps[i] * 2;
|
||||
}
|
||||
|
||||
sum += m_taps[i];
|
||||
|
||||
for (i = 0; i < (int)m_taps.size(); i++) {
|
||||
m_taps[i] /= sum;
|
||||
}
|
||||
}
|
||||
|
||||
Type filter(Type sample)
|
||||
{
|
||||
Type acc = 0;
|
||||
int a = m_ptr;
|
||||
int b = a - 1;
|
||||
int i, n_taps, size;
|
||||
|
||||
m_samples[m_ptr] = sample;
|
||||
size = m_samples.size(); // Valgrind optim (2)
|
||||
|
||||
while (b < 0) {
|
||||
b += size;
|
||||
}
|
||||
|
||||
n_taps = m_taps.size() - 1; // Valgrind optim
|
||||
|
||||
for (i = 0; i < n_taps; i++)
|
||||
{
|
||||
acc += (m_samples[a] + m_samples[b]) * m_taps[i];
|
||||
a++;
|
||||
|
||||
while (a >= size) {
|
||||
a -= size;
|
||||
}
|
||||
|
||||
b--;
|
||||
|
||||
while(b < 0) {
|
||||
b += size;
|
||||
}
|
||||
}
|
||||
|
||||
acc += m_samples[a] * m_taps[i];
|
||||
m_ptr++;
|
||||
|
||||
while(m_ptr >= size) {
|
||||
m_ptr -= size;
|
||||
}
|
||||
|
||||
return acc;
|
||||
}
|
||||
|
||||
private:
|
||||
std::vector<Real> m_taps;
|
||||
std::vector<Type> m_samples;
|
||||
int m_ptr;
|
||||
};
|
||||
|
||||
#endif // INCLUDE_LOWPASS_H
|
Loading…
Reference in New Issue
Block a user