From 41d4f083cd1275702969cdb63e79dd4e9e711a77 Mon Sep 17 00:00:00 2001 From: Lukas Winkler Date: Wed, 5 Aug 2020 19:43:45 +0200 Subject: [PATCH] initial commit --- .gitignore | 7 ++ main.py | 189 ++++++++++++++++++++++++++++++++++++++++++++++++ spectrograph.py | 12 +++ 3 files changed, 208 insertions(+) create mode 100644 .gitignore create mode 100644 main.py create mode 100644 spectrograph.py diff --git a/.gitignore b/.gitignore new file mode 100644 index 0000000..d96cc8d --- /dev/null +++ b/.gitignore @@ -0,0 +1,7 @@ +.idea/ +sound +*.mp3 +*.ogg +*.wav +*.png +tmp.py diff --git a/main.py b/main.py new file mode 100644 index 0000000..5edc66a --- /dev/null +++ b/main.py @@ -0,0 +1,189 @@ +from typing import List + +import numpy as np +from matplotlib import pyplot as plt +from numpy import ndarray +from scipy.interpolate import interp1d +from scipy.io import wavfile +from scipy.optimize import curve_fit + + +def gaus(x, a, x0, sigma): + return a * np.exp(-(x - x0) ** 2 / (2 * sigma ** 2)) + + +def peak_finder(signal: ndarray, samples_per_second: int, plot: bool = False, cut_borders: bool = False) -> float: + if cut_borders: + signal = signal[20:] + signal = signal[:-20] + + sp: ndarray = np.fft.rfft(signal) + real_spectrum: ndarray = np.abs(sp) + freq = np.fft.rfftfreq(len(signal), 1 / samples_per_second) + peak = freq[np.argmax(real_spectrum)] + # return peak + try: + popt, pcov = curve_fit(gaus, freq, real_spectrum, p0=[3.5e5, peak, 1.6e2], maxfev=2000) + except RuntimeError: + return 0 + # print(popt[1], peak) + if plot: + newx = np.linspace(0, np.max(freq), 1000) + plt.figure() + plt.plot(freq, real_spectrum, linewidth=0.3) + plt.plot(newx, gaus(newx, *popt), linewidth=0.3) + plt.show() + + return float(popt[1]) + + +def even_parity_check(data: List[bool], parity: bool) -> bool: + total = 0 + for bit in data: + if bit: + total += 1 + print(total) + print(total % 2) + even = total % 2 == 0 + parity_bit = not even + return parity_bit == parity + + +assert even_parity_check([False, False, True, True, False, True, False], True) + + +def lsb_first_binary(data: List[bool]) -> int: + i = 1 + total = 0 + for bit in data: + if bit: + total += i + i *= 2 + return total + + +# 4 + 8 + 32 = 44 +assert lsb_first_binary([False, False, True, True, False, True, False]) == 44 + + +def main(): + data: ndarray + sunset = False + samples_per_second, data = wavfile.read("sound/ambient/dinosaur3.wav") + # samples_per_second, data = wavfile.read("SSTV_sunset_audio.wav") + + spm = samples_per_second // 1000 + + # silence_len = 134 if sunset else 510 # ms + silence_len = 885 if sunset else 598 # ms + header_len = 300 # ms + # header_break_len = 10 # ms + header_break_len = 30 if sunset else 10 # ms + vis_bit_length = 30 # ms + + header1_start = silence_len + header1_end = header1_start + header_len + header_break_start = header1_end + header_break_end = header_break_start + header_break_len + header2_start = header_break_end + header2_end = header2_start + header_len + + headers = [ + data[header1_start * spm: header1_end * spm], + data[header2_start * spm: header2_end * spm] + ] + + for header in headers: + header_freq = peak_finder(header, samples_per_second) + print(header_freq) + assert 1895 < header_freq < 1905 # should be 1900 hz + + print(header_break_start) + header_break = data[header_break_start * spm: header_break_end * spm] + header_freq = peak_finder(header_break, samples_per_second) + print(header_freq) + assert 1195 < header_freq < 1206 # should be 1200 hz + + print("--- VIS start ---") + pos = header2_end + bits: List[bool] = [] + for i in range(10): + print(pos, pos + vis_bit_length) + bit_data = data[pos * spm: (pos + vis_bit_length) * spm] + pos += vis_bit_length + pixel_freq = peak_finder(bit_data, samples_per_second, cut_borders=True) + print(pixel_freq) + if i in [0, 9]: # start bit at 1200 Hz, stop bit at 1200 Hz + assert 1190 < pixel_freq < 1210 + continue + if 1090 < pixel_freq < 1110: + bits.append(True) + elif 1290 < pixel_freq < 1310: + bits.append(False) + else: + raise ValueError(f"{pixel_freq} hz is not an unique bit") + print("---") + print(pos) + print(bits) + print(bits[:7], bits[7]) + assert even_parity_check(bits[:7], bits[7]) + mode = lsb_first_binary(bits[:7]) + print(mode) + if mode == 44: + martin_m1 = True + robot36 = False + row_width = 447 + print("Martin M1 detected") + elif mode == 8: + martin_m1 = False + robot36 = True + row_width = 150 + print("Robot 36 detected") + else: + martin_m1 = False + robot36 = False + row_width = 0 + print(f"unknown mode: {mode}") + exit() + row_width_padded = row_width + 100 + + image = [] + overlap = 40 + row = [] + while True: + pixel_data = data[pos * spm - overlap: (pos + 1) * spm + overlap] + pixel_data = pixel_data * np.hamming(len(pixel_data)) + pos += 1 + if len(pixel_data) != 1 * spm + overlap * 2: + break + pixel_freq = peak_finder(pixel_data, samples_per_second, False) + if pixel_freq < 1300: + if len(row): + if len(row) < row_width_padded: + print(row_width_padded - len(row)) + row.extend([0] * (row_width_padded - len(row))) + if len(row) > row_width_padded: + row = row[:250] + image.append(row) + print(len(row), pos) + row = [] + else: + row.append(pixel_freq) + + image = np.asarray(image) + print(image.shape) + if martin_m1: + grid = 147 + image = np.array([image[::, (2 * grid):(3 * grid)], image[::, 0:grid], image[::, grid:2 * grid]]) + image -= 1500 + image /= 2300 - 1500 + image = np.moveaxis(image, 0, -1) + image = image.clip(0, 1) + plt.imshow(image, aspect="auto") + # plt.imsave("output.png", image) + # plt.colorbar() + plt.show() + + +if __name__ == '__main__': + main() diff --git a/spectrograph.py b/spectrograph.py new file mode 100644 index 0000000..ba94f0e --- /dev/null +++ b/spectrograph.py @@ -0,0 +1,12 @@ +from scipy.io import wavfile +from numpy import ndarray +from matplotlib import pyplot as plt + +data: ndarray +samples_per_second, data = wavfile.read("sound/ambient/dinosaur3.wav") + +samples_per_ms = samples_per_second // 1000 +print(samples_per_ms) +plt.specgram(data[0 * samples_per_ms:20 * samples_per_second], Fs=samples_per_second, scale_by_freq=False) +plt.colorbar() +plt.show()