From dbdf45ef059590e9030c6861cfed156b8a83e1df Mon Sep 17 00:00:00 2001 From: David Lenfesty Date: Thu, 5 Sep 2019 09:56:02 -0600 Subject: [PATCH] Trying to debug issues with carrier estimation --- channel.py | 55 ++++++++++++++++++++++++++++++++++++++++++------------ main.py | 7 +++---- modem.py | 9 +++++++++ qam.py | 3 ++- serpar.py | 3 +++ 5 files changed, 60 insertions(+), 17 deletions(-) create mode 100644 modem.py diff --git a/channel.py b/channel.py index b895a77..19cb34a 100644 --- a/channel.py +++ b/channel.py @@ -4,13 +4,13 @@ import scipy.interpolate import matplotlib.pyplot as plt # Dirac delta function response, no change -#channel_response = np.array([0, 0, 0, 1, 0, 0, 0]) -channel_response = np.array([-1 + 1j, 0, 1, 0.5, 0.25]) +channel_response = np.array([0, 0, 0, 1, 0, 0, 0]) +#channel_response = np.array([-1 + 1j, 0, 1, 0.5, 0.25]) # How do I sync this across two files? # figure it out later -pilot_value = 1 + 1j +pilot_value = 5 + 5j # 15dB seems to be around the minimum for error-free transmission snr_db = 300 @@ -65,6 +65,8 @@ def estimate(in_data, pilots=0): every time. This also has the advantage of being able to synchronise the symbols. Since I will be implementing some sort of protocol anyways, I think this will be a good idea. As well, we move slow enough that the channel will not likely change significantly over a single packet. + + Actually maybe not? Looked at a paper, check the drive. """ all_carriers = np.arange(len(in_data[0])) @@ -73,24 +75,50 @@ def estimate(in_data, pilots=0): pilot_carriers = all_carriers[::(len(all_carriers)) // pilots] pilot_carriers = np.delete(pilot_carriers, 0) + print(pilot_carriers) + # start averaging H_est = 0 - for i in range(len(in_data)): - # Obtain channel response at pilot carriers - H_est_pilots = in_data[i][pilot_carriers] / pilot_value + # for i in range(len(in_data)): + # # Obtain channel response at pilot carriers + # H_est_pilots = in_data[i][pilot_carriers] / pilot_value - # Interpolate estimates based on what we get from the few pilot values - H_est_abs = scipy.interpolate.interp1d(pilot_carriers, abs(H_est_pilots), kind='linear', fill_value='extrapolate')(all_carriers) - H_est_phase = scipy.interpolate.interp1d(pilot_carriers, np.angle(H_est_pilots), kind='linear', fill_value='extrapolate')(all_carriers) - # Take angular form and turn into rectangular form - H_est += H_est_abs * np.exp(1j*H_est_phase) + # # Interpolate estimates based on what we get from the few pilot values + # H_est_abs = scipy.interpolate.interp1d(pilot_carriers, abs(H_est_pilots), kind='linear', fill_value='extrapolate')(all_carriers) + # H_est_phase = scipy.interpolate.interp1d(pilot_carriers, np.angle(H_est_pilots), kind='linear', fill_value='extrapolate')(all_carriers) - H_est = H_est / len(in_data) + # # Take angular form and turn into rectangular form + # H_est += H_est_abs * np.exp(1j*H_est_phase) + + H_est_pilots = np.ndarray((len(pilot_carriers)), dtype=np.csingle) + + j = 0 + # Obtain channel response at pilot carriers + for i in pilot_carriers: + H_est_pilots[j] = in_data[0][i] / pilot_value + j += 1 + + + # Interpolate estimates based on what we get from the few pilot values + H_est_abs = scipy.interpolate.interp1d(pilot_carriers, abs(H_est_pilots), kind='linear', fill_value='extrapolate')(all_carriers) + H_est_phase = scipy.interpolate.interp1d(pilot_carriers, np.angle(H_est_pilots), kind='linear', fill_value='extrapolate')(all_carriers) + + # Take angular form and turn into rectangular form + H_est += H_est_abs * np.exp(1j*H_est_phase) + + # for an average channel estimate + # H_est = H_est / len(in_data) # Exact channel response H_exact = np.fft.fft(channel_response, len(in_data[0])) + + for i in range(3): + plt.plot(abs(in_data[i]), "b") + plt.plot(np.abs(in_data[i]), "r") + plt.show(block=True) + plt.plot(all_carriers, np.real(H_exact), "b") plt.plot(all_carriers, np.real(H_est), "r") plt.plot(pilot_carriers, np.real(H_est_pilots), "go") @@ -101,6 +129,9 @@ def estimate(in_data, pilots=0): plt.plot(pilot_carriers, np.imag(H_est_pilots), "go") plt.show(block=True) + print(H_est_pilots) + print(H_exact[pilot_carriers]) + return H_est diff --git a/main.py b/main.py index 2ee78bc..198011d 100755 --- a/main.py +++ b/main.py @@ -20,7 +20,6 @@ Right now though it's just proof of concept to see if I can get a reliable signa import numpy as np import matplotlib.pyplot as plt -plt.ion() import channel import qam @@ -69,7 +68,7 @@ if __name__ == '__main__': parallel = parallelise(64, bytes) # modulate data with a QAM scheme - modulated = qam.modulate(parallel, pilots=20) + modulated = qam.modulate(parallel, pilots=10) # Run IFFT to get a time-domain signal to send ofdm_time = np.fft.ifft(modulated) @@ -87,13 +86,13 @@ if __name__ == '__main__': to_equalize = np.fft.fft(ofdm_cp_removed) # Find an estimate for channel effect - H_est = channel.estimate(to_equalize, pilots=20) + H_est = channel.estimate(to_equalize, pilots=10) # Equalise based on estimated channel to_decode = channel.equalize(to_equalize, H_est) # Demodulate symbol into output data - to_serialise = qam.demodulate(to_decode, pilots=20) + to_serialise = qam.demodulate(to_decode, pilots=10) # Turn data back into string data = serialise(64, to_serialise) diff --git a/modem.py b/modem.py new file mode 100644 index 0000000..c65c691 --- /dev/null +++ b/modem.py @@ -0,0 +1,9 @@ +""" +The modem for communications. + +I need to define my protocol here first. + +Packet: + +Prefix +""" diff --git a/qam.py b/qam.py index 7ce7a9b..876ea88 100644 --- a/qam.py +++ b/qam.py @@ -1,7 +1,7 @@ import numpy as np from scipy.spatial.distance import euclidean -pilot_value = 1 + 1j +pilot_value = 5 + 5j qam_mapping_table = { 0 : 1 + 1j, @@ -36,6 +36,7 @@ def modulate(in_data, pilots=0): pilot_carriers = all_carriers[::(num_data_carriers + pilots)//pilots] pilot_carriers = np.delete(pilot_carriers, 0) # not sure how to not have this line data_carriers = np.delete(all_carriers, pilot_carriers) + print(pilot_carriers) else: data_carriers = all_carriers diff --git a/serpar.py b/serpar.py index 9a44537..ffd02d7 100644 --- a/serpar.py +++ b/serpar.py @@ -67,7 +67,10 @@ def serialise(n, in_data): new_byte = np.uint8(0) for k in range(4): + #try: new_byte |= in_data[i][(j*4) + k] << (k * 2) + #except: + # print(i) out_data.append(new_byte)