veeresht / CommPy

Digital Communication with Python
http://veeresht.github.com/CommPy
BSD 3-Clause "New" or "Revised" License
538 stars 176 forks source link

size mismatch in channel coding README #90

Closed londumas closed 3 years ago

londumas commented 3 years ago

size mismatch in channel coding README

ValueError                                Traceback (most recent call last)
<ipython-input-8-e1f7b72f48a6> in <module>
     36 
     37 
---> 38     NumErr, BER_soft[cntr] = BER_calc(message_bits, decoded_soft[:-(L-1)]) # bit-error ratio (soft decision)
     39     NumErr, BER_hard[cntr] = BER_calc(message_bits, decoded_hard[:-(L-1)]) # bit-error ratio (hard decision)
     40     NumErr, BER_uncoded[cntr] = BER_calc(message_bits, demodulated_uncoded) # bit-error ratio (uncoded case)

<ipython-input-1-a3bc5a317da7> in BER_calc(a, b)
      4 
      5 def BER_calc(a, b):
----> 6     num_ber = np.sum(np.abs(a - b))
      7     ber = np.mean(np.abs(a - b))
      8     return int(num_ber), ber

ValueError: operands could not be broadcast together with shapes (100000,) (99998,)
EbNo = 5 # energy per bit to noise power spectral density ratio (in dB)
snrdB = EbNo + 10*np.log10(k*rate) # Signal-to-Noise ratio (in dB)
noiseVar = 10**(-snrdB/10) # noise variance (power)

N_c = 10 # number of trials

BER_soft = np.empty((N_c,))
BER_hard = np.empty((N_c,))
BER_uncoded = np.empty((N_c,))

for cntr in range(N_c):

    message_bits = np.random.randint(0, 2, N) # message
    coded_bits = cc.conv_encode(message_bits, trellis) # encoding

    modulated = modem.modulate(coded_bits) # modulation
    modulated_uncoded = modem.modulate(message_bits) # modulation (uncoded case)

    Es = np.mean(np.abs(modulated)**2) # symbol energy
    No = Es/((10**(EbNo/10))*np.log2(M)) # noise spectrum density

    noisy = modulated + np.sqrt(No/2)*\
        (np.random.randn(modulated.shape[0])+\
         1j*np.random.randn(modulated.shape[0])) # AWGN

    noisy_uncoded = modulated_uncoded + np.sqrt(No/2)*\
        (np.random.randn(modulated_uncoded.shape[0])+\
         1j*np.random.randn(modulated_uncoded.shape[0])) # AWGN (uncoded case)

    demodulated_soft = modem.demodulate(noisy, demod_type='soft', noise_var=noiseVar) # demodulation (soft output)
    demodulated_hard = modem.demodulate(noisy, demod_type='hard') # demodulation (hard output)
    demodulated_uncoded = modem.demodulate(noisy_uncoded, demod_type='hard') # demodulation (uncoded case)

    decoded_soft = cc.viterbi_decode(demodulated_soft, trellis, tb_depth, decoding_type='unquantized') # decoding (soft decision)
    decoded_hard = cc.viterbi_decode(demodulated_hard, trellis, tb_depth, decoding_type='hard') # decoding (hard decision)

    NumErr, BER_soft[cntr] = BER_calc(message_bits, decoded_soft[:-(L-1)]) # bit-error ratio (soft decision)
    NumErr, BER_hard[cntr] = BER_calc(message_bits, decoded_hard[:-(L-1)]) # bit-error ratio (hard decision)
    NumErr, BER_uncoded[cntr] = BER_calc(message_bits, demodulated_uncoded) # bit-error ratio (uncoded case)

mean_BER_soft = np.mean(BER_soft) # averaged bit-error ratio (soft decision)
mean_BER_hard = np.mean(BER_hard) # averaged bit-error ratio (hard decision)
mean_BER_uncoded = np.mean(BER_uncoded) # averaged bit-error ratio (uncoded case)

print("Soft decision:\n"+str(mean_BER_soft)+"\n")
print("Hard decision:\n"+str(mean_BER_hard)+"\n")
print("Uncoded message:\n"+str(mean_BER_uncoded)+"\n")
londumas commented 3 years ago

The sizes are message_bits.shape = (100000,) and decoded_soft.shape = (100004,). This is very strange. Taking all but the last 4 bits gives a very bad result of ~ 0.5 error, indeed the bits after the first few are very different. I simply replaced by:

NumErr, BER_soft[cntr] = BER_calc(message_bits, decoded_soft[:-4]) # bit-error ratio (soft decision)
NumErr, BER_hard[cntr] = BER_calc(message_bits, decoded_hard[:-4]) # bit-error ratio (hard decision)

It gives very poor results:

Soft decision:
0.498325

Hard decision:
0.49850700000000003

Uncoded message:
0.008961

np.mean(np.abs(message_bits[:10]-decoded_soft[:10]))
0.0

np.mean(np.abs(message_bits[:100]-decoded_soft[:100]))
0.0

np.mean(np.abs(message_bits[:1000]-decoded_soft[:1000]))
0.368
londumas commented 3 years ago

Fixed in PR https://github.com/veeresht/CommPy/pull/92

londumas commented 3 years ago

I do not know where the hexadecimal matrix came from np.array([[0o171, 0o133]]), but it produces a wrong looking modulation. The bits are only encoded into the bottom symbols.

example-issue

import numpy as np
import commpy.channelcoding.convcode as cc
import commpy.modulation as modulation
import matplotlib.pyplot as plt

###
N = 124

###
M = 4
modem = modulation.PSKModem(M)

###
message_bits = np.random.randint(0, 2, N)
modulated_uncoded = modem.modulate(message_bits)
plt.errorbar(np.real(modulated_uncoded),np.imag(modulated_uncoded),fmt='o',label='modulated+uncoded')

###
generator_matrix = np.array([[5, 7]])
trellis = cc.Trellis(np.array([M]), generator_matrix)
coded_bits = cc.conv_encode(message_bits, trellis)
modulated = modem.modulate(coded_bits)
plt.errorbar(np.real(modulated)+0.05,np.imag(modulated)+0.05,fmt='o',label='modulated+coded matrix [5, 7]')

###
generator_matrix = np.array([[0o171, 0o133]])
trellis = cc.Trellis(np.array([M]), generator_matrix)
coded_bits = cc.conv_encode(message_bits, trellis)
modulated = modem.modulate(coded_bits)
plt.errorbar(np.real(modulated)-0.05,np.imag(modulated)-0.05,fmt='o',label='modulated+coded matrix [0o171, 0o133]')

###
plt.xlabel('I')
plt.ylabel('Q')
plt.legend()
plt.grid()
plt.show()