Updates pull request from Roman Tsypuk to convert the build system from

ant to gradle.  Removes linux and windows application launch scripts and
produces an executable jar with all resources included.
This commit is contained in:
Denny
2017-12-17 07:28:03 -05:00
parent a1535c4157
commit 9316599dec
981 changed files with 8368 additions and 9049 deletions

View File

@@ -0,0 +1,77 @@
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp;
import io.github.dsheirer.sample.Listener;
public class NRZDecoder implements Listener<Boolean>
{
public final static boolean MODE_NORMAL = true;
public final static boolean MODE_INVERTED = false;
private Listener<Boolean> mListener;
private boolean mMode = MODE_NORMAL;
private boolean mPrevious;
/**
* Non Return To Zero (NRZ) and NRZ Inverted (NRZ-I) Decoder.
*
* Use default constructor for NRZ output. Use the alternate constructor
* with the MODE_INVERTED argument for NRZ-I output.
*
* Performs XOR of incoming bit with previous output and sends the result
* to the registered listener.
*/
public NRZDecoder( boolean mode )
{
mMode = mode;
}
public NRZDecoder()
{
this( MODE_NORMAL );
}
public void dispose()
{
mListener = null;
}
@Override
public void receive( Boolean bit )
{
boolean result = mPrevious ^ bit;
if( mListener != null )
{
mListener.receive( mMode ? result : !result );
}
mPrevious = result;
}
public void setListener( Listener<Boolean> listener )
{
mListener = listener;
}
public void removeListener( Listener<Boolean> listener )
{
mListener = null;
}
}

View File

@@ -0,0 +1,22 @@
package io.github.dsheirer.dsp.am;
public abstract class AMDemodulator
{
private float mGain;
public AMDemodulator( float gain )
{
mGain = gain;
}
public float demodulate( float inphase, float quadrature )
{
return (float)Math.sqrt( ( inphase * inphase ) +
( quadrature * quadrature ) ) * mGain;
}
public void setGain( float gain )
{
mGain = gain;
}
}

View File

@@ -0,0 +1,69 @@
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014,2015 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp.am;
import io.github.dsheirer.sample.Listener;
import io.github.dsheirer.sample.complex.ComplexBuffer;
import io.github.dsheirer.sample.real.RealBuffer;
public class AMDemodulator_CB extends AMDemodulator implements Listener<ComplexBuffer>
{
private Listener<RealBuffer> mListener;
public AMDemodulator_CB( float gain )
{
super( gain );
}
@Override
public void receive( ComplexBuffer buffer )
{
if( mListener != null )
{
float[] samples = buffer.getSamples();
int half = samples.length / 2;
float[] demodulated = new float[ half ];
for( int x = 0; x < half; x++ )
{
int index = x * 2;
demodulated[ x ] = demodulate( samples[ index ], samples[ index + 1 ] );
}
mListener.receive( new RealBuffer( demodulated ) );
}
}
public void setListener( Listener<RealBuffer> listener )
{
mListener = listener;
}
public void removeListener()
{
mListener = null;
}
public void dispose()
{
mListener = null;
}
}

View File

@@ -0,0 +1,76 @@
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp.filter;
import io.github.dsheirer.sample.Listener;
import java.util.BitSet;
/**
* Averages (low pass filter) boolean values over the defined size
*/
public class BooleanAveragingFilter implements Listener<Boolean>
{
private BitSet mBuffer = new BitSet();
private int mBufferPointer;
private int mBufferLength;
private int mThreshold;
private Listener<Boolean> mListener;
public BooleanAveragingFilter( int length )
{
mBufferLength = length;
/* Round up */
mThreshold = (int)( length / 2) + ( length % 2 );
}
public void receive( Boolean value )
{
if( value )
{
mBuffer.set( mBufferPointer );
}
else
{
mBuffer.clear( mBufferPointer );
}
mBufferPointer++;
if( mBufferPointer >= mBufferLength )
{
mBufferPointer = 0;
}
if( mListener != null )
{
mListener.receive( mBuffer.cardinality() >= mThreshold );
}
}
public void setListener( Listener<Boolean> listener )
{
mListener = listener;
}
public void removeListener( Listener<Boolean> listener )
{
mListener = null;
}
}

View File

@@ -0,0 +1,55 @@
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp.filter;
import io.github.dsheirer.sample.Listener;
import io.github.dsheirer.sample.complex.Complex;
@Deprecated
public abstract class ComplexFilter implements Listener<Complex>
{
private Listener<Complex> mListener;
public void setListener( Listener<Complex> listener )
{
mListener = listener;
}
public void dispose()
{
mListener = null;
}
public Listener<Complex> getListener()
{
return mListener;
}
public boolean hasListener()
{
return mListener != null;
}
protected void send( Complex sample )
{
if( mListener != null )
{
mListener.receive( sample );
}
}
}

View File

@@ -0,0 +1,760 @@
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014,2015 Dennis Sheirer
*
* Root Raised Cosine filter designer:
* Copyright 2002,2007,2008,2012,2013 Free Software Foundation, Inc.
* http://gnuradio.org/redmine/projects/gnuradio/repository/changes/gr-filter
* /lib/firdes.cc?rev=435b1d166f0c7092bbd5e1f788e75dbb6ade3a4b
*
* 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp.filter;
import io.github.dsheirer.dsp.filter.cic.ComplexPrimeCICDecimate;
import io.github.dsheirer.dsp.filter.design.FilterDesignException;
import io.github.dsheirer.dsp.filter.fir.FIRFilterSpecification;
import io.github.dsheirer.dsp.filter.fir.remez.RemezFIRFilterDesigner;
import org.jtransforms.fft.FloatFFT_1D;
import org.slf4j.Logger;
import org.slf4j.LoggerFactory;
import java.util.Arrays;
import java.util.Set;
import java.util.TreeSet;
public class FilterFactory
{
private final static Logger mLog =
LoggerFactory.getLogger(FilterFactory.class);
/**
* Generates coefficients for a unity-gain, windowed low-pass filter
*
* @param sampleRate - hertz
* @param frequency - cutoff in hertz
* @param length - filter length
* @param window - to apply against the coefficients
* @return
*/
public static float[] getSinc(int sampleRate,
long frequency,
int length,
Window.WindowType window)
{
//Ensure we have an odd length
assert (length % 2 == 0);
//Get unity response array (one element longer to align with IDFT size)
float[] frequencyResponse = getUnityResponseArray(sampleRate, frequency, length + 1);
//Apply Inverse DFT against frequency response unity values, leaving the
//IDFT bin results in the frequency response array
FloatFFT_1D idft = new FloatFFT_1D(length + 1);
idft.realInverseFull(frequencyResponse, false);
//Transfer the IDFT results to the return array
float[] coefficients = new float[length];
int middleCoefficient = (int) (length / 2);
//Bin 0 of the idft is our center coefficient
coefficients[middleCoefficient] = frequencyResponse[0];
//The remaining idft bins from 1 to (middle - 1) are the mirror image
//coefficients
for(int x = 1; x < middleCoefficient; x++)
{
coefficients[middleCoefficient + x] = frequencyResponse[2 * x];
coefficients[middleCoefficient - x] = frequencyResponse[2 * x];
}
//Apply the window against the coefficients
coefficients = Window.apply(window, coefficients);
//Normalize to unity (1) gain
coefficients = normalize(coefficients);
return coefficients;
}
/**
* Normalizes all filter coefficients to achieve unity (1) gain, by ensuring
* that the sum of the absolute value of all coefficients adds up to 1.
*
* @param coefficients
* @return
*/
private static float[] normalize(float[] coefficients)
{
float accumulator = 0;
for(int x = 0; x < coefficients.length; x++)
{
accumulator += Math.abs(coefficients[x]);
}
for(int x = 0; x < coefficients.length; x++)
{
coefficients[x] = coefficients[x] / accumulator;
}
return coefficients;
}
/**
* Constructs an array of unity (1) response values representing the
* desired (pre-windowing) frequency response, used by the sync function
*
* Returns an array twice the length, filled with unity (1) response values
* in the desired pass-band with the positive frequency response starting at
* the lower end of the array, and the negative frequency response at the
* higher end of the first half of the array. The remaining zero-valued
* indexes in the second half will store the results of the JTransforms
* inverse DFT operation
*
* @param sampleRate
* @param frequency
* @param length
* @return
*/
public static float[] getUnityResponseArray(int sampleRate,
long frequency,
int length)
{
float[] unityArray = new float[length * 2];
int binCount = (int) ((Math.round(
(float) frequency / (float) sampleRate * (float) length)));
if(length % 2 == 0) //even length
{
for(int x = 0; x < binCount; x++)
{
unityArray[x] = 1.0f;
unityArray[length - 1 - x] = 1.0f;
}
}
else //odd length
{
unityArray[0] = 1.0f;
for(int x = 1; x <= binCount; x++)
{
unityArray[x] = 1.0f;
unityArray[length - x] = 1.0f;
}
}
return unityArray;
}
public static float[] getSine(float sampleRate, float frequency, int length)
{
float[] retVal = new float[length];
float radianFrequency = (float) (2.0d * Math.PI * (frequency / sampleRate));
int middle = (int) (length / 2);
for(int x = 0; x < middle; x++)
{
float val = (float) Math.sin(radianFrequency * x);
retVal[middle + x] = val;
retVal[middle - x] = -val;
}
return retVal;
}
/**
* Applies a repeating sequence of 1, -1 to the coefficients to invert
* the frequency response of the filter. Used in converting a low-pass
* filter to a high-pass filter.
*
* Inverts the sign of all odd index coefficients ( 1, 3, 5, etc.)
* returning:
* Index 0: same
* Index 1: inverted
* Index 2: same
* ...
* Index length - 1: same
*/
public static float[] invert(float[] coefficients)
{
for(int x = 1; x < coefficients.length; x += 2)
{
coefficients[x] = -coefficients[x];
}
return coefficients;
}
/**
* Generates filter coefficients for a unity-gain, odd-length, windowed,
* low pass filter with passband from 0-hertz to the cutoff frequency.
*
* @param sampleRate - hertz
* @param cutoff - frequency in hertz
* @param filterLength - odd filter length
* @param windowType - window to apply against the generated coefficients
* @return
*/
public static float[] getLowPass(int sampleRate,
long cutoff,
int filterLength,
Window.WindowType windowType)
{
if(filterLength % 2 == 0) //even length
{
float[] values = getSinc(sampleRate, cutoff, filterLength + 2, windowType);
//throw away the 0 index and the last index
return Arrays.copyOfRange(values, 1, values.length - 2);
}
else
{
float[] values = getSinc(sampleRate, cutoff, filterLength + 1, windowType);
//throw away the 0 index
return Arrays.copyOfRange(values, 1, values.length);
}
}
/**
* Creates a low-pass filter with ~ 0.1 db ripple in the pass band
*
* Note: stop frequency - pass frequency defines the transition band.
*/
/**
* Creates a low-pass filter with ~0.1 db ripple in the pass band. The
* transition region is defined by the stop frequency minus the pass
* frequency.
*
* Requires:
* - passFrequency < stopFrequency
* - stopFrequency <= sampleRate/2
*/
public static float[] getLowPass(int sampleRate,
int passFrequency,
int stopFrequency,
int attenuation,
Window.WindowType windowType,
boolean forceOddLength)
{
if(stopFrequency < passFrequency || stopFrequency > (sampleRate / 2))
{
throw new IllegalArgumentException("FilterFactory - low pass "
+ "filter pass frequency [" + passFrequency + "] must be "
+ "less than the stop frequency [" + stopFrequency + "] "
+ "and must be less than half [" + (int) (sampleRate / 2) +
"] of the sample rate [" + sampleRate + "]");
}
int tapCount = getTapCount(sampleRate,
passFrequency,
stopFrequency,
attenuation);
if(forceOddLength)
{
if(tapCount % 2 == 0)
{
tapCount--;
}
}
return getLowPass(sampleRate, passFrequency, tapCount, windowType);
}
/**
* Generates filter coefficients for a unity-gain, odd-length, windowed,
* high pass filter with passband from cutoff frequency to half the sample
* rate.
*
* @param sampleRate - hertz
* @param cutoff - frequency in hertz
* @param filterLength - odd filter length
* @param windowType - window to apply against the generated coefficients
* @return
*/
public static float[] getHighPass(int sampleRate,
long cutoff,
int filterLength,
Window.WindowType windowType)
{
//Convert the high frequency cutoff to its low frequency cutoff
//equivalent, so that when we generate the low pass filter, prior to
//inversion, its at the correct frequency
long convertedCutoff = sampleRate / 2 - cutoff;
return invert(getSinc(sampleRate,
convertedCutoff,
filterLength,
windowType));
}
public static float[] getHighPass(int sampleRate,
long stopFrequency,
long passFrequency,
int attenuation,
Window.WindowType windowType,
boolean forceOddLength)
{
/* reverse the stop and pass frequency to get the low pass variant */
int tapCount = getTapCount(sampleRate,
stopFrequency,
passFrequency,
attenuation);
if(forceOddLength)
{
if(tapCount % 2 == 0)
{
tapCount--;
}
}
return invert(getLowPass(sampleRate, stopFrequency, tapCount, windowType));
}
/**
* Utility to log the arrays of doubles with line breaks
*/
public static String arrayToString(float[] array, boolean breaks)
{
StringBuilder sb = new StringBuilder();
for(int x = 0; x < array.length; x++)
{
sb.append(x + ": " + array[x]);
if(breaks)
{
sb.append("\n");
}
else
{
if(x % 8 == 7)
{
sb.append("\n");
}
else
{
sb.append("\t");
}
}
}
return sb.toString();
}
/**
* Determines the number of fir filter taps required to produce the
* specified frequency response with passband ripple near .1dB.
*
* Implements the algorithm from Understanding Digital Signal Processing, 3e
* , Lyons, section 5.10.5.
*
* @param sampleRate in hertz
* @param pass pass frequency in hertz
* @param stop stop frequency in hertz
* @param attenuation in dB
* @return
*/
public static int getTapCount(int sampleRate,
long pass,
long stop,
int attenuation)
{
double frequency = ((double) stop - (double) pass) / (double) sampleRate;
return (int) (Math.round((double) attenuation / (22.0d * frequency)));
}
public static ComplexPrimeCICDecimate getDecimationFilter(int sampleRate,
int decimatedRate, int order, int passFrequency, int attenuation,
Window.WindowType windowType)
{
int decimationRate = (int) (sampleRate / decimatedRate);
return new ComplexPrimeCICDecimate(decimationRate, order,
passFrequency, attenuation, windowType);
}
/**
* Determines decimation rate(s) for a polyphase decimation filter.
*
* @param sampleRate - starting sample rate
* @param decimatedRate - final (decimated) output rate
* @return - set of integer decimation rates for a single or multi-stage
* polyphase filter decimation chain
* @throws - AssertionException if sample rate is not a multiple of 48 kHz
*/
public static int[] getPolyphaseDecimationRates(int sampleRate,
int decimatedRate,
long passFrequency,
long stopFrequency)
{
int[] rates;
if(sampleRate % decimatedRate != 0)
{
throw new IllegalArgumentException("Decimated rate must be an "
+ "integer multiple of sample rate");
}
int decimation = (int) (sampleRate / decimatedRate);
//Decimation rates below 20 will use single stage polyphase filter
if(decimation < 20)
{
rates = new int[1];
rates[0] = decimation;
return rates;
}
else
{
int optimalStage1 =
getOptimalStageOneRate(sampleRate, decimation, passFrequency, stopFrequency);
Set<Integer> factors = getFactors(decimation);
int stage1 = findClosest(optimalStage1, factors);
// mLog.info( "Decimation rate [" + decimation +
// "] stage1 optimal [" + optimalStage1 +
// "] stage1 actual [" + stage1 +
// "]");
if(stage1 == decimation || stage1 == 1)
{
rates = new int[1];
rates[0] = decimation;
return rates;
}
else
{
rates = new int[2];
int stage2 = (int) (decimation / stage1);
if(stage1 > stage2)
{
rates[0] = stage1;
rates[1] = stage2;
}
else
{
rates[0] = stage2;
rates[1] = stage1;
}
return rates;
}
}
}
/**
* Finds the factor that is closest to the desired factor, from an ordered
* list of factors.
*/
private static int findClosest(int desiredFactor, Set<Integer> factors)
{
int bestFactor = 1;
int bestDelta = desiredFactor;
for(Integer factor : factors)
{
int testDelta = Math.abs(desiredFactor - factor);
if(testDelta < bestDelta)
{
bestDelta = testDelta;
bestFactor = factor;
}
}
return bestFactor;
}
/**
* Determines the factors that make up an integer. Uses a brute force
* method to iterate all integers from 1 to value, determining which factors
* are evenly divisible into the value.
*
* @param value - integer decimation value
* @return - ordered set of factors for value
*/
private static Set<Integer> getFactors(int value)
{
Set<Integer> factors = new TreeSet<Integer>();
/* Brute force */
for(int x = 1; x <= value; x++)
{
int remainder = (int) (value / x);
if(remainder * x == value)
{
factors.add(x);
}
}
return factors;
}
/**
* Determines the optimal decimation rate for stage 1 of a two-stage
* poly-phase decimation filter chain, to produce a final sample rate of
* 48 kHz using a pass bandwidth of 25 kHz.
*
* Use for total decimation rates of 20 or higher.
*
* Implements the algorithm described in Lyons, Understanding Digital Signal
* Processing, 3e, section 10.2.1, page 511.
*
* @param sampleRate
* @param decimation
* @param passFrequency frequency of the pass band
* @return optimum integer decimation rate for the first stage decimation
* filter
*/
public static int getOptimalStageOneRate(int sampleRate,
int decimation,
long passFrequency,
long stopFrequency)
{
double ratio = getBandwidthRatio(passFrequency, stopFrequency);
double numerator = 1.0d - (
Math.sqrt((double) decimation * ratio / (2.0d - ratio)));
double denominator = 2.0d - (ratio * (decimation + 1.0d));
int retVal = (int) (2.0d * decimation * (numerator / denominator));
// mLog.info( "Optimal Stage 1 Decimation - rate [" + sampleRate +
// "] pass [" + passFrequency +
// "] bw ratio [" + ratio +
// "] optimal [" + retVal +
// "]" );
return retVal;
}
/**
* Determines the F ratio as described in Lyons, Understanding Digital
* Signal Processing, 3e, section 10.2.1, page 511
*
* Used in conjunction with the optimal stage one decimation rate method
* above.
*/
private static double getBandwidthRatio(long passFrequency, long stopFrequency)
{
assert (passFrequency < stopFrequency);
return ((double) (stopFrequency - passFrequency) /
(double) stopFrequency);
}
/**
* Assumes that the pass band is 1/4 of the output sample rate.
*
* Assumes the stop band is: pass + (pass * .25).
*
* @param outputSampleRate
* @return
*/
public static float[] getCICCleanupFilter(int outputSampleRate,
int passFrequency,
int attenuation,
Window.WindowType window)
{
int taps = getTapCount(outputSampleRate, passFrequency, passFrequency + 1500,
attenuation);
/* Make tap count odd */
if(taps % 2 == 0)
{
taps++;
}
float[] frequencyResponse =
getCICResponseArray(outputSampleRate, passFrequency, taps);
//Apply Inverse DFT against frequency response unity values, leaving the
//IDFT bin results in the frequency response array
FloatFFT_1D idft = new FloatFFT_1D(taps);
idft.realInverseFull(frequencyResponse, false);
//Transfer the IDFT results to the odd length return array
float[] coefficients = new float[taps];
int middleCoefficient = (int) (taps / 2);
//Bin 0 of the idft is our center coefficient
coefficients[middleCoefficient] = frequencyResponse[0];
//The remaining idft bins from 1 to (middle - 1) are the mirror image
//coefficients
for(int x = 1; x <= middleCoefficient; x++)
{
coefficients[middleCoefficient + x] = frequencyResponse[2 * x];
coefficients[middleCoefficient - x] = frequencyResponse[2 * x];
}
//Apply the window against the coefficients
// coefficients = Window.apply( window, coefficients );
normalize(coefficients);
return coefficients;
}
public static float[] getCICResponseArray(int sampleRate,
int frequency,
int length)
{
float[] cicArray = new float[length * 2];
int binCount = (int) ((Math.round(
(double) frequency / (double) sampleRate * 2.0d * (double) length)));
cicArray[0] = 1.0f;
float unityResponse = (float) (Math.sin(1.0d / (double) length) /
(1.0d / (double) length));
for(int x = 1; x <= binCount; x++)
{
/* Calculate unity response plus amplification for droop */
float compensated = 1.0f + (unityResponse -
(float) (Math.sin((double) x / (double) length) /
((double) x / (double) length)));
cicArray[x] = compensated;
cicArray[length - x] = compensated;
}
return cicArray;
}
/**
* Creates root raised cosine filter coefficients with a tap count equal
* to the symbols x samplesPerSymbol + 1.
*
* Symbol count should be an even integer.
*
* Ported to java from gnuradio/filter/firdes.cc
*
* For 40db attenuation, calculate the number of symbols based on the
* following formula:
*
* Symbols = -44 * alpha + 33
*
* @param samplesPerSymbol - number of samples per symbol
* @param symbols - number of symbols - must be even
* @param alpha - roll-off factor
* @return - filter coefficients
*/
public static float[] getRootRaisedCosine(int samplesPerSymbol,
int symbols,
float alpha)
{
int taps = samplesPerSymbol * symbols + 1;
float scale = 0;
float[] coefficients = new float[taps];
for(int x = 0; x < taps; x++)
{
float index = (float) x - ((float) taps / 2.0f);
float x1 = (float) Math.PI * index / (float) samplesPerSymbol;
float x2 = 4.0f * alpha * index / (float) samplesPerSymbol;
float x3 = x2 * x2 - 1.0f;
float numerator, denominator;
if(Math.abs(x3) >= 0.000001)
{
if(x != taps / 2)
{
numerator = (float) Math.cos((1.0 + alpha) * x1) +
(float) Math.sin((1.0f - alpha) * x1) /
(4.0f * alpha * index / (float) samplesPerSymbol);
}
else
{
numerator = (float) Math.cos((1.0f + alpha) * x1) +
(1.0f - alpha) * (float) Math.PI / (4.0f * alpha);
}
denominator = x3 * (float) Math.PI;
}
else
{
if(alpha == 1.0f)
{
coefficients[x] = -1.0f;
continue;
}
x3 = (1.0f - alpha) * x1;
x2 = (1.0f + alpha) * x1;
numerator = (float) (Math.sin(x2) * (1.0f + alpha) * Math.PI -
Math.cos(x3) * ((1.0 - alpha) * Math.PI *
(double) samplesPerSymbol) / (4 * alpha * index) +
Math.sin(x3) * (double) samplesPerSymbol *
(double) samplesPerSymbol / (4.0 * alpha * index * index));
denominator = (float) (-32.0 * Math.PI * alpha * alpha * index /
(double) samplesPerSymbol);
}
coefficients[x] = 4.0f * alpha * numerator / denominator;
scale += coefficients[x];
}
/**
* Normalize (scale) coefficients to 1.0 sum and apply gain
*/
for(int x = 0; x < taps; x++)
{
coefficients[x] = coefficients[x] / scale;
}
return coefficients;
}
/**
* Creates a filter from the filter specification
*
* @param specification
* @return
* @throws FilterDesignException if the filter cannot be designed
*/
public static float[] getTaps(FIRFilterSpecification specification) throws FilterDesignException
{
RemezFIRFilterDesigner designer = new RemezFIRFilterDesigner(specification);
if(designer.isValid())
{
return designer.getImpulseResponse();
}
return null;
}
}

File diff suppressed because it is too large Load Diff

View File

@@ -0,0 +1,178 @@
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp.filter;
import io.github.dsheirer.sample.real.RealSampleListener;
import org.apache.commons.lang3.Validate;
import java.util.ArrayList;
public class FloatFIRFilter implements RealSampleListener
{
private RealSampleListener mListener;
private ArrayList<Float> mBuffer;
private int mBufferSize = 1; //Temporary initial value
private int mBufferPointer = 0;
private float[] mCoefficients;
private int[][] mIndexMap;
private int mCenterCoefficient;
private int mCenterCoefficientMapIndex;
private float mGain;
public FloatFIRFilter( float[] coefficients, float gain )
{
mCoefficients = coefficients;
mBuffer = new ArrayList<Float>();
mBufferSize = mCoefficients.length;
mGain = gain;
//Fill the buffer with zero valued samples, so we don't have to check for null
for( int x = 0; x < mCoefficients.length; x++ )
{
mBuffer.add( 0.0f );
}
generateIndexMap( mCoefficients.length );
}
public void dispose()
{
mListener = null;
mBuffer.clear();
}
public int getTapCount()
{
return mCoefficients.length;
}
public void setListener( RealSampleListener listener )
{
mListener = listener;
}
public void receive( float newSample )
{
send( get( newSample ) );
}
public float get( float newSample )
{
//Add the new sample to the buffer
mBuffer.set( mBufferPointer, newSample );
//Increment & Adjust the buffer pointer for circular wrap around
mBufferPointer++;
if( mBufferPointer >= mBufferSize )
{
mBufferPointer = 0;
}
//Convolution - multiply filter coefficients by the circular buffer
//samples to calculate a new filtered value
float accumulator = 0.0f;
//Start with the center tap value
accumulator += mCoefficients[ mCenterCoefficient ] *
mBuffer.get( mIndexMap[ mBufferPointer ][ mCenterCoefficientMapIndex ] );
//For the remaining coefficients, add the symmetric samples, oldest and newest
//first, then multiply by the single coefficient
for( int x = 0; x < mCenterCoefficient; x++ )
{
accumulator += mCoefficients[ x ] *
( mBuffer.get( mIndexMap[ mBufferPointer ][ x ] ) +
mBuffer.get( mIndexMap[ mBufferPointer ][ x + mCenterCoefficient ] ) );
}
//We're almost finished ... apply gain, cast the doubles to shorts and
//return the value
return accumulator * mGain;
}
/**
* Dispatch the filtered sample to all registered listeners
*/
private void send( float sample )
{
if( mListener != null )
{
mListener.receive( sample );
}
}
public String printIndexMap()
{
StringBuilder sb = new StringBuilder();
for( int x = 0; x < mCoefficients.length; x++ )
{
for( int y = 0; y < mCoefficients.length; y++ )
{
sb.append( mIndexMap[ x ][ y ] + " " );
}
sb.append( "\n" );
}
return sb.toString();
}
/**
* @param odd-sized number of filter taps (ie coefficients) and buffer
*/
private void generateIndexMap( int size )
{
//Ensure we have an odd size
Validate.isTrue(size % 2 == 1);
mIndexMap = new int[ size ][ size ];
//Last column will be the center coefficient index value for each row
mCenterCoefficientMapIndex = size - 1;
mCenterCoefficient = (int)( size / 2 );
//Setup the first row. Offset is the first row's center coefficient value
//and will become the index value we place in the last column
mIndexMap[ 0 ][ mCenterCoefficient ] = mCenterCoefficient;
//Indexes 0 to 1/2 are their index value, and indexes 1/2 to end are
//an offset of the first half of the values
for( int x = 0; x < mCenterCoefficient; x++ )
{
mIndexMap[ 0 ][ x ] = x;
mIndexMap[ 0 ][ x + mCenterCoefficient ] = size - 1 - x;
}
//For each subsequent map row, increment the value of the preceding row
//same column by 1, wrapping to zero when we exceed the size value
for( int x = 1; x < size; x++ )
{
for( int y = 0; y < size; y++ )
{
mIndexMap[ x ][ y ] = mIndexMap[ x - 1 ][ y ] + 1;
if( mIndexMap[ x ][ y ] >= size )
{
mIndexMap[ x ][ y ] = 0;
}
}
}
}
}

View File

@@ -0,0 +1,204 @@
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp.filter;
import io.github.dsheirer.sample.real.RealSampleListener;
import org.apache.commons.lang3.Validate;
import java.util.ArrayList;
public class FloatHalfBandFilter implements RealSampleListener
{
private RealSampleListener mListener;
private ArrayList<Float> mBuffer;
private int mBufferSize = 1; //Temporary initial value
private int mBufferPointer = 0;
private float mGain;
private float[] mCoefficients;
private int[][] mIndexMap;
private int mCenterCoefficient;
private int mCenterCoefficientMapIndex;
private boolean mDispatchFlag = false;
/**
* Half-Band filter with decimation by 2 against real valued floats.
*
* Takes advantage of the symmetrical nature of FIR filter coefficients by
* adding oldest and newest sample first, then multiplying once by the
* corresponding coefficient
*
* Also, takes advantage of the 0-valued FIR half-band coefficents inherent
* in the half-band filter, and does not calculate those coefficients.
*
* This reduces the workload to (tap-size - 1) / 4 + 1 calculations per sample.
*
* @param filter - filter coefficients
* @param gain - gain multiplier. Use 1.0 for unity/no gain
*/
public FloatHalfBandFilter( Filters filter, float gain )
{
mCoefficients = filter.getCoefficients();
mBuffer = new ArrayList<Float>();
mBufferSize = mCoefficients.length;
//Fill the buffer with zero valued samples
for( int x = 0; x < mCoefficients.length; x++ )
{
mBuffer.add( 0.0f );
}
generateIndexMap( mCoefficients.length );
mGain = gain;
}
public void dispose()
{
mListener = null;
}
/**
* Calculate the filtered value by applying the coefficients against
* the complex samples in mBuffer
*/
public void receive( float newSample )
{
//Add the new sample to the buffer
mBuffer.set( mBufferPointer, newSample );
//Increment & Adjust the buffer pointer for circular wrap around
mBufferPointer++;
if( mBufferPointer >= mBufferSize )
{
mBufferPointer = 0;
}
//Toggle the flag every time, so that we can calculate and dispatch a
//sample every other time, when the flag is true (ie decimate by 2)
mDispatchFlag = !mDispatchFlag;
//Calculate a filtered sample when the flag is true
if( mDispatchFlag )
{
//Convolution - multiply filter coefficients by the circular buffer
//samples to calculate a new filtered value
float accumulator = 0;
//Start with the center tap value
accumulator += mCoefficients[ mCenterCoefficient ] *
mBuffer.get( mIndexMap[ mBufferPointer ][ mCenterCoefficientMapIndex ] );
//For the remaining coefficients, add the symmetric samples, oldest and newest
//first, then multiply by the single coefficient
for( int x = 0; x < mCenterCoefficientMapIndex; x += 2 )
{
accumulator += mCoefficients[ x ] *
( mBuffer.get( mIndexMap[ mBufferPointer ][ x ] ) +
mBuffer.get( mIndexMap[ mBufferPointer ][ x + 1 ] ) );
}
//We're almost finished ... apply gain, cast the floats to shorts and
//send it on it's merry way
if( mListener != null )
{
mListener.receive( (float)( accumulator * mGain ) );
}
}
}
/**
* Creates an n X (n + 1 / 2) index map enabling quick access to the
* circular buffer samples.
*
* As the buffer shifts right with each subsequent sample, we have to move
* the index pointers with it, for efficient access of the samples.
*
* The first array index value in the index map corresponds to the current
* buffer pointer location.
*
* The second array index value points to the samples that should be
* multiplied by the coefficients as follows:
*
* 0 = center tap sample, to be multiplied by center coefficient
*
* 0 = sample( 1 )
* 1 = sample( size - 1 )
*
* Indexes 0 and 1 will be multiplied by coefficient( 0 ).
*
* Subsequent indexes 3, 4, etc, point to the oldest and newest samples that
* correspond to the matching ( 3 ) coefficient index.
*
* @param odd-sized number of filter taps (ie coefficients) and buffer
*/
private void generateIndexMap( int size )
{
//Ensure we have an odd size
Validate.isTrue(size % 2 == 1);
int mapWidth = ( ( size + 1 ) / 2 ) + 1;
//Set center tap index for coefficients array
mCenterCoefficient = ( size - 1 ) / 2;
mCenterCoefficientMapIndex = mCenterCoefficient + 1;
mIndexMap = new int[ size ][ mapWidth ];
//Setup the first row, buffer pointer index 0, as a starting point
for( int x = 0; x < mapWidth - 2; x += 2 )
{
mIndexMap[ 0 ][ x ] = x;
mIndexMap[ 0 ][ x + 1 ] = size - 1 - x;
}
//Place center tap index in last element
mIndexMap[ 0 ][ mCenterCoefficientMapIndex ] = mCenterCoefficient;
//For each subsequent row, increment the previous row's value by 1,
//subtracting size as needed, to keep the values between 0 and size - 1
for( int x = 1; x < size; x++ )
{
for( int y = 0; y < mapWidth; y++ )
{
mIndexMap[ x ][ y ] = mIndexMap[ x - 1 ][ y ] + 1;
if( mIndexMap[ x ][ y ] >= size )
{
mIndexMap[ x ][ y ] -= size;
}
}
}
}
/**
* Registers a listener for filtered samples
*/
public void setListener( RealSampleListener listener )
{
mListener = listener;
}
/**
* Removes (if exists) a registered filtered sample listener
*/
public void clearListener()
{
mListener = null;
}
}

View File

@@ -0,0 +1,194 @@
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp.filter;
import io.github.dsheirer.sample.real.RealSampleListener;
import org.apache.commons.lang3.Validate;
import java.util.ArrayList;
public class FloatHalfBandNoDecimateFilter implements RealSampleListener
{
private RealSampleListener mListener;
private ArrayList<Float> mBuffer;
private int mBufferSize = 1; //Temporary initial value
private int mBufferPointer = 0;
private float mGain;
private float[] mCoefficients;
private int[][] mIndexMap;
private int mCenterCoefficient;
private int mCenterCoefficientMapIndex;
/**
* Half-Band filter with decimation by a real valued float.
*
* Takes advantage of the symmetrical nature of FIR filter coefficients by
* adding oldest and newest sample first, then multiplying once by the
* corresponding coefficient
*
* Also, takes advantage of the 0-valued FIR half-band coefficents inherent
* in the half-band filter, and does not calculate those coefficients.
*
* This reduces the workload to (tap-size - 1) / 4 + 1 calculations per sample.
*
* @param filter - filter coefficients
* @param gain - gain multiplier. Use 1.0 for unity/no gain
*/
public FloatHalfBandNoDecimateFilter( Filters filter, float gain )
{
mCoefficients = filter.getCoefficients();
mBuffer = new ArrayList<Float>();
mBufferSize = mCoefficients.length;
//Fill the buffer with zero valued samples
for( int x = 0; x < mCoefficients.length; x++ )
{
mBuffer.add( 0.0f );
}
generateIndexMap( mCoefficients.length );
mGain = gain;
}
public void dispose()
{
mListener = null;
}
/**
* Calculate the filtered value by applying the coefficients against
* the complex samples in mBuffer
*/
public void receive( float newSample )
{
//Add the new sample to the buffer
mBuffer.set( mBufferPointer, newSample );
//Increment & Adjust the buffer pointer for circular wrap around
mBufferPointer++;
if( mBufferPointer >= mBufferSize )
{
mBufferPointer = 0;
}
//Convolution - multiply filter coefficients by the circular buffer
//samples to calculate a new filtered value
float accumulator = 0;
//Start with the center tap value
accumulator += mCoefficients[ mCenterCoefficient ] *
mBuffer.get( mIndexMap[ mBufferPointer ][ mCenterCoefficientMapIndex ] );
//For the remaining coefficients, add the symmetric samples, oldest and newest
//first, then multiply by the single coefficient
for( int x = 0; x < mCenterCoefficientMapIndex; x += 2 )
{
accumulator += mCoefficients[ x ] *
( mBuffer.get( mIndexMap[ mBufferPointer ][ x ] ) +
mBuffer.get( mIndexMap[ mBufferPointer ][ x + 1 ] ) );
}
//We're almost finished ... apply gain, cast the floats to floats and
//send it on it's merry way
if( mListener != null )
{
mListener.receive( (float)( accumulator * mGain ) );
}
}
/**
* Creates an n X (n + 1 / 2) index map enabling quick access to the
* circular buffer samples.
*
* As the buffer shifts right with each subsequent sample, we have to move
* the index pointers with it, for efficient access of the samples.
*
* The first array index value in the index map corresponds to the current
* buffer pointer location.
*
* The second array index value points to the samples that should be
* multiplied by the coefficients as follows:
*
* 0 = center tap sample, to be multiplied by center coefficient
*
* 0 = sample( 1 )
* 1 = sample( size - 1 )
*
* Indexes 0 and 1 will be multiplied by coefficient( 0 ).
*
* Subsequent indexes 3, 4, etc, point to the oldest and newest samples that
* correspond to the matching ( 3 ) coefficient index.
*
* @param odd-sized number of filter taps (ie coefficients) and buffer
*/
private void generateIndexMap( int size )
{
//Ensure we have an odd size
Validate.isTrue(size % 2 == 1);
int mapWidth = ( ( size + 1 ) / 2 ) + 1;
//Set center tap index for coefficients array
mCenterCoefficient = ( size - 1 ) / 2;
mCenterCoefficientMapIndex = mCenterCoefficient + 1;
mIndexMap = new int[ size ][ mapWidth ];
//Setup the first row, buffer pointer index 0, as a starting point
for( int x = 0; x < mapWidth - 2; x += 2 )
{
mIndexMap[ 0 ][ x ] = x;
mIndexMap[ 0 ][ x + 1 ] = size - 1 - x;
}
//Place center tap index in last element
mIndexMap[ 0 ][ mCenterCoefficientMapIndex ] = mCenterCoefficient;
//For each subsequent row, increment the previous row's value by 1,
//subtracting size as needed, to keep the values between 0 and size - 1
for( int x = 1; x < size; x++ )
{
for( int y = 0; y < mapWidth; y++ )
{
mIndexMap[ x ][ y ] = mIndexMap[ x - 1 ][ y ] + 1;
if( mIndexMap[ x ][ y ] >= size )
{
mIndexMap[ x ][ y ] -= size;
}
}
}
}
/**
* Registers a listener for filtered samples
*/
public void setListener( RealSampleListener listener )
{
mListener = listener;
}
/**
* Removes (if exists) a registered filtered sample listener
*/
public void clearListener()
{
mListener = null;
}
}

View File

@@ -0,0 +1,161 @@
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp.filter;
import io.github.dsheirer.dsp.filter.Window.WindowType;
/**
* Implements the Goertzel Filter as described in
* http://en.wikipedia.org/wiki/Goertzel_algorithm
*
* Calculates the magnitude of a specific frequency component within a window of
* samples, using a Goertzel Filter (ie an optimized DFT).
*/
public class GoertzelFilter
{
private int mSampleRate;
private long mTargetFrequency;
private int mBlockSize;
private double[] mWindowCoefficients;
private double mCoefficient;
/**
* Constructs the Goertzel Filter to use for block processing of
* time-domain signal samples.
*
* @param sampleRate
* - sample rate (Hz)
* @param targetFrequency
* - target frequency (Hz)
* @param blockSize
* - number of samples per block of processing
* @param window
* - window to apply against the samples prior to calculating
* the magnitude
*/
public GoertzelFilter( int sampleRate,
long targetFrequency,
int blockSize,
WindowType window )
{
mSampleRate = sampleRate;
mTargetFrequency = targetFrequency;
mBlockSize = blockSize;
mWindowCoefficients = Window.getWindow( window, blockSize );
init();
}
/**
* Establish the pre-calculated values to use in the filter
*/
private void init()
{
double normalizedFrequency = (double) mTargetFrequency / mSampleRate;
mCoefficient = 2.0D * Math.cos( 2 * Math.PI * normalizedFrequency );
}
/**
* Returns the power (dB) of the target frequency within the block of
* samples as normalized against the bin(0) magnitude value
*
* @param samples
* - array of time-domain samples -- array size must be the same
* size as the windowSize parameter passed upon construction
*
* @return - power measurement in dB
* @throws IllegalArgumentException if the sample array size is not equal
* to the defined block size
*/
public int getPower( float[] samples ) throws IllegalArgumentException
{
// Verify size of samples array against block size
if ( samples.length != mBlockSize )
{
throw new IllegalArgumentException(
"Sample array size does not equal Block Size" );
}
// Apply the window against the samples
samples = Window.apply( mWindowCoefficients, samples );
// Process the samples
double s = 0.0D;
double s_prev = 0.0D;
double s_prev2 = 0.0D;
for ( int x = 0; x < samples.length; x++ )
{
s = samples[x] + ( mCoefficient * s_prev ) - s_prev2;
s_prev2 = s_prev;
s_prev = s;
}
//power = s_prev2 * s_prev2 + s_prev * s_prev - coeff * s_prev * s_prev2 ;
double magnitude = ( s_prev2 * s_prev2 ) + ( s_prev * s_prev) - (mCoefficient * s_prev * s_prev2);
int binZero = getBinZeroPower( samples );
int power = (int) (20 * Math.log10( magnitude / binZero ) );
return power;
}
/**
* Sums all of the sample values to derive the bin 0 power
*/
private int getBinZeroPower( float[] samples )
{
int retVal = 0;
for( int x = 0; x < samples.length; x++ )
{
retVal += samples[ x ];
}
return retVal;
}
/**
* @return the Sample Rate
*/
public int getSampleRate()
{
return mSampleRate;
}
/**
* @return the Target Frequency
*/
public long getTargetFrequency()
{
return mTargetFrequency;
}
/**
* @return the Block Size
*/
public int getBlockSize()
{
return mBlockSize;
}
}

View File

@@ -0,0 +1,205 @@
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp.filter;
/***************************************************************************
* Copyright (C) 2011 by Paul Lutus *
* lutusp@arachnoid.com *
* *
* 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 2 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, write to the *
* Free Software Foundation, Inc., *
* 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. *
***************************************************************************/
// http://en.wikipedia.org/wiki/Digital_biquad_filter
final public class IIRBiQuadraticFilter {
public enum Type {
BANDPASS, LOWPASS, HIGHPASS, NOTCH, PEAK, LOWSHELF, HIGHSHELF
};
double a0, a1, a2, b0, b1, b2;
double x1, x2, y, y1, y2;
double gain_abs;
Type type;
double center_freq, sample_rate, Q, gainDB;
public IIRBiQuadraticFilter() {
}
public IIRBiQuadraticFilter(Type type, double center_freq, double sample_rate, double Q, double gainDB) {
configure(type, center_freq, sample_rate, Q, gainDB);
}
// constructor without gain setting
public IIRBiQuadraticFilter(Type type, double center_freq, double sample_rate, double Q)
{
configure(type, center_freq, sample_rate, Q, 0);
}
public void reset()
{
x1 = x2 = y1 = y2 = 0;
}
public double frequency()
{
return center_freq;
}
public void configure(Type type, double center_freq, double sample_rate, double Q, double gainDB)
{
reset();
Q = (Q == 0) ? 1e-9 : Q;
this.type = type;
this.sample_rate = sample_rate;
this.Q = Q;
this.gainDB = gainDB;
reconfigure(center_freq);
}
public void configure(Type type, double center_freq, double sample_rate, double Q)
{
configure(type, center_freq, sample_rate, Q, 0);
}
// allow parameter change while running
public void reconfigure(double cf) {
center_freq = cf;
// only used for peaking and shelving filters
gain_abs = Math.pow(10, gainDB / 40);
double omega = 2 * Math.PI * cf / sample_rate;
double sn = Math.sin(omega);
double cs = Math.cos(omega);
double alpha = sn / (2 * Q);
double beta = Math.sqrt(gain_abs + gain_abs);
switch (type) {
case BANDPASS:
b0 = alpha;
b1 = 0;
b2 = -alpha;
a0 = 1 + alpha;
a1 = -2 * cs;
a2 = 1 - alpha;
break;
case LOWPASS:
b0 = (1 - cs) / 2;
b1 = 1 - cs;
b2 = (1 - cs) / 2;
a0 = 1 + alpha;
a1 = -2 * cs;
a2 = 1 - alpha;
break;
case HIGHPASS:
b0 = (1 + cs) / 2;
b1 = -(1 + cs);
b2 = (1 + cs) / 2;
a0 = 1 + alpha;
a1 = -2 * cs;
a2 = 1 - alpha;
break;
case NOTCH:
b0 = 1;
b1 = -2 * cs;
b2 = 1;
a0 = 1 + alpha;
a1 = -2 * cs;
a2 = 1 - alpha;
break;
case PEAK:
b0 = 1 + (alpha * gain_abs);
b1 = -2 * cs;
b2 = 1 - (alpha * gain_abs);
a0 = 1 + (alpha / gain_abs);
a1 = -2 * cs;
a2 = 1 - (alpha / gain_abs);
break;
case LOWSHELF:
b0 = gain_abs * ((gain_abs + 1) - (gain_abs - 1) * cs + beta * sn);
b1 = 2 * gain_abs * ((gain_abs - 1) - (gain_abs + 1) * cs);
b2 = gain_abs * ((gain_abs + 1) - (gain_abs - 1) * cs - beta * sn);
a0 = (gain_abs + 1) + (gain_abs - 1) * cs + beta * sn;
a1 = -2 * ((gain_abs - 1) + (gain_abs + 1) * cs);
a2 = (gain_abs + 1) + (gain_abs - 1) * cs - beta * sn;
break;
case HIGHSHELF:
b0 = gain_abs * ((gain_abs + 1) + (gain_abs - 1) * cs + beta * sn);
b1 = -2 * gain_abs * ((gain_abs - 1) + (gain_abs + 1) * cs);
b2 = gain_abs * ((gain_abs + 1) + (gain_abs - 1) * cs - beta * sn);
a0 = (gain_abs + 1) - (gain_abs - 1) * cs + beta * sn;
a1 = 2 * ((gain_abs - 1) - (gain_abs + 1) * cs);
a2 = (gain_abs + 1) - (gain_abs - 1) * cs - beta * sn;
break;
}
// prescale flter constants
b0 /= a0;
b1 /= a0;
b2 /= a0;
a1 /= a0;
a2 /= a0;
}
// provide a static amplitude result for testing
public double result(double f) {
double phi = Math.pow((Math.sin(2.0 * Math.PI * f / (2.0 * sample_rate))), 2.0);
return (Math.pow(b0 + b1 + b2, 2.0) - 4.0 * (b0 * b1 + 4.0 * b0 * b2 + b1 * b2) * phi + 16.0 * b0 * b2 * phi * phi) / (Math.pow(1.0 + a1 + a2, 2.0) - 4.0 * (a1 + 4.0 * a2 + a1 * a2) * phi + 16.0 * a2 * phi * phi);
}
// provide a static decibel result for testing
public double log_result(double f)
{
double r;
try {
r = 10 * Math.log10(result(f));
} catch (Exception e) {
r = -100;
}
if(Double.isInfinite(r) || Double.isNaN(r)) {
r = -100;
}
return r;
}
// return the constant set for this filter
public double[] constants()
{
return new double[]{b0, b1, b2, a1, a2};
}
// perform one filtering step
public double filter(double x)
{
y = b0 * x + b1 * x1 + b2 * x2 - a1 * y1 - a2 * y2;
x2 = x1;
x1 = x;
y2 = y1;
y1 = y;
return (y);
}
}

View File

@@ -0,0 +1,103 @@
package io.github.dsheirer.dsp.filter;
import io.github.dsheirer.sample.Listener;
import java.util.BitSet;
public class LTRPulseShapingFilter implements Listener<Boolean>
{
private int[][] mIndexes = { { 0, 1, 2, 3, 4 },
{ 1, 2, 3, 4, 0 },
{ 2, 3, 4, 0, 1 },
{ 3, 4, 0, 1, 2 },
{ 4, 0, 1, 2, 3 } };
private BitSet mBuffer;
private int mBufferPointer = 0;
private Listener<Boolean> mListener;
/**
* Specialized pulse shaping filter that detects a weak 0 or 1 pulse that
* is characterized by two consecutive values, true or false, within a
* sequence of five boolean values where the remainder are the opposite
* value, and elongates that weak pulse to be three consecutive values,
* which is enough to cause the slicer to render a correct value.
*/
public LTRPulseShapingFilter()
{
mBuffer = new BitSet();
mBuffer.clear();
}
@Override
public void receive( Boolean current )
{
boolean previous = mBuffer.get( mBufferPointer );
if( current )
{
mBuffer.set( mBufferPointer );
}
else
{
mBuffer.clear( mBufferPointer );
}
mBufferPointer++;
if( mBufferPointer >= 5 )
{
mBufferPointer = 0;
}
/* Check for weak positive pulse*/
if( mBuffer.cardinality() == 2 )
{
if( mBuffer.get( mIndexes[ mBufferPointer ][ 1 ] ) &&
mBuffer.get( mIndexes[ mBufferPointer ][ 2 ] ) )
{
mBuffer.set( mIndexes[ mBufferPointer][ 3 ] );
}
}
/* Or, a weak negative pulse */
else if( mBuffer.cardinality() == 3 )
{
if( !mBuffer.get( mIndexes[ mBufferPointer][ 1 ] ) &&
!mBuffer.get( mIndexes[ mBufferPointer][ 2 ] ) )
{
mBuffer.clear( mIndexes[ mBufferPointer][ 3 ] );
}
}
if( mListener != null )
{
mListener.receive( previous );
}
}
public String getIndexSet( int pointer )
{
StringBuilder sb = new StringBuilder();
sb.append( "[" );
for( int x = 0; x < 5; x++ )
{
sb.append( mIndexes[ pointer ][ x ] );
sb.append( "," );
}
sb.append( "[" );
return sb.toString();
}
public void setListener( Listener<Boolean> listener )
{
mListener = listener;
}
public void removeListener( Listener<Boolean> listener )
{
mListener = null;
}
}

View File

@@ -0,0 +1,52 @@
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp.filter;
import io.github.dsheirer.sample.Listener;
import io.github.dsheirer.sample.real.RealSampleListener;
public class SquaringFilter implements RealSampleListener
{
private Listener<Boolean> mListener;
/**
* Transforms a short sample waveform into a square wave
*/
public SquaringFilter()
{
}
@Override
public void receive( float sample )
{
if( mListener != null )
{
mListener.receive( sample >= 0 );
}
}
public void setListener( Listener<Boolean> listener )
{
mListener = listener;
}
public void clearListener( Listener<Boolean> listener )
{
mListener = null;
}
}

View File

@@ -0,0 +1,234 @@
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp.filter;
import java.util.Arrays;
/**
* Window creator
*
* Blackman, Hamming and Hanning windows
* based on windowing methods described at:
* http://www.mstarlabs.com/dsp/goertzel/goertzel.html
* and
* Understanding Digital Signal Processing 3e, Lyons, p.91
* and
* Digital Filters for Everyone, Allred, p. 132
*
*/
public class Window
{
public static double[] getWindow( WindowType type, int length )
{
switch( type )
{
case BLACKMAN:
return getBlackmanWindow( length );
case COSINE:
return getCosineWindow( length );
case HAMMING:
return getHammingWindow( length );
case HANNING:
return getHanningWindow( length );
case NONE:
default:
return getRectangularWindow( length );
}
}
public static double[] getRectangularWindow( int length )
{
double[] coefficients = new double[ length ];
Arrays.fill( coefficients, 1.0d );
return coefficients;
}
public static double[] getCosineWindow( int length )
{
double[] coefficients = new double[ length ];
if( length % 2 == 0 ) //Even length
{
int half = (int)( ( length - 1 ) / 2 );
for( int x = -half; x < length / 2 + 1; x++ )
{
coefficients[ x + half ] = Math.cos(
( (double)x * Math.PI ) / ( (double)length + 1.0d ) );
}
}
else //Odd length
{
int half = (int) length / 2;
for( int x = -half; x < half + 1; x++ )
{
coefficients[ x + half ] = Math.cos(
( (double)x * Math.PI ) / ( (double)length + 1.0d ) );
}
}
return coefficients;
}
public static double[] getBlackmanWindow( int length )
{
double[] coefficients = new double[ length ];
for( int x = 0; x < length; x++ )
{
coefficients[ x ] = .426591D -
( .496561D * Math.cos( ( Math.PI * 2.0D * (double)x ) / (double)( length - 1 ) ) +
( .076848D * Math.cos( ( Math.PI * 4.0D * (double)x ) / (double)( length - 1 ) ) ) );
}
return coefficients;
}
public static double[] getHammingWindow( int length )
{
double[] coefficients = new double[ length ];
if( length % 2 == 0 ) //Even length
{
for( int x = 0; x < length; x++ )
{
coefficients[ x ] = .54D - ( .46D *
Math.cos( ( Math.PI * ( 2.0D * (double)x + 1.0d ) ) / (double)( length - 1 ) ) );
}
}
else //Odd length
{
for( int x = 0; x < length; x++ )
{
coefficients[ x ] = .54D - ( .46D *
Math.cos( ( 2.0D * Math.PI * (double)x ) / (double)( length - 1 ) ) );
}
}
return coefficients;
}
public static double[] getHanningWindow( int length )
{
double[] coefficients = new double[ length ];
if( length % 2 == 0 ) //Even length
{
for( int x = 0; x < length; x++ )
{
coefficients[ x ] = .5D - ( .5D *
Math.cos( ( Math.PI * ( 2.0D * (double)x + 1 ) ) / (double)( length - 1 ) ) );
}
}
else //Odd length
{
for( int x = 0; x < length; x++ )
{
coefficients[ x ] = .5D - ( .5D *
Math.cos( ( 2.0D * Math.PI * (double)x ) / (double)( length - 1 ) ) );
}
}
return coefficients;
}
/**
* Apply the window against an array of float-type samples
*/
public static float[] apply( double[] coefficients, float[] samples )
{
for( int x = 0; x < coefficients.length; x++ )
{
samples[ x ] = (float)( samples[ x ] * coefficients[ x ] );
}
return samples;
}
public static float[] apply( WindowType type, float[] samples )
{
double[] coefficients = getWindow( type, samples.length );
return apply( coefficients, samples );
}
/**
* Apply the window against an array of float-type samples
*/
public static double[] apply( double[] coefficients, double[] samples )
{
for( int x = 0; x < coefficients.length; x++ )
{
samples[ x ] = samples[ x ] * coefficients[ x ];
}
return samples;
}
public static double[] apply( WindowType type, double[] samples )
{
double[] coefficients = getWindow( type, samples.length );
return apply( coefficients, samples );
}
/**
* Window types
*/
public enum WindowType
{
BLACKMAN, COSINE, HAMMING, HANNING, NONE;
}
/**
* Utility to log the double arrays
*/
public static String arrayToString( double[] array, boolean breaks )
{
StringBuilder sb = new StringBuilder();
for( int x = 0; x < array.length; x++ )
{
sb.append( x + ":" + array[ x ] );
if( breaks )
{
sb.append( "\n" );
}
else
{
if( x % 8 == 7 )
{
sb.append( "\n" );
}
else
{
sb.append( "\t" );
}
}
}
return sb.toString();
}
}

View File

@@ -0,0 +1,412 @@
/*******************************************************************************
* 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 <http://www.gnu.org/licenses/>
*
******************************************************************************/
package io.github.dsheirer.dsp.filter.cic;
import io.github.dsheirer.dsp.filter.FilterFactory;
import io.github.dsheirer.dsp.filter.Filters;
import io.github.dsheirer.dsp.filter.Window;
import io.github.dsheirer.dsp.filter.fir.complex.ComplexFIRFilter_CB_CB;
import io.github.dsheirer.dsp.filter.halfband.complex.HalfBandFilter_CB_CB;
import io.github.dsheirer.sample.Listener;
import io.github.dsheirer.sample.complex.ComplexBuffer;
import io.github.dsheirer.sample.complex.ComplexSampleListener;
import io.github.dsheirer.sample.complex.ComplexToComplexBufferAssembler;
import io.github.dsheirer.sample.decimator.ComplexDecimator;
import org.apache.commons.lang3.Validate;
import java.util.ArrayList;
import java.util.List;
public class ComplexPrimeCICDecimate implements Listener<ComplexBuffer>
{
public static int[] PRIMES = {2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53,
59, 61, 67, 71, 73, 79, 83, 89, 97, 101, 103, 107, 109, 113, 127, 131, 137, 139, 149, 151,
157, 163, 167, 173, 179, 181, 191, 193, 197, 199, 211, 223, 227, 229, 233, 239, 241, 251,
257, 263, 269, 271, 277, 281, 283, 293, 307, 311, 313, 317, 331, 337, 347, 349, 353, 359,
367, 373, 379, 383, 389, 397, 401, 409, 419, 421, 431, 433, 439, 443, 449, 457, 461, 463,
467, 479, 487, 491, 499, 503, 509, 521, 523, 541, 547, 557, 563, 569, 571, 577, 587, 593,
599, 601, 607, 613, 617, 619, 631, 641, 643, 647, 653, 659, 661, 673, 677, 683, 691};
private ArrayList<DecimatingStage> mDecimatingStages =
new ArrayList<DecimatingStage>();
private DecimatingStage mFirstDecimatingStage;
private Output mOutput;
/**
* Non-Recursive Prime-Factor CIC Filter with float sample array inputs and
* decimated, single paired i/q float sample output.
*
* Implements the CIC filter described in Understanding Digital Signal
* Processing, 3e, Lyons, on page 769. This filter is comprised of multiple
* decimating stages each with a prime factor decimation rate. Multiple
* stages are cascaded to achieve the overall decimation rate.
*
* This filter supports a maximum decimation rate of 700. This filter can
* be adapted to higher decimation rates by adding additional prime factors
* to the PRIMES array.
*
* @param decimation - overall decimation rate
* @param order - filter order
*/
public ComplexPrimeCICDecimate(int decimation, int order,
int passFrequency, int attenuation, Window.WindowType windowType)
{
Validate.isTrue(decimation <= 700);
List<Integer> stageSizes = getPrimeFactors(decimation);
for(int x = 0; x < stageSizes.size(); x++)
{
DecimatingStage stage = new DecimatingStage(stageSizes.get(x), order);
mDecimatingStages.add(stage);
if(x == 0)
{
/* Reference to first stage -- will receive all samples */
mFirstDecimatingStage = stage;
}
else
{
/* Wire the current stage to the previous stage */
mDecimatingStages.get(x - 1).setListener(stage);
}
}
mOutput = new Output(48000, passFrequency, attenuation, windowType);
mDecimatingStages.get(mDecimatingStages.size() - 1)
.setListener(mOutput);
}
public void dispose()
{
for(DecimatingStage stage : mDecimatingStages)
{
stage.dispose();
}
mDecimatingStages.clear();
mDecimatingStages = null;
mFirstDecimatingStage = null;
mOutput.dispose();
mOutput = null;
}
/**
* Adds a listener to receive the output of this CIC decimation filter
*/
public void setListener(Listener<ComplexBuffer> listener)
{
mOutput.setListener(listener);
}
/**
* Removes listener from output of this CIC decimation filter
*/
public void removeListener()
{
mOutput.removeListener();
}
/**
* Calculates the prime factors of the decimation rate up to a maximum
* decimation rate of 700. If you wish to have decimation rates higher
* than 700, then add additional prime factors to the PRIMES array.
*
* @param decimation - integral decimation rate
* @return - ordered list (smallest to largest) of prime factors
*/
public static List<Integer> getPrimeFactors(int decimation)
{
ArrayList<Integer> stages = new ArrayList<Integer>();
int pointer = 0;
while(decimation > 0 && pointer < PRIMES.length)
{
int prime = PRIMES[pointer];
if(decimation % prime == 0)
{
stages.add(prime);
decimation /= prime;
}
else
{
pointer++;
}
}
return stages;
}
/**
* Primary input method for receiving sample arrays composed as I,Q,I,Q, etc.
*/
@Override
public void receive(ComplexBuffer buffer)
{
if(mFirstDecimatingStage != null)
{
float[] samples = buffer.getSamples();
for(int x = 0; x < samples.length; x += 2)
{
mFirstDecimatingStage.receive(samples[x], samples[x + 1]);
}
}
}
/**
* Decimating stage combines multiple CIC stages with a decimator. The
* number of stages is indicated by the order value and the size indicates
* the decimation rate of this stage.
*/
public class DecimatingStage implements ComplexSampleListener
{
private ArrayList<Stage> mStages = new ArrayList<Stage>();
private Stage mFirstStage;
private ComplexDecimator mDecimator;
public DecimatingStage(int size, int order)
{
for(int x = 0; x < order; x++)
{
Stage stage;
if(size == 2)
{
stage = new TwoStage();
}
else
{
stage = new Stage(size);
}
mStages.add(stage);
if(x == 0)
{
mFirstStage = stage;
}
else
{
mStages.get(x - 1).setListener(stage);
}
}
mDecimator = new ComplexDecimator(size);
mStages.get(mStages.size() - 1).setListener(mDecimator);
}
public void dispose()
{
for(Stage stage : mStages)
{
stage.dispose();
}
mStages.clear();
mDecimator.dispose();
mDecimator = null;
mFirstStage = null;
mStages = null;
}
@Override
public void receive(float i, float q)
{
mFirstStage.receive(i, q);
}
public void setListener(ComplexSampleListener listener)
{
mDecimator.setListener(listener);
}
}
/**
* Single non-decimating CIC stage component. Uses a circular buffer and a
* running average internally to implement the stage so that stage size has
* essentially no impact on the computational requirements of the stage
*/
public class Stage implements ComplexSampleListener
{
protected ComplexSampleListener mListener;
private float[] mISamples;
private float[] mQSamples;
protected float mISum;
protected float mQSum;
private int mSamplePointer = 0;
private int mSize;
protected float mGain;
public Stage()
{
}
public Stage(int size)
{
mSize = size - 1;
mISamples = new float[mSize];
mQSamples = new float[mSize];
mGain = 1.0f / (float)size;
}
public void dispose()
{
mListener = null;
}
public void receive(float i, float q)
{
/* Subtract the oldest sample and add back in the newest */
mISum = mISum - mISamples[mSamplePointer] + i;
mQSum = mQSum - mQSamples[mSamplePointer] + q;
/* Overwrite the oldest sample with the newest */
mISamples[mSamplePointer] = i;
mQSamples[mSamplePointer] = q;
mSamplePointer++;
if(mSamplePointer >= mSize)
{
mSamplePointer = 0;
}
if(mListener != null)
{
mListener.receive((mISum * mGain), (mQSum * mGain));
}
}
public void setListener(ComplexSampleListener listener)
{
mListener = listener;
}
}
/**
* Size 2 stage that removes the unnecessary circular buffer management for
* a two-stage.
*/
public class TwoStage extends Stage
{
public TwoStage()
{
super();
mGain = 0.5f;
}
public void receive(float i, float q)
{
float iSum = mISum + i;
float qSum = mQSum + q;
mISum = i;
mQSum = q;
if(mListener != null)
{
mListener.receive((iSum * mGain), (qSum * mGain));
}
}
}
/**
* Output adapter - applies gain correction and cleanup filter and broadcast
* to registered listener.
*/
public class Output implements ComplexSampleListener
{
/* Decimated output buffers will contain 1024 samples */
private ComplexToComplexBufferAssembler mAssembler =
new ComplexToComplexBufferAssembler(2048);
private ComplexFIRFilter_CB_CB mCleanupFilter;
private HalfBandFilter_CB_CB mHalfBandFilter = new HalfBandFilter_CB_CB(
Filters.FIR_HALF_BAND_31T_ONE_EIGHTH_FCO.getCoefficients(), 0.4f, false);
public Output(int outputSampleRate, int passFrequency, int attenuation,
Window.WindowType windowType)
{
mCleanupFilter = new ComplexFIRFilter_CB_CB(FilterFactory
.getCICCleanupFilter(outputSampleRate,
passFrequency,
attenuation,
windowType), 0.4f);
//Bypassing the CIC cleanup filter for now
mAssembler.setListener(mHalfBandFilter);
// mAssembler.setListener( mCleanupFilter );
// mCleanupFilter.setListener( mHalfBandFilter );
}
public void dispose()
{
mAssembler.dispose();
mCleanupFilter.dispose();
mHalfBandFilter.dispose();
}
/**
* Receiver method for the output adapter to receive a filtered,
* decimated sample, apply gain correction, apply cleanup filtering
* and output the sample values as a complex sample.
*/
@Override
public void receive(float inphase, float quadrature)
{
mAssembler.receive(inphase, quadrature);
}
/**
* Adds a listener to receive output samples
*/
public void setListener(Listener<ComplexBuffer> listener)
{
mHalfBandFilter.setListener(listener);
}
/**
* Removes the listener from receiving output samples
*/
public void removeListener()
{
mHalfBandFilter.removeListener();
}
}
}

View File

@@ -0,0 +1,414 @@
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2015 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp.filter.cic;
import io.github.dsheirer.dsp.filter.FilterFactory;
import io.github.dsheirer.dsp.filter.Filters;
import io.github.dsheirer.dsp.filter.halfband.real.HalfBandFilter_RB_RB;
import io.github.dsheirer.sample.Listener;
import io.github.dsheirer.sample.decimator.RealDecimator;
import io.github.dsheirer.sample.real.RealBuffer;
import io.github.dsheirer.sample.real.RealSampleListener;
import io.github.dsheirer.sample.real.RealToRealBufferAssembler;
import org.apache.commons.lang3.Validate;
import org.slf4j.Logger;
import org.slf4j.LoggerFactory;
import io.github.dsheirer.dsp.filter.Window.WindowType;
import io.github.dsheirer.dsp.filter.fir.real.RealFIRFilter_RB_RB;
import java.util.ArrayList;
import java.util.List;
public class RealPrimeCICDecimate
{
private final static Logger mLog =
LoggerFactory.getLogger( RealPrimeCICDecimate.class );
public static int[] PRIMES = { 2,3,5,7,11,13,17,19,23,29,31,37,41,43,47,53,
59,61,67,71,73,79,83,89,97,101,103,107,109,113,127,131,137,139,149,151,
157,163,167,173,179,181,191,193,197,199,211,223,227,229,233,239,241,251,
257,263,269,271,277,281,283,293,307,311,313,317,331,337,347,349,353,359,
367,373,379,383,389,397,401,409,419,421,431,433,439,443,449,457,461,463,
467,479,487,491,499,503,509,521,523,541,547,557,563,569,571,577,587,593,
599,601,607,613,617,619,631,641,643,647,653,659,661,673,677,683,691 };
private ArrayList<DecimatingStage> mDecimatingStages =
new ArrayList<DecimatingStage>();
private DecimatingStage mFirstDecimatingStage;
private Output mOutput;
private int mOutputBufferSize = 2048;
/**
* Non-Recursive Prime-Factor CIC Filter with float sample array inputs and
* decimated float sample output.
*
* Implements the CIC filter described in Understanding Digital Signal
* Processing, 3e, Lyons, on page 769. This filter is comprised of multiple
* decimating stages each with a prime factor decimation rate. Multiple
* stages are cascaded to achieve the overall decimation rate.
*
* This filter supports a maximum decimation rate of 700. This filter can
* be adapted to higher decimation rates by adding additional prime factors
* to the PRIMES array.
*
* @param decimation - overall decimation rate
* @param order - filter order
*/
public RealPrimeCICDecimate( int decimation, int order, int passFrequency,
int attenuation, WindowType windowType )
{
this( decimation, order, passFrequency, attenuation, windowType, 2048 );
}
public RealPrimeCICDecimate( int decimation, int order, int passFrequency,
int attenuation, WindowType windowType, int outputBufferSize )
{
Validate.isTrue(decimation <= 700);
mOutputBufferSize = outputBufferSize;
List<Integer> stageSizes = getPrimeFactors( decimation );
for( int x = 0; x < stageSizes.size(); x++ )
{
DecimatingStage stage = new DecimatingStage( stageSizes.get( x ), order );
mDecimatingStages.add( stage );
if( x == 0 )
{
/* Reference to first stage -- will receive all samples */
mFirstDecimatingStage = stage;
}
else
{
/* Wire the current stage to the previous stage */
mDecimatingStages.get( x - 1 ).setListener( stage );
}
}
mOutput = new Output( 48000, passFrequency, attenuation, windowType,
outputBufferSize );
mDecimatingStages.get( mDecimatingStages.size() - 1 )
.setListener( mOutput );
}
public void dispose()
{
for( DecimatingStage stage: mDecimatingStages )
{
stage.dispose();
}
mDecimatingStages.clear();
mDecimatingStages = null;
mFirstDecimatingStage = null;
mOutput.dispose();
mOutput = null;
}
/**
* Adds a listener to receive the output of this CIC decimation filter
*/
public void setListener( Listener<RealBuffer> listener )
{
mOutput.setListener( listener );
}
/**
* Removes listener from output of this CIC decimation filter
*/
public void removeListener()
{
mOutput.removeListener();
}
/**
* Calculates the prime factors of the decimation rate up to a maximum
* decimation rate of 700. If you wish to have decimation rates higher
* than 700, then add additional prime factors to the PRIMES array.
*
* @param decimation - integral decimation rate
* @return - ordered list (smallest to largest) of prime factors
*/
public static List<Integer> getPrimeFactors( int decimation )
{
ArrayList<Integer> stages = new ArrayList<Integer>();
int pointer = 0;
while( decimation > 0 && pointer < PRIMES.length )
{
int prime = PRIMES[ pointer ];
if( decimation % prime == 0 )
{
stages.add( prime );
decimation /= prime;
}
else
{
pointer++;
}
}
return stages;
}
/**
* Primary input method for receiving real buffers.
*/
public void receive( RealBuffer buffer )
{
if( mFirstDecimatingStage != null )
{
for( float sample: buffer.getSamples() )
{
mFirstDecimatingStage.receive( sample );
}
}
}
/**
* Decimating stage combines multiple CIC stages with a decimator. The
* number of stages is indicated by the order value and the size indicates
* the decimation rate of this stage.
*/
public class DecimatingStage implements RealSampleListener
{
private ArrayList<Stage> mStages = new ArrayList<Stage>();
private Stage mFirstStage;
private RealDecimator mDecimator;
public DecimatingStage( int size, int order )
{
for( int x = 0; x < order; x++ )
{
Stage stage;
if( size == 2 )
{
stage = new TwoStage();
}
else
{
stage = new Stage( size );
}
mStages.add( stage );
if( x == 0 )
{
mFirstStage = stage;
}
else
{
mStages.get( x - 1 ).setListener( stage );
}
}
mDecimator = new RealDecimator( size );
mStages.get( mStages.size() - 1 ).setListener( mDecimator );
}
public void dispose()
{
for( Stage stage: mStages )
{
stage.dispose();
}
mStages.clear();
mDecimator.dispose();
mDecimator = null;
mFirstStage = null;
mStages = null;
}
@Override
public void receive( float sample )
{
mFirstStage.receive( sample );
}
public void setListener( RealSampleListener listener )
{
mDecimator.setListener( listener );
}
}
/**
* Single non-decimating CIC stage component. Uses a circular buffer and a
* running average internally to implement the stage so that stage size has
* essentially no impact on the computational requirements of the stage
*/
public class Stage implements RealSampleListener
{
protected RealSampleListener mListener;
private float[] mSamples;
protected float mSum;
private int mSamplePointer = 0;
private int mSize;
protected float mGain;
public Stage()
{
}
public Stage( int size )
{
mSize = size - 1;
mSamples = new float[ mSize ];
mGain = 1.0f / (float)size;
}
public void dispose()
{
mListener = null;
}
public void receive( float sample )
{
/* Subtract the oldest sample and add back in the newest */
mSum = mSum - mSamples[ mSamplePointer ] + sample;
/* Overwrite the oldest sample with the newest */
mSamples[ mSamplePointer ] = sample;
mSamplePointer++;
if( mSamplePointer >= mSize )
{
mSamplePointer = 0;
}
if( mListener != null )
{
mListener.receive( mSum * mGain );
}
}
public void setListener( RealSampleListener listener )
{
mListener = listener;
}
}
/**
* Size 2 stage that removes the unnecessary circular buffer management for
* a two-stage.
*/
public class TwoStage extends Stage
{
public TwoStage()
{
super();
mGain = 0.5f;
}
public void receive( float sample )
{
float sum = mSum + sample;
mSum = sample;
if( mListener != null )
{
mListener.receive( sum * mGain );
}
}
}
/**
* Output adapter - applies gain correction, applies cleanup filter, casts
* double-valued samples into float-valued complex sample and sends the
* result to the registered listener.
*/
public class Output implements RealSampleListener
{
/* Decimated output buffers will contain 1024 samples */
private RealToRealBufferAssembler mAssembler;
private RealFIRFilter_RB_RB mCleanupFilter;
private HalfBandFilter_RB_RB mHalfBandFilter = new HalfBandFilter_RB_RB(
Filters.FIR_HALF_BAND_31T_ONE_EIGHTH_FCO.getCoefficients(), 0.4f, false );
public Output( int outputSampleRate, int passFrequency, int attenuation,
WindowType windowType, int outputBufferSize )
{
mCleanupFilter = new RealFIRFilter_RB_RB( FilterFactory
.getCICCleanupFilter( outputSampleRate,
passFrequency,
attenuation,
windowType ), 0.4f );
mAssembler = new RealToRealBufferAssembler( outputBufferSize );
mAssembler.setListener( mCleanupFilter );
mCleanupFilter.setListener( mHalfBandFilter );
}
public void dispose()
{
mAssembler.dispose();
mCleanupFilter.dispose();
mHalfBandFilter.dispose();
}
/**
* Receiver method for the output adapter to receive a filtered,
* decimated sample, apply gain correction, apply cleanup filtering
* and output the sample values as a complex sample.
*/
@Override
public void receive( float sample )
{
mAssembler.receive( sample * 32.0f );
}
/**
* Adds a listener to receive output samples
*/
public void setListener( Listener<RealBuffer> listener )
{
mHalfBandFilter.setListener( listener );
}
/**
* Removes the listener from receiving output samples
*/
public void removeListener()
{
mHalfBandFilter.removeListener();
}
}
}

View File

@@ -0,0 +1,41 @@
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014,2015 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp.filter.dc;
public class DCRemovalFilter
{
protected float mAverage;
protected float mRatio;
public DCRemovalFilter( float ratio )
{
mRatio = ratio;
}
public void reset()
{
mAverage = 0.0f;
}
public float filter( float sample )
{
float filtered = sample - mAverage;
mAverage += mRatio * filtered;
return filtered;
}
}

View File

@@ -0,0 +1,65 @@
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp.filter.dc;
import io.github.dsheirer.sample.real.RealSampleListener;
/**
* Implements the real-time general DC removal algorithm described on page
* 762, figure 13-62d, in Understanding Digital Signal Processing 3e, Lyons.
*/
public class DCRemovalFilterLyons implements RealSampleListener
{
private float mAlpha;
private float mPrevious;
private RealSampleListener mListener;
public DCRemovalFilterLyons( float alpha )
{
mAlpha = alpha;
}
public void dispose()
{
mListener = null;
}
@Override
public void receive( float sample )
{
send( sample + ( mPrevious * mAlpha ) - mPrevious );
mPrevious = sample;
}
/**
* Sends the filtered sample to the listener
*/
private void send( float sample )
{
if( mListener != null )
{
mListener.receive( sample );
}
}
public void setListener( RealSampleListener listener )
{
mListener = listener;
}
}

View File

@@ -0,0 +1,71 @@
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp.filter.dc;
import io.github.dsheirer.sample.Listener;
import io.github.dsheirer.sample.complex.ComplexBuffer;
public class DCRemovalFilter_CB implements Listener<ComplexBuffer>
{
private Listener<ComplexBuffer> mListener;
private DCRemovalFilter mIFilter;
private DCRemovalFilter mQFilter;
public DCRemovalFilter_CB( float ratio )
{
mIFilter = new DCRemovalFilter( ratio );
mQFilter = new DCRemovalFilter( ratio );
}
public void dispose()
{
mListener = null;
}
public float[] filter( float[] samples )
{
for( int x = 0; x < samples.length; x += 2 )
{
samples[ x ] = mIFilter.filter( samples[ x ] );
samples[ x + 1 ] = mQFilter.filter( samples[ x + 1 ] );
}
return samples;
}
public ComplexBuffer filter( ComplexBuffer buffer )
{
filter( buffer.getSamples() );
return buffer;
}
@Override
public void receive( ComplexBuffer buffer )
{
if( mListener != null )
{
mListener.receive( filter( buffer ) );
}
}
public void setListener( Listener<ComplexBuffer> listener )
{
mListener = listener;
}
}

View File

@@ -0,0 +1,61 @@
package io.github.dsheirer.dsp.filter.dc;
import io.github.dsheirer.sample.Listener;
import io.github.dsheirer.sample.real.RealBuffer;
public class DCRemovalFilter_RB extends DCRemovalFilter
implements Listener<RealBuffer>
{
protected Listener<RealBuffer> mListener;
public DCRemovalFilter_RB( float ratio )
{
super( ratio );
}
public void setListener( Listener<RealBuffer> listener )
{
mListener = listener;
}
public void removeListener()
{
mListener = null;
}
public void dispose()
{
mListener = null;
}
public float[] filter( float[] samples )
{
for( int x = 0; x < samples.length; x++ )
{
samples[ x ] = filter( samples[ x ] );
}
return samples;
}
public RealBuffer filter( RealBuffer buffer )
{
float[] samples = buffer.getSamples();
for( int x = 0; x < samples.length; x++ )
{
samples[ x ] = filter( samples[ x ] );
}
return buffer;
}
@Override
public void receive( RealBuffer buffer )
{
if( mListener != null )
{
mListener.receive( filter( buffer ) );
}
}
}

View File

@@ -0,0 +1,54 @@
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp.filter.dc;
import io.github.dsheirer.sample.real.RealSampleListener;
public class DCRemovalFilter_R_R extends DCRemovalFilter implements RealSampleListener
{
RealSampleListener mListener;
public DCRemovalFilter_R_R( float ratio )
{
super( ratio );
}
public void dispose()
{
mListener = null;
}
@Override
public void receive( float sample )
{
if( mListener != null )
{
mListener.receive( filter( sample ) );
}
}
public void setListener( RealSampleListener listener )
{
mListener = listener;
}
public void removeListener()
{
mListener = null;
}
}

View File

@@ -0,0 +1,60 @@
package io.github.dsheirer.dsp.filter.dc;
import io.github.dsheirer.sample.Listener;
import io.github.dsheirer.sample.real.RealBuffer;
public class DCRemovalFilter_RealBuffer implements Listener<RealBuffer>
{
private Listener<RealBuffer> mListener;
/**
* Removes DC bias present across the samples in the buffer by calculating
* the mean DC component and subtracting that component from each of the
* samples.
*/
public DCRemovalFilter_RealBuffer()
{
}
public void dispose()
{
mListener = null;
}
public RealBuffer filter( RealBuffer buffer )
{
double accumulator = 0.0d;
for( float sample: buffer.getSamples() )
{
accumulator += sample;
}
float mean = (float)( accumulator / (double)buffer.getSamples().length );
for( float sample: buffer.getSamples() )
{
sample -= mean;
}
return buffer;
}
@Override
public void receive( RealBuffer buffer )
{
if( mListener != null )
{
mListener.receive( filter( buffer ) );
}
}
public void setListener( Listener<RealBuffer> listener )
{
mListener = listener;
}
public void removeListener()
{
mListener = null;
}
}

View File

@@ -0,0 +1,45 @@
package io.github.dsheirer.dsp.filter.dc;
import io.github.dsheirer.sample.real.RealSampleListener;
public class IIRSinglePoleDCRemovalFilter implements RealSampleListener
{
private float mAlpha;
private float mPreviousInput = 0.0f;
private float mPreviousOutput = 0.0f;
private RealSampleListener mListener;
/**
* IIR single-pole DC removal filter, as described by J M de Freitas in
* 29Jan2007 paper at:
*
* http://www.mathworks.com/matlabcentral/fileexchange/downloads/72134/download
*
* @param alpha 0.0 - 1.0 float - the closer alpha is to unity, the closer
* the cutoff frequency is to DC.
*/
public IIRSinglePoleDCRemovalFilter( float alpha )
{
mAlpha = alpha;
}
@Override
public void receive( float currentInput )
{
float currentOutput = ( currentInput - mPreviousInput ) +
( mAlpha * mPreviousOutput );
if( mListener != null )
{
mListener.receive( currentOutput );
}
mPreviousInput = currentInput;
mPreviousOutput = currentOutput;
}
public void setListener( RealSampleListener listener )
{
mListener = listener;
}
}

View File

@@ -0,0 +1,69 @@
package io.github.dsheirer.dsp.filter.dc;
import io.github.dsheirer.sample.Listener;
import io.github.dsheirer.sample.real.RealBuffer;
public class IIRSinglePoleDCRemovalFilter_RB implements Listener<RealBuffer>
{
private float mAlpha;
private float mPreviousInput = 0.0f;
private float mPreviousOutput = 0.0f;
private Listener<RealBuffer> mListener;
/**
* IIR single-pole DC removal filter, as described by J M de Freitas in
* 29Jan2007 paper at:
*
* http://www.mathworks.com/matlabcentral/fileexchange/downloads/72134/download
*
* @param alpha 0.0 - 1.0 float - the closer alpha is to unity, the closer
* the cutoff frequency is to DC.
*/
public IIRSinglePoleDCRemovalFilter_RB( float alpha )
{
mAlpha = alpha;
}
@Override
public void receive( RealBuffer buffer )
{
if( mListener != null )
{
mListener.receive( filter( buffer ) );
}
}
public RealBuffer filter( RealBuffer buffer )
{
float[] samples = buffer.getSamples();
for( int x = 0; x < samples.length; x++ )
{
samples[ x ] = filter( samples[ x ] );
}
return buffer;
}
public float filter( float sample )
{
float currentOutput = ( sample - mPreviousInput ) +
( mAlpha * mPreviousOutput );
mPreviousInput = sample;
mPreviousOutput = currentOutput;
return currentOutput;
}
public void reset()
{
mPreviousInput = 0.0f;
mPreviousOutput = 0.0f;
}
public void setListener( Listener<RealBuffer> listener )
{
mListener = listener;
}
}

View File

@@ -0,0 +1,16 @@
package io.github.dsheirer.dsp.filter.design;
public class FilterDesignException extends Exception
{
private static final long serialVersionUID = 1L;
public FilterDesignException( String text )
{
super( text );
}
public FilterDesignException( String text, Throwable throwable )
{
super( text, throwable );
}
}

View File

@@ -0,0 +1,39 @@
/* ....Show License.... */
.chart {
-fx-background-image: url("../../../Curve-fitted-background.png");
-fx-padding: 15 25 15 15;
}
.chart-plot-background {
-fx-background-color: red;
-fx-background-image: url("../../../Curve-fitted-chart-background.png"),
url("../../../Curve-fitted-graph-gridlines.png");
-fx-background-size: cover, auto;
-fx-background-repeat: no-repeat, repeat;
-fx-background-position: 0% 0%, 0% 100%;
-fx-border-color: black black transparent transparent;
}
.chart-area-symbol {
-fx-background-color: white;
}
.chart-series-area-line {
-fx-stroke: white;
-fx-stroke-width: 2px;
}
.chart-series-area-fill {
-fx-fill: linear-gradient(to right, white, rgba(255,255,255,0));
-fx-blend-mode: OVERLAY;
}
.axis {
-fx-tick-mark-visible: false;
-fx-minor-tick-visible: false;
-fx-tick-length: 3;
-fx-minor-tick-length: 0;
-fx-border-color: transparent;
-fx-tick-label-fill: white;
}
.axis:bottom {
-fx-border-color: black transparent transparent transparent;
}
.axis:left {
-fx-border-color: transparent black transparent transparent;
}

View File

@@ -0,0 +1,166 @@
package io.github.dsheirer.dsp.filter.design;
import io.github.dsheirer.gui.control.CurveFittedAreaChart;
import javafx.collections.ObservableList;
import javafx.scene.chart.NumberAxis;
import javafx.scene.chart.XYChart;
import javafx.scene.chart.XYChart.Series;
import javafx.scene.layout.BorderPane;
import javafx.scene.layout.HBox;
import javafx.scene.layout.Priority;
import org.jtransforms.fft.FloatFFT_1D;
import org.slf4j.Logger;
import org.slf4j.LoggerFactory;
public class FilterView extends BorderPane
{
private final static Logger mLog = LoggerFactory.getLogger( FilterView.class );
private CurveFittedAreaChart mImpulseChart;
private NumberAxis mImpulseXAxis;
private NumberAxis mImpulseYAxis;
private CurveFittedAreaChart mDFTChart;
private NumberAxis mDFTXAxis;
private NumberAxis mDFTYAxis;
private static final int FFT_SIZE = 4096;
private FloatFFT_1D mFFT = new FloatFFT_1D( FFT_SIZE );
public FilterView( float[] taps )
{
init();
if( taps != null )
{
setImpulseTaps( taps );
float[] dft = new float[ FFT_SIZE ];
System.arraycopy( taps, 0, dft, 0, taps.length );
mFFT.realForward( dft );
setDFTTaps( convertDFTToDecibel( dft, taps.length ) );
}
}
private void init()
{
HBox hbox = new HBox();
mImpulseXAxis = new NumberAxis( 0.0, 0.5, 1.0 );
mImpulseYAxis = new NumberAxis( -0.2, 1.2, .1 );
mImpulseChart = new CurveFittedAreaChart( mImpulseXAxis, mImpulseYAxis );
mImpulseChart.setLegendVisible( false );
mImpulseChart.setHorizontalGridLinesVisible( false );
mImpulseChart.setVerticalGridLinesVisible(false);
mImpulseChart.setAlternativeColumnFillVisible(false);
mImpulseChart.setAlternativeRowFillVisible(false);
String css = FilterView.class.getResource( "FilterView.css" ).toExternalForm();
mImpulseChart.getStylesheets().add( css );
hbox.getChildren().add( mImpulseChart );
HBox.setHgrow( mImpulseChart, Priority.ALWAYS );
mDFTXAxis = new NumberAxis( 0.0, FFT_SIZE / 2, FFT_SIZE / 8 );
mDFTYAxis = new NumberAxis( -200.0, 10.0, 10.0 );
mDFTChart = new CurveFittedAreaChart( mDFTXAxis, mDFTYAxis );
mDFTChart.setLegendVisible( false );
mDFTChart.setHorizontalGridLinesVisible( false );
mDFTChart.setVerticalGridLinesVisible(false);
mDFTChart.setAlternativeColumnFillVisible(false);
mDFTChart.setAlternativeRowFillVisible(false);
mDFTChart.getStylesheets().add( css );
hbox.getChildren().add( mDFTChart );
HBox.setHgrow( mDFTChart, Priority.ALWAYS );
setCenter( hbox );
}
public void setImpulseTaps( double[] taps )
{
float[] convertedTaps = new float[ taps.length ];
for( int x = 0; x < taps.length; x++ )
{
convertedTaps[ x ] = (float)taps[ x ];
}
setImpulseTaps( convertedTaps );
}
public void setImpulseTaps( float[] taps )
{
double smallest = 0.0;
double largest = 0.0;
final XYChart.Series<Number, Number> series = new XYChart.Series<>();
for( int x = 0; x < taps.length; x++ )
{
if( taps[ x ] < smallest )
{
smallest = taps[ x ];
}
if( taps[ x ] > largest )
{
largest = taps[ x ];
}
series.getData().add( new XYChart.Data<Number, Number>( x, taps[ x ] ) );
}
mImpulseXAxis.setLowerBound( 0.0 );
mImpulseXAxis.setUpperBound( taps.length - 1 );
mImpulseXAxis.setTickUnit( 1.0 );
double padding = Math.abs( largest - smallest ) * .1;
mImpulseYAxis.setLowerBound( smallest - padding );
mImpulseYAxis.setUpperBound( largest + padding );
mImpulseYAxis.setTickUnit( padding );
ObservableList<Series<Number,Number>> list = mImpulseChart.getData();
list.add( series );
}
public void setDFTTaps( float[] taps )
{
double smallest = 0.0;
final XYChart.Series<Number, Number> series = new XYChart.Series<>();
for( int x = 0; x < taps.length; x++ )
{
if( taps[ x ] < smallest )
{
smallest = taps[ x ];
}
series.getData().add( new XYChart.Data<Number, Number>( x, taps[ x ] ) );
}
ObservableList<Series<Number,Number>> list = mDFTChart.getData();
list.add( series );
}
public float[] convertDFTToDecibel( float[] dft, int tapLength )
{
float[] decibels = new float[ dft.length / 2 ];
int index = 0;
for( int x = 0; x < decibels.length; x ++ )
{
index = x * 2;
decibels[ x ] = 20.0f * (float)Math.log10( ( ( dft[ index ] * dft[ index ] ) +
( dft[ index + 1 ] * dft[ index + 1 ] ) ) );
}
return decibels;
}
}

View File

@@ -0,0 +1,94 @@
package io.github.dsheirer.dsp.filter.design;
import javafx.application.Application;
import javafx.scene.Scene;
import javafx.stage.Stage;
import org.slf4j.Logger;
import org.slf4j.LoggerFactory;
import io.github.dsheirer.dsp.filter.fir.FIRFilterSpecification;
import io.github.dsheirer.dsp.filter.fir.remez.RemezFIRFilterDesigner;
public class FilterViewer extends Application
{
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);
}
/**
* 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(2500)
.passBandAmplitude(1.0)
.passBandRipple(0.01)
.stopBandStart(4000)
.stopBandAmplitude(0.0)
.stopBandRipple(0.008)
.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 )
// .stopFrequency1( 200 )
// .passFrequencyBegin( 400 )
// .passFrequencyEnd( 3200 )
// .stopFrequency2( 3400 )
// .stopRipple( 0.0003 )
// .passRipple( 0.008)
// .build();
float[] taps = null;
// taps = FilterFactory.getLowPass(48000, 6750, 7500, 60,
// Window.WindowType.HAMMING, true);
try
{
RemezFIRFilterDesigner designer = new RemezFIRFilterDesigner(specification);
if(designer.isValid())
{
taps = designer.getImpulseResponse();
}
}
catch(FilterDesignException fde)
{
mLog.error("Filter design error", fde);
}
return taps;
}
}

View File

@@ -0,0 +1,102 @@
package io.github.dsheirer.dsp.filter.equalizer;
import io.github.dsheirer.sample.Listener;
import io.github.dsheirer.sample.complex.Complex;
import org.slf4j.Logger;
import org.slf4j.LoggerFactory;
public class CMAEqualizer implements Listener<Complex>
{
private final static Logger mLog = LoggerFactory.getLogger( CMAEqualizer.class );
public static final int TAP_COUNT = 11;
private float mModulus;
private float mMu;
private Complex mError = new Complex( 0.0f, 0.0f );
private Complex[] mTaps = new Complex[ TAP_COUNT ];
private Complex[] mBuffer = new Complex[ TAP_COUNT ];
private int mBufferPointer = 0;
private Listener<Complex> mListener;
public CMAEqualizer( float modulus, float mu )
{
mModulus = modulus;
mMu = mu;
/* Set first tap to 1,0 and all other taps to 0,0 as a starting point.
* This means that the output is solely determined by the last sample
* entering the filter, initially. */
mTaps[ 0 ] = new Complex( 1.0f, 0.0f );
for( int x = 1; x < TAP_COUNT; x++ )
{
mTaps[ x ] = new Complex( 0.0f, 0.0f );
mBuffer[ x ] = new Complex( 0.0f, 0.0f );
}
}
private void updateTaps( Complex sample )
{
/* Calculate the error. The ideal sample magnitude is unity (1.0).
* When the sample is not at unity, the error vector points in the
* opposite direction with a magnitude equal to 1 - sample magnitude.*/
mError = Complex.multiply( sample, sample.norm() - mModulus );
/* Ensure we don't exceed unity */
mError.clip( 1.0f );
for( int x = 0; x < TAP_COUNT; x++ )
{
// tap -= mu * conj(sample) * error
mTaps[ x ] = Complex.subtract( mTaps[ x ],
Complex.multiply( Complex.multiply(
mBuffer[ x ].conjugate(), mError ), mMu ) );
}
}
@Override
public void receive( Complex sample )
{
mBuffer[ mBufferPointer ] = sample;
mBufferPointer++;
mBufferPointer = mBufferPointer % TAP_COUNT;
Complex equalized = filter();
updateTaps( equalized );
if( mListener != null )
{
mListener.receive( equalized );
}
mLog.debug( "Tap 0:" + mTaps[ 0 ].toString() + " error:" + mError.toString() );
}
private Complex filter()
{
Complex sum = new Complex( 0.0f, 0.0f );
for( int x = 0; x < TAP_COUNT; x++ )
{
sum = Complex.add( sum, Complex.multiply(
mBuffer[ ( x + mBufferPointer ) % TAP_COUNT ],
mTaps[ TAP_COUNT - x - 1 ] ) );
}
return sum;
}
public void setListener( Listener<Complex> listener )
{
mListener = listener;
}
}

View File

@@ -0,0 +1,6 @@
package io.github.dsheirer.dsp.filter.fir;
public abstract class FIRFilter
{
public abstract void dispose();
}

File diff suppressed because it is too large Load Diff

View File

@@ -0,0 +1,54 @@
package io.github.dsheirer.dsp.filter.fir.complex;
import io.github.dsheirer.dsp.filter.fir.FIRFilter;
import io.github.dsheirer.dsp.filter.fir.real.RealFIRFilter;
import io.github.dsheirer.sample.complex.Complex;
public class ComplexFIRFilter extends FIRFilter
{
private RealFIRFilter mIFilter;
private RealFIRFilter mQFilter;
/**
* Complex FIR Filter for processing complex sample pairs. Wraps two real
* FIR filters for processing each of the inphase and quadrature samples.
*
* @param coefficients - filter taps
* @param gain - gain to apply to filtered outputs - use 1.0f for no gain
*/
public ComplexFIRFilter( float[] coefficients, float gain )
{
mIFilter = new RealFIRFilter( coefficients, gain );
mQFilter = new RealFIRFilter( coefficients, gain );
}
public float[] getCoefficients()
{
return mIFilter.getCoefficients();
}
public float filterInphase( float sample )
{
return mIFilter.filter( sample );
}
public float filterQuadrature( float sample )
{
return mQFilter.filter( sample );
}
public Complex filter( Complex sample )
{
float i = filterInphase( sample.inphase() );
float q = filterQuadrature( sample.quadrature() );
return new Complex( i, q );
}
@Override
public void dispose()
{
mIFilter.dispose();
mQFilter.dispose();
}
}

View File

@@ -0,0 +1,54 @@
package io.github.dsheirer.dsp.filter.fir.complex;
import io.github.dsheirer.sample.Listener;
import io.github.dsheirer.sample.complex.ComplexBuffer;
public class ComplexFIRFilter_CB_CB extends ComplexFIRFilter
implements Listener<ComplexBuffer>
{
private Listener<ComplexBuffer> mListener;
/**
* Complex FIR Filter for processing complex buffers. Wraps two real FIR
* filters for processing each of the inphase and quadrature samples.
*
* @param coefficients - odd length symmetric filter taps
* @param gain - gain value to apply to each of the filtered samples
*/
public ComplexFIRFilter_CB_CB( float[] coefficients, float gain )
{
super( coefficients, gain );
}
public void dispose()
{
mListener = null;
}
@Override
public void receive( ComplexBuffer buffer )
{
if( mListener != null )
{
float[] samples = buffer.getSamples();
for( int x = 0; x < samples.length; x += 2 )
{
samples[ x ] = filterInphase( samples[ x ] );
samples[ x + 1 ] = filterQuadrature( samples[ x + 1 ] );
}
mListener.receive( buffer );
}
}
public void setListener( Listener<ComplexBuffer> listener )
{
mListener = listener;
}
public void removeListener()
{
mListener = null;
}
}

View File

@@ -0,0 +1,105 @@
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2015 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp.filter.fir.real;
import io.github.dsheirer.dsp.filter.fir.FIRFilter;
public class RealFIRFilter extends FIRFilter
{
private float[] mBuffer;
private int mBufferSize = 1;
private int mBufferPointer = 0;
private int[][] mIndexMap;
private float[] mCoefficients;
private float mGain;
/**
* Float sample FIR filter base class.
*
* @param coefficients - filter coefficients
* @param gain value to apply to the filtered output - use 1.0f for no gain
*/
public RealFIRFilter( float[] coefficients, float gain )
{
mCoefficients = coefficients;
mGain = gain;
mBufferSize = mCoefficients.length;
mBuffer = new float[ mBufferSize ];
mBufferPointer = mBufferSize - 1;
generateIndexMap( mBufferSize );
}
public float[] getCoefficients()
{
return mCoefficients;
}
@Override
public void dispose()
{
mCoefficients = null;
mIndexMap = null;
mBuffer = null;
}
public float filter( float sample )
{
mBuffer[ mBufferPointer ] = sample;
float accumulator = 0.0f;
for( int x = 0; x < mBufferSize; x++ )
{
accumulator += mCoefficients[ x ] *
mBuffer[ mIndexMap[ mBufferPointer][ x ] ];
}
mBufferPointer--;
if( mBufferPointer < 0 )
{
mBufferPointer += mBufferSize;
}
/* Apply gain and return the filtered value */
return accumulator * mGain;
}
/**
* Generates a circular buffer index map to support lookup of the translated
* index based on the current buffer pointer and the desired sample index.
*/
private void generateIndexMap( int size )
{
mIndexMap = new int[ size ][ size ];
for( int x = 0; x < size; x++ )
{
for( int y = 0; y < size; y++ )
{
int z = x + y;
mIndexMap[ x ][ y ] = z < size ? z : z - size;
}
}
}
}

View File

@@ -0,0 +1,67 @@
package io.github.dsheirer.dsp.filter.fir.real;
import io.github.dsheirer.sample.Listener;
import io.github.dsheirer.sample.real.RealBuffer;
public class RealFIRFilter_RB_RB extends RealFIRFilter
implements Listener<RealBuffer>
{
private Listener<RealBuffer> mListener;
/**
* Real Buffer Float FIR filter with streaming float buffer provider and
* listener interfaces.
*
* @param coefficients
* @param gain
*/
public RealFIRFilter_RB_RB( float[] coefficients, float gain )
{
super( coefficients, gain );
}
@Override
public void dispose()
{
super.dispose();
mListener = null;
}
@Override
public void receive( RealBuffer buffer )
{
if( mListener != null )
{
mListener.receive( filter( buffer ) );
}
}
/**
* Filters the buffer samples
*
* @param buffer
* @return
*/
public RealBuffer filter( RealBuffer buffer )
{
float[] samples = buffer.getSamples();
for( int x = 0; x < samples.length; x++ )
{
samples[ x ] = filter( samples[ x ] );
}
return buffer;
}
public void setListener( Listener<RealBuffer> listener )
{
mListener = listener;
}
public void removeListener()
{
mListener = null;
}
}

View File

@@ -0,0 +1,56 @@
package io.github.dsheirer.dsp.filter.fir.real;
import io.github.dsheirer.sample.real.RealSampleListener;
public class RealFIRFilter_R_R extends RealFIRFilter implements RealSampleListener
{
private RealSampleListener mListener;
/**
* Float FIR filter with streaming float provider and listener interfaces.
*
* @param coefficients
* @param gain
*/
public RealFIRFilter_R_R( float[] coefficients, float gain )
{
super( coefficients, gain );
}
@Override
public void dispose()
{
super.dispose();
mListener = null;
}
@Override
public void receive( float sample )
{
if( mListener != null )
{
mListener.receive( filter( sample ) );
}
}
public float[] filter( float[] samples )
{
for( int x = 0; x < samples.length; x++ )
{
samples[ x ] = filter( samples[ x ] );
}
return samples;
}
public void setListener( RealSampleListener listener )
{
mListener = listener;
}
public void removeListener()
{
mListener = null;
}
}

View File

@@ -0,0 +1,56 @@
package io.github.dsheirer.dsp.filter.fir.remez;
public enum FIRLinearPhaseFilterType
{
/**
* TYPE 1
* Length: ODD
* Order: EVEN
* Symmetry: SYMMETRICAL
*
* Valid for: LOW PASS, HIGH PASS
*/
TYPE_1_ODD_LENGTH_EVEN_ORDER_SYMMETRICAL,
/**
* TYPE 2
* Length: EVEN
* Order: ODD
* Symmetry: SYMMETRICAL
* Valid for: LOW PASS
*
* Notes:
* -Frequency response at FS/2 (Pi) must be 0 -- not suitable for High Pass.
*/
TYPE_2_EVEN_LENGTH_ODD_ORDER_SYMMETRICAL,
/**
* TYPE 3
* Length: ODD
* Order: EVEN
* Symmetry: ANTI-SYMMETRICAL
* Valid for:
*
* Notes:
* -Introduces a 90 degree phase shift such that h[N-n] = -h[n] for all n
* -Frequency response at 0 and FS/2 (Pi) must be 0 -- not suitable for Low or High Pass.
*/
TYPE_3_ODD_LENGTH_EVEN_ORDER_ANTI_SYMMETRICAL,
/**
* TYPE 4
* Length: EVEN
* Order: ODD
* Symmetry: ANTI-SYMMETRICAL
* Valid for: HIGH PASS
*
* Notes:
* -Introduces a 90 degree phase shift such that h[N-n] = -h[n] for all n
* -Frequency response at 0 must be zero -- not suitable for Low Pass
*/
TYPE_4_EVEN_LENGTH_ODD_ORDER_ANTI_SYMMETRICAL,
BANDPASS,
DIFFERENTIATOR,
HILBERT;
}

View File

@@ -0,0 +1,111 @@
package io.github.dsheirer.dsp.filter.fir.remez;
import io.github.dsheirer.dsp.filter.fir.FIRFilterSpecification;
public class Grid
{
public static final double TWO_PI = 2.0 * Math.PI;
private double[] mFrequencyGrid;
private double[] mCosineOfFrequencyGrid;
private double[] mDesiredResponse;
private double[] mWeight;
/**
* Creates a dense frequency grid from the filter specification for those frequency bands
* where the weighted value is non-zero. This grid is appropriate for use in the Remez
* exchange algorithm.
*/
public Grid( FIRFilterSpecification specification )
{
create( specification );
}
private void create( FIRFilterSpecification specification )
{
int gridSize = specification.getGridSize();
mFrequencyGrid = new double[ gridSize ];
mCosineOfFrequencyGrid = new double[ gridSize ];
mDesiredResponse = new double[ gridSize ];
mWeight = new double[ gridSize ];
double gridFrequencyInterval = specification.getGridFrequencyInterval();
double bandZeroEdge = specification.getFrequencyBands().get( 0 ).getStart();
// For differentiator, hilbert, symmetry is odd and Grid[0] = max(delf, bands[0])
double grid0 = specification.isAntiSymmetrical() &&
( gridFrequencyInterval > bandZeroEdge ) ? gridFrequencyInterval : bandZeroEdge;
int j = 0;
double maxRipple = specification.getMaxBandAmplitude();
for( int x = 0; x < specification.getFrequencyBands().size(); x++ )
{
FIRFilterSpecification.FrequencyBand band = specification.getFrequencyBands().get( x );
double lowFrequency = ( x == 0 ? grid0 : band.getStart() );
for( int i = 0; i < band.getGridSize(); i++ )
{
mDesiredResponse[ j ] = band.getAmplitude();
//Set the weighting function equally across this band's dense grids
mWeight[ j ] = band.getWeight( maxRipple );
//Frequency grid is linear taper from start frequency to stop frequency in
//gridFrequencyInterval increments
mFrequencyGrid[ j ] = lowFrequency;
mCosineOfFrequencyGrid[ j ] = Math.cos( TWO_PI * mFrequencyGrid[ j ] );
lowFrequency += gridFrequencyInterval;
j++;
}
int lastIndex = j - 1;
mFrequencyGrid[ lastIndex ] = band.getEnd();
}
}
/**
* Size of the dense frequency grid
*/
public int getSize()
{
return mFrequencyGrid.length;
}
/**
* The dense frequency grid
*/
public double[] getFrequencyGrid()
{
return mFrequencyGrid;
}
public double[] getCosineFrequencyGrid()
{
return mCosineOfFrequencyGrid;
}
/**
* Desired frequency response
*/
public double[] getDesiredResponse()
{
return mDesiredResponse;
}
/**
* Weighting function
*/
public double[] getWeight()
{
return mWeight;
}
}

View File

@@ -0,0 +1,674 @@
package io.github.dsheirer.dsp.filter.fir.remez;
import io.github.dsheirer.dsp.filter.design.FilterDesignException;
import io.github.dsheirer.dsp.filter.fir.FIRFilterSpecification;
import org.slf4j.Logger;
import org.slf4j.LoggerFactory;
import java.text.DecimalFormat;
import java.util.ArrayList;
import java.util.Iterator;
import java.util.List;
public class RemezFIRFilterDesigner
{
private final static Logger mLog = LoggerFactory.getLogger( RemezFIRFilterDesigner.class );
private static final double CONVERGENCE_THRESHOLD = 0.0001;
public static final int MAXIMUM_ITERATION_COUNT = 40;
public static final double TWO_PI = 2.0 * Math.PI;
private static final DecimalFormat DECIMAL_FORMATTER = new DecimalFormat( "0.00000000" );
private FIRFilterSpecification mSpecification;
private Grid mGrid;
private List<Integer> mExtremalIndices;
private double[] mD;
private double[] mGridErrors;
private double[] mGridFrequencyResponse;
private double[] mIdealFrequencyResponse;
private double mDelta;
private boolean mConverged;
/**
* Linear FIR filter designer using the Parks-McClelland Remez exchange algorithm.
*
* Construct a filter specification using the builder methods in the FIRFilterSpecification
* class and construct an instance of this class with the specification. Check the isValid()
* method after construction to see if the filter design was successful. If successful, use the
* getFrequencyResponse() and getImpulseResponse() methods to access the filter taps as the
* specified length.
*
* If the filter cannot be designed from the specification the isValid() method will return
* false and any attempts to access the getFrequencyResponse() or getImpulseResponse() methods
* will throw a FilterDesignException.
*
* @param specification that defines a Type 1-4 linear phase FIR Filter
*/
public RemezFIRFilterDesigner( FIRFilterSpecification specification )
{
mSpecification = specification;
design();
}
/**
* Designs the filter according to the specification
*/
private void design()
{
mGrid = new Grid( mSpecification );
mExtremalIndices = getInitialExtremalIndices();
int iterationCount = 0;
try
{
do
{
calculateGridFrequencyResponse();
calculateGridError();
findExtremalIndices();
checkConvergence();
iterationCount++;
}
while( !mConverged && iterationCount < MAXIMUM_ITERATION_COUNT );
}
catch( FilterDesignException fde )
{
mLog.error( "Filter design error - couldn't find extremal indices at count [" +
iterationCount + "] try changing filter order up/down from [" +
mSpecification.getOrder() + "]" );
mConverged = false;
}
if( mConverged )
{
//Update frequency response using the final set of extremal indices
calculateGridFrequencyResponse();
}
}
/**
* Indicates if the filter was successfully created from the filter specification. Check this
* method before accessing frequency response or filter response, otherwise a FilterDesignException
* is thrown.
*
* @return true if the filter was successfully designed.
*/
public boolean isValid()
{
return mConverged;
}
/**
* Frequency response of the designed filter with double precision
*
* @return half ( length / 2 ) of filter frequency response sampled at 0.0 - PI radians
* @throws FilterDesignException if the specified filter cannot be designed
*/
public double[] getFrequencyResponse() throws FilterDesignException
{
if( !mConverged )
{
throw new FilterDesignException( "Can't create filter from specification - failed "
+ "to converge" );
}
//Resample the final polynomial at the desired filter length
double[] resampled = resample();
//Apply Parks/McClelland frequency response corrections according to the filter type
return correctFrequencyResponse( resampled );
}
/**
* Impulse response of the designed filter with floating point precision.
*
* @return filter impulse response
* @throws FilterDesignException if the specified filter cannot be designed
*/
public float[] getImpulseResponse() throws FilterDesignException
{
return convertToFloatArray( getImpulseResponseDoubles() );
}
/**
* Impulse response of the designed filter with double precision
*/
public double[] getImpulseResponseDoubles() throws FilterDesignException
{
double[] frequencyResponse = getFrequencyResponse();
int length = mSpecification.getFilterLength();
double[] impulseResponse = new double[ length ];
double M;
switch( mSpecification.getFilterType() )
{
case TYPE_1_ODD_LENGTH_EVEN_ORDER_SYMMETRICAL:
M = ( (double)length - 1.0 ) / 2.0;
for( int n = 0; n < length; n++ )
{
double accumulator = frequencyResponse[ 0 ];
double frequency = TWO_PI * ( n - M ) / length;
for( int k = 1; k <= M; k++ )
{
accumulator += 2.0 * frequencyResponse[ k ] * Math.cos( frequency * (double)k );
}
impulseResponse[ n ] = accumulator / (double)length;
}
break;
case TYPE_2_EVEN_LENGTH_ODD_ORDER_SYMMETRICAL:
double offset = (double)( length - 1 ) / 2.0;
for( int n = 0; n < length; n++ )
{
double accumulator = frequencyResponse[ 0 ];
double frequency = TWO_PI * ( (double)n - offset ) / (double)length;
for( int k = 1; k < frequencyResponse.length; k++ )
{
accumulator += 2.0 * frequencyResponse[ k ] * Math.cos( frequency * (double)k );
}
impulseResponse[ n ] = accumulator / (double)length;
}
break;
}
return impulseResponse;
}
/**
* Calculates the frequency response of this filter at the specified cosine of the frequency
*
* Implements Oppenheim/Schafer Discrete Time Signal Processing, 3e, 2016, equation 116a to
* determine the actual response of the specified frequency cosine using Lagrange
* interpolation of the polynomial formed at the extremal index set.
*
* @param cosineOfFrequency to calculate the frequency response
*
* @return actual frequency response for the specified frequency with the current set of
* extreaml indices
*/
private double getFrequencyResponse( double cosineOfFrequency )
{
double numerator = 0.0;
double denominator = 0.0;
for( int k = 0; k < mExtremalIndices.size() - 1; k++ )
{
double cosineDelta = cosineOfFrequency -
mGrid.getCosineFrequencyGrid()[ mExtremalIndices.get( k ) ];
//If this frequency is or is close to one of the polynomial points, use c[k] for that point
if( Math.abs( cosineDelta ) < 1.0e-7 )
{
numerator = mIdealFrequencyResponse[ k ];
denominator = 1.0;
break;
}
else
{
double dkOverCosineDelta = mD[ k ] / cosineDelta;
numerator += dkOverCosineDelta * mIdealFrequencyResponse[ k ];
denominator += dkOverCosineDelta;
}
}
return numerator / denominator;
}
/**
* Calculated Delta value
*
* @return delta value
*/
public double getDelta()
{
return mDelta;
}
/**
* Creates initial array of extremal frequency indices spaced at intervals of the grid frequency
* set for all frequencies where a frequency band weighting is specified.
*/
private List<Integer> getInitialExtremalIndices()
{
int count = mSpecification.getExtremaCount();
int density = mSpecification.getGridDensity();
List<Integer> extremalIndices = new ArrayList<>();
for( int i = 0; i < count; i++ )
{
extremalIndices.add( i * density );
}
return extremalIndices;
}
/**
* Calculates the actual frequency response of a dense frequency grid using the polynomial
* represented by the extremal indices.
*/
private void calculateGridFrequencyResponse()
{
double[] b = calculateB();
calculateDelta( b );
calculateC();
calculateD( b );
updateGridFrequencyResponse();
}
/**
* Updates the frequency response of frequencies specified in the dense frequency grid
*/
private void updateGridFrequencyResponse()
{
mGridFrequencyResponse = new double[ mGrid.getSize() ];
double[] gridFrequencyCosines = mGrid.getCosineFrequencyGrid();
for( int i = 0; i < mGridFrequencyResponse.length; i++ )
{
mGridFrequencyResponse[ i ] = getFrequencyResponse( gridFrequencyCosines[ i ] );
}
}
/**
* Calculates the value of b across the set of extremal indices.
*
* Implements Oppenheim/Schafer Discrete Time Signal Processing, 3e, 2016, equation 115
*
* @return the array of b values of length L+2
*/
private double[] calculateB()
{
int length = mExtremalIndices.size();
double[] b = new double[ length ];
for( int k = 0; k < length; k++ )
{
b[ k ] = 1.0;
double xk = mGrid.getCosineFrequencyGrid()[ mExtremalIndices.get( k ) ];
for( int i = 0; i < length; i++ )
{
if( i != k )
{
double xi = mGrid.getCosineFrequencyGrid()[ mExtremalIndices.get( i ) ];
double denominator = xk - xi;
if( Math.abs( denominator ) < 0.00001 )
{
denominator = 0.00001;
}
b[ k ] *= 1.0 / denominator;
}
}
}
return b;
}
/**
* Calculates the value of delta that represents the maximum ripple for the current set of
* extremal indices (L+2)
*
* Implements Oppenheim/Schafer Discrete Time Signal Processing, 3e, 2016, equation 114
*
* @param b calculated for the grid from equation 115
*/
private void calculateDelta( double[] b )
{
double numerator = 0.0;
double denominator = 0.0;
double sign = 1.0;
for( int k = 0; k < b.length; k++ )
{
if( k < mExtremalIndices.size() )
{
int extremalIndex = mExtremalIndices.get( k );
numerator += ( b[ k ] * mGrid.getDesiredResponse()[ extremalIndex ] );
denominator += b[ k ] * sign / mGrid.getWeight()[ extremalIndex ];
sign = -sign;
}
else
{
mLog.error( "Something went wrong -- the length of b exceeds the set of extremal indices" );
}
}
mDelta = numerator / denominator;
}
/**
* Calculates the ideal frequency response (C) of the current extremal index set. This reference
* is used to determine the level of error in the set of actual frequency responses calculated
* from this set of extremal indices to determine convergence.
*
* Implements Oppenheim/Schafer Discrete Time Signal Processing, 3e, 2016, equation 116b
*/
private void calculateC()
{
int length = mSpecification.getExtremaCount() - 1;
mIdealFrequencyResponse = new double[ length ];
double sign = 1.0;
for( int k = 0; k < length; k++ )
{
if( k < mExtremalIndices.size() )
{
int index = mExtremalIndices.get( k );
mIdealFrequencyResponse[ k ] = mGrid.getDesiredResponse()[ index ] -
( sign * mDelta / mGrid.getWeight()[ index ] );
sign = -sign;
}
}
}
/**
* Calculates the set of D values for the current extremal index set.
*
* Implements Oppenheim/Schafer Discrete Time Signal Processing, 3e, 2016, equation 116c
*/
private void calculateD( double[] b )
{
int length = mExtremalIndices.size() - 1;
mD = new double[ length ];
for( int k = 0; k < length; k++ )
{
//Note: extremalCosines array is one index longer than d array, so we use length as
//a pointer to the final (L+2) extremalCosine index
mD[k] = b[k] * ( mGrid.getCosineFrequencyGrid()[ mExtremalIndices.get( k ) ] -
mGrid.getCosineFrequencyGrid()[ mExtremalIndices.get( length ) ] );
}
}
/**
* Calculates the weighted error between the specification desired frequency response and the
* actual frequency response of the current extremal index set across the frequencies in the
* dense frequency grid.
*
* Implements Oppenheim/Schafer Discrete Time Signal Processing, 3e, 2016, equation 112
*/
public void calculateGridError()
{
int length = mGridFrequencyResponse.length;
mGridErrors = new double[ length ];
for( int i = 0; i < length; i++ )
{
mGridErrors[ i ] = mGrid.getWeight()[ i ] *
( mGrid.getDesiredResponse()[ i ] - mGridFrequencyResponse[ i ] );
}
}
/**
* Identifies the indexes of the extremal error values using the Alternation theorem
*
* @throws FilterDesignException if the full set of extremal indices can't be identified
*/
private void findExtremalIndices() throws FilterDesignException
{
mExtremalIndices.clear();
//Check for extremum at error index 0
if( ( ( mGridErrors[ 0 ] > 0.0 && mGridErrors[ 0 ] > mGridErrors[ 1 ] ) ||
( mGridErrors[ 0 ] < 0.0 && mGridErrors[ 0 ] < mGridErrors[ 1 ] ) ) &&
isGTEDelta( mGridErrors[ 0 ] ) )
{
mExtremalIndices.add( 0 );
}
//Check for extrema in middle error indices )
for( int x = 1; x < mGridErrors.length - 1; x++ )
{
if( ( ( mGridErrors[ x ] > 0.0 && ( mGridErrors[ x - 1 ] <= mGridErrors[ x ] && mGridErrors[ x ] > mGridErrors[ x + 1 ] ) ) ||
( mGridErrors[ x ] < 0.0 && ( mGridErrors[ x - 1 ] >= mGridErrors[ x ] && mGridErrors[ x ] < mGridErrors[ x + 1 ] ) ) ) &&
isGTEDelta( mGridErrors[ x ] ) )
{
mExtremalIndices.add( x );
}
}
//Check for extremum at final error index
int last = mGridErrors.length - 1;
if( ( ( mGridErrors[ last ] > 0.0 && ( mGridErrors[ last ] > mGridErrors[ last - 1 ] ) ) ||
( mGridErrors[ last ] < 0.0 && ( mGridErrors[ last ] < mGridErrors[ last - 1 ] ) ) ) &&
isGTEDelta( mGridErrors[ last ] ) )
{
mExtremalIndices.add( last );
}
//Ensure we have the minimum of extremals before we continue
if( mExtremalIndices.size() < mSpecification.getExtremaCount() )
{
throw new FilterDesignException( "Couldn't find the minimum extremal frequencies in "
+ "error set" );
}
//Enforce alternation theory -- only one extremal for each transition about the zero axis
List<Integer> indicesToRemove = new ArrayList<>();
Iterator<Integer> it = mExtremalIndices.iterator();
Integer current = it.next();
Integer next;
boolean positiveAxis = mGridErrors[ current ] > 0.0;
while( it.hasNext() )
{
next = it.next();
//Check for consecutive (redundant) indices on same side of zero axis - retain largest
if( !( positiveAxis ^ ( mGridErrors[ next ] > 0.0 ) ) )
{
if( Math.abs( mGridErrors[ next ] ) <= Math.abs( mGridErrors[ current ] ))
{
//Remove next error index that is less than current index and on
it.remove();
next = current;
}
else
{
//Since we can't remove the current index with the iterator, queue for later removal
indicesToRemove.add( current );
}
}
else
{
positiveAxis = !positiveAxis;
}
current = next;
}
//Remove redundant extremal indices that couldn't be removed via the iterator
mExtremalIndices.removeAll( indicesToRemove );
//Truncate the list to one larger than needed by removing excess tailing indices
while( mExtremalIndices.size() > mSpecification.getExtremaCount() )
{
mExtremalIndices.remove( mExtremalIndices.size() - 1 );
}
//If we have one too many indices, delete the smaller of the first or last
if( mExtremalIndices.size() > mSpecification.getExtremaCount() )
{
int lastIndex = mExtremalIndices.size() - 1;
if( Math.abs( mGridErrors[ mExtremalIndices.get( 0 )] ) >
Math.abs( mGridErrors[ mExtremalIndices.get( lastIndex )] ) )
{
mExtremalIndices.remove( lastIndex );
}
else
{
mExtremalIndices.remove( 0 );
}
}
//Detect if we have too few indices
if( mExtremalIndices.size() < mSpecification.getExtremaCount() )
{
throw new FilterDesignException( "Couldn't find the minimum extremal frequencies in "
+ "error set" );
}
}
/**
* Indicates if the absolute value of the argument is greater than or equal to the delta
* value with accuracy to 14 digits of precision. This avoids rounding errors at 15 digits or
* precision or greater.
*
* @param value to evaluate
*
* @return true if the absolute value is greater than or equal to the absolute delta value
*/
private boolean isGTEDelta( double value )
{
return Math.abs( value ) - Math.abs( mDelta ) > -1.0e-5;
}
/**
* Checks for convergence of the frequency response of the current set of extremal indices to
* the filter specification by comparing the maximum absolute error value against the delta
* value to determine if these two values are within the convergence threshold.
*/
public void checkConvergence()
{
double maximum = Math.abs( mGridErrors[ mExtremalIndices.get( 0 ) ] );
for( int i = 1; i < mExtremalIndices.size(); i++ )
{
double current = Math.abs( mGridErrors[ mExtremalIndices.get( i ) ] );
if( current > maximum )
{
maximum = current;
}
}
double convergence = maximum - Math.abs( mDelta );
mConverged = convergence < CONVERGENCE_THRESHOLD;
}
/**
* Coverts/casts the double array to a float array
*/
private static float[] convertToFloatArray( double[] samples )
{
float[] converted = new float[ samples.length ];
for( int x = 0; x < samples.length; x++ )
{
converted[ x ] = (float)samples[ x ];
}
return converted;
}
/**
* Resamples the designed filter at half of the filter length at evenly spaced intervals of the
* frequency spectrum from 0 to PI radians.
*
* @return frequency response of the filter sampled the desired filter length and order
*/
private double[] resample()
{
int length = mSpecification.getFilterLength();
if( length % 2 == 0 )
{
length--;
}
double half = (double)length / 2.0;
double[] resampled = new double[ (int)Math.ceil( half ) ] ;
for( int x = 0; x < resampled.length; x++ )
{
resampled[ x ] = getFrequencyResponse( Math.cos( Math.PI * (double)x / half ) );
}
return resampled;
}
/**
* Applies Parks/McClelland filter type correction to the final frequency response array
*
* @param frequencyResponse of the resampled filter design from 0 <> ~ PI radians
*
* @return corrected frequency response
*/
private double[] correctFrequencyResponse( double[] frequencyResponse )
{
double filterLength = mSpecification.getFilterLength();
switch( mSpecification.getFilterType() )
{
case TYPE_1_ODD_LENGTH_EVEN_ORDER_SYMMETRICAL:
//No correction needed
break;
case TYPE_2_EVEN_LENGTH_ODD_ORDER_SYMMETRICAL:
//No correction needed
break;
case TYPE_3_ODD_LENGTH_EVEN_ORDER_ANTI_SYMMETRICAL:
for( int x = 0; x < frequencyResponse.length; x++ )
{
double frequencyRadians = Math.PI * ( ( (double)x ) / filterLength );
frequencyResponse[ x ] *= Math.sin( TWO_PI * frequencyRadians );
}
break;
case TYPE_4_EVEN_LENGTH_ODD_ORDER_ANTI_SYMMETRICAL:
for( int x = 0; x < frequencyResponse.length; x++ )
{
double frequencyRadians = Math.PI * ( ( (double)x ) / filterLength );
frequencyResponse[ x ] *= Math.sin( Math.PI * frequencyRadians );
}
break;
default:
break;
}
return frequencyResponse;
}
}

View File

@@ -0,0 +1,117 @@
package io.github.dsheirer.dsp.filter.halfband.complex;
import io.github.dsheirer.sample.Listener;
import io.github.dsheirer.sample.complex.ComplexBuffer;
import io.github.dsheirer.dsp.filter.halfband.real.HalfBandFilter;
public class HalfBandFilter_CB_CB implements Listener<ComplexBuffer>
{
private HalfBandFilter mIFilter;
private HalfBandFilter mQFilter;
private boolean mDecimate;
private boolean mDecimateFlag;
private Listener<ComplexBuffer> mListener;
/**
*
* @param coefficients
* @param gain
* @param decimate
*/
public HalfBandFilter_CB_CB( float[] coefficients, float gain, boolean decimate )
{
mIFilter = new HalfBandFilter( coefficients, gain );
mQFilter = new HalfBandFilter( coefficients, gain );
mDecimate = decimate;
}
public void dispose()
{
mIFilter.dispose();
mQFilter.dispose();
mListener = null;
}
@Override
public void receive( ComplexBuffer buffer )
{
if( mListener != null )
{
if( mDecimate )
{
float[] samples = buffer.getSamples();
int half = samples.length / 2;
float[] decimated;
if( half % 2 == 0 || mDecimateFlag )
{
decimated = new float[ half ];
}
else
{
/* If inbound buffer size is odd-length, then we have to
* adjust when the first operation is non-decimation, since
* that will produce an outbound buffer 1 sample larger */
decimated = new float[ half + 1 ];
}
int decimatedPointer = 0;
for( int x = 0; x < samples.length; x += 2 )
{
/* Insert the sample but don't filter */
if( mDecimateFlag )
{
mIFilter.insert( samples[ x ] );
mQFilter.insert( samples[ x + 1 ] );
}
else
{
decimated[ decimatedPointer++ ] = mIFilter.filter( samples[ x ] );
decimated[ decimatedPointer++ ] = mQFilter.filter( samples[ x + 1 ] );
}
/* Toggle the decimation flag for every sample */
mDecimateFlag = !mDecimateFlag;
}
if( mListener != null )
{
mListener.receive( new ComplexBuffer( decimated ) );
}
}
else
{
/* Non Decimate - filter each of the values and return the buffer */
float[] samples = buffer.getSamples();
for( int x = 0; x < samples.length; x += 2 )
{
samples[ x ] = mIFilter.filter( samples[ x ] );
samples[ x + 1 ] = mQFilter.filter( samples[ x + 1 ] );
}
if( mListener != null )
{
mListener.receive( buffer );
}
}
}
}
public void setListener( Listener<ComplexBuffer> listener )
{
mListener = listener;
}
public void removeListener()
{
mListener = null;
}
}

View File

@@ -0,0 +1,162 @@
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014,2015 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp.filter.halfband.real;
import org.apache.commons.lang3.Validate;
public class HalfBandFilter
{
private float[] mBuffer;
private int mBufferSize = 1; //Temporary initial value
private int mBufferPointer = 0;
private float mGain;
private float[] mCoefficients;
private int[][] mIndexMap;
private int mCenterCoefficient;
private int mCenterCoefficientMapIndex;
/**
* Half-Band filter for processing real float samples.
*
* Takes advantage of the 0-valued FIR half-band coefficents inherent
* in the half-band filter and does not calculate those coefficients.
*
* This reduces the workload to (tap-size - 1) / 4 + 1 calculations per sample.
*
* @param filter - filter coefficients
* @param gain - gain multiplier. Use 1.0 for unity/no gain
*/
public HalfBandFilter( float[] coefficients, float gain )
{
mCoefficients = coefficients;
mGain = gain;
mBufferSize = mCoefficients.length;
mBuffer = new float[ mBufferSize ];
generateIndexMap( mCoefficients.length );
}
public void dispose()
{
}
/**
* Inserts the sample into the buffer without calculating a filtered value
*/
public void insert( float sample )
{
mBuffer[ mBufferPointer ] = sample;
mBufferPointer++;
mBufferPointer = mBufferPointer % mBufferSize;
}
/**
* Calculate the filtered value by applying the coefficients against
* the complex samples in mBuffer
*/
public float filter( float sample )
{
insert( sample );
float accumulator = 0;
//Start with the center tap value
accumulator += mCoefficients[ mCenterCoefficient ] *
mBuffer[ mIndexMap[ mBufferPointer ][ mCenterCoefficientMapIndex ] ];
//For the remaining coefficients, add the symmetric samples, oldest and newest
//first, then multiply by the single coefficient
for( int x = 0; x < mCenterCoefficientMapIndex; x += 2 )
{
accumulator += mCoefficients[ x ] *
( mBuffer[ mIndexMap[ mBufferPointer ][ x ] ] +
mBuffer[ mIndexMap[ mBufferPointer ][ x + 1 ] ] );
}
return accumulator * mGain;
}
/**
* Creates an n X (n + 1 / 2) index map enabling quick access to the
* circular buffer samples.
*
* As the buffer shifts right with each subsequent sample, we have to move
* the index pointers with it, for efficient access of the samples.
*
* The first array index value in the index map corresponds to the current
* buffer pointer location.
*
* The second array index value points to the samples that should be
* multiplied by the coefficients as follows:
*
* 0 = center tap sample, to be multiplied by center coefficient
*
* 0 = sample( 1 )
* 1 = sample( size - 1 )
*
* Indexes 0 and 1 will be multiplied by coefficient( 0 ).
*
* Subsequent indexes 3, 4, etc, point to the oldest and newest samples that
* correspond to the matching ( 3 ) coefficient index.
*
* @param odd-sized number of filter taps (ie coefficients) and buffer
*/
private void generateIndexMap( int size )
{
//Ensure we have an odd size
Validate.isTrue(size % 2 == 1);
int mapWidth = ( ( size + 1 ) / 2 ) + 1;
//Set center tap index for coefficients array
mCenterCoefficient = ( size - 1 ) / 2;
mCenterCoefficientMapIndex = mCenterCoefficient + 1;
mIndexMap = new int[ size ][ mapWidth ];
//Setup the first row, buffer pointer index 0, as a starting point
for( int x = 0; x < mapWidth - 2; x += 2 )
{
mIndexMap[ 0 ][ x ] = x;
mIndexMap[ 0 ][ x + 1 ] = size - 1 - x;
}
//Place center tap index in last element
mIndexMap[ 0 ][ mCenterCoefficientMapIndex ] = mCenterCoefficient;
//For each subsequent row, increment the previous row's value by 1,
//subtracting size as needed, to keep the values between 0 and size - 1
for( int x = 1; x < size; x++ )
{
for( int y = 0; y < mapWidth; y++ )
{
mIndexMap[ x ][ y ] = mIndexMap[ x - 1 ][ y ] + 1;
if( mIndexMap[ x ][ y ] >= size )
{
mIndexMap[ x ][ y ] -= size;
}
}
}
}
}

View File

@@ -0,0 +1,117 @@
package io.github.dsheirer.dsp.filter.halfband.real;
import io.github.dsheirer.sample.Listener;
import io.github.dsheirer.sample.real.RealBuffer;
public class HalfBandFilter_RB_RB extends HalfBandFilter implements Listener<RealBuffer>
{
private boolean mDecimate;
private boolean mDecimateFlag;
private Listener<RealBuffer> mListener;
/**
* Halfband filter for real buffers with optional decimation by two.
*
* If this filter is used in non-decimating mode, the filter will overwrite
* the inbound sample buffer with the filtered samples. This will cause
* issues if the buffers are being used across several processes are sharing
* the buffer, since this filter will modify the contents of the object, so
* consider making a copy of the buffer prior to feeding the buffer to this
* filter.
*
* If used in decimating mode, a new (half-length) buffer will be created
* for the filtered samples.
*
* @param coefficients - filter kernel
* @param gain - gain to apply to the outputs
* @param decimate - true for decimate by 2, or false for no decimation
*/
public HalfBandFilter_RB_RB( float[] coefficients, float gain, boolean decimate )
{
super( coefficients, gain );
mDecimate = decimate;
}
public void dispose()
{
super.dispose();
mListener = null;
}
@Override
public void receive( RealBuffer buffer )
{
if( mListener != null )
{
mListener.receive( filter( buffer ) );
}
}
public RealBuffer filter( RealBuffer buffer )
{
if( mDecimate )
{
float[] samples = buffer.getSamples();
int half = samples.length / 2;
float[] decimated;
if( half % 2 == 0 || mDecimateFlag )
{
decimated = new float[ half ];
}
else
{
/* If inbound buffer size is odd-length, then we have to
* adjust when the first operation is non-decimation, since
* that will produce an outbound buffer 1 sample larger */
decimated = new float[ half + 1 ];
}
int decimatedPointer = 0;
for( float sample: samples )
{
/* Insert the sample but don't filter */
if( mDecimateFlag )
{
insert( sample );
}
else
{
decimated[ decimatedPointer++ ] = filter( sample );
}
/* Toggle the decimation flag for every sample */
mDecimateFlag = !mDecimateFlag;
}
return new RealBuffer( decimated );
}
else
{
/* Non Decimate - filter each of the values and return the buffer */
float[] samples = buffer.getSamples();
for( int x = 0; x < samples.length; x++ )
{
samples[ x ] = filter( samples[ x ] );
}
return buffer;
}
}
public void setListener( Listener<RealBuffer> listener )
{
mListener = listener;
}
public void removeListener()
{
mListener = null;
}
}

View File

@@ -0,0 +1,281 @@
package io.github.dsheirer.dsp.filter.hilbert;
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2015 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 <http://www.gnu.org/licenses/>
******************************************************************************/
import io.github.dsheirer.dsp.filter.Filters;
import io.github.dsheirer.sample.Listener;
import io.github.dsheirer.sample.complex.ComplexBuffer;
import io.github.dsheirer.sample.real.RealBuffer;
import org.slf4j.Logger;
import org.slf4j.LoggerFactory;
import java.util.Arrays;
public class HilbertTransform implements Listener<RealBuffer>
{
private static final Logger mLog = LoggerFactory.getLogger( HilbertTransform.class );
private float[] mHilbertFilter;
private Listener<ComplexBuffer> mListener;
private float[] mBuffer;
private int mBufferSize;
private int mBufferPointer;
private int[][] mIndexMap;
private int mMapHeight;
private int mCenterTapIndex;
private float mCenterCoefficent;
private boolean mInvertFlag = false;
/**
* Hilbert transform filter used for converting real-valued samples into
* complex valued samples using frequency translation (FS/4) and a half-band
* filter.
*
* This filter uses a circular sample buffer and a pre-calculated
* index map to correctly map the filter coefficients to the circular delay
* buffer contents as each new sample is added to the buffer and processed.
*
* Half-band filter coefficients used in this filter must be of length N
* where (N + 1) is a multiple of 4. This filter uses a pre-defined half
* band filter with 47 coefficients.
*
* This transform process is described in Understanding Digital Signal
* Processing, Lyons, 3e, 2011, sections 13.1.2 and 13.1.3 (p 674-678) and
* implemented as described in Section 13.37.1 and 13.37.2 (p 802-804)
*/
public HilbertTransform()
{
convertHalfBandToHilbert( Filters.HALF_BAND_FILTER_47T.getCoefficients() );
mBufferSize = mHilbertFilter.length + 1;
mBuffer = new float[ mBufferSize ];
mCenterCoefficent = mHilbertFilter[ mHilbertFilter.length / 2 ];
generateIndexMap( mHilbertFilter.length );
}
@Override
public void receive( RealBuffer buffer )
{
if( mListener != null )
{
mListener.receive( filter( buffer ) );
}
}
/**
* Sets the listener to receive filtered complex buffers.
*/
public void setListener( Listener<ComplexBuffer> listener )
{
mListener = listener;
}
/**
* Applies Hilbert transform to the real sample buffer and returns a
* complex sample buffer
*/
public ComplexBuffer filter( RealBuffer buffer )
{
return new ComplexBuffer( filter( buffer.getSamples() ) );
}
/**
* Filters the continuous real-sample stream by treating the stream as a
* set of interleaved values and applying two hilbert filters against each
* interleaved set (I & Q). Frequency translation has already been
* pre-applied to the hilbert filter coefficients.
*
* Filtering of the inphase component is a simple multiply by the center
* coefficient (0.5) since all of the other coefficients are zero-valued.
*
* Filtering of the quadrature component uses a folded FIR structure and
* expects the hilbert coefficients to be symmetric with all negative valued
* coefficients below the center coefficient and all positive-valued coefficients
* above the center coefficient. So, convolution only applies the filter
* coefficient against the difference of the two samples.
*
* ( coefficient * value1 ) + (-coefficient * value2 ) = coefficient * ( value1 - value2 )
*
* Performs frequency translation by FS/2 on the final filtered values by
* applying a sequence of 1,-1 (sign change) to each I/Q output sample.
*/
public float[] filter( float[] samples )
{
float accumulator;
for( int y = 0; y < samples.length; y +=2 )
{
insert( samples[ y ] );
insert( samples[ y + 1 ] );
accumulator = 0.0f;
int index = mBufferPointer / 2;
for( int x = 0; x < mHilbertFilter.length / 2; x += 2 )
{
accumulator += mHilbertFilter[ x ] *
( mBuffer[ mIndexMap[ index ][ x + 1] ] -
mBuffer[ mIndexMap[ index ][ x ] ] );
}
//Perform FS/2 frequency translation on the final filtered values
if( mInvertFlag )
{
//inphase
samples[ y ] = -( mBuffer[ mIndexMap[ index ][ mCenterTapIndex ] ] );
//quadrature
samples[ y + 1 ] = -accumulator;
}
else
{
//inphase
samples[ y ] = mBuffer[ mIndexMap[ index ][ mCenterTapIndex ] ];
//quadrature
samples[ y + 1 ] = accumulator;
}
mInvertFlag = !mInvertFlag;
}
return samples;
}
/**
* Inserts the sample into the circular buffer, overwriting the oldest value
*/
private void insert( float sample )
{
mBuffer[ mBufferPointer++ ] = sample;
mBufferPointer = mBufferPointer % mBufferSize;
}
/**
* Creates an index map to support efficient lookup of sample indexes from
* the circular buffer.
*
* The first index of the map corresponds to the current buffer pointer
* setting.
*
* The second index is structured in sample pairs to align with
* the hilbert filter coefficients. For example, the first iteration of
* the filter would use filter coefficient 0 and the index map values 0 and
* 1 would be pointing to samples 46 and 0.
*
* The final element of each index set is the index of the center sample.
*/
private void generateIndexMap( int size )
{
mMapHeight = size / 2 + 1;
int mapWidth = mMapHeight + 1;
mIndexMap = new int[ mMapHeight ][ mapWidth ];
//Setup the first row with buffer pointer index 0 as a starting point
for( int x = 0; x < mapWidth - 1; x += 2 )
{
mIndexMap[ 0 ][ x ] = size - 1 - x;
mIndexMap[ 0 ][ x + 1 ] = x;
}
//Place the center index at the end of the array
mCenterTapIndex = mapWidth - 1;
mIndexMap[ 0 ][ mCenterTapIndex ] = size / 2;
//For each subsequent row, increment the previous row's value by 2,
//wrapping as needed, to keep the values between 0 and size - 1
for( int x = 1; x < mMapHeight; x++ )
{
for( int y = 0; y < mapWidth; y++ )
{
mIndexMap[ x ][ y ] = mIndexMap[ x - 1 ][ y ] + 2;
if( mIndexMap[ x ][ y ] >= size )
{
mIndexMap[ x ][ y ] -= size + 1;
//Special handling for center tap wrap around
if( y == mCenterTapIndex && mIndexMap[ x ][ y ] < 0 )
{
mIndexMap[ x ][ y ] = size;
}
}
}
}
}
/**
* Logs the index map contents
*/
public void logIndexMap()
{
for( int[] indexSet: mIndexMap )
{
mLog.debug( "Row:" + Arrays.toString( indexSet ) );
}
}
/**
* Converts the half-band filter coefficients for use as hilbert transform
* filter coefficients. Sets all even-numbered coefficients left of center
* to negative and all even-numbered coefficients right of center to positive
* which applies a FS/4 frequency translation to the coefficients.
*
* Note: even though the coefficients above the center coefficient are being
* defined, they are unused during the convolution process since we're using
* a folded FIR structure to exploit the symmetric nature of the half-band
* filter coefficients. The polarity difference between the upper and lower
* halves is accounted for in the filter() method.
*
* Note: we apply a 2.0 gain to the coefficients to compensate for the loss
* of splitting the signal via two filters.
*/
private void convertHalfBandToHilbert( float[] coefficients )
{
mHilbertFilter = new float[ coefficients.length ];
int middle = coefficients.length / 2;
for( int x = 0; x < coefficients.length; x++ )
{
if( x < middle )
{
mHilbertFilter[ x ] = 2.0f * -Math.abs( coefficients[ x ] );
}
else if( x > middle )
{
mHilbertFilter[ x ] = 2.0f * Math.abs( coefficients[ x ] );
}
else
{
mHilbertFilter[ x ] = 2.0f * coefficients[ x ];
}
}
}
}

View File

@@ -0,0 +1,280 @@
package io.github.dsheirer.dsp.filter.hilbert;
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2015 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 <http://www.gnu.org/licenses/>
*
* ----------------------------------------------------------------------------
* Half-band filter coefficients retrieved October 2015 from:
* https://github.com/airspy/host/libairspy/src/filters.h
*
* Copyright (C) 2014, Youssef Touil <youssef@airspy.com>
*
* Permission is hereby granted, free of charge, to any person obtaining a copy
* of this software and associated documentation files (the "Software"), to deal
* in the Software without restriction, including without limitation the rights
* to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
* copies of the Software, and to permit persons to whom the Software is
* furnished to do so, subject to the following conditions:
*
* The above copyright notice and this permission notice shall be included in
* all copies or substantial portions of the Software.
*
* THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
* IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
* FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
* AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
* LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
* OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
* THE SOFTWARE.
******************************************************************************/
import io.github.dsheirer.sample.Listener;
import io.github.dsheirer.sample.complex.ComplexBuffer;
import io.github.dsheirer.sample.real.RealBuffer;
import org.slf4j.Logger;
import org.slf4j.LoggerFactory;
import java.util.Arrays;
public class InterpolatingHilbertTransform implements Listener<RealBuffer>
{
private static final Logger mLog = LoggerFactory.getLogger( InterpolatingHilbertTransform.class );
private static final float[] HALF_BAND_FILTER =
{
-0.000998606272947510f, 0.0f, 0.001695637278417295f, 0.0f,
-0.003054430179754289f, 0.0f, 0.005055504379767936f, 0.0f,
-0.007901319195893647f, 0.0f, 0.011873357051047719f, 0.0f,
-0.017411159379930066f, 0.0f, 0.025304817427568772f, 0.0f,
-0.037225225204559217f, 0.0f, 0.057533286997004301f, 0.0f,
-0.102327462004259350f, 0.0f, 0.317034472508947400f, 0.5f,
0.317034472508947400f, 0.0f, -0.102327462004259350f, 0.0f,
0.057533286997004301f, 0.0f, -0.037225225204559217f, 0.0f,
0.025304817427568772f, 0.0f, -0.017411159379930066f, 0.0f,
0.011873357051047719f, 0.0f, -0.007901319195893647f, 0.0f,
0.005055504379767936f, 0.0f, -0.003054430179754289f, 0.0f,
0.001695637278417295f, 0.0f, -0.000998606272947510f
};
private static final float GAIN = 2.0f;
private Listener<ComplexBuffer> mListener;
private float[] mBuffer;
private int mBufferSize;
private int mBufferPointer;
private int mGroupDelayedInphaseTapPointIndex;
private float[] mCoefficients;
private int[][] mIndexMap;
/**
* Hilbert transform filter used for converting real-valued samples into
* complex valued samples using frequency translation (FS/4) and a half-band
* filter.
*
* This filter uses a circular sample buffer and a pre-calculated
* index map to correctly map the filter coefficients to the circular delay
* buffer contents as each new sample is added to the buffer and processed.
*
* Half-band filter coefficients used in this filter must be of length N
* where (N + 1) is a multiple of 4. This filter uses a pre-defined half
* band filter with 47 coefficients.
*
* This transform process is described in Understanding Digital Signal
* Processing, Lyons, 3e, 2011, sections 13.1.2 and 13.1.3 (p 674-678) and
* implemented as described in Section 13.37.1 and 13.37.2 (p 802-804)
*/
public InterpolatingHilbertTransform()
{
mBufferSize = HALF_BAND_FILTER.length;
mBuffer = new float[ mBufferSize ];
createHilbertCoefficients();
generateIndexMap( mBufferSize );
}
@Override
public void receive( RealBuffer buffer )
{
if( mListener != null )
{
mListener.receive( filter( buffer ) );
}
}
/**
* Sets the listener to receive filtered complex buffers.
*/
public void setListener( Listener<ComplexBuffer> listener )
{
mListener = listener;
}
/**
* Applies Hilbert transform to the real sample buffer and returns a
* complex sample buffer
*/
public ComplexBuffer filter( RealBuffer buffer )
{
return new ComplexBuffer( filter( buffer.getSamples() ) );
}
/**
* Applies a hilbert transform to an array of real samples and returns an
* array twice as long containing I,Q,I,Q ... complex samples
*/
public float[] filter( float[] samples )
{
float[] filtered = new float[ samples.length * 2 ];
float accumulator;
for( int y = 0; y < samples.length; y++ )
{
insert( samples[ y ] );
accumulator = 0.0f;
for( int x = 0; x < mCoefficients.length; x++ )
{
accumulator += mCoefficients[ x ] *
mBuffer[ mIndexMap[ mBufferPointer ][ x ] ] ;
}
int index = y * 2;
/* We use the buffered sample at the half group delay point as the
* Inphase value and the accumulator as the quadrature value */
filtered[ index ] = mBuffer[ mIndexMap[ mBufferPointer ]
[ mGroupDelayedInphaseTapPointIndex ] ];
filtered[ index + 1 ] = accumulator;
}
return filtered;
}
/**
* Inserts the sample into the circular buffer, overwriting the oldest value
*/
private void insert( float sample )
{
mBuffer[ mBufferPointer ] = sample;
mBufferPointer++;
mBufferPointer = mBufferPointer % mBufferSize;
}
/**
* Creates an N x (N + 1 / 2) index map enabling quick access to the
* circular buffer samples to support convolution.
*
* As the buffer pointer shifts right, as each sample is added, we use the
* index map pointed to by the buffer pointer, so that the filter coefficients
* are pointing to the correct circular buffer samples.
*
* We add an extra index to each index map with a value that points to the
* current index of the middle sample in the buffer, so that we can retrieve
* that value as the half-group delayed sample that forms the inphase
* component.
*/
private void generateIndexMap( int size )
{
int mapWidth = ( ( size + 1 ) / 2 ) + 1;
int center = ( size - 1 ) / 2;
//Set center tap index for coefficients array
mGroupDelayedInphaseTapPointIndex = center + 1;
mIndexMap = new int[ size ][ mapWidth ];
//Setup the first row with buffer pointer index 0 as a starting point
for( int x = 0; x < mapWidth - 1; x++ )
{
mIndexMap[ 0 ][ x ] = size - 1 - ( x * 2 );
}
//Place the group delayed in-phase sample tap point in last element of
//the first row
mIndexMap[ 0 ][ mGroupDelayedInphaseTapPointIndex ] = center;
//For each subsequent row, increment the previous row's value by 1,
//subtracting size as needed, to keep the values between 0 and size - 1
for( int x = 1; x < size; x++ )
{
for( int y = 0; y < mapWidth; y++ )
{
mIndexMap[ x ][ y ] = mIndexMap[ x - 1 ][ y ] + 1;
if( mIndexMap[ x ][ y ] >= size )
{
mIndexMap[ x ][ y ] -= size;
}
}
}
}
public void logIndexMap()
{
for( int[] indexSet: mIndexMap )
{
mLog.debug( "Row:" + Arrays.toString( indexSet ) );
}
}
/**
* Converts the half-band filter coefficients into a packed set of hilbert
* coefficients by removing all of the odd index, zero-valued coefficients
* and the center 0.5 coefficient. Then, all coefficients left of the
* original center are made negative and all coefficients right of the
* original center are made positive. These coefficient sign modifications
* effectively incorporate a frequency translation by FS/4 into the hilbert
* coefficients. A gain of 2.0 is applied to the coefficients to offset the
* gain of 0.5 induced by the filter.
*/
private void createHilbertCoefficients()
{
int half = ( HALF_BAND_FILTER.length + 1 ) / 2;
int middle = half / 2;
mCoefficients = new float[ half ];
/* Pack all of the even numbered HB coefficients into our coefficient
* array and pre-apply a gain of two to offset the loss in the filter */
for( int x = 0; x < half; x++ )
{
mCoefficients[ x ] = HALF_BAND_FILTER[ x * 2 ] * GAIN;
/* Make all of the left of center coefficients negative */
if( x < middle && mCoefficients[ x ] > 0.0f )
{
mCoefficients[ x ] = -mCoefficients[ x ];
}
/* Make all of the right of center coefficients positive */
else if( x >= middle && mCoefficients[ x ] < 0.0f )
{
mCoefficients[ x ] = -mCoefficients[ x ];
}
}
}
}

View File

@@ -0,0 +1,60 @@
package io.github.dsheirer.dsp.filter.iir;
import io.github.dsheirer.sample.real.RealBuffer;
public class DeemphasisFilter
{
private static final float MAX_SAMPLE_VALUE = 0.95f;
private float mAlpha;
private float mGain;
private float mPrevious = 0.0f;
public DeemphasisFilter( float sampleRate, float cutoff, float gain )
{
mAlpha = (float)Math.exp( -2.0 * Math.PI * cutoff * ( 1.0 / sampleRate ) );
mGain = gain;
}
public float filter( float sample )
{
mPrevious = sample + ( mAlpha * mPrevious );
return declip( mPrevious * mGain );
}
private float declip( float value )
{
if( value > MAX_SAMPLE_VALUE )
{
return MAX_SAMPLE_VALUE;
}
else if( value < -MAX_SAMPLE_VALUE )
{
return -MAX_SAMPLE_VALUE;
}
else return value;
}
public float[] filter( float[] samples )
{
for( int x = 0; x < samples.length; x++ )
{
samples[ x ] = filter( samples[ x ] );
}
return samples;
}
public RealBuffer filter( RealBuffer buffer )
{
float[] samples = buffer.getSamples();
for( int x = 0; x < samples.length; x++ )
{
samples[ x ] = filter( samples[ x ] );
}
return buffer;
}
}

View File

@@ -0,0 +1,161 @@
package io.github.dsheirer.dsp.filter.interpolator;
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2015 Dennis Sheirer
*
* Copyright (C) 2002,2012 Free Software Foundation, Inc.
* Ported to java from gnuradio - interpolator_taps.h
*
* 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 <http://www.gnu.org/licenses/>
******************************************************************************/
public abstract class Interpolator
{
public static final int NTAPS = 8;
public static final int NSTEPS = 128;
public static final float TAPS[][] = new float[][]
{
// -4 -3 -2 -1 0 1 2 3 mu
{ 0.00000e+00f, 0.00000e+00f, 0.00000e+00f, 0.00000e+00f, 1.00000e+00f, 0.00000e+00f, 0.00000e+00f, 0.00000e+00f }, // 0/128
{ -1.54700e-04f, 8.53777e-04f, -2.76968e-03f, 7.89295e-03f, 9.98534e-01f, -5.41054e-03f, 1.24642e-03f, -1.98993e-04f }, // 1/128
{ -3.09412e-04f, 1.70888e-03f, -5.55134e-03f, 1.58840e-02f, 9.96891e-01f, -1.07209e-02f, 2.47942e-03f, -3.96391e-04f }, // 2/128
{ -4.64053e-04f, 2.56486e-03f, -8.34364e-03f, 2.39714e-02f, 9.95074e-01f, -1.59305e-02f, 3.69852e-03f, -5.92100e-04f }, // 3/128
{ -6.18544e-04f, 3.42130e-03f, -1.11453e-02f, 3.21531e-02f, 9.93082e-01f, -2.10389e-02f, 4.90322e-03f, -7.86031e-04f }, // 4/128
{ -7.72802e-04f, 4.27773e-03f, -1.39548e-02f, 4.04274e-02f, 9.90917e-01f, -2.60456e-02f, 6.09305e-03f, -9.78093e-04f }, // 5/128
{ -9.26747e-04f, 5.13372e-03f, -1.67710e-02f, 4.87921e-02f, 9.88580e-01f, -3.09503e-02f, 7.26755e-03f, -1.16820e-03f }, // 6/128
{ -1.08030e-03f, 5.98883e-03f, -1.95925e-02f, 5.72454e-02f, 9.86071e-01f, -3.57525e-02f, 8.42626e-03f, -1.35627e-03f }, // 7/128
{ -1.23337e-03f, 6.84261e-03f, -2.24178e-02f, 6.57852e-02f, 9.83392e-01f, -4.04519e-02f, 9.56876e-03f, -1.54221e-03f }, // 8/128
{ -1.38589e-03f, 7.69462e-03f, -2.52457e-02f, 7.44095e-02f, 9.80543e-01f, -4.50483e-02f, 1.06946e-02f, -1.72594e-03f }, // 9/128
{ -1.53777e-03f, 8.54441e-03f, -2.80746e-02f, 8.31162e-02f, 9.77526e-01f, -4.95412e-02f, 1.18034e-02f, -1.90738e-03f }, // 10/128
{ -1.68894e-03f, 9.39154e-03f, -3.09033e-02f, 9.19033e-02f, 9.74342e-01f, -5.39305e-02f, 1.28947e-02f, -2.08645e-03f }, // 11/128
{ -1.83931e-03f, 1.02356e-02f, -3.37303e-02f, 1.00769e-01f, 9.70992e-01f, -5.82159e-02f, 1.39681e-02f, -2.26307e-03f }, // 12/128
{ -1.98880e-03f, 1.10760e-02f, -3.65541e-02f, 1.09710e-01f, 9.67477e-01f, -6.23972e-02f, 1.50233e-02f, -2.43718e-03f }, // 13/128
{ -2.13733e-03f, 1.19125e-02f, -3.93735e-02f, 1.18725e-01f, 9.63798e-01f, -6.64743e-02f, 1.60599e-02f, -2.60868e-03f }, // 14/128
{ -2.28483e-03f, 1.27445e-02f, -4.21869e-02f, 1.27812e-01f, 9.59958e-01f, -7.04471e-02f, 1.70776e-02f, -2.77751e-03f }, // 15/128
{ -2.43121e-03f, 1.35716e-02f, -4.49929e-02f, 1.36968e-01f, 9.55956e-01f, -7.43154e-02f, 1.80759e-02f, -2.94361e-03f }, // 16/128
{ -2.57640e-03f, 1.43934e-02f, -4.77900e-02f, 1.46192e-01f, 9.51795e-01f, -7.80792e-02f, 1.90545e-02f, -3.10689e-03f }, // 17/128
{ -2.72032e-03f, 1.52095e-02f, -5.05770e-02f, 1.55480e-01f, 9.47477e-01f, -8.17385e-02f, 2.00132e-02f, -3.26730e-03f }, // 18/128
{ -2.86289e-03f, 1.60193e-02f, -5.33522e-02f, 1.64831e-01f, 9.43001e-01f, -8.52933e-02f, 2.09516e-02f, -3.42477e-03f }, // 19/128
{ -3.00403e-03f, 1.68225e-02f, -5.61142e-02f, 1.74242e-01f, 9.38371e-01f, -8.87435e-02f, 2.18695e-02f, -3.57923e-03f }, // 20/128
{ -3.14367e-03f, 1.76185e-02f, -5.88617e-02f, 1.83711e-01f, 9.33586e-01f, -9.20893e-02f, 2.27664e-02f, -3.73062e-03f }, // 21/128
{ -3.28174e-03f, 1.84071e-02f, -6.15931e-02f, 1.93236e-01f, 9.28650e-01f, -9.53307e-02f, 2.36423e-02f, -3.87888e-03f }, // 22/128
{ -3.41815e-03f, 1.91877e-02f, -6.43069e-02f, 2.02814e-01f, 9.23564e-01f, -9.84679e-02f, 2.44967e-02f, -4.02397e-03f }, // 23/128
{ -3.55283e-03f, 1.99599e-02f, -6.70018e-02f, 2.12443e-01f, 9.18329e-01f, -1.01501e-01f, 2.53295e-02f, -4.16581e-03f }, // 24/128
{ -3.68570e-03f, 2.07233e-02f, -6.96762e-02f, 2.22120e-01f, 9.12947e-01f, -1.04430e-01f, 2.61404e-02f, -4.30435e-03f }, // 25/128
{ -3.81671e-03f, 2.14774e-02f, -7.23286e-02f, 2.31843e-01f, 9.07420e-01f, -1.07256e-01f, 2.69293e-02f, -4.43955e-03f }, // 26/128
{ -3.94576e-03f, 2.22218e-02f, -7.49577e-02f, 2.41609e-01f, 9.01749e-01f, -1.09978e-01f, 2.76957e-02f, -4.57135e-03f }, // 27/128
{ -4.07279e-03f, 2.29562e-02f, -7.75620e-02f, 2.51417e-01f, 8.95936e-01f, -1.12597e-01f, 2.84397e-02f, -4.69970e-03f }, // 28/128
{ -4.19774e-03f, 2.36801e-02f, -8.01399e-02f, 2.61263e-01f, 8.89984e-01f, -1.15113e-01f, 2.91609e-02f, -4.82456e-03f }, // 29/128
{ -4.32052e-03f, 2.43930e-02f, -8.26900e-02f, 2.71144e-01f, 8.83893e-01f, -1.17526e-01f, 2.98593e-02f, -4.94589e-03f }, // 30/128
{ -4.44107e-03f, 2.50946e-02f, -8.52109e-02f, 2.81060e-01f, 8.77666e-01f, -1.19837e-01f, 3.05345e-02f, -5.06363e-03f }, // 31/128
{ -4.55932e-03f, 2.57844e-02f, -8.77011e-02f, 2.91006e-01f, 8.71305e-01f, -1.22047e-01f, 3.11866e-02f, -5.17776e-03f }, // 32/128
{ -4.67520e-03f, 2.64621e-02f, -9.01591e-02f, 3.00980e-01f, 8.64812e-01f, -1.24154e-01f, 3.18153e-02f, -5.28823e-03f }, // 33/128
{ -4.78866e-03f, 2.71272e-02f, -9.25834e-02f, 3.10980e-01f, 8.58189e-01f, -1.26161e-01f, 3.24205e-02f, -5.39500e-03f }, // 34/128
{ -4.89961e-03f, 2.77794e-02f, -9.49727e-02f, 3.21004e-01f, 8.51437e-01f, -1.28068e-01f, 3.30021e-02f, -5.49804e-03f }, // 35/128
{ -5.00800e-03f, 2.84182e-02f, -9.73254e-02f, 3.31048e-01f, 8.44559e-01f, -1.29874e-01f, 3.35600e-02f, -5.59731e-03f }, // 36/128
{ -5.11376e-03f, 2.90433e-02f, -9.96402e-02f, 3.41109e-01f, 8.37557e-01f, -1.31581e-01f, 3.40940e-02f, -5.69280e-03f }, // 37/128
{ -5.21683e-03f, 2.96543e-02f, -1.01915e-01f, 3.51186e-01f, 8.30432e-01f, -1.33189e-01f, 3.46042e-02f, -5.78446e-03f }, // 38/128
{ -5.31716e-03f, 3.02507e-02f, -1.04150e-01f, 3.61276e-01f, 8.23188e-01f, -1.34699e-01f, 3.50903e-02f, -5.87227e-03f }, // 39/128
{ -5.41467e-03f, 3.08323e-02f, -1.06342e-01f, 3.71376e-01f, 8.15826e-01f, -1.36111e-01f, 3.55525e-02f, -5.95620e-03f }, // 40/128
{ -5.50931e-03f, 3.13987e-02f, -1.08490e-01f, 3.81484e-01f, 8.08348e-01f, -1.37426e-01f, 3.59905e-02f, -6.03624e-03f }, // 41/128
{ -5.60103e-03f, 3.19495e-02f, -1.10593e-01f, 3.91596e-01f, 8.00757e-01f, -1.38644e-01f, 3.64044e-02f, -6.11236e-03f }, // 42/128
{ -5.68976e-03f, 3.24843e-02f, -1.12650e-01f, 4.01710e-01f, 7.93055e-01f, -1.39767e-01f, 3.67941e-02f, -6.18454e-03f }, // 43/128
{ -5.77544e-03f, 3.30027e-02f, -1.14659e-01f, 4.11823e-01f, 7.85244e-01f, -1.40794e-01f, 3.71596e-02f, -6.25277e-03f }, // 44/128
{ -5.85804e-03f, 3.35046e-02f, -1.16618e-01f, 4.21934e-01f, 7.77327e-01f, -1.41727e-01f, 3.75010e-02f, -6.31703e-03f }, // 45/128
{ -5.93749e-03f, 3.39894e-02f, -1.18526e-01f, 4.32038e-01f, 7.69305e-01f, -1.42566e-01f, 3.78182e-02f, -6.37730e-03f }, // 46/128
{ -6.01374e-03f, 3.44568e-02f, -1.20382e-01f, 4.42134e-01f, 7.61181e-01f, -1.43313e-01f, 3.81111e-02f, -6.43358e-03f }, // 47/128
{ -6.08674e-03f, 3.49066e-02f, -1.22185e-01f, 4.52218e-01f, 7.52958e-01f, -1.43968e-01f, 3.83800e-02f, -6.48585e-03f }, // 48/128
{ -6.15644e-03f, 3.53384e-02f, -1.23933e-01f, 4.62289e-01f, 7.44637e-01f, -1.44531e-01f, 3.86247e-02f, -6.53412e-03f }, // 49/128
{ -6.22280e-03f, 3.57519e-02f, -1.25624e-01f, 4.72342e-01f, 7.36222e-01f, -1.45004e-01f, 3.88454e-02f, -6.57836e-03f }, // 50/128
{ -6.28577e-03f, 3.61468e-02f, -1.27258e-01f, 4.82377e-01f, 7.27714e-01f, -1.45387e-01f, 3.90420e-02f, -6.61859e-03f }, // 51/128
{ -6.34530e-03f, 3.65227e-02f, -1.28832e-01f, 4.92389e-01f, 7.19116e-01f, -1.45682e-01f, 3.92147e-02f, -6.65479e-03f }, // 52/128
{ -6.40135e-03f, 3.68795e-02f, -1.30347e-01f, 5.02377e-01f, 7.10431e-01f, -1.45889e-01f, 3.93636e-02f, -6.68698e-03f }, // 53/128
{ -6.45388e-03f, 3.72167e-02f, -1.31800e-01f, 5.12337e-01f, 7.01661e-01f, -1.46009e-01f, 3.94886e-02f, -6.71514e-03f }, // 54/128
{ -6.50285e-03f, 3.75341e-02f, -1.33190e-01f, 5.22267e-01f, 6.92808e-01f, -1.46043e-01f, 3.95900e-02f, -6.73929e-03f }, // 55/128
{ -6.54823e-03f, 3.78315e-02f, -1.34515e-01f, 5.32164e-01f, 6.83875e-01f, -1.45993e-01f, 3.96678e-02f, -6.75943e-03f }, // 56/128
{ -6.58996e-03f, 3.81085e-02f, -1.35775e-01f, 5.42025e-01f, 6.74865e-01f, -1.45859e-01f, 3.97222e-02f, -6.77557e-03f }, // 57/128
{ -6.62802e-03f, 3.83650e-02f, -1.36969e-01f, 5.51849e-01f, 6.65779e-01f, -1.45641e-01f, 3.97532e-02f, -6.78771e-03f }, // 58/128
{ -6.66238e-03f, 3.86006e-02f, -1.38094e-01f, 5.61631e-01f, 6.56621e-01f, -1.45343e-01f, 3.97610e-02f, -6.79588e-03f }, // 59/128
{ -6.69300e-03f, 3.88151e-02f, -1.39150e-01f, 5.71370e-01f, 6.47394e-01f, -1.44963e-01f, 3.97458e-02f, -6.80007e-03f }, // 60/128
{ -6.71985e-03f, 3.90083e-02f, -1.40136e-01f, 5.81063e-01f, 6.38099e-01f, -1.44503e-01f, 3.97077e-02f, -6.80032e-03f }, // 61/128
{ -6.74291e-03f, 3.91800e-02f, -1.41050e-01f, 5.90706e-01f, 6.28739e-01f, -1.43965e-01f, 3.96469e-02f, -6.79662e-03f }, // 62/128
{ -6.76214e-03f, 3.93299e-02f, -1.41891e-01f, 6.00298e-01f, 6.19318e-01f, -1.43350e-01f, 3.95635e-02f, -6.78902e-03f }, // 63/128
{ -6.77751e-03f, 3.94578e-02f, -1.42658e-01f, 6.09836e-01f, 6.09836e-01f, -1.42658e-01f, 3.94578e-02f, -6.77751e-03f }, // 64/128
{ -6.78902e-03f, 3.95635e-02f, -1.43350e-01f, 6.19318e-01f, 6.00298e-01f, -1.41891e-01f, 3.93299e-02f, -6.76214e-03f }, // 65/128
{ -6.79662e-03f, 3.96469e-02f, -1.43965e-01f, 6.28739e-01f, 5.90706e-01f, -1.41050e-01f, 3.91800e-02f, -6.74291e-03f }, // 66/128
{ -6.80032e-03f, 3.97077e-02f, -1.44503e-01f, 6.38099e-01f, 5.81063e-01f, -1.40136e-01f, 3.90083e-02f, -6.71985e-03f }, // 67/128
{ -6.80007e-03f, 3.97458e-02f, -1.44963e-01f, 6.47394e-01f, 5.71370e-01f, -1.39150e-01f, 3.88151e-02f, -6.69300e-03f }, // 68/128
{ -6.79588e-03f, 3.97610e-02f, -1.45343e-01f, 6.56621e-01f, 5.61631e-01f, -1.38094e-01f, 3.86006e-02f, -6.66238e-03f }, // 69/128
{ -6.78771e-03f, 3.97532e-02f, -1.45641e-01f, 6.65779e-01f, 5.51849e-01f, -1.36969e-01f, 3.83650e-02f, -6.62802e-03f }, // 70/128
{ -6.77557e-03f, 3.97222e-02f, -1.45859e-01f, 6.74865e-01f, 5.42025e-01f, -1.35775e-01f, 3.81085e-02f, -6.58996e-03f }, // 71/128
{ -6.75943e-03f, 3.96678e-02f, -1.45993e-01f, 6.83875e-01f, 5.32164e-01f, -1.34515e-01f, 3.78315e-02f, -6.54823e-03f }, // 72/128
{ -6.73929e-03f, 3.95900e-02f, -1.46043e-01f, 6.92808e-01f, 5.22267e-01f, -1.33190e-01f, 3.75341e-02f, -6.50285e-03f }, // 73/128
{ -6.71514e-03f, 3.94886e-02f, -1.46009e-01f, 7.01661e-01f, 5.12337e-01f, -1.31800e-01f, 3.72167e-02f, -6.45388e-03f }, // 74/128
{ -6.68698e-03f, 3.93636e-02f, -1.45889e-01f, 7.10431e-01f, 5.02377e-01f, -1.30347e-01f, 3.68795e-02f, -6.40135e-03f }, // 75/128
{ -6.65479e-03f, 3.92147e-02f, -1.45682e-01f, 7.19116e-01f, 4.92389e-01f, -1.28832e-01f, 3.65227e-02f, -6.34530e-03f }, // 76/128
{ -6.61859e-03f, 3.90420e-02f, -1.45387e-01f, 7.27714e-01f, 4.82377e-01f, -1.27258e-01f, 3.61468e-02f, -6.28577e-03f }, // 77/128
{ -6.57836e-03f, 3.88454e-02f, -1.45004e-01f, 7.36222e-01f, 4.72342e-01f, -1.25624e-01f, 3.57519e-02f, -6.22280e-03f }, // 78/128
{ -6.53412e-03f, 3.86247e-02f, -1.44531e-01f, 7.44637e-01f, 4.62289e-01f, -1.23933e-01f, 3.53384e-02f, -6.15644e-03f }, // 79/128
{ -6.48585e-03f, 3.83800e-02f, -1.43968e-01f, 7.52958e-01f, 4.52218e-01f, -1.22185e-01f, 3.49066e-02f, -6.08674e-03f }, // 80/128
{ -6.43358e-03f, 3.81111e-02f, -1.43313e-01f, 7.61181e-01f, 4.42134e-01f, -1.20382e-01f, 3.44568e-02f, -6.01374e-03f }, // 81/128
{ -6.37730e-03f, 3.78182e-02f, -1.42566e-01f, 7.69305e-01f, 4.32038e-01f, -1.18526e-01f, 3.39894e-02f, -5.93749e-03f }, // 82/128
{ -6.31703e-03f, 3.75010e-02f, -1.41727e-01f, 7.77327e-01f, 4.21934e-01f, -1.16618e-01f, 3.35046e-02f, -5.85804e-03f }, // 83/128
{ -6.25277e-03f, 3.71596e-02f, -1.40794e-01f, 7.85244e-01f, 4.11823e-01f, -1.14659e-01f, 3.30027e-02f, -5.77544e-03f }, // 84/128
{ -6.18454e-03f, 3.67941e-02f, -1.39767e-01f, 7.93055e-01f, 4.01710e-01f, -1.12650e-01f, 3.24843e-02f, -5.68976e-03f }, // 85/128
{ -6.11236e-03f, 3.64044e-02f, -1.38644e-01f, 8.00757e-01f, 3.91596e-01f, -1.10593e-01f, 3.19495e-02f, -5.60103e-03f }, // 86/128
{ -6.03624e-03f, 3.59905e-02f, -1.37426e-01f, 8.08348e-01f, 3.81484e-01f, -1.08490e-01f, 3.13987e-02f, -5.50931e-03f }, // 87/128
{ -5.95620e-03f, 3.55525e-02f, -1.36111e-01f, 8.15826e-01f, 3.71376e-01f, -1.06342e-01f, 3.08323e-02f, -5.41467e-03f }, // 88/128
{ -5.87227e-03f, 3.50903e-02f, -1.34699e-01f, 8.23188e-01f, 3.61276e-01f, -1.04150e-01f, 3.02507e-02f, -5.31716e-03f }, // 89/128
{ -5.78446e-03f, 3.46042e-02f, -1.33189e-01f, 8.30432e-01f, 3.51186e-01f, -1.01915e-01f, 2.96543e-02f, -5.21683e-03f }, // 90/128
{ -5.69280e-03f, 3.40940e-02f, -1.31581e-01f, 8.37557e-01f, 3.41109e-01f, -9.96402e-02f, 2.90433e-02f, -5.11376e-03f }, // 91/128
{ -5.59731e-03f, 3.35600e-02f, -1.29874e-01f, 8.44559e-01f, 3.31048e-01f, -9.73254e-02f, 2.84182e-02f, -5.00800e-03f }, // 92/128
{ -5.49804e-03f, 3.30021e-02f, -1.28068e-01f, 8.51437e-01f, 3.21004e-01f, -9.49727e-02f, 2.77794e-02f, -4.89961e-03f }, // 93/128
{ -5.39500e-03f, 3.24205e-02f, -1.26161e-01f, 8.58189e-01f, 3.10980e-01f, -9.25834e-02f, 2.71272e-02f, -4.78866e-03f }, // 94/128
{ -5.28823e-03f, 3.18153e-02f, -1.24154e-01f, 8.64812e-01f, 3.00980e-01f, -9.01591e-02f, 2.64621e-02f, -4.67520e-03f }, // 95/128
{ -5.17776e-03f, 3.11866e-02f, -1.22047e-01f, 8.71305e-01f, 2.91006e-01f, -8.77011e-02f, 2.57844e-02f, -4.55932e-03f }, // 96/128
{ -5.06363e-03f, 3.05345e-02f, -1.19837e-01f, 8.77666e-01f, 2.81060e-01f, -8.52109e-02f, 2.50946e-02f, -4.44107e-03f }, // 97/128
{ -4.94589e-03f, 2.98593e-02f, -1.17526e-01f, 8.83893e-01f, 2.71144e-01f, -8.26900e-02f, 2.43930e-02f, -4.32052e-03f }, // 98/128
{ -4.82456e-03f, 2.91609e-02f, -1.15113e-01f, 8.89984e-01f, 2.61263e-01f, -8.01399e-02f, 2.36801e-02f, -4.19774e-03f }, // 99/128
{ -4.69970e-03f, 2.84397e-02f, -1.12597e-01f, 8.95936e-01f, 2.51417e-01f, -7.75620e-02f, 2.29562e-02f, -4.07279e-03f }, // 100/128
{ -4.57135e-03f, 2.76957e-02f, -1.09978e-01f, 9.01749e-01f, 2.41609e-01f, -7.49577e-02f, 2.22218e-02f, -3.94576e-03f }, // 101/128
{ -4.43955e-03f, 2.69293e-02f, -1.07256e-01f, 9.07420e-01f, 2.31843e-01f, -7.23286e-02f, 2.14774e-02f, -3.81671e-03f }, // 102/128
{ -4.30435e-03f, 2.61404e-02f, -1.04430e-01f, 9.12947e-01f, 2.22120e-01f, -6.96762e-02f, 2.07233e-02f, -3.68570e-03f }, // 103/128
{ -4.16581e-03f, 2.53295e-02f, -1.01501e-01f, 9.18329e-01f, 2.12443e-01f, -6.70018e-02f, 1.99599e-02f, -3.55283e-03f }, // 104/128
{ -4.02397e-03f, 2.44967e-02f, -9.84679e-02f, 9.23564e-01f, 2.02814e-01f, -6.43069e-02f, 1.91877e-02f, -3.41815e-03f }, // 105/128
{ -3.87888e-03f, 2.36423e-02f, -9.53307e-02f, 9.28650e-01f, 1.93236e-01f, -6.15931e-02f, 1.84071e-02f, -3.28174e-03f }, // 106/128
{ -3.73062e-03f, 2.27664e-02f, -9.20893e-02f, 9.33586e-01f, 1.83711e-01f, -5.88617e-02f, 1.76185e-02f, -3.14367e-03f }, // 107/128
{ -3.57923e-03f, 2.18695e-02f, -8.87435e-02f, 9.38371e-01f, 1.74242e-01f, -5.61142e-02f, 1.68225e-02f, -3.00403e-03f }, // 108/128
{ -3.42477e-03f, 2.09516e-02f, -8.52933e-02f, 9.43001e-01f, 1.64831e-01f, -5.33522e-02f, 1.60193e-02f, -2.86289e-03f }, // 109/128
{ -3.26730e-03f, 2.00132e-02f, -8.17385e-02f, 9.47477e-01f, 1.55480e-01f, -5.05770e-02f, 1.52095e-02f, -2.72032e-03f }, // 110/128
{ -3.10689e-03f, 1.90545e-02f, -7.80792e-02f, 9.51795e-01f, 1.46192e-01f, -4.77900e-02f, 1.43934e-02f, -2.57640e-03f }, // 111/128
{ -2.94361e-03f, 1.80759e-02f, -7.43154e-02f, 9.55956e-01f, 1.36968e-01f, -4.49929e-02f, 1.35716e-02f, -2.43121e-03f }, // 112/128
{ -2.77751e-03f, 1.70776e-02f, -7.04471e-02f, 9.59958e-01f, 1.27812e-01f, -4.21869e-02f, 1.27445e-02f, -2.28483e-03f }, // 113/128
{ -2.60868e-03f, 1.60599e-02f, -6.64743e-02f, 9.63798e-01f, 1.18725e-01f, -3.93735e-02f, 1.19125e-02f, -2.13733e-03f }, // 114/128
{ -2.43718e-03f, 1.50233e-02f, -6.23972e-02f, 9.67477e-01f, 1.09710e-01f, -3.65541e-02f, 1.10760e-02f, -1.98880e-03f }, // 115/128
{ -2.26307e-03f, 1.39681e-02f, -5.82159e-02f, 9.70992e-01f, 1.00769e-01f, -3.37303e-02f, 1.02356e-02f, -1.83931e-03f }, // 116/128
{ -2.08645e-03f, 1.28947e-02f, -5.39305e-02f, 9.74342e-01f, 9.19033e-02f, -3.09033e-02f, 9.39154e-03f, -1.68894e-03f }, // 117/128
{ -1.90738e-03f, 1.18034e-02f, -4.95412e-02f, 9.77526e-01f, 8.31162e-02f, -2.80746e-02f, 8.54441e-03f, -1.53777e-03f }, // 118/128
{ -1.72594e-03f, 1.06946e-02f, -4.50483e-02f, 9.80543e-01f, 7.44095e-02f, -2.52457e-02f, 7.69462e-03f, -1.38589e-03f }, // 119/128
{ -1.54221e-03f, 9.56876e-03f, -4.04519e-02f, 9.83392e-01f, 6.57852e-02f, -2.24178e-02f, 6.84261e-03f, -1.23337e-03f }, // 120/128
{ -1.35627e-03f, 8.42626e-03f, -3.57525e-02f, 9.86071e-01f, 5.72454e-02f, -1.95925e-02f, 5.98883e-03f, -1.08030e-03f }, // 121/128
{ -1.16820e-03f, 7.26755e-03f, -3.09503e-02f, 9.88580e-01f, 4.87921e-02f, -1.67710e-02f, 5.13372e-03f, -9.26747e-04f }, // 122/128
{ -9.78093e-04f, 6.09305e-03f, -2.60456e-02f, 9.90917e-01f, 4.04274e-02f, -1.39548e-02f, 4.27773e-03f, -7.72802e-04f }, // 123/128
{ -7.86031e-04f, 4.90322e-03f, -2.10389e-02f, 9.93082e-01f, 3.21531e-02f, -1.11453e-02f, 3.42130e-03f, -6.18544e-04f }, // 124/128
{ -5.92100e-04f, 3.69852e-03f, -1.59305e-02f, 9.95074e-01f, 2.39714e-02f, -8.34364e-03f, 2.56486e-03f, -4.64053e-04f }, // 125/128
{ -3.96391e-04f, 2.47942e-03f, -1.07209e-02f, 9.96891e-01f, 1.58840e-02f, -5.55134e-03f, 1.70888e-03f, -3.09412e-04f }, // 126/128
{ -1.98993e-04f, 1.24642e-03f, -5.41054e-03f, 9.98534e-01f, 7.89295e-03f, -2.76968e-03f, 8.53777e-04f, -1.54700e-04f }, // 127/128
{ 0.00000e+00f, 0.00000e+00f, 0.00000e+00f, 1.00000e+00f, 0.00000e+00f, 0.00000e+00f, 0.00000e+00f, 0.00000e+00f }, // 128/128
};
}

View File

@@ -0,0 +1,58 @@
package io.github.dsheirer.dsp.filter.interpolator;
import io.github.dsheirer.sample.complex.Complex;
import org.apache.commons.lang3.Validate;
import org.slf4j.Logger;
import org.slf4j.LoggerFactory;
public class RealInterpolator extends Interpolator
{
private final static Logger mLog = LoggerFactory.getLogger( RealInterpolator.class );
private float mGain;
public RealInterpolator( float gain )
{
mGain = gain;
}
/**
* Calculates an interpolated value from the 8 samples beginning at the
* offset index. Uses the mu value to determine which of 128 filter kernels
* to use in order to closely approximate the ideal frequency response for
* the value of mu, ranging from 0.0 to 1.0, indicating the interpolated
* position in the sample set ranging from offset to offset + 7.
*
* @param samples - sample array of length at least offset + 7
*
* @param mu - interpolated sample position between 0 and 1.0
*
* @return - interpolated sample value
*/
public float filter( float[] samples, int offset, float mu )
{
/* Ensure we have enough samples in the array */
Validate.isTrue(samples.length >= offset + 7);
/* Identify the filter bank that corresponds to mu */
int index = (int)( NSTEPS * mu );
float accumulator = ( TAPS[ index ][ 7 ] * samples[ offset ] );
accumulator += ( TAPS[ index ][ 6 ] * samples[ offset + 1 ] );
accumulator += ( TAPS[ index ][ 5 ] * samples[ offset + 2 ] );
accumulator += ( TAPS[ index ][ 4 ] * samples[ offset + 3 ] );
accumulator += ( TAPS[ index ][ 3 ] * samples[ offset + 4 ] );
accumulator += ( TAPS[ index ][ 2 ] * samples[ offset + 5 ] );
accumulator += ( TAPS[ index ][ 1 ] * samples[ offset + 6 ] );
accumulator += ( TAPS[ index ][ 0 ] * samples[ offset + 7 ] );
return accumulator * mGain;
}
public Complex filter( float[] iSamples, float[] qSamples, int offset, float mu )
{
float i = filter( iSamples, offset, mu );
float q = filter( qSamples, offset, mu );
return new Complex( i, q );
}
}

View File

@@ -0,0 +1,116 @@
/*******************************************************************************
* sdrtrunk
* Copyright (C) 2014-2016 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp.filter.polyphase;
public class FilterStage
{
private float[] mBuffer;
private int mBufferSize;
private int mBufferPointer;
private int[][] mIndexMap;
private float[] mCoefficients;
/**
* Single stage of a polyphase filter implemented with a circular sample buffer that uses an index map to ensure
* correct alignment of samples to filter taps.
*
* @param coefficients
*/
public FilterStage(float[] coefficients)
{
mCoefficients = coefficients;
mBufferSize = mCoefficients.length;
mBuffer = new float[mBufferSize];
//Note: buffer pointer starts at an illegal index but is decremented before first usage
mBufferPointer = mBufferSize;
generateIndexMap(mBufferSize);
}
/**
* Coefficients used by this filter stage.
*/
public float[] getCoefficients()
{
return mCoefficients;
}
/**
* Cleanup prior to disposal
*/
public void dispose()
{
mCoefficients = null;
mIndexMap = null;
mBuffer = null;
}
/**
* Loads the sample into the sample buffer
*/
private void load(float sample)
{
mBufferPointer--;
if (mBufferPointer < 0)
{
mBufferPointer += mBufferSize;
}
mBuffer[mBufferPointer] = sample;
}
/**
* Applies convolution against the currently loaded buffer samples.
*/
public float filter(float sample)
{
load(sample);
float accumulator = 0.0f;
for (int x = 0; x < mBufferSize; x++)
{
accumulator += mCoefficients[x] * mBuffer[mIndexMap[mBufferPointer][x]];
}
return accumulator;
}
/**
* Generates a circular buffer index map to support lookup of the translated index based on the current buffer
* pointer and the desired sample index.
*/
private void generateIndexMap(int size)
{
mIndexMap = new int[size][size];
for (int x = 0; x < size; x++)
{
for (int y = 0; y < size; y++)
{
int z = x + y;
mIndexMap[x][y] = z < size ? z : z - size;
}
}
}
}

View File

@@ -0,0 +1,211 @@
/*******************************************************************************
* sdrtrunk
* Copyright (C) 2014-2016 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp.filter.polyphase;
import io.github.dsheirer.dsp.filter.design.FilterDesignException;
import io.github.dsheirer.dsp.filter.fir.FIRFilter;
import io.github.dsheirer.dsp.filter.fir.FIRFilterSpecification;
import io.github.dsheirer.dsp.filter.fir.remez.RemezFIRFilterDesigner;
import io.github.dsheirer.sample.Listener;
import io.github.dsheirer.sample.real.RealBuffer;
import io.github.dsheirer.sample.real.RealToRealBufferAssembler;
import org.slf4j.Logger;
import org.slf4j.LoggerFactory;
import io.github.dsheirer.dsp.mixer.Oscillator;
import java.util.ArrayList;
import java.util.Arrays;
import java.util.List;
public class PolyphaseFIRDecimatingFilter_RB extends FIRFilter implements Listener<RealBuffer>
{
private final static Logger mLog = LoggerFactory.getLogger( PolyphaseFIRDecimatingFilter_RB.class );
private List<FilterStage> mFilterStages = new ArrayList<>();
private int mFilterStagePointer;
private Listener<RealBuffer> mListener;
private RealToRealBufferAssembler mAssembler;
private int mDecimationRatio;
private float mAccumulator = 0.0f;
private float mGain = 1.0f;
/**
* Polyphase filter stage for decimating real-valued sample buffers by the decimation ratio. Creates a number of
* polyphase filter arms/stages equal to the decimation ratio.
*
* Note: output real-valued buffer sizes are determined by the first input buffer that is presented to the filter
* where the output buffer size is equal to the floor of the input buffer size divided by the decimation ratio.
*
* @param coefficients of a low pass filter designed for the input sample rate with a cutoff frequency that conforms
* to the nyquist frequency for an output sample rate equal to the input sample rate divided by the decimation ratio.
*
* @param decimationRatio to decimate the input sample rate by.
* @param gain to apply to the filtered output. Use a value of 1.0 to apply no gain.
*/
public PolyphaseFIRDecimatingFilter_RB(float[] coefficients, int decimationRatio, float gain)
{
mDecimationRatio = decimationRatio;
mGain = gain;
createFilterStages(coefficients, mDecimationRatio);
mFilterStagePointer = mFilterStages.size() - 1;
}
/**
* Creates polyphase filter stages (arms) from the set of coefficients.
*
* @param coefficients for the low-pass filter
* @param decimationRatio to determine the number of polyphase filter stages.
*/
private void createFilterStages(float[] coefficients, int decimationRatio)
{
int tapSize = (int)Math.ceil((double)coefficients.length / (double)decimationRatio);
float[][] coefficientSets = new float[decimationRatio][tapSize];
int stagePointer = 0;
int coefficientPointer = 0;
//Split the coefficients up into polyphase stage coefficient sets
for(int x = 0; x < coefficients.length; x++)
{
coefficientSets[stagePointer++][coefficientPointer] = coefficients[x];
if(stagePointer >= decimationRatio)
{
stagePointer = 0;
coefficientPointer++;
}
}
//Create the stages
for(int x = 0; x < decimationRatio; x++)
{
mFilterStages.add(new FilterStage(coefficientSets[x]));
}
}
public void dispose()
{
mListener = null;
}
@Override
public void receive(RealBuffer buffer)
{
if(mAssembler == null)
{
mAssembler = new RealToRealBufferAssembler(buffer.getSamples().length / mDecimationRatio);
mAssembler.setListener(mListener);
}
for(float sample: buffer.getSamples())
{
filter(sample);
}
int a = 0;
}
private void filter(float sample)
{
FilterStage stage = mFilterStages.get(mFilterStagePointer--);
mAccumulator += stage.filter(sample);
if(mFilterStagePointer < 0)
{
mAssembler.receive(mAccumulator * mGain);
mAccumulator = 0.0f;
mFilterStagePointer = mFilterStages.size() - 1;
}
}
public void setListener(Listener<RealBuffer> listener)
{
mListener = listener;
if(mAssembler != null)
{
mAssembler.setListener(mListener);
}
}
public void removeListener()
{
mListener = null;
if(mAssembler != null)
{
mAssembler.setListener(null);
}
}
public static void main(String[] args)
{
FIRFilterSpecification specification = FIRFilterSpecification.lowPassBuilder()
.sampleRate( 48000 )
.gridDensity( 16 )
.passBandCutoff( 3000 )
.passBandAmplitude( 1.0 )
.passBandRipple( 0.01 )
.stopBandStart( 4000 )
.stopBandAmplitude( 0.0 )
.stopBandRipple( 0.027 )
.build();
RemezFIRFilterDesigner designer = new RemezFIRFilterDesigner(specification);
try
{
if(!designer.isValid())
{
throw new FilterDesignException("Couldn't design the filter");
}
float[] coefficients = designer.getImpulseResponse();
PolyphaseFIRDecimatingFilter_RB filter = new PolyphaseFIRDecimatingFilter_RB(coefficients, 6, 1.0f);
filter.setListener(new Listener<RealBuffer>()
{
@Override
public void receive(RealBuffer realBuffer)
{
mLog.info("Filtered: " + Arrays.toString(realBuffer.getSamples()));
}
});
Oscillator oscillator = new Oscillator(3400, 48000);
float[] samples = new float[500];
for(int x = 0; x < 500; x++)
{
samples[x] = oscillator.getFloat();
oscillator.rotate();
}
RealBuffer buffer = new RealBuffer(samples);
filter.receive(buffer);
}
catch(FilterDesignException e)
{
mLog.debug("Error designing filter", e);
}
}
}

View File

@@ -0,0 +1,203 @@
/*******************************************************************************
* sdrtrunk
* Copyright (C) 2014-2016 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 <http://www.gnu.org/licenses/>
*
******************************************************************************/
package io.github.dsheirer.dsp.filter.polyphase;
import org.apache.commons.lang3.Validate;
import org.slf4j.Logger;
import org.slf4j.LoggerFactory;
/**
* Polyphase interpolating FIR filter with an internal circular buffer as
* described by the dspGuru at:
*
* http://www.dspguru.com/dsp/faqs/multirate/interpolation
* See sections 3.4.2 and 3.4.3 )
*/
public class PolyphaseFIRInterpolatingFilter
{
private final static Logger mLog =
LoggerFactory.getLogger( PolyphaseFIRInterpolatingFilter.class );
private double[] mBuffer;
private int mBufferPointer = 0;
private double[] mTaps;
private Phase[] mPhases;
/**
* Constructs the interpolating filter. Coefficients length must be an
* integer multiple of the interpolation value. Each sample processed by
* the filter will generate a quantity of samples equal to the interpolation
* value.
*
* @param coefficients - low pass filter coefficients sized for the desired
* interpolated sample rate with a cutoff frequency of no more than half of
* the original sample rate.
*
* @param interpolation - integral interpolation quantity
*/
public PolyphaseFIRInterpolatingFilter( double[] coefficients, int interpolation )
{
/* Ensure coefficients are an integral multiple of the interpolation factor */
Validate.isTrue(coefficients.length % interpolation == 0);
mTaps = coefficients;
mPhases = new Phase[ interpolation ];
int phaseSize = coefficients.length / interpolation;
mBuffer = new double[ phaseSize ];
for( int x = 0; x < interpolation; x++ )
{
int coefficientIndex = 0;
int[] indexes = new int[ phaseSize ];
for( int y = x; y < coefficients.length; y += interpolation )
{
indexes[ coefficientIndex++ ] = y;
}
/* Create a phase using the tap indexes and the interpolation value
* for the gain value, to offset the loss induced in upsampling from
* one sample to (interpolation) samples. */
mPhases[ interpolation - x - 1 ] = new Phase( indexes, interpolation );
}
}
/**
* Inserts the sample into the filter and returns (interpolation) number
* of interpolated samples in the return array.
*
* @param sample - float sample
* @return - array of interpolated samples, size = interpolation quantity
*/
public float[] interpolate( float sample )
{
mBuffer[ mBufferPointer ] = (double)sample;
mBufferPointer++;
mBufferPointer %= mBuffer.length;
float[] samples = new float[ mPhases.length ];
for( int x = 0; x < mPhases.length; x++ )
{
samples[ x ] = mPhases[ x ].filter();
}
return samples;
}
/**
* Upsamples and interpolates the array of samples, returning an array that
* is upsampled and low-pass filtered by the interpolation quantity.
*
* @param samples - array of float samples
* @return - array of interpolated samples, size = interpolation * samples size
*/
public float[] interpolate( float[] samples )
{
int pointer = 0;
float[] interpolated = new float[ samples.length * mPhases.length ];
for( float sample: samples )
{
mBuffer[ mBufferPointer ] = (double)sample;
mBufferPointer++;
mBufferPointer %= mBuffer.length;
for( int x = 0; x < mPhases.length; x++ )
{
interpolated[ pointer ] = mPhases[ x ].filter();
pointer++;
}
}
return interpolated;
}
/**
* Single phase element of a polyphase interpolating FIR filter. Each phase
* element shares a common buffer. The tap indexes argument defines the
* filter taps used by this phase element to generate an interpolated value
* from the shared circular buffer of samples.
*/
public class Phase
{
private int[] mIndexes;
private double mGain;
public Phase( int[] tapIndexes, int gain )
{
mIndexes = tapIndexes;
mGain = gain;
}
public float filter()
{
double accumulator = 0;
for( int x = 0; x < mIndexes.length; x++ )
{
accumulator += mTaps[ mIndexes[ x ] ] *
mBuffer[ ( x + mBufferPointer ) % mIndexes.length ];
}
return (float)( accumulator * mGain );
}
}
// public static void main( String[] args )
// {
// int interpolation = 6;
//
// PolyphaseFIRInterpolatingFilter p =
// new PolyphaseFIRInterpolatingFilter( IMBESynthesizer.INTERPOLATION_TAPS, interpolation );
//
// Oscillator o = new Oscillator( 500, 8000 );
//
// int iterations = 40;
//
// float[] in = new float[ iterations ];
// float[] out = new float[ iterations * interpolation ];
//
// int inIndex = 0;
// int outIndex = 0;
//
// for( int x = 0; x < iterations; x++ )
// {
// float next = o.nextFloat();
//
// in[ inIndex++ ] = next;
//
// float[] samples = p.interpolate( next );
//
// System.arraycopy( samples, 0, out, outIndex, samples.length );
//
// outIndex += interpolation;
// }
//
// mLog.debug( "In: " + Arrays.toString( in ) );
// mLog.debug( "Out: " + Arrays.toString( out ) );
// }
}

View File

@@ -0,0 +1,39 @@
package io.github.dsheirer.dsp.filter.smoothing;
import java.util.HashMap;
import java.util.Map;
public class GaussianSmoothingFilter extends SmoothingFilter
{
private static Map<Integer,float[]> mMap = new HashMap<>();
static
{
mMap.put( 3, new float[] { 0.1842153217f,0.6315693566f,0.1842153217f } );
mMap.put( 5, new float[] { 0.0646709835f,0.2446878423f,0.3812823483f,0.2446878423f,0.0646709835f } );
mMap.put( 7, new float[] { 0.0355966134f,0.1103628062f,0.2176061504f,0.27286886f,0.2176061504f,0.1103628062f,0.0355966134f } );
mMap.put( 9, new float[] { 0.0237625774f,0.061955345f,0.1228439993f,0.185233293f,0.2124095706f,0.185233293f,0.1228439993f,0.061955345f,0.0237625774f } );
mMap.put( 11, new float[] { 0.0175873706f,0.040124091f,0.0762094006f,0.120506627f,0.1586398759f,0.1738652697f,0.1586398759f,0.120506627f,0.0762094006f,0.040124091f,0.0175873706f } );
mMap.put( 13, new float[] { 0.0138646514f,0.0285342951f,0.0515030495f,0.0815279159f,0.1131845993f,0.137808478f,0.1471540217f,0.137808478f,0.1131845993f,0.0815279159f,0.0515030495f,0.0285342951f,0.0138646514f } );
mMap.put( 15, new float[] { 0.011400054f,0.0216349783f,0.0372047506f,0.0579739252f,0.0818576722f,0.1047317717f,0.1214199044f,0.1275538873f,0.1214199044f,0.1047317717f,0.0818576722f,0.0579739252f,0.0372047506f,0.0216349783f,0.011400054f } );
mMap.put( 17, new float[] { 0.0096581767f,0.0171732764f,0.0282802625f,0.0431306476f,0.0609201251f,0.07969074f,0.0965444243f,0.1083225312f,0.1125596326f,0.1083225312f,0.0965444243f,0.07969074f,0.0609201251f,0.0431306476f,0.0282802625f,0.0171732764f,0.0096581767f } );
mMap.put( 19, new float[] { 0.0083665107f,0.0141035977f,0.0223580953f,0.0333318211f,0.0467307368f,0.0616120364f,0.0763919939f,0.0890736966f,0.097672084f,0.1007188551f,0.097672084f,0.0890736966f,0.0763919939f,0.0616120364f,0.0467307368f,0.0333318211f,0.0223580953f,0.0141035977f,0.0083665107f } );
mMap.put( 21, new float[] { 0.0073728906f,0.0118884085f,0.0182292596f,0.0265811467f,0.0368585098f,0.0486027908f,0.0609458331f,0.0726751947f,0.0824114898f,0.0888686744f,0.091131604f,0.0888686744f,0.0824114898f,0.0726751947f,0.0609458331f,0.0486027908f,0.0368585098f,0.0265811467f,0.0182292596f,0.0118884085f,0.0073728906f } );
mMap.put( 23, new float[] { 0.0065861356f,0.010228416f,0.0152327635f,0.0217541246f,0.0297918261f,0.0391241787f,0.0492703835f,0.0595003129f,0.0689041144f,0.0765180086f,0.0814844543f,0.0832105633f,0.0814844543f,0.0765180086f,0.0689041144f,0.0595003129f,0.0492703835f,0.0391241787f,0.0297918261f,0.0217541246f,0.0152327635f,0.010228416f,0.0065861356f } );
mMap.put( 25, new float[] { 0.0059484873f,0.0089460248f,0.0129850288f,0.0181905149f,0.0245944051f,0.0320934871f,0.0404191102f,0.0491298959f,0.0576360481f,0.0652577058f,0.0713113508f,0.07520986f,0.076556162f,0.07520986f,0.0713113508f,0.0652577058f,0.0576360481f,0.0491298959f,0.0404191102f,0.0320934871f,0.0245944051f,0.0181905149f,0.0129850288f,0.0089460248f,0.0059484873f } );
mMap.put( 27, new float[] { 0.0054216802f,0.007930248f,0.0112519401f,0.0154865908f,0.0206762609f,0.0267778655f,0.0336409056f,0.0409965368f,0.0484634619f,0.0555737199f,0.0618176187f,0.0667026128f,0.0698170011f,0.0708871156f,0.0698170011f,0.0667026128f,0.0618176187f,0.0555737199f,0.0484634619f,0.0409965368f,0.0336409056f,0.0267778655f,0.0206762609f,0.0154865908f,0.0112519401f,0.007930248f,0.0054216802f } );
mMap.put( 29, new float[] { 0.0049794093f,0.0071086836f,0.0098843442f,0.013386093f,0.0176566017f,0.022683377f,0.0283828243f,0.0345900198f,0.0410575774f,0.0474660547f,0.0534466215f,0.0586144488f,0.0626089499f,0.0651351564f,0.0659996766f,0.0651351564f,0.0626089499f,0.0586144488f,0.0534466215f,0.0474660547f,0.0410575774f,0.0345900198f,0.0283828243f,0.022683377f,0.0176566017f,0.013386093f,0.0098843442f,0.0071086836f,0.0049794093f } );
mMap.put( 31, new float[] { 0.0046030224f,0.0064323659f,0.0087836656f,0.0117208269f,0.0152833344f,0.0194740097f,0.0242476721f,0.0295027244f,0.0350777394f,0.0407547727f,0.0462703471f,0.0513339215f,0.0556523529f,0.0589576307f,0.0610342929f,0.0617426429f,0.0610342929f,0.0589576307f,0.0556523529f,0.0513339215f,0.0462703471f,0.0407547727f,0.0350777394f,0.0295027244f,0.0242476721f,0.0194740097f,0.0152833344f,0.0117208269f,0.0087836656f,0.0064323659f,0.0046030224f } );
}
public GaussianSmoothingFilter()
{
super( mMap, 3 );
}
@Override
public SmoothingType getSmoothingType()
{
return SmoothingType.GAUSSIAN;
}
}

View File

@@ -0,0 +1,31 @@
package io.github.dsheirer.dsp.filter.smoothing;
public class NoSmoothingFilter extends SmoothingFilter
{
public NoSmoothingFilter()
{
super( null, 3 );
}
public float[] filter( float[] data )
{
return data;
}
public static float[] getCoefficients( int points )
{
return null;
}
public int getPointSize()
{
return SMOOTHING_DEFAULT;
}
@Override
public SmoothingType getSmoothingType()
{
return SmoothingType.NONE;
}
}

View File

@@ -0,0 +1,49 @@
package io.github.dsheirer.dsp.filter.smoothing;
import java.util.Arrays;
import java.util.HashMap;
import java.util.Map;
public class RectangularSmoothingFilter extends SmoothingFilter
{
private static Map<Integer,float[]> mMap = new HashMap<>();
static
{
mMap.put( 3, getCoefficients( 3 ) );
mMap.put( 5, getCoefficients( 5 ) );
mMap.put( 7, getCoefficients( 7 ) );
mMap.put( 9, getCoefficients( 9 ) );
mMap.put( 11, getCoefficients( 11 ) );
mMap.put( 13, getCoefficients( 13 ) );
mMap.put( 15, getCoefficients( 15 ) );
mMap.put( 17, getCoefficients( 17 ) );
mMap.put( 19, getCoefficients( 19 ) );
mMap.put( 21, getCoefficients( 21 ) );
mMap.put( 23, getCoefficients( 23 ) );
mMap.put( 25, getCoefficients( 25 ) );
mMap.put( 27, getCoefficients( 27 ) );
mMap.put( 29, getCoefficients( 29 ) );
mMap.put( 31, getCoefficients( 31 ) );
}
public RectangularSmoothingFilter()
{
super( mMap, 3 );
}
public static float[] getCoefficients( int points )
{
float[] coefficients = new float[ points ];
Arrays.fill( coefficients, 1.0f / (float)points );
return coefficients;
}
@Override
public SmoothingType getSmoothingType()
{
return SmoothingType.RECTANGLE;
}
}

View File

@@ -0,0 +1,97 @@
package io.github.dsheirer.dsp.filter.smoothing;
import java.util.Map;
import java.util.Set;
public abstract class SmoothingFilter
{
public static final int SMOOTHING_MINIMUM = 3;
public static final int SMOOTHING_MAXIMUM = 29;
public static final int SMOOTHING_DEFAULT = 9;
public enum SmoothingType
{
NONE,
RECTANGLE,
TRIANGLE,
GAUSSIAN;
}
private Map<Integer,float[]> mCoefficientMap;
private float[] mCoefficients;
private float[] mNewCoefficients;
private boolean mNewCoefficientsAvailable;
public SmoothingFilter( Map<Integer,float[]> coefficients, int index )
{
mCoefficientMap = coefficients;
if( mCoefficientMap != null )
{
mCoefficients = mCoefficientMap.get( index );
}
}
public abstract SmoothingType getSmoothingType();
public Set<Integer> getPointSizeList()
{
return mCoefficientMap.keySet();
}
public int getPointSize()
{
return mCoefficients.length;
}
public void setPointSize( int pointSize )
{
if( mCoefficientMap != null )
{
if( mCoefficientMap.containsKey( pointSize ) )
{
mNewCoefficients = mCoefficientMap.get( pointSize );
mNewCoefficientsAvailable = true;
}
else
{
throw new IllegalArgumentException( "Point size [" + pointSize +
"] is not valid in the current coefficient map" );
}
}
}
public float[] filter( float[] data )
{
if( mNewCoefficientsAvailable )
{
mCoefficients = mNewCoefficients;
mNewCoefficientsAvailable = false;
}
int middle = mCoefficients.length / 2;
float[] filtered = new float[ data.length ];
int toCopy = middle;
System.arraycopy( data, 0, filtered, 0, toCopy );
System.arraycopy( data, data.length - toCopy, filtered, filtered.length - toCopy, toCopy );
float accumulator;
for( int x = 0; x < data.length - mCoefficients.length + 1; x++ )
{
accumulator = 0.0f;
for( int y = 0; y < mCoefficients.length; y++ )
{
accumulator += data[ x + y ] * mCoefficients[ y ];
}
filtered[ x + middle ] = accumulator;
}
return filtered;
}
}

View File

@@ -0,0 +1,38 @@
package io.github.dsheirer.dsp.filter.smoothing;
import java.util.HashMap;
import java.util.Map;
public class TriangularSmoothingFilter extends SmoothingFilter
{
private static Map<Integer,float[]> mMap = new HashMap<>();
static
{
mMap.put( 3, new float[]{ 0.25f,0.5f,0.25f } );
mMap.put( 5, new float[]{ 0.1111111111f,0.2222222222f,0.3333333333f,0.2222222222f,0.1111111111f } );
mMap.put( 7, new float[]{ 0.0625f,0.125f,0.1875f,0.25f,0.1875f,0.125f,0.0625f } );
mMap.put( 9, new float[]{ 0.04f,0.08f,0.12f,0.16f,0.2f,0.16f,0.12f,0.08f,0.04f } );
mMap.put( 11, new float[]{ 0.0277777778f,0.0555555556f,0.0833333333f,0.1111111111f,0.1388888889f,0.1666666667f,0.1388888889f,0.1111111111f,0.0833333333f,0.0555555556f,0.0277777778f } );
mMap.put( 13, new float[]{ 0.0204081633f,0.0408163265f,0.0612244898f,0.0816326531f,0.1020408163f,0.1224489796f,0.1428571429f,0.1224489796f,0.1020408163f,0.0816326531f,0.0612244898f,0.0408163265f,0.0204081633f } );
mMap.put( 15, new float[]{ 0.015625f,0.03125f,0.046875f,0.0625f,0.078125f,0.09375f,0.109375f,0.125f,0.109375f,0.09375f,0.078125f,0.0625f,0.046875f,0.03125f,0.015625f } );
mMap.put( 17, new float[]{ 0.012345679f,0.024691358f,0.037037037f,0.049382716f,0.0617283951f,0.0740740741f,0.0864197531f,0.0987654321f,0.1111111111f,0.0987654321f,0.0864197531f,0.0740740741f,0.0617283951f,0.049382716f,0.037037037f,0.024691358f,0.012345679f } );
mMap.put( 19, new float[]{ 0.01f,0.02f,0.03f,0.04f,0.05f,0.06f,0.07f,0.08f,0.09f,0.1f,0.09f,0.08f,0.07f,0.06f,0.05f,0.04f,0.03f,0.02f,0.01f } );
mMap.put( 21, new float[]{ 0.0082644628f,0.0165289256f,0.0247933884f,0.0330578512f,0.041322314f,0.0495867769f,0.0578512397f,0.0661157025f,0.0743801653f,0.0826446281f,0.0909090909f,0.0826446281f,0.0743801653f,0.0661157025f,0.0578512397f,0.0495867769f,0.041322314f,0.0330578512f,0.0247933884f,0.0165289256f,0.0082644628f } );
mMap.put( 23, new float[]{ 0.0069444444f,0.0138888889f,0.0208333333f,0.0277777778f,0.0347222222f,0.0416666667f,0.0486111111f,0.0555555556f,0.0625f,0.0694444444f,0.0763888889f,0.0833333333f,0.0763888889f,0.0694444444f,0.0625f,0.0555555556f,0.0486111111f,0.0416666667f,0.0347222222f,0.0277777778f,0.0208333333f,0.0138888889f,0.0069444444f } );
mMap.put( 25, new float[]{ 0.0059171598f,0.0118343195f,0.0177514793f,0.0236686391f,0.0295857988f,0.0355029586f,0.0414201183f,0.0473372781f,0.0532544379f,0.0591715976f,0.0650887574f,0.0710059172f,0.0769230769f,0.0710059172f,0.0650887574f,0.0591715976f,0.0532544379f,0.0473372781f,0.0414201183f,0.0355029586f,0.0295857988f,0.0236686391f,0.0177514793f,0.0118343195f,0.0059171598f } );
mMap.put( 27, new float[]{ 0.0051020408f,0.0102040816f,0.0153061224f,0.0204081633f,0.0255102041f,0.0306122449f,0.0357142857f,0.0408163265f,0.0459183673f,0.0510204082f,0.056122449f,0.0612244898f,0.0663265306f,0.0714285714f,0.0663265306f,0.0612244898f,0.056122449f,0.0510204082f,0.0459183673f,0.0408163265f,0.0357142857f,0.0306122449f,0.0255102041f,0.0204081633f,0.0153061224f,0.0102040816f,0.0051020408f } );
mMap.put( 29, new float[]{ 0.0044444444f,0.0088888889f,0.0133333333f,0.0177777778f,0.0222222222f,0.0266666667f,0.0311111111f,0.0355555556f,0.04f,0.0444444444f,0.0488888889f,0.0533333333f,0.0577777778f,0.0622222222f,0.0666666667f,0.0622222222f,0.0577777778f,0.0533333333f,0.0488888889f,0.0444444444f,0.04f,0.0355555556f,0.0311111111f,0.0266666667f,0.0222222222f,0.0177777778f,0.0133333333f,0.0088888889f,0.0044444444f } );
}
public TriangularSmoothingFilter()
{
super( mMap, 3 );
}
@Override
public SmoothingType getSmoothingType()
{
return SmoothingType.TRIANGLE;
}
}

View File

@@ -0,0 +1,78 @@
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014,2015 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp.fm;
public abstract class FMDemodulator
{
private float mPreviousI = 0.0f;
private float mPreviousQ = 0.0f;
protected float mGain;
public FMDemodulator( float gain )
{
mGain = gain;
}
public float demodulate( float currentI, float currentQ )
{
/**
* Multiply the current sample against the complex conjugate of the
* previous sample to derive the phase delta between the two samples
*
* Negating the previous sample quadrature produces the conjugate
*/
double inphase = ( currentI * mPreviousI ) - ( currentQ * -mPreviousQ );
double quadrature = ( currentQ * mPreviousI ) + ( currentI * -mPreviousQ );
double angle = 0.0f;
//Check for divide by zero
if( inphase != 0 )
{
/**
* Use the arc-tangent of quadrature divided by inphase to
* get the phase angle (+/-) which was directly manipulated by the
* original message waveform during the modulation. This value now
* serves as the instantaneous amplitude of the demodulated signal
*/
double denominator = 1.0d / inphase;
angle = Math.atan( (double)quadrature * denominator );
}
/**
* Store the current sample to use during the next iteration
*/
mPreviousI = currentI;
mPreviousQ = currentQ;
return (float)( angle * mGain );
}
public abstract void dispose();
public void reset()
{
mPreviousI = 0.0f;
mPreviousQ = 0.0f;
}
public void setGain( float gain )
{
mGain = gain;
}
}

View File

@@ -0,0 +1,79 @@
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014,2015 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp.fm;
import io.github.dsheirer.sample.Listener;
import io.github.dsheirer.sample.complex.ComplexBuffer;
import io.github.dsheirer.sample.real.RateCounter_RB;
import io.github.dsheirer.sample.real.RealBuffer;
import org.slf4j.Logger;
import org.slf4j.LoggerFactory;
public class FMDemodulator_CB extends FMDemodulator implements Listener<ComplexBuffer>
{
private final static Logger mLog =
LoggerFactory.getLogger( FMDemodulator_CB.class );
private Listener<RealBuffer> mListener;
private RateCounter_RB mRateCounter = new RateCounter_RB( "FM DEMOD OUTPUT" );
public FMDemodulator_CB( float gain )
{
super( gain );
}
@Override
public void receive( ComplexBuffer buffer )
{
if( mListener != null )
{
float[] samples = buffer.getSamples();
int half = samples.length / 2;
float[] demodulated = new float[ half ];
for( int x = 0; x < half; x++ )
{
int index = x * 2;
demodulated[ x ] = demodulate( samples[ index ],
samples[ index + 1 ] );
}
mListener.receive( new RealBuffer( demodulated ) );
}
}
public void setListener( Listener<RealBuffer> listener )
{
mListener = listener;
}
public void removeListener()
{
mListener = null;
}
@Override
public void dispose()
{
mListener = null;
}
}

View File

@@ -0,0 +1,102 @@
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014,2015 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp.fm;
import io.github.dsheirer.sample.Listener;
import io.github.dsheirer.sample.complex.Complex;
import io.github.dsheirer.sample.real.RealSampleListener;
@Deprecated
public class FMDiscriminator implements Listener<Complex>
{
private RealSampleListener mListener;
private Complex mPreviousSample = new Complex( 0.0f, 0.0f );
private double mGain;
/**
* Implements a polar discriminator with ArcTangent angle estimator to
* demodulate complex sampled frequency modulated signals.
*
* @param gain - gain to be applied to the demodulated output - can be
* dynamically applied to continuously adjust audio output
*/
public FMDiscriminator( double gain )
{
mGain = gain;
}
public void dispose()
{
mListener = null;
mPreviousSample = null;
}
public synchronized void setGain( double gain )
{
mGain = gain;
}
@Override
public void receive( Complex currentSample )
{
/**
* Multiply the current sample against the complex conjugate of the
* previous sample to derive the phase delta between the two samples
*
* Negating the previous sample quadrature produces the conjugate
*/
double i = ( currentSample.inphase() * mPreviousSample.inphase() ) -
( currentSample.quadrature() * -mPreviousSample.quadrature() );
double q = ( currentSample.quadrature() * mPreviousSample.inphase() ) +
( currentSample.inphase() * -mPreviousSample.quadrature() );
double angle;
//Check for divide by zero
if( i == 0 )
{
angle = 0.0;
}
else
{
/**
* Use the arc-tangent of imaginary (q) divided by real (i) to
* get the phase angle (+/-) which was directly manipulated by the
* original message waveform during the modulation. This value now
* serves as the instantaneous amplitude of the demodulated signal
*/
double denominator = 1.0d / i;
angle = Math.atan( (double)q * denominator );
}
if( mListener != null )
{
mListener.receive( (float)( angle * mGain ) );
}
/**
* Store the current sample to use during the next iteration
*/
mPreviousSample = currentSample;
}
public void setListener( RealSampleListener listener )
{
mListener = listener;
}
}

View File

@@ -0,0 +1,302 @@
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp.fsk;
import io.github.dsheirer.buffer.BooleanAveragingBuffer;
import io.github.dsheirer.instrument.Instrumentable;
import io.github.dsheirer.instrument.tap.Tap;
import io.github.dsheirer.instrument.tap.TapGroup;
import io.github.dsheirer.instrument.tap.stream.SymbolEventTap;
import io.github.dsheirer.sample.Listener;
import io.github.dsheirer.sample.real.RealBuffer;
import org.slf4j.Logger;
import org.slf4j.LoggerFactory;
import io.github.dsheirer.dsp.symbol.SymbolEvent;
import io.github.dsheirer.dsp.symbol.SymbolEvent.Shift;
import java.util.ArrayList;
import java.util.BitSet;
import java.util.List;
/**
* Binary Frequency Shift Keying (FSK) decoder. Implements a BFSK correlation
* decoder. Provides normal or inverted decoded output. Requires the symbol
* rate to be an integer multiple of the sample rate.
*
* Implements a correlation decoder that correlates a modulated FSK signal upon
* against a one-baud delayed version of itself, converting all float samples to
* binary (0 or 1), and uses XOR to produce the correlation value. Low pass
* filtering smoothes the correlated output.
*
* Automatically aligns to symbol timing during symbol transitions (0 to 1, or
* 1 to 0) by inspecting the samples at the symbol edges and advancing or
* retarding symbol window to maintain continuous symbol alignment.
*
* Use a DC-removal filter prior to this decoder to ensure samples don't have
* a DC component.
*
* Implements instrumentable interface, so that slice events can be received
* externally to analyze decoder performance.
*/
public class FSK2Decoder implements Instrumentable, Listener<RealBuffer>
{
private final static Logger mLog = LoggerFactory.getLogger(FSK2Decoder.class);
public enum Output
{
NORMAL, INVERTED
}
/* Instrumentation taps */
private static final String INSTRUMENT_DECISION =
"Tap Point: FSK2 Symbol Decision";
private List<TapGroup> mAvailableTaps;
private ArrayList<SymbolEventTap> mTaps = new ArrayList<SymbolEventTap>();
private Listener<Boolean> mListener;
private BooleanAveragingBuffer mDelayBuffer;
private BooleanAveragingBuffer mLowPassFilter;
private Slicer mSlicer;
private boolean mNormalOutput;
private int mSamplesPerSymbol;
private int mSymbolRate;
public FSK2Decoder(int sampleRate, int symbolRate, Output output)
{
/* Ensure we're using an integral symbolRate of the sampleRate */
assert (sampleRate % symbolRate == 0);
mSamplesPerSymbol = (int) (sampleRate / symbolRate);
mNormalOutput = (output == Output.NORMAL);
mSymbolRate = symbolRate;
mDelayBuffer = new BooleanAveragingBuffer(mSamplesPerSymbol);
mLowPassFilter = new BooleanAveragingBuffer(mSamplesPerSymbol);
mSlicer = new Slicer(mSamplesPerSymbol);
}
/**
* Disposes of all references to prepare for garbage collection
*/
public void dispose()
{
mListener = null;
}
/**
* Primary sample input
*/
@Override
public void receive(RealBuffer buffer)
{
for(float sample : buffer.getSamples())
{
/* Square the sample. Greater than zero is a 1 (true) and less than
* zero is a 0 (false) */
boolean bitSample = (sample >= 0.0f);
/* Feed the delay buffer and fetch the one-baud delayed sample */
boolean delayedBitSample = mDelayBuffer.get(bitSample);
/* Correlation: xor current bit with delayed bit */
boolean softBit = bitSample ^ delayedBitSample;
/* Low pass filter to smooth the correlated values */
boolean filteredSoftBit = mLowPassFilter.getAverage(softBit);
/* Send the filtered correlated bit to the slicer */
mSlicer.receive(filteredSoftBit);
}
}
/**
* Registers a listener to receive the decoded FSK bits
*/
public void setListener(Listener<Boolean> listener)
{
mListener = listener;
}
/**
* Removes the listener
*/
public void removeListener(Listener<Boolean> listener)
{
mListener = null;
}
/**
* Symbol slicer with auto-aligning baud timing
*/
public class Slicer
{
private BitSet mBitSet = new BitSet();
private int mSymbolLength;
private int mDecisionThreshold;
private int mSampleCounter;
public Slicer(int samplesPerSymbol)
{
mSymbolLength = samplesPerSymbol;
mDecisionThreshold = (int) (mSymbolLength / 2);
/* Adjust for an odd number of samples per baud */
if(mSymbolLength % 2 == 1)
{
mDecisionThreshold++;
}
}
public void receive(boolean softBit)
{
if(mSampleCounter >= 0)
{
if(softBit)
{
mBitSet.set(mSampleCounter);
}
else
{
mBitSet.clear(mSampleCounter);
}
}
mSampleCounter++;
if(mSampleCounter >= mSymbolLength)
{
boolean decision = mBitSet.cardinality() >= mDecisionThreshold;
send(decision);
/* Shift timing left if the left bit in the bitset is opposite
* the decision and the right bit is the same */
if((mBitSet.get(0) ^ decision) &&
(!(mBitSet.get(mSymbolLength - 1) ^ decision)))
{
sendTapEvent(mBitSet, Shift.LEFT, decision);
reset();
mSampleCounter--;
}
/* Shift timing right if the left bit is the same as the
* decision and the right bit is opposite */
else if((!(mBitSet.get(0) ^ decision)) &&
(mBitSet.get(mSymbolLength - 1) ^ decision))
{
sendTapEvent(mBitSet, Shift.RIGHT, decision);
/* Last bit from previous symbol to pre-fill next symbol */
boolean previousSoftBit = mBitSet.get(mSymbolLength - 1);
reset();
if(previousSoftBit)
{
mBitSet.set(0);
}
mSampleCounter++;
}
/* No shift */
else
{
sendTapEvent(mBitSet, Shift.NONE, decision);
reset();
}
}
}
/**
* Sends the bit decision to the listener
*/
private void send(boolean decision)
{
if(mListener != null)
{
mListener.receive(mNormalOutput ? decision : !decision);
}
}
private void reset()
{
mBitSet.clear();
mSampleCounter = 0;
}
/**
* Sends instrumentation tap event to all registered listeners
*/
private void sendTapEvent(BitSet bitset, Shift shift, boolean decision)
{
for(SymbolEventTap tap : mTaps)
{
SymbolEvent event =
new SymbolEvent(bitset.get(0, mSymbolLength),
mSymbolLength,
decision,
shift);
tap.receive(event);
}
}
}
/**
* Get instrumentation taps
*/
@Override
public List<TapGroup> getTapGroups()
{
if(mAvailableTaps == null)
{
mAvailableTaps = new ArrayList<TapGroup>();
TapGroup tapGroup = new TapGroup("FSK2 Decoder");
tapGroup.add(new SymbolEventTap(INSTRUMENT_DECISION, 0, .025f));
mAvailableTaps.add(tapGroup);
}
return mAvailableTaps;
}
/**
* Add instrumentation tap
*/
@Override
public void registerTap(Tap tap)
{
if(tap instanceof SymbolEventTap && !mTaps.contains(tap))
{
mTaps.add((SymbolEventTap) tap);
}
}
/**
* Remove instrumentation tap
*/
@Override
public void unregisterTap(Tap tap)
{
mTaps.remove(tap);
}
}

View File

@@ -0,0 +1,253 @@
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp.fsk;
import io.github.dsheirer.dsp.filter.Filters;
import io.github.dsheirer.dsp.filter.FloatHalfBandNoDecimateFilter;
import io.github.dsheirer.dsp.filter.LTRPulseShapingFilter;
import io.github.dsheirer.dsp.filter.SquaringFilter;
import io.github.dsheirer.instrument.Instrumentable;
import io.github.dsheirer.instrument.tap.Tap;
import io.github.dsheirer.instrument.tap.TapGroup;
import io.github.dsheirer.instrument.tap.stream.BinaryTap;
import io.github.dsheirer.instrument.tap.stream.FloatBufferTap;
import io.github.dsheirer.instrument.tap.stream.FloatTap;
import io.github.dsheirer.instrument.tap.stream.SymbolEventTap;
import io.github.dsheirer.sample.Listener;
import io.github.dsheirer.sample.real.RealBuffer;
import io.github.dsheirer.dsp.filter.FloatHalfBandFilter;
import io.github.dsheirer.dsp.filter.Window.WindowType;
import io.github.dsheirer.dsp.filter.cic.RealPrimeCICDecimate;
import io.github.dsheirer.dsp.filter.dc.IIRSinglePoleDCRemovalFilter;
import io.github.dsheirer.dsp.symbol.Slicer;
import io.github.dsheirer.dsp.symbol.Slicer.Output;
import java.util.ArrayList;
import java.util.List;
public class LTRFSKDecoder implements Listener<RealBuffer>, Instrumentable
{
private RealPrimeCICDecimate mCICDecimator;
private FloatHalfBandFilter mHBFilter1;
private FloatHalfBandFilter mHBFilter2;
private FloatHalfBandFilter mHBFilter3;
private FloatHalfBandFilter mHBFilter4;
private FloatHalfBandFilter mHBFilter5;
private FloatHalfBandNoDecimateFilter mHBFilter6;
private IIRSinglePoleDCRemovalFilter mDCFilter;
private SquaringFilter mSquaringFilter;
private LTRPulseShapingFilter mPulseShaper;
private Slicer mSlicer;
private List<TapGroup> mAvailableTaps;
private final String TAP_F1_F2 = "LTR FSK Demod Filter1 >< Filter2";
private final String TAP_F2_F3 = "LTR FSK Demod Filter2 >< Filter3";
private final String TAP_F3_F4 = "LTR FSK Demod Filter3 >< Filter4";
private final String TAP_F4_F5 = "LTR FSK Demod Filter4 >< Filter5";
private final String TAP_F5_F6 = "LTR FSK Demod Filter5 >< Filter6";
private final String TAP_F6_DT = "LTR FSK Demod Filter6 >< DC Filter";
private final String TAP_CIC = "LTR FSK Demod CIC Decimator";
private final String TAP_DT_SF1 = "LTR FSK Demod DC Filter >< Squaring Filter";
private final String TAP_SF1_PS = "LTR FSK Demod Squaring Filter >< Pulse Shaper";
private final String TAP_PS_SLICER = "LTR FSK Demod Pulse Shaper >< Slicer";
private final String TAP_SYMBOL_EVENT = "LTR FSK Demod Symbol Event";
/**
* Implements a Logic Trunked Radio sub-audible 300 baud FSK signaling
* decoder. Expects a 48000 sample rate input.
*/
public LTRFSKDecoder()
{
mCICDecimator = new RealPrimeCICDecimate( 32, 1, 375, 60, WindowType.HAMMING, 2 );
mHBFilter1 = new FloatHalfBandFilter(
Filters.FIR_HALF_BAND_31T_ONE_EIGHTH_FCO, 1.1002f );
mHBFilter2 = new FloatHalfBandFilter(
Filters.FIR_HALF_BAND_31T_ONE_EIGHTH_FCO, 1.1002f );
mHBFilter1.setListener( mHBFilter2 );
mHBFilter3 = new FloatHalfBandFilter(
Filters.FIR_HALF_BAND_31T_ONE_EIGHTH_FCO, 1.1002f );
mHBFilter2.setListener( mHBFilter3 );
mHBFilter4 = new FloatHalfBandFilter(
Filters.FIR_HALF_BAND_31T_ONE_EIGHTH_FCO, 1.1002f );
mHBFilter3.setListener( mHBFilter4 );
mHBFilter5 = new FloatHalfBandFilter(
Filters.FIR_HALF_BAND_31T_ONE_EIGHTH_FCO, 1.1002f );
mHBFilter4.setListener( mHBFilter5 );
mHBFilter6 = new FloatHalfBandNoDecimateFilter(
Filters.FIR_HALF_BAND_31T_ONE_EIGHTH_FCO, 1.1002f );
mHBFilter5.setListener( mHBFilter6 );
mDCFilter = new IIRSinglePoleDCRemovalFilter( 0.9946f );
mHBFilter6.setListener( mDCFilter );
mSquaringFilter = new SquaringFilter();
mDCFilter.setListener( mSquaringFilter );
mPulseShaper = new LTRPulseShapingFilter();
mSquaringFilter.setListener( mPulseShaper );
mSlicer = new Slicer( Output.NORMAL, 5 );
mPulseShaper.setListener( mSlicer );
}
@Override
public void receive( RealBuffer buffer )
{
for( float sample: buffer.getSamples() )
{
mHBFilter1.receive( sample );
}
mCICDecimator.receive(buffer);
}
public void addListener( Listener<Boolean> listener )
{
mSlicer.setListener( listener );
}
public void removeListener( Listener<Boolean> listener )
{
mSlicer.removeListener( listener );
}
@Override
public List<TapGroup> getTapGroups()
{
if( mAvailableTaps == null )
{
mAvailableTaps = new ArrayList<TapGroup>();
TapGroup group = new TapGroup( "LTR FSK Decoder" );
group.add( new FloatTap( TAP_F1_F2, 0, 1.0f ) );
group.add( new FloatTap( TAP_F2_F3, 8, 0.5f ) );
group.add( new FloatTap( TAP_F3_F4, 12, 0.25f ) );
group.add( new FloatTap( TAP_F4_F5, 14, 0.125f ) );
group.add( new FloatTap( TAP_F5_F6, 15, 0.0625f ) );
group.add( new FloatTap( TAP_F6_DT, 31, 0.0625f ) );
group.add( new FloatBufferTap( TAP_CIC, 31, 0.0625f ) );
group.add( new FloatTap( TAP_DT_SF1, 31, 0.0625f ) );
group.add( new BinaryTap( TAP_SF1_PS, 33, 0.0625f ) );
group.add( new BinaryTap( TAP_PS_SLICER, 33, 0.0625f ) );
group.add( new SymbolEventTap( TAP_SYMBOL_EVENT, 7, 0.0125f ) );
mAvailableTaps.add( group );
}
return mAvailableTaps;
}
@Override
public void registerTap( Tap tap )
{
switch( tap.getName() )
{
case TAP_CIC:
mCICDecimator.setListener( (FloatBufferTap)tap );
break;
case TAP_F1_F2:
mHBFilter1.setListener( (FloatTap)tap );
((FloatTap)tap).setListener( mHBFilter2 );
break;
case TAP_F2_F3:
mHBFilter2.setListener( (FloatTap)tap );
((FloatTap)tap).setListener( mHBFilter3 );
break;
case TAP_F3_F4:
mHBFilter3.setListener( (FloatTap)tap );
((FloatTap)tap).setListener( mHBFilter4 );
break;
case TAP_F4_F5:
mHBFilter4.setListener( (FloatTap)tap );
((FloatTap)tap).setListener( mHBFilter5 );
break;
case TAP_F5_F6:
mHBFilter5.setListener( (FloatTap)tap );
((FloatTap)tap).setListener( mHBFilter6 );
break;
case TAP_F6_DT:
mHBFilter6.setListener( (FloatTap)tap );
((FloatTap)tap).setListener( mDCFilter );
break;
case TAP_DT_SF1:
mDCFilter.setListener( (FloatTap)tap );
((FloatTap)tap).setListener( mSquaringFilter );
break;
case TAP_SF1_PS:
mSquaringFilter.setListener( (BinaryTap)tap );
((BinaryTap)tap).setListener( mPulseShaper );
break;
case TAP_PS_SLICER:
mPulseShaper.setListener( (BinaryTap)tap );
((BinaryTap)tap).setListener( mSlicer );
break;
case TAP_SYMBOL_EVENT:
mSlicer.addTap( (SymbolEventTap)tap );
break;
}
}
@Override
public void unregisterTap( Tap tap )
{
switch( tap.getName() )
{
case TAP_CIC:
mCICDecimator.removeListener();
break;
case TAP_F1_F2:
mHBFilter1.setListener( mHBFilter2 );
break;
case TAP_F2_F3:
mHBFilter2.setListener( mHBFilter3 );
break;
case TAP_F3_F4:
mHBFilter3.setListener( mHBFilter4 );
break;
case TAP_F4_F5:
mHBFilter4.setListener( mHBFilter5 );
break;
case TAP_F5_F6:
mHBFilter5.setListener( mHBFilter6 );
break;
case TAP_F6_DT:
mHBFilter6.setListener( mDCFilter );
break;
case TAP_DT_SF1:
mDCFilter.setListener( mSquaringFilter );
break;
case TAP_SF1_PS:
mSquaringFilter.setListener( mPulseShaper );
break;
case TAP_PS_SLICER:
mPulseShaper.setListener( mSlicer );
break;
case TAP_SYMBOL_EVENT:
mSlicer.removeTap();
break;
}
}
}

View File

@@ -0,0 +1,224 @@
package io.github.dsheirer.dsp.gain;
import io.github.dsheirer.buffer.DoubleCircularBuffer;
import io.github.dsheirer.buffer.FloatCircularBuffer;
import java.util.concurrent.atomic.AtomicBoolean;
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014 Dennis Sheirer
* Copyright (C) 2011 Alex Csete
* Copyright (C) 2010 Moe Wheatley
*
* 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 <http://www.gnu.org/licenses/>
*
* ---------------------------------------------------
* Ported from Alex Csete's gqrx source code at:
* https://github.com/csete/gqrx/blob/74ba2742c89b0dfb66854fc9f72d19e5c6e355b6/dsp/agc_impl.cpp
*
* This Software is released under the "Simplified BSD License" + + +
* Copyright 2010 Moe Wheatley. All rights reserved.
*
* Redistribution and use in source and binary forms, with or without
* modification, are permitted provided that the following conditions are
* met:
*
* 1. Redistributions of source code must retain the above copyright notice,
* this list of conditions and the following disclaimer.
*
* 2. Redistributions in binary form must reproduce the above copyright
* notice, this list of conditions and the following disclaimer in the
* documentation and/or other materials provided with the distribution.
*
* THIS SOFTWARE IS PROVIDED BY Moe Wheatley ``AS IS'' AND ANY EXPRESS OR
* IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
* OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
* IN NO EVENT SHALL Moe Wheatley OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT,
* INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
* (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
* SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
* HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT,
* STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
* ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
* POSSIBILITY OF SUCH DAMAGE.
*
* The views and conclusions contained in the software and documentation are
* those of the authors and should not be interpreted as representing
* official policies, either expressed or implied, of Moe Wheatley.
******************************************************************************/
public class AutomaticGainControl
{
private static final double SAMPLE_RATE = 48000;
/* Signal delay line - time delay in seconds */
// private static final double DELAY_TIME_CONSTANT = 0.015;
private static final double DELAY_TIME_CONSTANT = 0.01;
/* Peak detector window - time delay in seconds */
private static final double WINDOW_TIME_CONSTANT = 0.018;
/* Attack time constants in seconds */
private static final double ATTACK_RISE_TIME_CONSTANT = 0.002;
private static final double ATTACK_FALL_TIME_CONSTANT = 0.005;
private static final double ATTACK_RISE_ALPHA = 1.0 -
Math.exp( -1.0 / SAMPLE_RATE * ATTACK_RISE_TIME_CONSTANT );
private static final double ATTACK_FALL_ALPHA = 1.0 -
Math.exp( -1.0 / SAMPLE_RATE * ATTACK_FALL_TIME_CONSTANT );
/* AGC decay value in milliseconds (20 to 5000) */
private static final double DECAY = 200;
/* Ratio between rise and fall times of decay time constants - adjust for
* best action with SSB */
private static final double DECAY_RISEFALL_RATIO = 0.3;
private static final double DECAY_RISE_ALPHA = 1.0 -
Math.exp( -1.0 / ( SAMPLE_RATE * DECAY * .001 * DECAY_RISEFALL_RATIO ) );
private static final double DECAY_FALL_ALPHA = 1.0 -
Math.exp( -1.0 / ( SAMPLE_RATE * DECAY * .001 ) );
/* Hang timer release decay time constant in seconds */
@SuppressWarnings( "unused" )
private static final double RELEASE_TIME_CONSTANT = 0.05;
/* Specifies the AGC Knee in dB if AGC is active (nominal range -160 to 0 dB) */
private static final double THRESHOLD = -100;
/* Limit output to about 3db of maximum */
private static final double AGC_OUT_SCALE = 0.7;
/* Keep max input and output the same */
private static final double MAX_AMPLITUDE = 32767.0; //1.0;
private static final double MAX_MANUAL_AMPLITUDE = 32767.0; //1.0;
/* Specifies AGC manual gain in dB if AGC is not active ( 0 to 100 dB) */
private static final double MANUAL_GAIN = 0.0;
private static final double MANUAL_AGC_GAIN = MAX_MANUAL_AMPLITUDE *
Math.pow( 10.0, MANUAL_GAIN / 20.0 );
/* Specifies dB reduction in output at knee from max output level (0 - 10dB) */
private static final double SLOPE_FACTOR = 2.0;
private static final double KNEE = THRESHOLD / 20.0;
private static final double GAIN_SLOPE = SLOPE_FACTOR / 100.0;
private static final double FIXED_GAIN = AGC_OUT_SCALE *
Math.pow( 10.0, KNEE * ( GAIN_SLOPE - 1.0 ) );
/* Constant for calc log() so that a value of 0 magnitude = -8 */
private static final double MIN_CONSTANT = 3.2767E-4;
private AtomicBoolean mAGCEnabled = new AtomicBoolean( true );
private double mPeakMagnitude = 0.0;
private double mAttackAverage = 0.0;
private double mDecayAverage = 0.0;
private FloatCircularBuffer mDelayBuffer =
new FloatCircularBuffer( (int)( SAMPLE_RATE * DELAY_TIME_CONSTANT ) );
private DoubleCircularBuffer mMagnitudeBuffer =
new DoubleCircularBuffer( (int)( SAMPLE_RATE * WINDOW_TIME_CONSTANT ) );
public AutomaticGainControl()
{
}
public float process( float currentSample )
{
float delayedSample = mDelayBuffer.get( currentSample );
double gain = MANUAL_AGC_GAIN;
if( mAGCEnabled.get() )
{
double currentMagnitude = Math.log10( Math.abs( currentSample ) +
MIN_CONSTANT ) - Math.log10( MAX_AMPLITUDE );
double delayedMagnitude = mMagnitudeBuffer.get( currentMagnitude );
if( currentMagnitude > mPeakMagnitude )
{
/* Use current magnitude as peak if it's larger */
mPeakMagnitude = currentMagnitude;
}
else if( delayedMagnitude == mPeakMagnitude )
{
/* If delayed magnitude is the current peak, then find a new peak */
mPeakMagnitude = mMagnitudeBuffer.max();
}
/* Exponential decay mode */
if( mPeakMagnitude > mAttackAverage )
{
mAttackAverage = ( ( 1.0 - ATTACK_RISE_ALPHA ) * mAttackAverage ) +
( ATTACK_RISE_ALPHA * mPeakMagnitude );
}
else
{
mAttackAverage = ( ( 1.0 - ATTACK_FALL_ALPHA ) * mAttackAverage ) +
( ATTACK_FALL_ALPHA * mPeakMagnitude );
}
if( mPeakMagnitude > mDecayAverage )
{
mDecayAverage = ( ( 1.0 - DECAY_RISE_ALPHA ) * mDecayAverage ) +
( DECAY_RISE_ALPHA * mPeakMagnitude );
}
else
{
mDecayAverage = ( ( 1.0 - DECAY_FALL_ALPHA ) * mDecayAverage ) +
( DECAY_RISE_ALPHA * mPeakMagnitude );
}
double magnitude = ( mAttackAverage > mDecayAverage ) ?
mAttackAverage : mDecayAverage;
if( magnitude < KNEE )
{
gain = FIXED_GAIN;
}
else
{
gain = AGC_OUT_SCALE * Math.pow( 10.0,
magnitude * ( GAIN_SLOPE - 1.0 ) );
}
}
return (float)( delayedSample * gain );
}
/**
* Enables or disables Automatic Gain Control (AGC).
*/
public void setAGCEnabled( boolean enabled )
{
mAGCEnabled.set( enabled );
}
/**
* Indicates if AGC is enabled
* @return true=AGC, false=MANUAL GAIN
*/
public boolean isAGCEnabled()
{
return mAGCEnabled.get();
}
}

View File

@@ -0,0 +1,53 @@
package io.github.dsheirer.dsp.gain;
import io.github.dsheirer.sample.Listener;
import io.github.dsheirer.sample.real.RealBuffer;
public class AutomaticGainControl_RB extends AutomaticGainControl
implements Listener<RealBuffer>
{
private Listener<RealBuffer> mListener;
public AutomaticGainControl_RB()
{
}
public void dispose()
{
mListener = null;
}
/**
* Applies AGC to a buffer of real (float) samples
*/
public RealBuffer process( RealBuffer buffer )
{
float[] samples = buffer.getSamples();
for( int x = 0; x < samples.length; x++ )
{
samples[ x ] = process( samples[ x ] );
}
return buffer;
}
@Override
public void receive( RealBuffer buffer )
{
if( mListener != null )
{
mListener.receive( process( buffer ) );
}
}
public void setListener( Listener<RealBuffer> listener )
{
mListener = listener;
}
public void removeListener()
{
mListener = null;
}
}

View File

@@ -0,0 +1,249 @@
package io.github.dsheirer.dsp.gain;
import io.github.dsheirer.buffer.ComplexCircularBuffer;
import io.github.dsheirer.buffer.DoubleCircularBuffer;
import io.github.dsheirer.sample.Listener;
import io.github.dsheirer.sample.Provider;
import io.github.dsheirer.sample.complex.Complex;
import java.util.concurrent.atomic.AtomicBoolean;
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014 Dennis Sheirer
* Copyright (C) 2011 Alex Csete
* Copyright (C) 2010 Moe Wheatley
*
* 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 <http://www.gnu.org/licenses/>
*
* ---------------------------------------------------
* Ported from Alex Csete's gqrx source code at:
* https://github.com/csete/gqrx/blob/74ba2742c89b0dfb66854fc9f72d19e5c6e355b6/dsp/agc_impl.cpp
*
* This Software is released under the "Simplified BSD License" + + +
* Copyright 2010 Moe Wheatley. All rights reserved.
*
* Redistribution and use in source and binary forms, with or without
* modification, are permitted provided that the following conditions are
* met:
*
* 1. Redistributions of source code must retain the above copyright notice,
* this list of conditions and the following disclaimer.
*
* 2. Redistributions in binary form must reproduce the above copyright
* notice, this list of conditions and the following disclaimer in the
* documentation and/or other materials provided with the distribution.
*
* THIS SOFTWARE IS PROVIDED BY Moe Wheatley ``AS IS'' AND ANY EXPRESS OR
* IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
* OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
* IN NO EVENT SHALL Moe Wheatley OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT,
* INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
* (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
* SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
* HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT,
* STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
* ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
* POSSIBILITY OF SUCH DAMAGE.
*
* The views and conclusions contained in the software and documentation are
* those of the authors and should not be interpreted as representing
* official policies, either expressed or implied, of Moe Wheatley.
******************************************************************************/
public class ComplexAutomaticGainControl implements Listener<Complex>,
Provider<Complex>
{
private static final double SAMPLE_RATE = 48000;
/* Signal delay line - time delay in seconds */
private static final double DELAY_TIME_CONSTANT = 0.015;
/* Peak detector window - time delay in seconds */
private static final double WINDOW_TIME_CONSTANT = 0.018;
/* Attack time constants in seconds */
private static final double ATTACK_RISE_TIME_CONSTANT = 0.002;
private static final double ATTACK_FALL_TIME_CONSTANT = 0.005;
private static final double ATTACK_RISE_ALPHA = 1.0 -
Math.exp( -1.0 / SAMPLE_RATE * ATTACK_RISE_TIME_CONSTANT );
private static final double ATTACK_FALL_ALPHA = 1.0 -
Math.exp( -1.0 / SAMPLE_RATE * ATTACK_FALL_TIME_CONSTANT );
/* AGC decay value in milliseconds (20 to 5000) */
private static final double DECAY = 200;
/* Ratio between rise and fall times of decay time constants - adjust for
* best action with SSB */
private static final double DECAY_RISEFALL_RATIO = 0.3;
private static final double DECAY_RISE_ALPHA = 1.0 -
Math.exp( -1.0 / ( SAMPLE_RATE * DECAY * .001 * DECAY_RISEFALL_RATIO ) );
private static final double DECAY_FALL_ALPHA = 1.0 -
Math.exp( -1.0 / ( SAMPLE_RATE * DECAY * .001 ) );
/* Hang timer release decay time constant in seconds */
@SuppressWarnings( "unused" )
private static final double RELEASE_TIME_CONSTANT = 0.05;
/* Specifies the AGC Knee in dB if AGC is active (nominal range -160 to 0 dB) */
private static final double THRESHOLD = -100;
/* Limit output to about 3db of maximum */
private static final double AGC_OUT_SCALE = 0.7;
/* Keep max input and output the same */
private static final double MAX_AMPLITUDE = 1.0;
private static final double MAX_MANUAL_AMPLITUDE = 1.0;
/* Specifies AGC manual gain in dB if AGC is not active ( 0 to 100 dB) */
private static final double MANUAL_GAIN = 0.0;
private static final double MANUAL_AGC_GAIN = MAX_MANUAL_AMPLITUDE *
Math.pow( 10.0, MANUAL_GAIN / 20.0 );
/* Specifies dB reduction in output at knee from max output level (0 - 10dB) */
private static final double SLOPE_FACTOR = 2.0;
private static final double KNEE = THRESHOLD / 20.0;
private static final double GAIN_SLOPE = SLOPE_FACTOR / 100.0;
private static final double FIXED_GAIN = AGC_OUT_SCALE *
Math.pow( 10.0, KNEE * ( GAIN_SLOPE - 1.0 ) );
/* Constant for calc log() so that a value of 0 magnitude = -8 */
private static final double MIN_CONSTANT = 3.2767E-4;
private AtomicBoolean mAGCEnabled = new AtomicBoolean( true );
private double mPeakMagnitude = 0.0;
private double mAttackAverage = 0.0;
private double mDecayAverage = 0.0;
private ComplexCircularBuffer mDelayBuffer =
new ComplexCircularBuffer( (int)( SAMPLE_RATE * DELAY_TIME_CONSTANT ) );
private DoubleCircularBuffer mMagnitudeBuffer =
new DoubleCircularBuffer( (int)( SAMPLE_RATE * WINDOW_TIME_CONSTANT ) );
private Listener<Complex> mListener;
public ComplexAutomaticGainControl()
{
}
@Override
public void receive( Complex currentSample )
{
Complex delayedSample = mDelayBuffer.get( currentSample );
double gain = MANUAL_AGC_GAIN;
if( mAGCEnabled.get() )
{
float max = currentSample.maximumAbsolute();
double currentMagnitude = Math.log10( max + MIN_CONSTANT ) -
Math.log10( MAX_AMPLITUDE );
double delayedMagnitude = mMagnitudeBuffer.get( currentMagnitude );
if( currentMagnitude > mPeakMagnitude )
{
/* Use current magnitude as peak if it's larger */
mPeakMagnitude = currentMagnitude;
}
else if( delayedMagnitude == mPeakMagnitude )
{
/* If delayed magnitude is the current peak, then find a new peak */
mPeakMagnitude = mMagnitudeBuffer.max();
}
/* Exponential decay mode */
if( mPeakMagnitude > mAttackAverage )
{
mAttackAverage = ( ( 1.0 - ATTACK_RISE_ALPHA ) * mAttackAverage ) +
( ATTACK_RISE_ALPHA * mPeakMagnitude );
}
else
{
mAttackAverage = ( ( 1.0 - ATTACK_FALL_ALPHA ) * mAttackAverage ) +
( ATTACK_FALL_ALPHA * mPeakMagnitude );
}
if( mPeakMagnitude > mDecayAverage )
{
mDecayAverage = ( ( 1.0 - DECAY_RISE_ALPHA ) * mDecayAverage ) +
( DECAY_RISE_ALPHA * mPeakMagnitude );
}
else
{
mDecayAverage = ( ( 1.0 - DECAY_FALL_ALPHA ) * mDecayAverage ) +
( DECAY_RISE_ALPHA * mPeakMagnitude );
}
double magnitude = ( mAttackAverage > mDecayAverage ) ?
mAttackAverage : mDecayAverage;
if( magnitude < KNEE )
{
gain = FIXED_GAIN;
}
else
{
gain = AGC_OUT_SCALE * Math.pow( 10.0,
magnitude * ( GAIN_SLOPE - 1.0 ) );
}
}
delayedSample.multiply( (float)gain );
if( mListener != null )
{
mListener.receive( delayedSample );
}
}
/**
* Enables or disables Automatic Gain Control (AGC).
*/
public void setAGCEnabled( boolean enabled )
{
mAGCEnabled.set( enabled );
}
/**
* Indicates if AGC is enabled
* @return true=AGC, false=MANUAL GAIN
*/
public boolean isAGCEnabled()
{
return mAGCEnabled.get();
}
@Override
public void setListener( Listener<Complex> listener )
{
mListener = listener;
}
@Override
public void removeListener( Listener<Complex> listener )
{
mListener = null;
}
}

View File

@@ -0,0 +1,112 @@
package io.github.dsheirer.dsp.gain;
import io.github.dsheirer.buffer.FloatCircularBuffer;
import io.github.dsheirer.sample.complex.Complex;
import io.github.dsheirer.sample.complex.ComplexSampleListener;
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2015 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 <http://www.gnu.org/licenses/>
*
* This class is modeled using concepts detailed in the GNURadio
* implementation of feed forward AGC class located at:
* https://github.com/gnuradio/gnuradio/blob/master/gr-analog/lib/
* feedforward_agc_cc_impl.cc
*
******************************************************************************/
public class ComplexFeedForwardGainControl implements ComplexSampleListener
{
public static final float OBJECTIVE_ENVELOPE = 1.0f;
public static final float MINIMUM_ENVELOPE = 0.0001f;
private ComplexSampleListener mListener;
private FloatCircularBuffer mEnvelopeHistory;
private float mMaxEnvelope = 0.0f;
private float mGain = 1.0f;
/**
* Dynamic gain control for incoming sample stream to amplify or attenuate
* all samples toward an objective unity)gain, using the maximum envelope
* value detected in the stream history window.
*
* Uses the specified damping factor to limit gain swings. Damping factor
* is applied against the delta between current gain value and a recalculated
* gain value to limit how quickly the gain value will increase or decrease.
*
* @param window - history size to use in detecting maximum envelope value
* @param damping - a value between 0 < damping <= 1.0;
*/
public ComplexFeedForwardGainControl( int window )
{
mEnvelopeHistory = new FloatCircularBuffer( window );
}
public void dispose()
{
mListener = null;
}
@Override
public void receive( float inphase, float quadrature )
{
float envelope = Complex.envelope( inphase, quadrature );
if( envelope > mMaxEnvelope )
{
mMaxEnvelope = envelope;
adjustGain();
}
/* Replace oldest envelope value with current envelope value */
float oldestEnvelope = mEnvelopeHistory.get( envelope );
/* If the oldest envelope value was the max envelope value, then we
* have to rediscover the max value from the envelope history */
if( mMaxEnvelope == oldestEnvelope && mMaxEnvelope != envelope )
{
mMaxEnvelope = MINIMUM_ENVELOPE;
for( float value: mEnvelopeHistory.getBuffer() )
{
if( value > mMaxEnvelope )
{
mMaxEnvelope = value;
}
}
adjustGain();
}
/* Apply current gain value to the sample and send to the listener */
if( mListener != null )
{
mListener.receive( inphase *= mGain, quadrature *= mGain );
}
}
private void adjustGain()
{
mGain = OBJECTIVE_ENVELOPE / mMaxEnvelope;
}
public void setListener( ComplexSampleListener listener )
{
mListener = listener;
}
}

View File

@@ -0,0 +1,31 @@
package io.github.dsheirer.dsp.gain;
import io.github.dsheirer.sample.Listener;
import io.github.dsheirer.sample.complex.Complex;
public class ComplexGain implements Listener<Complex>
{
private float mGain;
private Listener<Complex> mListener;
public ComplexGain( float gain )
{
mGain = gain;
}
@Override
public void receive( Complex sample )
{
if( mListener != null )
{
sample.multiply( mGain );
mListener.receive( sample );
}
}
public void setListener( Listener<Complex> listener )
{
mListener = listener;
}
}

View File

@@ -0,0 +1,109 @@
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp.gain;
import io.github.dsheirer.sample.real.RealSampleListener;
import org.apache.commons.lang3.Validate;
import org.slf4j.Logger;
import org.slf4j.LoggerFactory;
public class DirectGainControl implements GainController, RealSampleListener
{
private final static Logger mLog =
LoggerFactory.getLogger( DirectGainControl.class );
private float mDefaultGain;
private float mGain;
private float mMinimum;
private float mMaximum;
private float mIncrement;
private RealSampleListener mListener;
/**
* Gain Control.
*
* Note: does not perform null checking on output listener in order to
* streamline sample processing.
*
* @param gain - initial value
* @param minimum - minimum value greater than 0
* @param maximum - maximum permitted value
* @param increment - amount gain is adjusted when increase() or decrease()
* is invoked
*/
public DirectGainControl( float gain, float minimum, float maximum, float increment )
{
Validate.isTrue(minimum > 0);
mDefaultGain = gain;
mGain = gain;
mMinimum = minimum;
mMaximum = maximum;
mIncrement = increment;
}
public void dispose()
{
mListener = null;
}
public float correct( float sample )
{
return sample * mGain;
}
@Override
public void receive( float sample )
{
mListener.receive( correct( sample ) );
}
@Override
public void reset()
{
mGain = mDefaultGain;
}
@Override
public void increase()
{
mGain += mIncrement;
if( mGain > mMaximum )
{
mGain = mMaximum;
}
}
@Override
public void decrease()
{
mGain -= mIncrement;
if( mGain < mMinimum )
{
mGain = mMinimum;
}
}
public void setListener( RealSampleListener listener )
{
mListener = listener;
}
}

View File

@@ -0,0 +1,10 @@
package io.github.dsheirer.dsp.gain;
public interface GainController
{
public abstract void increase();
public abstract void decrease();
public abstract void reset();
}

View File

@@ -0,0 +1,65 @@
package io.github.dsheirer.dsp.gain;
import io.github.dsheirer.sample.real.RealBuffer;
public class NonClippingGain
{
private float mGain;
private float mMaxValue;
/**
* Applies a fixed gain value to the float samples.
*
* @param gain value
* @param maxValue max absolute value to prevent clipping ( 0.0 to 1.0 )
*/
public NonClippingGain( float gain, float maxValue )
{
mGain = gain;
mMaxValue = maxValue;
}
/**
* Applies gain to the sample
*/
public float apply( float sample )
{
float adjusted = sample * mGain;
if( adjusted > mMaxValue )
{
adjusted = mMaxValue;
}
if( adjusted < -mMaxValue )
{
adjusted = -mMaxValue;
}
return adjusted;
}
/**
* Applies gain to the sample array
*/
public float[] apply( float[] samples )
{
for( int x = 0; x < samples.length; x++ )
{
samples[ x ] = apply( samples[ x ]);
}
return samples;
}
/**
* Applies gain to the real buffer
*/
public RealBuffer apply( RealBuffer buffer )
{
float[] samples = buffer.getSamples();
apply( samples );
return buffer;
}
}

View File

@@ -0,0 +1,112 @@
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp.mixer;
import io.github.dsheirer.sample.complex.Complex;
import org.slf4j.Logger;
import org.slf4j.LoggerFactory;
public class Oscillator
{
private final static Logger mLog = LoggerFactory.getLogger( Oscillator.class );
private double mFrequency;
private double mSampleRate;
private Complex mAnglePerSample;
private Complex mCurrentAngle = new Complex( 0.0f, -1.0f );
/**
* Oscillator produces complex or float samples corresponding to a sine wave
* oscillating at the specified frequency and sample rate
*
* @param frequency - positive or negative frequency in hertz
* @param sampleRate - in hertz
*/
public Oscillator( long frequency, int sampleRate )
{
mSampleRate = (double)sampleRate;
mFrequency = (double)frequency;
update();
}
/**
* Updates the internal values after a frequency or sample rate change
*/
private void update()
{
float anglePerSample =
(float)( 2.0d * Math.PI * mFrequency / mSampleRate );
mAnglePerSample = Complex.fromAngle( anglePerSample );
}
/**
* Sets or changes the frequency of this oscillator
*/
public void setFrequency( long frequency )
{
mFrequency = (double)frequency;
update();
}
/**
* Sets or changes the sample rate of this oscillator
*/
public void setSampleRate( int sampleRate )
{
mSampleRate = (double)sampleRate;
update();
}
/**
* Steps the current angle by the angle per sample amount
*/
public void rotate()
{
mCurrentAngle.multiply( mAnglePerSample );
mCurrentAngle.fastNormalize();
}
public float inphase()
{
return mCurrentAngle.inphase();
}
public float quadrature()
{
return mCurrentAngle.quadrature();
}
/**
* Get next complex sample
*/
public Complex getComplex()
{
return mCurrentAngle.copy();
}
/**
* Get the next float sample
*/
public float getFloat()
{
return mCurrentAngle.real();
}
}

View File

@@ -0,0 +1,588 @@
package io.github.dsheirer.dsp.psk;
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014, 2015 Dennis Sheirer
*
* Gardner detector, Costas loop and Interpolator derived from:
*
* OP25 - gardner_costas_cc_impl.cc
* Copyright 2010, 2011, 2012, 2013 KA1RBI (gardner_costas_cc_impl.cc)
*
* OP25 - scope.py
* Copyright 2008-2011 Steve Glass
* Copyright 2011, 2012, 2013 KA1RBI
*
* GNURadio - control_loop and mpsk_receiver_cc
* Copyright 2011,2013 Free Software Foundation, Inc.
*
* 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 <http://www.gnu.org/licenses/>
******************************************************************************/
import io.github.dsheirer.buffer.FloatAveragingBuffer;
import io.github.dsheirer.instrument.Instrumentable;
import io.github.dsheirer.instrument.tap.Tap;
import io.github.dsheirer.instrument.tap.TapGroup;
import io.github.dsheirer.instrument.tap.stream.EyeDiagramData;
import io.github.dsheirer.instrument.tap.stream.EyeDiagramDataTap;
import io.github.dsheirer.sample.Listener;
import io.github.dsheirer.sample.complex.Complex;
import io.github.dsheirer.sample.complex.ComplexSampleListener;
import io.github.dsheirer.source.tuner.frequency.IFrequencyChangeListener;
import io.github.dsheirer.dsp.filter.interpolator.RealInterpolator;
import io.github.dsheirer.dsp.symbol.Dibit;
import java.util.ArrayList;
import java.util.Arrays;
import java.util.List;
/**
* Implements a LSM (Pi/4) demodulator using a Gardner Detector to determine
* optimal sampling timing and a Costas Loop as a phase locked loop synchronized
* with the incoming signal carrier frequency.
*
* Sample Rate: 48000
* Symbol Rate: 4800
*/
public class LSMDemodulator implements Instrumentable, ComplexSampleListener
{
private List<TapGroup> mAvailableTaps;
/* 45 degree rotation angle */
public static final float THETA = (float)( Math.PI / 4.0d );
/* 45 degree point */
public static final Complex POINT_45_DEGREES =
new Complex( (float)Math.sin( THETA ), (float)Math.cos( THETA ) );
private Listener<Complex> mSymbolListener;
private GardnerSymbolTiming mGardnerDetector = new GardnerSymbolTiming();
private CostasLoop mCostasLoop = new CostasLoop();
private IFrequencyChangeListener mFrequencyChangeListener;
private EyeDiagramDataTap mEyeDiagramDataTap;
public LSMDemodulator()
{
}
public void dispose()
{
mSymbolListener = null;
mFrequencyChangeListener = null;
}
public void setSymbolListener( Listener<Complex> listener )
{
mSymbolListener = listener;
}
public void removeSymbolListener()
{
mSymbolListener = null;
}
@Override
public void receive( float inphase, float quadrature )
{
mGardnerDetector.receive( inphase, quadrature );
}
/**
* Applies a phase correction value to the costas loop to correct when a
* phase lock error is detected in the binary output stream.
*
* @param correction - value in radians
*/
public void correctPhaseError( double correction )
{
mCostasLoop.correctPhaseError( correction );
}
/**
* Constrains value to the range of ( -maximum <> maximum )
*/
public static float clip( float value, float maximum )
{
if( value > maximum )
{
return maximum;
}
else if( value < -maximum )
{
return -maximum;
}
return value;
}
/**
* Constrains timing error to +/- the maximum value and corrects any
* floating point invalid numbers
*/
private float normalize( float error, float maximum )
{
if( Float.isNaN( error ) )
{
return 0.0f;
}
else
{
return clip( error, maximum );
}
}
public void addListener( IFrequencyChangeListener listener )
{
mFrequencyChangeListener = listener;
}
public void removeListener( IFrequencyChangeListener listener )
{
mFrequencyChangeListener = null;
}
/**
* Processes the current loop frequency of the costas loop and broadcasts
* tuner frequency offset adjustments 10 times a second as needed to keep
* the frequency offset below half of the maximum frequency
*/
public class FrequencyControl
{
/* Set the trigger threshold at half of the maximum costas loop control
* frequency */
private static final float CORRECTION_THRESHOLD_FREQUENCY =
( 2.0f * (float)Math.PI * 600.0f ) / 48000.0f;
private FloatAveragingBuffer mBuffer = new FloatAveragingBuffer( 40 );
private long mFrequencyError = 0;
private int mCounter = 0;
public void receive( float frequency )
{
// mCounter++;
//
// if( mCounter >= 480 )
// {
// float average = mBuffer.get( frequency );
// mLog.debug( "Avg:" + average + " Threshold:" + CORRECTION_THRESHOLD_FREQUENCY);
//
// boolean correctionNeeded = false;
//
// if( average > CORRECTION_THRESHOLD_FREQUENCY )
// {
// mFrequencyError--;
//
// correctionNeeded = true;
// }
// else if( average < -CORRECTION_THRESHOLD_FREQUENCY )
// {
// mFrequencyError++;
//
// correctionNeeded = true;
// }
//
// if( correctionNeeded && mFrequencyChangeListener != null )
// {
// mLog.debug( "Issuing Frequency Correction: " + mFrequencyError );
// mFrequencyChangeListener.frequencyChanged(
// new FrequencyChangeEvent( Attribute.FREQUENCY_ERROR,
// mFrequencyError ) );
// }
//
// mCounter = 0;
// }
// else
// {
// mBuffer.get( frequency );
// }
}
}
/**
* Gardner symbol timing detector
*/
public class GardnerSymbolTiming implements ComplexSampleListener
{
public static final int HALF_SAMPLES_PER_SYMBOL = 5;
public static final int SAMPLES_PER_SYMBOL = 10;
public static final int TWICE_SAMPLES_PER_SYMBOL = 20;
private float[] mDelayLineInphase = new float[ 2 * TWICE_SAMPLES_PER_SYMBOL ];
private float[] mDelayLineQuadrature = new float[ 2 * TWICE_SAMPLES_PER_SYMBOL ];
private int mDelayLinePointer = 0;
/* Sampling point */
private float mMu = 10.0f;
private float mGainMu = 0.05f;
/* Samples per symbol */
private float mOmega = 10.0f;
private float mGainOmega = 0.1f * mGainMu * mGainMu;
private float mOmegaRel = 0.005f;
private float mOmegaMid = 10.0f;
private RealInterpolator mInterpolator = new RealInterpolator( 1.0f );
private Complex mPreviousSample = new Complex( 0, 0 );
private Complex mPreviousMiddleSample = new Complex( 0, 0 );
private Complex mPreviousSymbol = new Complex( 0, 0 );
/**
* Provides symbol sampling timing control
*/
public GardnerSymbolTiming()
{
}
@Override
public void receive( float inphase, float quadrature )
{
/* Count down samples per symbol until we calculate the symbol */
mMu--;
mCostasLoop.increment();
/* Mix incoming sample with costas loop to remove any rotation
* that is present from a mis-tuned carrier frequency */
Complex costas = mCostasLoop.getCurrentVector();
float derotatedInphase = Complex.multiplyInphase( inphase, quadrature, costas.inphase(), costas.quadrature() );
float derotatedQuadrature = Complex.multiplyQuadrature( inphase, quadrature, costas.inphase(), costas.quadrature() );
/* Fill up the delay line to use with the interpolator */
mDelayLineInphase[ mDelayLinePointer ] = derotatedInphase;
mDelayLineInphase[ mDelayLinePointer + TWICE_SAMPLES_PER_SYMBOL ] = derotatedInphase;
mDelayLineQuadrature[ mDelayLinePointer ] = derotatedQuadrature;
mDelayLineQuadrature[ mDelayLinePointer + TWICE_SAMPLES_PER_SYMBOL ] = derotatedQuadrature;
/* Increment pointer and keep pointer in bounds */
mDelayLinePointer = ( mDelayLinePointer + 1 ) % TWICE_SAMPLES_PER_SYMBOL;
/* Calculate the symbol once we've stored enough samples */
if( mMu <= 1.0f )
{
float half_omega = mOmega / 2.0f;
int half_sps = (int)Math.floor( half_omega );
float half_mu = mMu + half_omega - (float)half_sps;
if( half_mu > 1.0 )
{
half_mu -= 1.0;
half_sps += 1;
}
/* Calculate interpolated middle sample and current sample */
float middleSampleInphase = mInterpolator.filter(
mDelayLineInphase, mDelayLinePointer, mMu );
float middleSampleQuadrature = mInterpolator.filter(
mDelayLineQuadrature, mDelayLinePointer, mMu );
int index = mDelayLinePointer + half_sps;
float currentSampleInphase = mInterpolator.filter(
mDelayLineInphase, index, half_mu );
float currentSampleQuadrature = mInterpolator.filter(
mDelayLineQuadrature, index, half_mu );
/* Multiply current sample and conjugate (negated quadrature) of
* previous sample to get symbols to use for gardner error feedback */
Complex middleSymbol = Complex.multiply( middleSampleInphase,
middleSampleQuadrature, mPreviousMiddleSample.conjugate() );
Complex currentSymbol = Complex.multiply( currentSampleInphase,
currentSampleQuadrature, mPreviousSample.conjugate() );
/* Set gain to unity */
middleSymbol.normalize();
currentSymbol.normalize();
/* Gardner timing error calculations */
float errorInphase = ( mPreviousSymbol.inphase() -
currentSymbol.inphase() ) * middleSymbol.inphase();
float errorQuadrature = ( mPreviousSymbol.quadrature() -
currentSymbol.quadrature() ) * middleSymbol.quadrature();
float gardnerError = normalize( errorInphase + errorQuadrature, 1.0f );
if( mEyeDiagramDataTap != null )
{
mEyeDiagramDataTap.receive(
new EyeDiagramData(
Arrays.copyOfRange( mDelayLineInphase, 0, 20 ),
Arrays.copyOfRange( mDelayLineQuadrature, 0, 20 ),
mMu, (float)half_sps + half_mu, gardnerError ) );
}
/* mOmega is samples per symbol and is constrained to floating
* between +/- .005 of the nominal 10.0 samples per symbol */
mOmega = mOmega + mGainOmega * gardnerError;
mOmega = mOmegaMid + clip( mOmega - mOmegaMid, mOmegaRel );
/* Adjust sample timing based on error of current sample */
mMu += mOmega + ( mGainMu * gardnerError );
/* Store current samples/symbols to use for the next period */
mPreviousSample.setValues( currentSampleInphase, currentSampleQuadrature );
mPreviousMiddleSample.setValues( middleSampleInphase, middleSampleQuadrature );
mPreviousSymbol = currentSymbol;
/* Update costas loop using phase error present in current
* symbol. The symbol is rotated from star orientation to polar
* orientation to simplify error calculation */
mCostasLoop.receive( currentSymbol );
// mCostasLoop.receive(
// ComplexSample.multiply( currentSymbol, POINT_45_DEGREES ) );
/* Dispatch the differentiated symbol to the registered listener */
if( mSymbolListener != null )
{
mSymbolListener.receive( currentSymbol );
}
}
}
}
/**
* Costas Loop - phase locked loop synchronized to the frequency offset of
* the incoming signal to enable the signal to be mixed down to zero offset
* for proper synchronization. The mFrequency value indicates the detected
* frequency offset. We attempt to keep that value close to zero by issuing
* frequency adjustments to the tuner channel source.
*
* Most of the costas loop code was ported from gnuradio/control_loop and
* gnuradio/mpsk_receiver_cc using initialization values from KA1RBI's
* OP25/cqpsk.py
*/
public class CostasLoop implements Listener<Complex>
{
public static final double TWO_PI = 2.0 * Math.PI;
private static final float MAXIMUM_FREQUENCY =
( 2400.0f * (float)TWO_PI ) / 48000.0f;
/* http://www.trondeau.com/blog/2011/8/13/control-loop-gain-values.html */
private float mDamping = (float)Math.sqrt( 2.0 ) / 2.0f;
/* Use denominator between 100 (most) and 400 (least) to adjust costas
* loop control level */
private float mLoopBandwidth = (float)( TWO_PI / 400.0d );
private float mAlphaGain = ( 4.0f * mDamping * mLoopBandwidth ) /
( 1.0f + ( 2.0f * mDamping * mLoopBandwidth ) +
( mLoopBandwidth * mLoopBandwidth ) );
private float mBetaGain = ( 4.0f * mLoopBandwidth * mLoopBandwidth ) /
( 1.0f + ( 2.0f * mDamping * mLoopBandwidth ) +
( mLoopBandwidth * mLoopBandwidth ) );
private float mLoopPhase = 0.0f;
private float mLoopFrequency = 0.0f;
private FrequencyControl mFrequencyControl = new FrequencyControl();
public CostasLoop()
{
}
/**
* Applies a phase correction value to the current phasor. Use this
* method to apply correction when a +/- 90, or 180 degree phase error
* is detected in the binary output stream.
*
* If the supplied correction value places the loop frequency outside
* of the max frequency, then the frequency will be corrected 360
* degrees in the opposite direction to maintain within the max
* frequency bounds.
*
* @param correction - correction value in radians
*/
public void correctPhaseError( double correction )
{
mLoopFrequency += correction;
if( mLoopFrequency > MAXIMUM_FREQUENCY )
{
mLoopFrequency -= 2.0d * MAXIMUM_FREQUENCY;
}
if( mLoopFrequency < -MAXIMUM_FREQUENCY )
{
mLoopFrequency += 2.0d * MAXIMUM_FREQUENCY;
}
}
/**
* Increments the phase of the loop
*/
public void increment()
{
mLoopPhase += mLoopFrequency;
/* Keep the loop phase in bounds */
unwrapPhase();
}
private void unwrapPhase()
{
while( mLoopPhase > TWO_PI )
{
mLoopPhase -= TWO_PI;
}
while( mLoopPhase < -TWO_PI )
{
mLoopPhase += TWO_PI;
}
}
/**
* Current vector of the loop rotated by THETA degrees.
*/
public Complex getCurrentVector()
{
// return ComplexSample.fromAngle( mLoopPhase + THETA );
return Complex.fromAngle( mLoopPhase );
}
@Override
public void receive( Complex complex )
{
/* Calculate phase error */
float phaseError = getPhaseError( complex );
Dibit d = QPSKPolarSlicer.decide( complex );
// mLog.debug( "Dibit: " + d.name() + " Error: " + phaseError );
/* Adjust for phase error */
adjust( phaseError );
}
/**
* Updates the costas loop frequency and phase to adjust for the phase
* error value
*
* @param phase_error - (-)= late and (+)= early
*/
private void adjust( float phase_error )
{
mLoopFrequency += mBetaGain * phase_error;
mLoopPhase += mLoopFrequency + mAlphaGain * phase_error;
/* Maintain phase between +/- 2 * PI */
unwrapPhase();
/* Limit frequency to +/- maximum loop frequency */
limitFrequency();
}
/**
* Constrains the frequency within the bounds of +/- loop frequency
*/
private void limitFrequency()
{
/* Check for and issue tuner offset correction */
mFrequencyControl.receive( mLoopFrequency );
if( mLoopFrequency > MAXIMUM_FREQUENCY )
{
mLoopFrequency = MAXIMUM_FREQUENCY;
}
if( mLoopFrequency < -MAXIMUM_FREQUENCY )
{
mLoopFrequency = -MAXIMUM_FREQUENCY;
}
}
/**
* Detects rotational error present when the costas loop is not aligned
* with the carrier frequency. Provides error feedback to adjust
* mixer frequency.
*/
public float getPhaseError( Complex complex )
{
float phase_error = 0;
if( Math.abs( complex.inphase() ) > Math.abs( complex.quadrature() ) )
{
if( complex.inphase() > 0 )
{
phase_error = -complex.quadrature();
}
else
{
phase_error = complex.quadrature();
}
}
else
{
if( complex.quadrature() > 0 )
{
phase_error = complex.inphase();
}
else
{
phase_error = -complex.inphase();
}
}
return phase_error;
}
}
@Override
public List<TapGroup> getTapGroups()
{
if( mAvailableTaps == null )
{
mAvailableTaps = new ArrayList<>();
TapGroup group = new TapGroup( "LSM Demodulator" );
group.add( new EyeDiagramDataTap( "Eye Diagram", 0, 4800 ) );
mAvailableTaps.add( group );
}
return mAvailableTaps;
}
@Override
public void registerTap( Tap tap )
{
if( tap instanceof EyeDiagramDataTap )
{
mEyeDiagramDataTap = (EyeDiagramDataTap)tap;
}
}
@Override
public void unregisterTap( Tap tap )
{
if( tap instanceof EyeDiagramDataTap )
{
mEyeDiagramDataTap = null;
}
}
}

View File

@@ -0,0 +1,92 @@
package io.github.dsheirer.dsp.psk;
import io.github.dsheirer.sample.Broadcaster;
import io.github.dsheirer.sample.Listener;
import io.github.dsheirer.sample.complex.Complex;
import io.github.dsheirer.dsp.symbol.Dibit;
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014,2015 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 <http://www.gnu.org/licenses/>
******************************************************************************/
public class QPSKPolarSlicer implements Listener<Complex>
{
private Broadcaster<Dibit> mBroadcaster = new Broadcaster<Dibit>();
/**
* Slices a ComplexSample representing a phase shifted symbol according to
* the following constellation pattern:
*
* \ 00 /
* \ /
* 01 \/ 10
* /\
* / \
* / 11 \
*/
public QPSKPolarSlicer()
{
}
public void dispose()
{
mBroadcaster.dispose();
mBroadcaster = null;
}
public void addListener( Listener<Dibit> listener )
{
mBroadcaster.addListener( listener );
}
public void removeListener( Listener<Dibit> listener )
{
mBroadcaster.removeListener( listener );
}
@Override
public void receive( Complex complex )
{
mBroadcaster.broadcast( decide( complex ) );
}
public static Dibit decide( Complex complex )
{
if( Math.abs( complex.inphase() ) > Math.abs( complex.quadrature() ) )
{
if( complex.inphase() > 0 )
{
return Dibit.D10_MINUS_1;
}
else
{
return Dibit.D01_PLUS_3;
}
}
else
{
if( complex.quadrature() > 0 )
{
return Dibit.D00_PLUS_1;
}
else
{
return Dibit.D11_MINUS_3;
}
}
}
}

View File

@@ -0,0 +1,58 @@
package io.github.dsheirer.dsp.psk;
import io.github.dsheirer.sample.Broadcaster;
import io.github.dsheirer.sample.Listener;
import io.github.dsheirer.sample.complex.Complex;
import io.github.dsheirer.dsp.symbol.Dibit;
public class QPSKStarSlicer implements Listener<Complex>
{
private Broadcaster<Dibit> mBroadcaster = new Broadcaster<Dibit>();
/**
* Slices a ComplexSample representing a phase shifted symbol according to
* the following constellation pattern:
*
* 00 | 01
* ---|---
* 10 | 11
*/
public QPSKStarSlicer()
{
}
public void dispose()
{
mBroadcaster.dispose();
mBroadcaster = null;
}
public void addListener( Listener<Dibit> listener )
{
mBroadcaster.addListener( listener );
}
public void removeListener( Listener<Dibit> listener )
{
mBroadcaster.removeListener( listener );
}
@Override
public void receive( Complex sample )
{
mBroadcaster.broadcast( decide( sample ) );
}
public static Dibit decide( Complex sample )
{
if( sample.inphase() > 0 )
{
return sample.quadrature() > 0 ? Dibit.D00_PLUS_1 : Dibit.D10_MINUS_1;
}
else
{
return sample.quadrature() > 0 ? Dibit.D01_PLUS_3 : Dibit.D11_MINUS_3;
}
}
}

View File

@@ -0,0 +1,112 @@
package io.github.dsheirer.dsp.symbol;
public enum Dibit
{
D01_PLUS_3( false, true, 1, 4 ),
D00_PLUS_1( false, false, 0, 0 ),
D10_MINUS_1( true, false, 2, 8 ),
D11_MINUS_3( true, true, 3, 12 );
private boolean mBit1;
private boolean mBit2;
private int mLowValue;
private int mHighValue;
private Dibit( boolean bit1, boolean bit2, int lowValue, int highValue )
{
mBit1 = bit1;
mBit2 = bit2;
mLowValue = lowValue;
mHighValue = highValue;
}
/**
* Allowable transition state 1 for trellis coded modulation. Indicates
* the dibit that differs from this dibit in the LSB bit position.
*/
public Dibit getAllowableTransition1()
{
switch( this )
{
case D00_PLUS_1:
return D01_PLUS_3;
case D01_PLUS_3:
return D00_PLUS_1;
case D10_MINUS_1:
return D11_MINUS_3;
case D11_MINUS_3:
return D10_MINUS_1;
default:
/* We'll never get here */
return D00_PLUS_1;
}
}
/**
* Allowable transition state 2 for trellis coded modulation. Indicates
* the dibit that differs from this dibit in the MSB bit position.
*/
public Dibit getAllowableTransition2()
{
switch( this )
{
case D00_PLUS_1:
return D10_MINUS_1;
case D01_PLUS_3:
return D11_MINUS_3;
case D10_MINUS_1:
return D00_PLUS_1;
case D11_MINUS_3:
return D01_PLUS_3;
default:
/* We'll never get here */
return D00_PLUS_1;
}
}
/**
* Indicates if the trellis coded modulation dibit is allowed to transition
* to the dibit argument.
*/
public boolean isAllowableTransition( Dibit dibit )
{
return dibit == getAllowableTransition1() ||
dibit == getAllowableTransition2();
}
public boolean getBit1()
{
return mBit1;
}
public boolean getBit2()
{
return mBit2;
}
public int getLowValue()
{
return mLowValue;
}
public int getHighValue()
{
return mHighValue;
}
public static Dibit inverted( Dibit symbol )
{
switch( symbol )
{
case D10_MINUS_1:
return D11_MINUS_3;
case D11_MINUS_3:
return D01_PLUS_3;
case D00_PLUS_1:
return D10_MINUS_1;
case D01_PLUS_3:
default:
return D00_PLUS_1;
}
}
}

View File

@@ -0,0 +1,21 @@
package io.github.dsheirer.dsp.symbol;
public enum FrameSync
{
P25_PHASE1_ERROR_90_CCW( 0xFFEFAFAAEEAAl ),
P25_PHASE1_NORMAL( 0x5575F5FF77FFl ), // +33333 -3 +33 -33 +33 -3333 +3 -3 +3 -33333
P25_PHASE1_ERROR_90_CW( 0x001050551155l ),
P25_PHASE1_ERROR_180( 0xAA8A0A008800l );
private long mSync;
private FrameSync( long sync )
{
mSync = sync;
}
public long getSync()
{
return mSync;
}
}

View File

@@ -0,0 +1,149 @@
package io.github.dsheirer.dsp.symbol;
import io.github.dsheirer.instrument.tap.stream.SymbolEventTap;
import io.github.dsheirer.sample.Listener;
import java.util.BitSet;
/**
* Symbol slicer with auto-aligning baud timing
*/
public class Slicer implements Listener<Boolean>
{
public enum Output{ NORMAL, INVERTED };
private BitSet mBitSet = new BitSet();
private int mSymbolLength;
private int mDecisionThreshold;
private int mSampleCounter;
private SymbolEventTap mSymbolEventTap;
private boolean mNormalOutput = true;
private Listener<Boolean> mListener;
public Slicer( Output output, int samplesPerSymbol )
{
mNormalOutput = ( output == Output.NORMAL );
mSymbolLength = samplesPerSymbol;
/* Round up */
mDecisionThreshold = (int)( mSymbolLength / 2 ) + ( mSymbolLength % 2 );
}
public void receive( Boolean sample )
{
if( mSampleCounter >= 0 )
{
if( sample )
{
mBitSet.set( mSampleCounter );
}
else
{
mBitSet.clear( mSampleCounter );
}
}
mSampleCounter++;
if( mSampleCounter >= mSymbolLength )
{
boolean decision = mBitSet.cardinality() >= mDecisionThreshold;
send( decision );
/* Shift timing left if the left bit in the bitset is opposite
* the decision and the right bit is the same */
if( ( mBitSet.get( 0 ) ^ decision ) &&
( !( mBitSet.get( mSymbolLength - 1 ) ^ decision ) ) )
{
sendTapEvent( mBitSet, SymbolEvent.Shift.LEFT, decision );
reset();
mSampleCounter--;
}
/* Shift timing right if the left bit is the same as the
* decision and the right bit is opposite */
else if( ( !( mBitSet.get( 0 ) ^ decision ) ) &&
( mBitSet.get( mSymbolLength - 1 ) ^ decision ) )
{
sendTapEvent( mBitSet, SymbolEvent.Shift.RIGHT, decision );
/* Last bit from previous symbol to pre-fill next symbol */
boolean previousSoftBit = mBitSet.get( mSymbolLength - 1 );
reset();
if( previousSoftBit )
{
mBitSet.set( 0 );
}
mSampleCounter++;
}
/* No shift */
else
{
sendTapEvent( mBitSet, SymbolEvent.Shift.NONE, decision );
reset();
}
}
}
/**
* Sends the bit decision to the listener
*/
private void send( boolean decision )
{
if( mListener != null )
{
mListener.receive( mNormalOutput ? decision : !decision );
}
}
private void reset()
{
mBitSet.clear();
mSampleCounter = 0;
}
public void setListener( Listener<Boolean> listener)
{
mListener = listener;
}
public void removeListener( Listener<Boolean> listener )
{
mListener = null;
}
/**
* Sends instrumentation tap event to all registered listeners
*/
private void sendTapEvent(BitSet bitset, SymbolEvent.Shift shift, boolean decision )
{
if( mSymbolEventTap != null )
{
SymbolEvent event =
new SymbolEvent( bitset.get( 0, mSymbolLength ),
mSymbolLength,
decision,
shift );
mSymbolEventTap.receive( event );
}
}
public void addTap( SymbolEventTap tap )
{
mSymbolEventTap = tap;
}
public void removeTap()
{
mSymbolEventTap = null;
}
}

View File

@@ -0,0 +1,77 @@
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp.symbol;
import java.util.BitSet;
public class SymbolEvent
{
private BitSet mBitset;
private int mSamplesPerSymbol;
private boolean mDecision;
private Shift mShift;
public SymbolEvent( BitSet bitset, int samplesPerSymbol, boolean decision, Shift shift )
{
mBitset = bitset;
mSamplesPerSymbol = samplesPerSymbol;
mDecision = decision;
mShift = shift;
}
public BitSet getBitSet()
{
return mBitset;
}
public int getSamplesPerSymbol()
{
return mSamplesPerSymbol;
}
public boolean getDecision()
{
return mDecision;
}
public Shift getShift()
{
return mShift;
}
public enum Shift
{
LEFT( "<" ),
AGGRESSIVE_LEFT( "<<" ),
RIGHT( ">" ),
AGGRESSIVE_RIGHT( ">>" ),
NONE( "=" );
private String mLabel;
private Shift( String label )
{
mLabel = label;
}
public String getLabel()
{
return mLabel;
}
}
}

View File

@@ -0,0 +1,27 @@
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp.symbol;
/**
* Interface to allow an external process indicate that it has detected a
* sync pattern in the bit stream
*/
public interface SyncDetectListener
{
public void syncDetected();
}

View File

@@ -0,0 +1,24 @@
/*******************************************************************************
* SDR Trunk
* Copyright (C) 2014 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 <http://www.gnu.org/licenses/>
******************************************************************************/
package io.github.dsheirer.dsp.symbol;
public interface SyncDetectProvider
{
public void setSyncDetectListener( SyncDetectListener listener );
}