330 WIP on redesigning the P25 C4FM decoder. Creates a new decoder

based on the P25 LSM decoder using a costas loop and a gardner detector.
Updates the Demodulator Viewer FX application to add in a symbol chart
and a phase error chart.  The decoder won't lock on yet ... still a wip.

src/main/java/io/github/dsheirer/dsp/psk/DQPSKSymbolPhaseErrorCalculator.java

@@ -0,0 +1,30 @@
+/*******************************************************************************
+ * sdr-trunk
+ * Copyright (C) 2014-2018 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.psk;
+
+import io.github.dsheirer.sample.complex.Complex;
+
+public class DQPSKSymbolPhaseErrorCalculator extends QPSKSymbolPhaseErrorCalculator
+{
+    //Negative 45 degrees rotation to offset for differential additive 45 degrees of rotation
+    public static final Complex DIFFERENTIAL_OFFSET = Complex.fromAngle(Math.PI / -4.0d);
+
+    @Override
+    public void adjust(Complex symbol)
+    {
+        symbol.multiply(DIFFERENTIAL_OFFSET);
+    }
+}

src/main/java/io/github/dsheirer/dsp/psk/GardnerDetector.java

@@ -0,0 +1,82 @@
+/*******************************************************************************
+ * sdr-trunk
+ * Copyright (C) 2014-2018 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.psk;
+
+import io.github.dsheirer.sample.complex.Complex;
+
+public class GardnerDetector
+{
+    private float mErrorInphase;
+    private float mErrorQuadrature;
+    private float mError;
+
+    /**
+     * Gardner detector for recovering symbol timing error.
+     */
+    public GardnerDetector()
+    {
+    }
+
+
+    /**
+     * Calculates the symbol timing error by evaluating the previous, middle and current symbol instants.
+     *
+     * @param previous symbol
+     * @param middle symbol (ie midway between the previous and current)
+     * @param current symbol
+     * @return symbol timing error
+     */
+    public float getError(Complex previous, Complex middle, Complex current)
+    {
+        mErrorInphase = (previous.inphase() - current.inphase()) * middle.inphase();
+        mErrorQuadrature = (previous.quadrature() - current.quadrature()) * middle.quadrature();
+        mError = normalize(mErrorInphase + mErrorQuadrature, 1.0f);
+
+        return mError;
+    }
+
+    /**
+     * Constrains timing error to +/- the maximum value and corrects any floating point invalid numbers
+     */
+    public static float normalize(float error, float maximum)
+    {
+        if(Float.isNaN(error))
+        {
+            return 0.0f;
+        }
+        else
+        {
+            return clip(error, maximum);
+        }
+    }
+
+    /**
+     * 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;
+    }
+}

src/main/java/io/github/dsheirer/dsp/psk/IQPSKSymbolDecoder.java

@@ -0,0 +1,29 @@
+/*******************************************************************************
+ * sdr-trunk
+ * Copyright (C) 2014-2018 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.psk;
+
+import io.github.dsheirer.dsp.symbol.Dibit;
+import io.github.dsheirer.sample.complex.Complex;
+
+public interface IQPSKSymbolDecoder
+{
+    /**
+     * Decodes the symbol and converts it into a Dibit
+     * @param symbol to decode
+     * @return decoded dibit
+     */
+    Dibit decode(Complex symbol);
+}

src/main/java/io/github/dsheirer/dsp/psk/ISymbolPhaseErrorCalculator.java

@@ -0,0 +1,35 @@
+/*******************************************************************************
+ * sdr-trunk
+ * Copyright (C) 2014-2018 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.psk;
+
+import io.github.dsheirer.sample.complex.Complex;
+
+public interface ISymbolPhaseErrorCalculator
+{
+    /**
+     * Adjusts the calculated symbol as necessary before calculating any phase error.  This might be useful for
+     * differential encoding to remove a 45 degree rotational spin.
+     * @param symbol to adjust
+     */
+    void adjust(Complex symbol);
+
+    /**
+     * Calculates the phase error of the symbol relative to an alignment for the constellation
+     * @param symbol to calculate error
+     * @return error value
+     */
+    float getPhaseError(Complex symbol);
+}

src/main/java/io/github/dsheirer/dsp/psk/InterpolatingSymbolBuffer.java

@@ -0,0 +1,187 @@
+/*******************************************************************************
+ * sdr-trunk
+ * Copyright (C) 2014-2018 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.psk;
+
+import io.github.dsheirer.dsp.filter.interpolator.RealInterpolator;
+import io.github.dsheirer.sample.complex.Complex;
+
+public class InterpolatingSymbolBuffer
+{
+    private Complex mMiddleSample = new Complex(0,0);
+    private Complex mCurrentSample = new Complex(0,0);
+
+    private float[] mDelayLineInphase;
+    private float[] mDelayLineQuadrature;
+    private float mMu; //Sample Counter
+    private float mGainMu = 0.05f; //Sample Counter Gain
+    private float mOmega;
+    private float mGainOmega = 0.1f * mGainMu * mGainMu;
+    private float mOmegaRel = 0.005f;
+    private float mOmegaMid;
+
+    private int mDelayLinePointer = 0;
+    private int mTwiceSamplesPerSymbol;
+    private RealInterpolator mInterpolator = new RealInterpolator(1.0f);
+
+    public InterpolatingSymbolBuffer(float samplesPerSymbol)
+    {
+        mMu = samplesPerSymbol;
+        mTwiceSamplesPerSymbol = (int)Math.floor(2.0 * samplesPerSymbol);
+        mDelayLineInphase = new float[2 * mTwiceSamplesPerSymbol];
+        mDelayLineQuadrature = new float[2 * mTwiceSamplesPerSymbol];
+        mOmega = samplesPerSymbol;
+        mOmegaMid = samplesPerSymbol;
+
+    }
+
+    public void receive(Complex currentSample)
+    {
+        mMu--;
+
+        //Fill up the delay line to use with the interpolator
+        mDelayLineInphase[mDelayLinePointer] = currentSample.inphase();
+        mDelayLineInphase[mDelayLinePointer + mTwiceSamplesPerSymbol] = currentSample.inphase();
+        mDelayLineQuadrature[mDelayLinePointer] = currentSample.quadrature();
+        mDelayLineQuadrature[mDelayLinePointer + mTwiceSamplesPerSymbol] = currentSample.quadrature();
+
+        //Increment pointer and keep pointer in bounds
+        mDelayLinePointer = (mDelayLinePointer + 1) % mTwiceSamplesPerSymbol;
+    }
+
+    /**
+     * Updates the internal sample counter with additional samples to add to the symbol sample counter.
+     * @param samplesToAdd to the sample counter
+     */
+    public void increaseSampleCounter(float samplesToAdd)
+    {
+        mMu += samplesToAdd;
+    }
+
+    /**
+     * Adjusts samples per symbol and symbol timing counters and increments the sample counter to collect another
+     * symbol.
+     *
+     * @param symbolTimingError from a symbol timing error detector
+     */
+    public void resetAndAdjust(float symbolTimingError)
+    {
+        //mOmega is samples per symbol and is constrained to floating between +/- .005 of the nominal samples per
+        //symbol
+        mOmega = mOmega + mGainOmega * symbolTimingError;
+        mOmega = mOmegaMid + clip(mOmega - mOmegaMid, mOmegaRel);
+
+        //Add another symbol's worth of samples to the counter and adjust timing based on gardner error
+        increaseSampleCounter((mOmega + (mGainMu * symbolTimingError)));
+    }
+
+    /**
+     * Indicates if this buffer has accumulated a enough samples to represent a full symbol
+     */
+    public boolean hasSymbol()
+    {
+        return mMu < 1.0f;
+    }
+
+    /**
+     * Interpolated middle sample for the symbol.  Note: this method should only be invoked after testing for a
+     * complete symbol with the hasSymbol() method.
+     */
+    public Complex getMiddleSample()
+    {
+        /* Calculate interpolated middle sample and current sample */
+        mMiddleSample.setValues(getInphase(mMu), getQuadrature(mMu));
+
+        return mMiddleSample;
+    }
+
+    /**
+     * Interpolated current sample for the symbol.  This sample represents the ideal symbol sampling point.
+     * Note: this method should only be invoked after testing for a complete symbol with the hasSymbol() method.
+     */
+    public Complex getCurrentSample()
+    {
+        float half_omega = mOmega / 2.0f;
+        int half_sps = (int)Math.floor(half_omega);
+        float half_mu = mMu + half_omega - half_sps;
+
+        if(half_mu > 1.0)
+        {
+            half_mu -= 1.0;
+            half_sps += 1;
+        }
+
+        mCurrentSample.setValues(getInphase(half_mu, half_sps), getQuadrature(half_mu, half_sps));
+
+        return mCurrentSample;
+    }
+
+    /**
+     * Returns the interpolated inphase value for the specified offset
+     * @param interpolation into the buffer to calculate the interpolated sample
+     * @return inphase value of the interpolated sample
+     */
+    public float getInphase(float interpolation)
+    {
+        return mInterpolator.filter(mDelayLineInphase, mDelayLinePointer, interpolation);
+    }
+
+    /**
+     * Returns the interpolated inphase value for the specified offset
+     * @param interpolation into the buffer to calculate the interpolated sample
+     * @return inphase value of the interpolated sample
+     */
+    public float getInphase(float interpolation, int indexOffet)
+    {
+        return mInterpolator.filter(mDelayLineInphase, mDelayLinePointer + indexOffet, interpolation);
+    }
+
+    /**
+     * Returns the interpolated quadrature value for the specified offset
+     * @param interpolation into the buffer to calculate the interpolated sample
+     * @return quadrature value of the interpolated sample
+     */
+    public float getQuadrature(float interpolation)
+    {
+        return mInterpolator.filter(mDelayLineQuadrature, mDelayLinePointer, interpolation);
+    }
+
+    /**
+     * Returns the interpolated quadrature value for the specified offset
+     * @param interpolation into the buffer to calculate the interpolated sample
+     * @return quadrature value of the interpolated sample
+     */
+    public float getQuadrature(float interpolation, int indexOffset)
+    {
+        return mInterpolator.filter(mDelayLineQuadrature, mDelayLinePointer + indexOffset, interpolation);
+    }
+
+    /**
+     * 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;
+    }
+}

src/main/java/io/github/dsheirer/dsp/psk/LSMConstellationDecoder.java

@@ -0,0 +1,28 @@
+/*******************************************************************************
+ * sdr-trunk
+ * Copyright (C) 2014-2018 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.psk;
+
+import io.github.dsheirer.dsp.symbol.Dibit;
+import io.github.dsheirer.sample.complex.Complex;
+
+public class LSMConstellationDecoder implements IQPSKSymbolDecoder
+{
+    @Override
+    public Dibit decode(Complex symbol)
+    {
+        return null;
+    }
+}

src/main/java/io/github/dsheirer/dsp/psk/QPSKDemodulator.java

@@ -0,0 +1,182 @@
+/*******************************************************************************
+ * sdr-trunk
+ * Copyright (C) 2014-2018 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.psk;
+
+import io.github.dsheirer.dsp.psk.pll.IPhaseLockedLoop;
+import io.github.dsheirer.dsp.symbol.Dibit;
+import io.github.dsheirer.sample.Listener;
+import io.github.dsheirer.sample.complex.Complex;
+import io.github.dsheirer.sample.complex.ComplexBuffer;
+import io.github.dsheirer.sample.complex.ComplexSampleListener;
+
+public class QPSKDemodulator implements ComplexSampleListener
+{
+    private InterpolatingSymbolBuffer mInterpolatingSymbolBuffer;
+    private IPhaseLockedLoop mPLL;
+    private ISymbolPhaseErrorCalculator mSymbolPhaseErrorCalculator;
+    private IQPSKSymbolDecoder mQPSKSymbolDecoder;
+    private GardnerDetector mGardnerDetector = new GardnerDetector();
+    private Listener<Complex> mSymbolListener;
+    private Listener<Double> mPLLErrorListener;
+    private Listener<Dibit> mDibitListener;
+
+    private Complex mPreviousSample = new Complex(0, 0);
+    private Complex mPreviousMiddleSample = new Complex(0, 0);
+    private Complex mCurrentSample = new Complex(0, 0);
+    private Complex mMiddleSample = new Complex(0, 0);
+    private Complex mPreviousSymbol = new Complex(0, 0);
+    private Complex mMiddleSymbol = new Complex(0, 0);
+    private Complex mCurrentSymbol = new Complex(0, 0);
+    private float mPhaseError;
+    private float mSymbolTimingError;
+
+    /**
+     * Decoder for Quaternary Phase Shift Keying (QPSK) and Differential QPSK (DQPSK).  This decoder uses both a
+     * phase-locked loop and a gardner symbol timing error detector to automatically align to the incoming carrier
+     * frequency and to adjust for any changes in symbol timing.
+     *
+     * @param phaseLockedLoop for tracking carrier frequency error
+     * @param symbolPhaseErrorCalculator to calculate symbol phase errors and optionally adjust for differential encoding
+     * @param symbolDecoder to decode demodulated symbols into dibits
+     * @param sampleRate of the incoming complex sample stream
+     * @param symbolRate of the decoded QPSK symbols.
+     */
+    public QPSKDemodulator(IPhaseLockedLoop phaseLockedLoop, ISymbolPhaseErrorCalculator symbolPhaseErrorCalculator,
+                           IQPSKSymbolDecoder symbolDecoder, double sampleRate, double symbolRate)
+    {
+        if(sampleRate < (symbolRate * 2))
+        {
+            throw new IllegalArgumentException("Sample rate [" + sampleRate +
+                "] must be at least 2 x symbol rate [" + symbolRate + "]");
+        }
+
+        mPLL = phaseLockedLoop;
+        mSymbolPhaseErrorCalculator = symbolPhaseErrorCalculator;
+        mQPSKSymbolDecoder = symbolDecoder;
+
+        float samplesPerSymbol = (float)(sampleRate / symbolRate);
+        mInterpolatingSymbolBuffer = new InterpolatingSymbolBuffer(samplesPerSymbol);
+    }
+
+    /**
+     * Submits a buffer containing complex samples for decoding
+     * @param complexBuffer with complex samples
+     */
+    public void receive(ComplexBuffer complexBuffer)
+    {
+        float[] samples = complexBuffer.getSamples();
+
+        for(int x = 0; x < samples.length; x += 2)
+        {
+            receive(samples[x], samples[x + 1]);
+        }
+    }
+
+    /**
+     * Submits a complex sample for decoding.
+     * @param inphase value for the sample
+     * @param quadrature value for the sample
+     */
+    @Override
+    public void receive(float inphase, float quadrature)
+    {
+        //Update current sample with values
+        mCurrentSample.setValues(inphase, quadrature);
+
+        //Mix current sample with costas loop to remove any rotation that is present from a mis-tuned carrier frequency
+        mCurrentSample.multiply(mPLL.incrementAndGetCurrentVector());
+
+        mInterpolatingSymbolBuffer.receive(mCurrentSample);
+
+        //Calculate the symbol once we've stored enough samples
+        if(mInterpolatingSymbolBuffer.hasSymbol())
+        {
+            //Get middle and current samples from the interpolating buffer
+            mMiddleSample = mInterpolatingSymbolBuffer.getMiddleSample();
+            mCurrentSample = mInterpolatingSymbolBuffer.getCurrentSample();
+
+            //Calculate middle and current symbols as the delta between the previous and current samples
+            mMiddleSymbol.setInphase(Complex.multiplyInphase(mMiddleSample.inphase(), mMiddleSample.quadrature(), mPreviousMiddleSample.inphase(), -mPreviousMiddleSample.quadrature()));
+            mMiddleSymbol.setQuadrature(Complex.multiplyQuadrature(mMiddleSample.inphase(), mMiddleSample.quadrature(), mPreviousMiddleSample.inphase(), -mPreviousMiddleSample.quadrature()));
+            mCurrentSymbol.setInphase(Complex.multiplyInphase(mCurrentSample.inphase(), mCurrentSample.quadrature(), mPreviousSample.inphase(), -mPreviousSample.quadrature()));
+            mCurrentSymbol.setQuadrature(Complex.multiplyQuadrature(mCurrentSample.inphase(), mCurrentSample.quadrature(), mPreviousSample.inphase(), -mPreviousSample.quadrature()));
+
+            //Set gain to unity before we calculate the error value
+            mMiddleSymbol.normalize();
+            mCurrentSymbol.normalize();
+
+            //Symbol timing error calculations
+            mSymbolTimingError = mGardnerDetector.getError(mPreviousSymbol, mMiddleSymbol, mCurrentSymbol);
+            mInterpolatingSymbolBuffer.resetAndAdjust(mSymbolTimingError);
+
+            //Store current samples/symbols to use for the next period
+            mPreviousSample.setValues(mCurrentSample);
+            mPreviousMiddleSample.setValues(mMiddleSample);
+            mPreviousSymbol.setValues(mCurrentSymbol);
+
+            //Adjust the current symbol as necessary (e.g. for differential encoding)
+            mSymbolPhaseErrorCalculator.adjust(mCurrentSymbol);
+
+            //Calculate the phase error of the current symbol relative to the expected constellation and provide
+            //feedback to the PLL
+            mPhaseError = mSymbolPhaseErrorCalculator.getPhaseError(mCurrentSymbol);
+            mPLL.adjust(mPhaseError);
+
+            if(mPLLErrorListener != null)
+            {
+                mPLLErrorListener.receive((double)mPhaseError);
+            }
+
+            //Decode the dibit from the symbol and send to the listener
+            if(mDibitListener != null)
+            {
+                mDibitListener.receive(mQPSKSymbolDecoder.decode(mCurrentSymbol));
+            }
+
+            //Send to an external constellation symbol listener when registered
+            if(mSymbolListener != null)
+            {
+                mSymbolListener.receive(mCurrentSymbol);
+            }
+
+            //TODO: Assemble dibits here and broadcast
+        }
+    }
+
+    /**
+     * Registers the listener to receive decoded QPSK symbols
+     */
+    public void setSymbolListener(Listener<Complex> listener)
+    {
+        mSymbolListener = listener;
+    }
+
+    /**
+     * Registers the listener to receive decoded QPSK dibits
+     */
+    public void setDibitListener(Listener<Dibit> listener)
+    {
+        mDibitListener = listener;
+    }
+
+    /**
+     * Registers the listener to receive PLL error values
+     */
+    public void setPLLErrorListener(Listener<Double> listener)
+    {
+        mPLLErrorListener = listener;
+    }
+}

src/main/java/io/github/dsheirer/dsp/psk/QPSKSymbolDecoder.java

@@ -1,3 +1,18 @@
+/*******************************************************************************
+ * sdr-trunk
+ * Copyright (C) 2014-2018 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.psk;
 
 import io.github.dsheirer.dsp.symbol.Dibit;
@@ -5,24 +20,7 @@ import io.github.dsheirer.sample.Broadcaster;
 import io.github.dsheirer.sample.Listener;
 import io.github.dsheirer.sample.complex.Complex;
 
-/*******************************************************************************
- *     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>
+public class QPSKSymbolDecoder implements Listener<Complex>, IQPSKSymbolDecoder
 {
 	private Broadcaster<Dibit> mBroadcaster = new Broadcaster<Dibit>();
 	
@@ -37,7 +35,7 @@ public class QPSKPolarSlicer implements Listener<Complex>
 	 *     /  \
 	 *    / 11 \
 	 */
-	public QPSKPolarSlicer()
+	public QPSKSymbolDecoder()
 	{
 	}
 	
@@ -60,10 +58,10 @@ public class QPSKPolarSlicer implements Listener<Complex>
 	@Override
 	public void receive( Complex complex )
 	{
-		mBroadcaster.broadcast( decide( complex ) );
+		mBroadcaster.broadcast( decode( complex ) );
 	}
 	
-	public static Dibit decide( Complex complex )
+	public Dibit decode(Complex complex )
 	{
 		if( Math.abs( complex.inphase() ) > Math.abs( complex.quadrature() ) )
 		{

src/main/java/io/github/dsheirer/dsp/psk/QPSKSymbolPhaseErrorCalculator.java

@@ -0,0 +1,63 @@
+/*******************************************************************************
+ * sdr-trunk
+ * Copyright (C) 2014-2018 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.psk;
+
+import io.github.dsheirer.sample.complex.Complex;
+
+/**
+ * Calculates the phase error for a QPSK polar-orientation constellation
+ */
+public class QPSKSymbolPhaseErrorCalculator implements ISymbolPhaseErrorCalculator
+{
+    private float mPhaseError;
+
+    @Override
+    public void adjust(Complex symbol)
+    {
+        //No-op - since there is no differential, we don't have to adjust
+    }
+
+    @Override
+    public float getPhaseError(Complex symbol)
+    {
+        mPhaseError = 0.0f;
+
+        if(Math.abs(symbol.inphase()) > Math.abs(symbol.quadrature()))
+        {
+            if(symbol.inphase() > 0)
+            {
+                mPhaseError = -symbol.quadrature();
+            }
+            else
+            {
+                mPhaseError = symbol.quadrature();
+            }
+        }
+        else
+        {
+            if(symbol.quadrature() > 0)
+            {
+                mPhaseError = symbol.inphase();
+            }
+            else
+            {
+                mPhaseError = -symbol.inphase();
+            }
+        }
+
+        return mPhaseError;
+    }
+}

src/main/java/io/github/dsheirer/dsp/psk/SymbolInversion.java

@@ -0,0 +1,35 @@
+/*******************************************************************************
+ * sdr-trunk
+ * Copyright (C) 2014-2018 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.psk;
+
+public enum SymbolInversion
+{
+    DEGREES_90_PLUS(Math.PI / 2.0),
+    DEGREES_90_MINUS(Math.PI / -2.0),
+    DEGREES_180(Math.PI);
+
+    private double mAngleRadians;
+
+    SymbolInversion(double angleRadians)
+    {
+        mAngleRadians = angleRadians;
+    }
+
+    public double getError()
+    {
+        return mAngleRadians;
+    }
+}

src/main/java/io/github/dsheirer/dsp/psk/pll/CostasLoop.java

@@ -0,0 +1,337 @@
+/*******************************************************************************
+ * sdr-trunk
+ * Copyright (C) 2014-2018 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.psk.pll;
+
+import io.github.dsheirer.buffer.DoubleCircularBuffer;
+import io.github.dsheirer.dsp.mixer.Oscillator;
+import io.github.dsheirer.sample.complex.Complex;
+import org.slf4j.Logger;
+import org.slf4j.LoggerFactory;
+
+import java.text.DecimalFormat;
+
+/**
+ * Costas Loop - phase locked loop designed to automatically synchronize to the incoming carrier frequency in order to
+ * zeroize any frequency offset inherent in the signal due to either mis-tuning or carrier frequency drift.
+ *
+ * Costas loop structure was modeled from KA1RBI's OP25/cqpsk.py
+ *
+ * Loop bandwidth and alpha/beta gain formulas derived from:
+ * http://www.trondeau.com/blog/2011/8/13/control-loop-gain-values.html
+ */
+public class CostasLoop implements IPhaseLockedLoop
+{
+    private final static Logger mLog = LoggerFactory.getLogger(CostasLoop.class);
+    public static final double TWO_PI = 2.0 * Math.PI;
+    private DoubleCircularBuffer mErrorBuffer = new DoubleCircularBuffer(10);
+    private double mDamping = Math.sqrt(2.0) / 2.0;
+    private double mAlphaGain;
+    private double mBetaGain;
+    private double mLoopPhase = 0.0;
+    private double mLoopFrequency = 0.0;
+    private double mMaximumLoopFrequency;
+    private Complex mCurrentVector = new Complex(0, 0);
+
+    private Tracking mTracking = Tracking.SEARCHING;
+    private ITrackingStateListener mTrackingStateListener;
+    private boolean mAutomaticTrackingEnabled = true;
+
+    /**
+     * Costas Loop for tracking and correcting frequency error in a received carrier signal.
+     *
+     * @param sampleRate of the incoming samples
+     * @param symbolRate of the digital signal
+     */
+    public CostasLoop(double sampleRate, double symbolRate)
+    {
+        mMaximumLoopFrequency = TWO_PI * (symbolRate / 2.0) / sampleRate;
+        updateLoopBandwidth();
+    }
+
+    /**
+     * Corrects the current phase tracker when an inverted output is detected.  The inversion will take the form of
+     * +/- 90 degrees or +/- 180 degrees, although the latter correction is equal regardless of the sign.
+     *
+     * 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 as measured in radians
+     */
+    public void correctInversion(double correction)
+    {
+        mLoopFrequency += correction;
+
+        if(mLoopFrequency > mMaximumLoopFrequency)
+        {
+            mLoopFrequency -= 2.0 * mMaximumLoopFrequency;
+        }
+
+        if(mLoopFrequency < -mMaximumLoopFrequency)
+        {
+            mLoopFrequency += 2.0 * mMaximumLoopFrequency;
+        }
+    }
+
+    /**
+     * Updates the loop bandwidth and alpha/beta gains according to the current loop synchronization state.
+     */
+    private void updateLoopBandwidth()
+    {
+        double bandwidth = TWO_PI / mTracking.getLoopBandwidth();
+
+        mAlphaGain = (4.0 * mDamping * bandwidth) / (1.0 + (2.0 * mDamping * bandwidth) + (bandwidth * bandwidth));
+        mBetaGain = (4.0 * bandwidth * bandwidth) / (1.0 + (2.0 * mDamping * bandwidth) + (bandwidth * bandwidth));
+    }
+
+    /**
+     * Sets the tracking state for this costas loop.  The tracking state affects the aggressiveness of the alpha/beta
+     * gain values in synchronizing with the signal carrier.
+     *
+     * @param tracking
+     */
+    public void setTracking(Tracking tracking)
+    {
+        mTracking = tracking;
+        updateLoopBandwidth();
+
+        mLog.debug("Tracking: " + tracking.name());
+
+        if(mTrackingStateListener != null)
+        {
+            mTrackingStateListener.trackingStateChanged(tracking);
+        }
+    }
+
+    /**
+     * Sets the listener to receive notifications when the tracking state changes
+     */
+    public void setTrackingStateListener(ITrackingStateListener listener)
+    {
+        mTrackingStateListener = listener;
+    }
+
+    /**
+     * Enables or disables automatic tracking state updates.  When enabled, the loop will automatically increase or
+     * decrease the loop alpha/beta gain values based on the variance of the error tracking values submitted to the
+     * loop control.  As the variance decreases, so will the gain values, and vice-versa.  Setting this to disabled
+     * will force the tracker to retain the gain values associated with either the default or user-specified tracking
+     * value.
+     *
+     * @param enabled true (default) to enable automatic gain tracking
+     */
+    public void setAutomaticTracking(boolean enabled)
+    {
+        mAutomaticTrackingEnabled = enabled;
+    }
+
+    /**
+     * Increments the phase of the loop for each sample received at the sample rate.
+     */
+    public void increment()
+    {
+        mLoopPhase += mLoopFrequency;
+
+        /* Keep the loop phase in bounds */
+        normalizePhase();
+    }
+
+    /**
+     * Current vector of the loop.  Note: this value is updated for the current angle in radians each time this method
+     * is invoked.
+     */
+    public Complex getCurrentVector()
+    {
+        mCurrentVector.setAngle(mLoopPhase);
+        return mCurrentVector;
+    }
+
+    @Override
+    public Complex incrementAndGetCurrentVector()
+    {
+        increment();
+        return getCurrentVector();
+    }
+
+    public double getLoopFrequency()
+    {
+        return mLoopFrequency;
+    }
+
+    public double getMaximumLoopFrequency()
+    {
+        return mMaximumLoopFrequency;
+    }
+
+    /**
+     * Updates the costas loop frequency and phase to adjust for the phase
+     * error value
+     *
+     * @param phaseError - (-)= late and (+)= early
+     */
+    public void adjust(double phaseError)
+    {
+        mLoopFrequency += (mBetaGain * phaseError);
+        mLoopPhase += mLoopFrequency + (mAlphaGain * phaseError);
+
+        /* Maintain phase between +/- 2 * PI */
+        normalizePhase();
+
+        /* Limit frequency to +/- maximum loop frequency */
+        limitFrequency();
+
+        mErrorBuffer.put(phaseError);
+
+        double variance = mErrorBuffer.variance();
+
+        switch(mTracking)
+        {
+            case SEARCHING:
+                if(Tracking.COARSE.meetsThreshold(variance))
+                {
+                    setTracking(Tracking.COARSE);
+                }
+                break;
+            case COARSE:
+                if(Tracking.COARSE.meetsThreshold(variance))
+                {
+                    if(Tracking.FINE.meetsThreshold(variance))
+                    {
+                        setTracking(Tracking.FINE);
+                    }
+                }
+                else
+                {
+                    setTracking(Tracking.SEARCHING);
+                }
+                break;
+            case FINE:
+                if(Tracking.FINE.meetsThreshold(variance))
+                {
+                    if(Tracking.LOCKED.meetsThreshold(variance))
+                    {
+                        setTracking(Tracking.LOCKED);
+                    }
+                }
+                else
+                {
+                    setTracking(Tracking.COARSE);
+                }
+                break;
+            case LOCKED:
+                if(!Tracking.LOCKED.meetsThreshold(variance))
+                {
+                    setTracking(Tracking.FINE);
+                }
+                break;
+        }
+    }
+
+    public double getErrorMean()
+    {
+        return mErrorBuffer.mean();
+    }
+
+    public double getErrorVariance()
+    {
+        return mErrorBuffer.variance();
+    }
+
+    /**
+     * Normalizes the phase tracker to maintain within +/- 2 Pi
+     */
+    private void normalizePhase()
+    {
+        while(mLoopPhase > TWO_PI)
+        {
+            mLoopPhase -= TWO_PI;
+        }
+
+        while(mLoopPhase < -TWO_PI)
+        {
+            mLoopPhase += TWO_PI;
+        }
+    }
+
+    /**
+     * Constrains the frequency within the bounds of +/- loop frequency
+     */
+    private void limitFrequency()
+    {
+        if(mLoopFrequency > mMaximumLoopFrequency)
+        {
+            mLoopFrequency = mMaximumLoopFrequency;
+        }
+
+        if(mLoopFrequency < -mMaximumLoopFrequency)
+        {
+            mLoopFrequency = -mMaximumLoopFrequency;
+        }
+    }
+
+    /**
+     * Tracking state listener interface for receiving updates when the tracking state of the costas loop changes
+     */
+    public interface ITrackingStateListener
+    {
+        void trackingStateChanged(Tracking tracking);
+    }
+
+    public static void main(String[] arguments)
+    {
+        DecimalFormat decimalFormat = new DecimalFormat("0.000000");
+        double sampleRate = 25000.0;
+        int errorFrequency = 411;
+        double loopBandwidth = 50.0;
+
+        Oscillator oscillator = new Oscillator(errorFrequency, (int)sampleRate);
+        CostasLoop costasLoop = new CostasLoop(sampleRate, 4800.0);
+        costasLoop.setTrackingStateListener(new ITrackingStateListener()
+        {
+            @Override
+            public void trackingStateChanged(Tracking tracking)
+            {
+                mLog.debug("Tracking State Changed: " + tracking.name());
+            }
+        });
+
+        for(int x = 0; x < 1000; x++)
+        {
+            Complex complex = oscillator.getComplex();
+
+            oscillator.rotate();
+            costasLoop.increment();
+
+            if(x % 5 == 0)
+            {
+                complex.multiply(costasLoop.getCurrentVector().conjugate());
+
+                double phaseError = complex.angle();
+                costasLoop.adjust(phaseError);
+
+                double loopFrequency = costasLoop.getLoopFrequency() * sampleRate / TWO_PI;
+                double maxFrequency = costasLoop.getMaximumLoopFrequency() * sampleRate / TWO_PI;
+
+                mLog.debug(x + " Frequency: " + errorFrequency +
+                    " Loop:" + decimalFormat.format(loopFrequency) +
+                    " Max:" + maxFrequency +
+                    " Detected Error:" + decimalFormat.format(phaseError) +
+                    " Avg:" + decimalFormat.format(costasLoop.getErrorMean()) +
+                    " Var:" + decimalFormat.format(costasLoop.getErrorVariance()));
+            }
+        }
+    }
+
+}

src/main/java/io/github/dsheirer/dsp/psk/pll/IPhaseLockedLoop.java

@@ -0,0 +1,56 @@
+/*******************************************************************************
+ * sdr-trunk
+ * Copyright (C) 2014-2018 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.psk.pll;
+
+import io.github.dsheirer.sample.complex.Complex;
+
+/**
+ * Phase-Locked Loop (PLL) interface
+ */
+public interface IPhaseLockedLoop
+{
+    /**
+     * Increments the PLL by one sample period
+     */
+    void increment();
+
+    /**
+     * Applies an error value to the PLL.
+     *
+     * @param error value in the range of -PI <> +PI
+     */
+    void adjust(double error);
+
+    /**
+     * Adjusts the current phase of the loop to account for +/-90 or 180 degree inverted phase locks when detected
+     * @param correction to apply to the phase of the loop
+     */
+    void correctInversion(double correction);
+
+    /**
+     * Current vector of the PLL to use in de-spinning the incoming samples.
+     * @return complex vector
+     */
+    Complex getCurrentVector();
+
+    /**
+     * Increments the phase rotation and gets the current vector of the PLL to use in de-spinning the incoming samples.
+     * @return complex vector
+     */
+    Complex incrementAndGetCurrentVector();
+
+
+}

src/main/java/io/github/dsheirer/dsp/psk/pll/Tracking.java

@@ -0,0 +1,48 @@
+/*******************************************************************************
+ * sdr-trunk
+ * Copyright (C) 2014-2018 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.psk.pll;
+
+public enum Tracking
+{
+    SEARCHING(50.0, Math.PI),
+    COARSE(100.0, 0.01),
+    FINE(150.0, 0.001),
+    LOCKED(200.0, 0.0001);
+
+    private double mLoopBandwidth;
+    private double mVarianceThreshold;
+
+    Tracking(double loopBandwidth, double varianceThreshold)
+    {
+        mLoopBandwidth = loopBandwidth;
+        mVarianceThreshold = varianceThreshold;
+    }
+
+    public double getLoopBandwidth()
+    {
+        return mLoopBandwidth;
+    }
+
+    public double getVarianceThreshold()
+    {
+        return mVarianceThreshold;
+    }
+
+    public boolean meetsThreshold(double variance)
+    {
+        return variance < getVarianceThreshold();
+    }
+}