diff --git a/src/dsp/filter/channelizer/PolyphaseChannelizer.java b/src/dsp/filter/channelizer/PolyphaseChannelizer.java new file mode 100644 index 00000000..daf257a6 --- /dev/null +++ b/src/dsp/filter/channelizer/PolyphaseChannelizer.java @@ -0,0 +1,223 @@ +/******************************************************************************* + * sdrtrunk + * Copyright (C) 2014-2017 Dennis Sheirer + * + * 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, either 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 for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see + * + ******************************************************************************/ +package dsp.filter.channelizer; + +import dsp.filter.design.FilterDesignException; +import dsp.filter.fir.FIRFilterSpecification; +import dsp.filter.fir.real.RealFIRFilter; +import dsp.filter.fir.remez.RemezFIRFilterDesigner; +import dsp.mixer.Oscillator; +import org.jtransforms.fft.FloatFFT_1D; +import org.slf4j.Logger; +import org.slf4j.LoggerFactory; +import sample.complex.ComplexBuffer; + +import java.text.DecimalFormat; + +public class PolyphaseChannelizer +{ + private final static Logger mLog = LoggerFactory.getLogger(PolyphaseChannelizer.class); + + private RealFIRFilter[] mFilters; + private int mFilterPointer = 0; + private int mChannelCount; + private int mBlockSize; + private float[] mFilteredSamples; + private boolean mToggle = false; + private FloatFFT_1D mFFT; + + private DecimalFormat DECIMAL_FORMAT = new DecimalFormat("0.0000"); + + /** + * Creates a polyphase filter bank that divides the input frequency band into evenly spaced channels that + * are oversampled by 2x for output. + * + * @param taps of a low-pass filter designed for the inbound sample rate with a cutoff frequency + * equal to the channel bandwidth (sample rate / filters). If you need to synthesize (combine two or more + * filte outputs) a new bandwidth signal from the outputs of this filter, then the filter should be designed + * as a nyquist filter with -6 dB attenuation at the channel bandwidth cutoff frequency + * + * @param channels - number of filters/channels to output. Since this filter bank oversamples each filter + * output, this number must be even (divisible by the oversample rate). + */ + public PolyphaseChannelizer(float[] taps, int channels) + { + assert(channels % 2 == 0); + + mChannelCount = channels; + + initFilters(taps); + } + + public void filter(float inphase, float quadrature) + { + int iIndex = 2 * mFilterPointer; + int qIndex = iIndex + 1; + + mFilteredSamples[iIndex] = mFilters[iIndex].filter(inphase); + mFilteredSamples[qIndex] = mFilters[qIndex].filter(quadrature); + + if(mFilterPointer == 0) + { + calculate(); + mFilterPointer = mChannelCount - 1; + mToggle = !mToggle; + } + else if(mFilterPointer == mBlockSize) + { + calculate(); + mFilterPointer--; + } + else + { + mFilterPointer--; + } + } + + public void filter(ComplexBuffer buffer) + { + float[] samples = buffer.getSamples(); + + for(int x = 0; x < samples.length; x += 2) + { + filter(samples[x], samples[x + 1]); + } + } + + private void calculate() + { + float[] samples = new float[mFilters.length]; + + if(mToggle) + { + System.arraycopy(mFilteredSamples, 0, samples, mChannelCount, mChannelCount); + System.arraycopy(mFilteredSamples, mChannelCount, samples, 0, mChannelCount); + } + else + { + System.arraycopy(mFilteredSamples, 0, samples, 0, mFilteredSamples.length); + } + + //FFT is executed in-place where the output overwrites the input + mFFT.complexInverse(samples, true); + + dispatch(samples); + } + + private void dispatch(float[] channels) + { + StringBuilder sb = new StringBuilder(); + + for(int x = 0; x < channels.length; x += 2) + { + sb.append((x/2)).append("-") + .append(DECIMAL_FORMAT.format(magnitude(channels[x], channels[x + 1], mChannelCount))).append(" "); + } + + mLog.debug("Out: " + sb.toString()); + } + + private double magnitude(float real, float imaginary, int fftSize) + { + return Math.sqrt(Math.pow(real, 2.0) + Math.pow(imaginary, 2.0)); + } + + /** + * Distributes the filter taps to each of the polyphase filters + * + * @param taps + */ + private void initFilters(float[] taps) + { + int tapCount = (int)Math.ceil((double)taps.length / (double)mChannelCount); + + float[][] coefficients = new float[mChannelCount][tapCount]; + + for(int filter = 0; filter < mChannelCount; filter++) + { + for(int tap = 0; tap < tapCount; tap++) + { + int index = filter * mChannelCount + tap; + + if(index < taps.length) + { + coefficients[filter][tap] = taps[index]; + } + } + } + + mFilters = new RealFIRFilter[mChannelCount * 2]; + + for(int filter = 0; filter < mChannelCount; filter++) + { + //Inphase filter + mFilters[2 * filter] = new RealFIRFilter(coefficients[filter], 1.0f); + + //Quadrature filter + mFilters[2 * filter + 1] = new RealFIRFilter(coefficients[filter], 1.0f); + } + + mFilteredSamples = new float[mFilters.length]; + mFFT = new FloatFFT_1D(mChannelCount); + mBlockSize = mChannelCount / 2; + mFilterPointer = mChannelCount -1; + } + + public static void main(String[] args) + { + mLog.debug("Starting ..."); + + FIRFilterSpecification specification = FIRFilterSpecification.lowPassBuilder() + .sampleRate(4000) + .gridDensity(16) + .passBandCutoff(950) + .passBandAmplitude(1.0) + .passBandRipple(0.01) + .stopBandStart(1050) + .stopBandAmplitude(0.0) + .stopBandRipple(0.005) + .build(); + + RemezFIRFilterDesigner filterDesigner = new RemezFIRFilterDesigner(specification); + + if(filterDesigner.isValid()) + { + try + { + float[] taps = filterDesigner.getImpulseResponse(); + + PolyphaseChannelizer channelizer = new PolyphaseChannelizer(taps, 8); + + Oscillator oscillator = new Oscillator(-1000, 8000); + + for(int x = 0; x < 2000; x++) + { + channelizer.filter(oscillator.getComplex().inphase(), oscillator.getComplex().quadrature()); + oscillator.rotate(); + } + } + catch(FilterDesignException fde) + { + mLog.error("Couldn't design filter to specification"); + } + } + + mLog.debug("Finished!"); + } +} diff --git a/src/dsp/filter/design/FilterViewer.java b/src/dsp/filter/design/FilterViewer.java index 3044ae3c..5526d824 100644 --- a/src/dsp/filter/design/FilterViewer.java +++ b/src/dsp/filter/design/FilterViewer.java @@ -1,65 +1,64 @@ package dsp.filter.design; -import org.slf4j.Logger; -import org.slf4j.LoggerFactory; - import dsp.filter.fir.FIRFilterSpecification; import dsp.filter.fir.remez.RemezFIRFilterDesigner; import javafx.application.Application; import javafx.scene.Scene; import javafx.stage.Stage; +import org.slf4j.Logger; +import org.slf4j.LoggerFactory; public class FilterViewer extends Application { - private final static Logger mLog = LoggerFactory.getLogger( FilterViewer.class ); + private final static Logger mLog = LoggerFactory.getLogger(FilterViewer.class); - /** - * Developer tool to visualize filter designs - */ - public FilterViewer() - { - } - - @Override - public void start( Stage primaryStage ) throws Exception - { - Scene scene = new Scene( new FilterView( getFilter() ) ); - - primaryStage.setTitle( "Filter Viewer" ); - primaryStage.setScene( scene ); - primaryStage.show(); - } - - public static void main( String[] args ) - { - launch( args ); - } + /** + * Developer tool to visualize filter designs + */ + public FilterViewer() + { + } - /** - * Provides the filter to visualize. Modify this method to visualize your filter design. - */ - private float[] getFilter() - { -// FIRFilterSpecification specification = FIRFilterSpecification.lowPassBuilder() -// .sampleRate( 48000 ) -// .gridDensity( 16 ) -// .passBandCutoff( 3200 ) -// .passBandAmplitude( 1.0 ) -// .passBandRipple( 0.02 ) -// .stopBandStart( 4000 ) + @Override + public void start(Stage primaryStage) throws Exception + { + Scene scene = new Scene(new FilterView(getFilter())); + + primaryStage.setTitle("Filter Viewer"); + primaryStage.setScene(scene); + primaryStage.show(); + } + + public static void main(String[] args) + { + launch(args); + } + + /** + * Provides the filter to visualize. Modify this method to visualize your filter design. + */ + private float[] getFilter() + { + FIRFilterSpecification specification = FIRFilterSpecification.lowPassBuilder() + .sampleRate(4000) + .gridDensity(16) + .passBandCutoff(950) + .passBandAmplitude(1.0) + .passBandRipple(0.01) + .stopBandStart(1050) + .stopBandAmplitude(0.0) + .stopBandRipple(0.005) + .build(); + +// FIRFilterSpecification specification = FIRFilterSpecification.highPassBuilder() +// .sampleRate( 24000 ) +// .stopBandCutoff( 800 ) // .stopBandAmplitude( 0.0 ) // .stopBandRipple( 0.03 ) +// .passBandStart( 1000 ) +// .passBandAmplitude( 1.0 ) +// .passBandRipple( 0.08 ) // .build(); - - FIRFilterSpecification specification = FIRFilterSpecification.highPassBuilder() - .sampleRate( 24000 ) - .stopBandCutoff( 800 ) - .stopBandAmplitude( 0.0 ) - .stopBandRipple( 0.03 ) - .passBandStart( 1000 ) - .passBandAmplitude( 1.0 ) - .passBandRipple( 0.08 ) - .build(); // FIRFilterSpecification specification = FIRFilterSpecification.bandPassBuilder() // .sampleRate( 48000 ) @@ -71,24 +70,24 @@ public class FilterViewer extends Application // .passRipple( 0.008) // .build(); - float[] taps = null; + float[] taps = null; // taps = FilterFactory.getLowPass( 16000, 2400, 73, WindowType.HANNING ); - try - { - RemezFIRFilterDesigner designer = new RemezFIRFilterDesigner( specification ); - - if( designer.isValid() ) - { - taps = designer.getImpulseResponse(); - } - } - catch( FilterDesignException fde ) - { - mLog.error( "Filter design error", fde ); - } - - return taps; - } + try + { + RemezFIRFilterDesigner designer = new RemezFIRFilterDesigner(specification); + + if(designer.isValid()) + { + taps = designer.getImpulseResponse(); + } + } + catch(FilterDesignException fde) + { + mLog.error("Filter design error", fde); + } + + return taps; + } }