Mercurial > hg > qm-dsp
changeset 156:edb86e0d850c
Lower filter cutoff to below target Nyquist when downsampling
| author | Chris Cannam |
|---|---|
| date | Fri, 01 Nov 2013 12:07:08 +0000 |
| parents | 529f42b50d81 |
| children | 152abaf17c62 |
| files | dsp/rateconversion/Resampler.cpp dsp/rateconversion/Resampler.h tests/TestResampler.cpp |
| diffstat | 3 files changed, 16 insertions(+), 14 deletions(-) [+] |
line wrap: on
line diff
--- a/dsp/rateconversion/Resampler.cpp Tue Oct 22 08:58:28 2013 +0100 +++ b/dsp/rateconversion/Resampler.cpp Fri Nov 01 12:07:08 2013 +0000 @@ -50,7 +50,7 @@ } // peakToPole -> length -> beta -> window -static map<int, map<int, map<double, vector<double> > > > +static map<double, map<int, map<double, vector<double> > > > knownFilters; static Mutex @@ -63,11 +63,15 @@ int lower = std::min(m_sourceRate, m_targetRate); m_gcd = MathUtilities::gcd(lower, higher); + m_peakToPole = higher / m_gcd; - int peakToPole = higher / m_gcd; + if (m_targetRate < m_sourceRate) { + // antialiasing filter, should be slightly below nyquist + m_peakToPole = m_peakToPole / (1.0 - bandwidth/2.0); + } KaiserWindow::Parameters params = - KaiserWindow::parametersForBandwidth(snr, bandwidth, peakToPole); + KaiserWindow::parametersForBandwidth(snr, bandwidth, higher / m_gcd); params.length = (params.length % 2 == 0 ? params.length + 1 : params.length); @@ -80,21 +84,21 @@ vector<double> filter; knownFilterMutex.lock(); - if (knownFilters[peakToPole][m_filterLength].find(params.beta) == - knownFilters[peakToPole][m_filterLength].end()) { + if (knownFilters[m_peakToPole][m_filterLength].find(params.beta) == + knownFilters[m_peakToPole][m_filterLength].end()) { KaiserWindow kw(params); - SincWindow sw(m_filterLength, peakToPole * 2); + SincWindow sw(m_filterLength, m_peakToPole * 2); filter = vector<double>(m_filterLength, 0.0); for (int i = 0; i < m_filterLength; ++i) filter[i] = 1.0; sw.cut(filter.data()); kw.cut(filter.data()); - knownFilters[peakToPole][m_filterLength][params.beta] = filter; + knownFilters[m_peakToPole][m_filterLength][params.beta] = filter; } - filter = knownFilters[peakToPole][m_filterLength][params.beta]; + filter = knownFilters[m_peakToPole][m_filterLength][params.beta]; knownFilterMutex.unlock(); int inputSpacing = m_targetRate / m_gcd; @@ -309,10 +313,7 @@ std::cerr << "process: buf siz " << m_buffer.size() << " filt siz for phase " << m_phase << " " << m_phaseData[m_phase].filter.size() << std::endl; #endif - double scaleFactor = 1.0; - if (m_targetRate < m_sourceRate) { - scaleFactor = double(m_targetRate) / double(m_sourceRate); - } + double scaleFactor = (double(m_targetRate) / m_gcd) / m_peakToPole; while (outidx < maxout && m_buffer.size() >= m_phaseData[m_phase].filter.size() + m_bufferOrigin) {
--- a/dsp/rateconversion/Resampler.h Tue Oct 22 08:58:28 2013 +0100 +++ b/dsp/rateconversion/Resampler.h Fri Nov 01 12:07:08 2013 +0000 @@ -75,6 +75,7 @@ int m_filterLength; int m_bufferLength; int m_latency; + double m_peakToPole; struct Phase { int nextPhase;
--- a/tests/TestResampler.cpp Tue Oct 22 08:58:28 2013 +0100 +++ b/tests/TestResampler.cpp Fri Nov 01 12:07:08 2013 +0000 @@ -197,9 +197,9 @@ outSpectrum[i] /= outrate; } - // Don't compare bins any higher than 99% of Nyquist freq of lower sr + // Don't compare bins any higher than 96% of Nyquist freq of lower sr int lengthOfInterest = (inrate < outrate ? inrate : outrate) / 2; - lengthOfInterest = lengthOfInterest - (lengthOfInterest / 100); + lengthOfInterest = lengthOfInterest - (lengthOfInterest / 25); for (int i = 0; i < lengthOfInterest; ++i) { BOOST_CHECK_SMALL(inSpectrum[i] - outSpectrum[i], 1e-7);