PySDR icon indicating copy to clipboard operation
PySDR copied to clipboard

End-to-End example doesn't seem to work properly

Open Fabudi opened this issue 1 year ago • 6 comments

I've tried to use examples from Textbook to process my own signal, after I tried full code example and provided file, I'm still getting wrong results as I can assume from constellation diagrams. Any suggestions where I go wrong way or misunderstood something? Figure_1 Code I've used:


import matplotlib.pyplot as plt

import numpy as np
from scipy.signal import resample_poly, firwin


def main():
    x = np.fromfile('waves/fm_rds_250k_1Msamples.iq', dtype=np.complex64)
    sample_rate = 250e3
    center_freq = 99.5e6
    plt.scatter(x.real, x.imag, marker='.', label="Original")

    x = 0.5 * np.angle(x[0:-1] * np.conj(x[1:]))

    N = len(x)
    f_o = -57e3
    t = np.arange(N) / sample_rate
    x = x * np.exp(2j * np.pi * f_o * t)

    taps = firwin(numtaps=101, cutoff=7.5e3, fs=sample_rate)
    x = np.convolve(x, taps, 'valid')
    plt.scatter(x.real, x.imag, marker='.', label="After filtering")

    x = x[::10]
    sample_rate = 25e3

    x = resample_poly(x, 19, 25)
    sample_rate = 19e3

    samples = x
    print(len(x))
    samples_interpolated = resample_poly(samples, 32, 1)
    sps = 16
    mu = 0.01
    out = np.zeros(len(samples) + 10, dtype=np.complex64)
    out_rail = np.zeros(len(samples) + 10, dtype=np.complex64)
    i_in = 0
    i_out = 2
    while i_out < len(samples) and i_in + 32 < len(samples):
        out[i_out] = samples_interpolated[i_in * 32 + int(mu * 32)]
        out_rail[i_out] = int(np.real(out[i_out]) > 0) + 1j * int(np.imag(out[i_out]) > 0)
        x = (out_rail[i_out] - out_rail[i_out - 2]) * np.conj(out[i_out - 1])
        y = (out[i_out] - out[i_out - 2]) * np.conj(out_rail[i_out - 1])
        mm_val = np.real(y - x)
        mu += sps + 0.01 * mm_val
        i_in += int(np.floor(mu))
        mu = mu - np.floor(mu)
        i_out += 1
    x = out[2:i_out]

    print(len(x))
    sample_rate /= 16

    samples = x
    N = len(samples)
    plt.scatter(x.real, x.imag, marker='.', label="After time sync")
    phase = 0
    freq = 0

    alpha = 8.0
    beta = 0.002
    out = np.zeros(N, dtype=np.complex64)
    freq_log = []
    for i in range(N):
        out[i] = samples[i] * np.exp(-1j * phase)
        error = np.real(out[i]) * np.imag(out[i])

        freq += (beta * error)
        freq_log.append(freq * sample_rate / (2 * np.pi))
        phase += freq + (alpha * error)

        while phase >= 2 * np.pi:
            phase -= 2 * np.pi
        while phase < 0:
            phase += 2 * np.pi
    x = out
    plt.scatter(x.real, x.imag, marker='.', label="After frequency sync")
    plt.legend(loc="upper right")
    plt.show()

Fabudi avatar May 14 '24 17:05 Fabudi

Did you try using the exact code given in the expandable part of https://pysdr.org/content/rds.html#wrap-up-and-final-code all you should have to change is the filepath. If you want to plot the constellation you'll do it right before the # Demod BPSK

777arc avatar May 14 '24 17:05 777arc

Yes, i've done it, and still get such constellation Figure_1

Fabudi avatar May 14 '24 17:05 Fabudi

I assume I find out where I could have done a mistake. In your example of Fine Frequency Tuning, you attached plot, where there is at least 50k samples, so loop have enough samples to settle: image But in real example there are about 5k samples, so with alpha=8 and beta=0.002 freq_log looks like that: Figure_1 If I set alpha=100 and beta=0.5, I get such freq_log and constellation plot seems to get right but not quite: Figure_8 Figure_2

Fabudi avatar May 14 '24 22:05 Fabudi

In this real example, are you talking about the one I provided a link to, or your own recording?

That last BPSK constellation looks really good, there should be nearly zero errors after you demod.

777arc avatar May 14 '24 23:05 777arc

Yeah, I'm talking about your example. If I continue with it, I get such output after decoding and parsing:

Sync State Detected
Still Sync-ed (Got  1  bad blocks on  50  total)
Still Sync-ed (Got  0  bad blocks on  50  total)
Still Sync-ed (Got  0  bad blocks on  50  total)
PTY: Top 40
program: 29
coverage_area: Regional 4
            ing.                                                 
            ing. Upb                                             
            ing. Upbeat.                                         
            ing. Upbeat. Rea

It seems to be much shorter when in your example, so I tried to tweak alpha and beta a little bit more. After I've set them to original alpha=8.0 and beta=0.02:

Sync State Detected
Still Sync-ed (Got  12  bad blocks on  50  total)
Still Sync-ed (Got  2  bad blocks on  50  total)
Still Sync-ed (Got  0  bad blocks on  50  total)
PTY: Top 40
program: 29
coverage_area: Regional 4
            ing.                                                 
            ing.    eat.                                         
            ing.    eat. Rea                                     
WayF        ing.    eat. Rea                                     
WayFM Up    ing.    eat. Rea                                     
WayFM Uplifting.    eat. Rea                                     
WayFM Uplifting.    eat. Rea                                     
WayFM Uplifting. Upbeat. Rea 

But in that case I once again getting that constellation 2110812384. So now I'm just not really understand, what I'm doing wrong.

Fabudi avatar May 15 '24 09:05 Fabudi

Ahhh I figured out what happened, when I was doing the initial testing/tweaking I had a much longer recording, so stuff like the freq_log was able to level out, but then the file was too large to distribute via github so I shortened it to just 1M samples, and I made sure it still outputted enough final messages but I didnt verify that the freq_log and constellation still looked good. So thanks for bringing this to my attention!

I just spent some time messing with alpha/beta and its interesting that alpha = 8.0 beta = 0.02 gives you the extra characters but its really sensitive to those values..

As far as plotting constellation, by the time you have done symbol sync you're left with about 5k samples. Not only is that a lot to plot together, but the first half or so are going to be before the freq sync finished so they will give you a circle. But if you plot the final 1k samples, after freq sync, it should look pretty good, but you do have to increase beta like you said, since its a short recording. And I'll update the code to reflect this!

... x = out plt.plot(np.real(x[-1000:]), np.imag(x[-1000:]), '.') plt.show()

image

777arc avatar May 16 '24 04:05 777arc

https://github.com/777arc/PySDR/commit/d6fe5c489586c4bbed7346e13eb80555a85a7eaa

777arc avatar May 31 '24 07:05 777arc