Filter out CTCSS tones in NFMDemod

Readme.md

@@ -105,14 +105,15 @@ Done since the fork
   - New plugin for BladeRF interfacing libbladeRF directly
   - Corrected the nasty audio band resampling bug preventing use of sample rates that are not power of 2 multiples of 48kHz. This was because the resampling ratio was calculated with an integer division instead of a float division. 
   - As a consequence of the above added more interesting values for the available sampling rates of the BladeRF plugin
+  - Variable span for the SSB demod down to 1.5 kHz
+  - Filter out CTCSS tones in NFMDemod
     
 =====
 To Do
 =====
 
   - Enhance WFM (stereo, RDS?)
-  - Possibility to completely undock the receiver in a separate window. Useful when there are many receivers
-  - Even larger decimation capability for narrowband and very narrowband work? (64, ...)
   - Even more demods ... 
   - Triggering capability like on expensive spectrum analyzers to trap burst signals
   - recording capability
+  - Tx channels for Rx/Tx boards like BladeRF

include-gpl/dsp/bandpass.h

@@ -0,0 +1,120 @@
+#ifndef INCLUDE_BANDPASS_H
+#define INCLUDE_BANDPASS_H
+
+#define _USE_MATH_DEFINES
+#include <math.h>
+#include "dsp/dsptypes.h"
+
+template <class Type> class Bandpass {
+public:
+	Bandpass() { }
+
+	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];
+
+		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;
+
+		m_samples[m_ptr] = sample;
+
+		while(b < 0) {
+			b += m_samples.size();
+		}
+
+		for(i = 0; i < (int)m_taps.size() - 1; i++)
+		{
+			acc +=  (m_samples[a] + m_samples[b]) * m_taps[i];
+			a++;
+
+			while(a >= (int)m_samples.size()) {
+				a -= m_samples.size();
+			}
+
+			b--;
+
+			while(b < 0) {
+				b += m_samples.size();
+			}
+		}
+
+		acc += m_samples[a] * m_taps[i];
+
+		m_ptr++;
+
+		while(m_ptr >= (int)m_samples.size()) {
+			m_ptr -= m_samples.size();
+		}
+
+		return acc;
+	}
+
+private:
+	std::vector<Real> m_taps;
+	std::vector<Type> m_samples;
+	int m_ptr;
+};
+
+#endif // INCLUDE_BANDPASS_H

include-gpl/dsp/highpass.h

@@ -0,0 +1,90 @@
+#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];
+		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;
+
+		m_samples[m_ptr] = sample;
+
+		while(b < 0)
+			b += m_samples.size();
+
+		for(i = 0; i < (int)m_taps.size() - 1; i++) {
+			acc +=  (m_samples[a] + m_samples[b]) * m_taps[i];
+			a++;
+			while(a >= (int)m_samples.size())
+				a -= m_samples.size();
+			b--;
+			while(b < 0)
+				b += m_samples.size();
+		}
+		acc += m_samples[a] * m_taps[i];
+
+		m_ptr++;
+		while(m_ptr >= (int)m_samples.size())
+			m_ptr -= m_samples.size();
+
+		return acc;
+	}
+
+private:
+	std::vector<Real> m_taps;
+	std::vector<Type> m_samples;
+	int m_ptr;
+};
+
+#endif // INCLUDE_HIGHPASS_H

plugins/channel/nfm/nfmdemod.cpp

@@ -147,17 +147,11 @@ void NFMDemod::feed(SampleVector::const_iterator begin, SampleVector::const_iter
 
 					// AF processing
 
-					demod = m_lowpass.filter(demod);
+					//demod = m_lowpass.filter(demod);
-
+					demod = m_bandpass.filter(demod);
-					/*
-					if(demod < -1)
-						demod = -1;
-					else if(demod > 1)
-						demod = 1;
-					*/
-
 					demod *= m_running.m_volume;
-					sample = demod * 32700;
+					//sample = demod * 32700;
+					sample = demod * ((1<<16)/301); // denominator = bandpass filter number of taps
 
 				} else {
 					m_AGC.close();
@@ -245,7 +239,8 @@ void NFMDemod::apply()
 
 	if((m_config.m_afBandwidth != m_running.m_afBandwidth) ||
 		(m_config.m_audioSampleRate != m_running.m_audioSampleRate)) {
-		m_lowpass.create(21, m_config.m_audioSampleRate, m_config.m_afBandwidth);
+		//m_lowpass.create(21, m_config.m_audioSampleRate, m_config.m_afBandwidth);
+		m_bandpass.create(301, m_config.m_audioSampleRate, 300.0, m_config.m_afBandwidth);
 	}
 
 	if(m_config.m_squelch != m_running.m_squelch) {

plugins/channel/nfm/nfmdemod.h

@@ -23,6 +23,7 @@
 #include "dsp/nco.h"
 #include "dsp/interpolator.h"
 #include "dsp/lowpass.h"
+#include "dsp/bandpass.h"
 #include "dsp/movingaverage.h"
 #include "dsp/agc.h"
 #include "audio/audiofifo.h"
@@ -111,7 +112,8 @@ private:
 	Interpolator m_interpolator;
 	Real m_interpolatorDistance;
 	Real m_interpolatorDistanceRemain;
-	Lowpass<Real> m_lowpass;
+	//Lowpass<Real> m_lowpass;
+	Bandpass<Real> m_bandpass;
 
 	Real m_squelchLevel;
 	int m_squelchState;