<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 FSK Modulation
Now that we've got some ideas for demodulating fsk, let's do some freakin' modulation!
First, how do we increase the frequency of a wave?
For starters, let's just try doubling the frequency of the carrier wave. What mathematical operation can we use to double a wave's frequency? Well, going back to high school trigonometric identities, sin(2x) = 2*sin(x)*cos(x). So if our starting wave is sin(x) and we want to double it's frequency, we can just multiply it by 2*cos(x).
Step1:  Interesting, so just multiplying by 2*cos(x) doubles the frequency. Or, to think of it another way, just multiplying by a 90-degree phase-shifted wave of the same frequency (cos(x)... or as it could alternatively be called, sin((pi/2)-x)) doubles the frequency of a wave, but when the frequency doubles, the amplitude halves, and we have to do a scalar doubling in order to restore the original amplitude.
So let's now try to use this to modulate a carrier wave...
Step2:  Excellent, so that worked well.
An easier approach
Another simpler way to do FSK modulation could be to just calculate arrays for each of the different frequencies, and then just concatenate those arrays together. In this case, we wouldn't really be modulating a "carrier" wave so much as building the modulated wave, one bit-worth at a time.
Let's try that approach...
Step3:  As can be seen, the results are practically equivalent, and this later approach is much less computationally intensive (as it mostly just entails memcpy's).
So let's now build a general FSK algorithm using this approach.
General FSK modulation algorithm
Let's start with BFSK (only 2 tones)...
Step4:  Let's see how to write the wave data out to a .wav file
Step5:  General FSK algorithm
Step6:  Now let's listen to one of audible frequencies (by just multiplying each of the frequencies by 10), and speed it up (by increasing the baud rate).
Step7:  FSK Demodulation
I want to be able to write a software FSK (Frequency Shift Key) Modem. First, I'll try to figure out the Demodulation step...
FSK Demodulation with Zero-crossing detection
The simplest approach is zero-crossing detection - just see how often the sign changes to determine the frequency in a particular point.
Step8:  Excellent! So using zero-crossing detection, we accurately determined the frequency to be 20Hz
FSK Demodulation with Multiplying by prospective frequencies
To try to determine the frequency, we can just take the dot product of a number of prospective frequencies, and see which has the largest dot product, which should correspond to the closest approximation of the frequency.
Step9:  Success. The 20Hz dot product is clearly much larger than the 10, 20, or 40Hz dot products, indicating the frequency is 20Hz. 
It's worth noting, however, that if the frequency is just a bit off, the dot product is still near zero, so there's not a lot of room for error in terms of the frequency with this approach. Only 0.5Hz off has a near-zero dot product.
Step10:  Relationship between frequency difference and dot product?
Out of curiosity, let's graph the relationship between the difference in the frequencies and the dot product...
Step11:  Interesting! So it looks like a reflected damped sine wave. Now, does the variation in the dot product change based on absolute deviation from the center frequency, or is it based on percentage? Let's try a higher frequency with the same deviations...
Step12:  So it looks like the dot product varies as you move away from the center frequency purely based on absolute difference in Hz, rather than based on percentage of deviation (relative to the frequency).
FSK Demodulation with FFT
The above method is kind of like a manual FFT. Let's see how the normal FFT fares.
Step13:  Just as with the manual wave dot products above, you can see that the FFT output corresponding to 20Hz is the largest, indicating the wave is 20Hz. 
Now, it's worth noting that the fft output length defaults to the num_samples/2. But for efficiency, we can constrain the fft to some lower number of output bins...
Step14:  And it works, but since the sample size didn't change, we need to divide the frequencies we're looking for by 10 as well (so below, fft_output100[2] corresponds to the 20Hz component.
Step15:  Let's do a full demodulation algorithm now
First, let's make a nice testing modulated signal
Step16:  How to do full demodulation algorithm
First, we have to locate the beginning of the frame. We'll postpone that for now, and just assume the beginning of the signal is the beginning of the frame.
Next, we need to split the the signal into the windows to represent each bit.
Finally, we'll determine the bit represented by each window.
<ASSISTANT_TASK:>
Python Code:
samp_rate = 1000
len_in_sec = 1
t = np.linspace(0, 1, samp_rate * len_in_sec)
hz_4 = 1*np.sin(4 * 2 * np.pi * t)
hz_8 = hz_4 * (2 * np.cos(4 * 2 * np.pi * t))
plt.plot(t, hz_4)
plt.show()
plt.plot(t, hz_8)
plt.show()
Explanation: FSK Modulation
Now that we've got some ideas for demodulating fsk, let's do some freakin' modulation!
First, how do we increase the frequency of a wave?
For starters, let's just try doubling the frequency of the carrier wave. What mathematical operation can we use to double a wave's frequency? Well, going back to high school trigonometric identities, sin(2x) = 2*sin(x)*cos(x). So if our starting wave is sin(x) and we want to double it's frequency, we can just multiply it by 2*cos(x).
End of explanation
samp_rate = 1000  # samples/second
len_in_sec = 1
space_freq = 8  # Hz
t = np.linspace(0, 1, samp_rate * len_in_sec)
carrier = np.sin(space_freq * 2 * np.pi * t)
# Note: in FSK venaculare, "space" represents 0 and "mark" 1.
mark_multiplier_array = 2*np.cos(space_freq * 2 * np.pi * t)
modulation_array = np.array([1]*500 + [i for i in mark_multiplier_array[500:]])
setup_graph(title='modulation array for "01"', fig_size=(10,4))
plt.plot(modulation_array)
plt.margins(0.1)
setup_graph(title='modulated "01"', fig_size=(10,4))
plt.plot(carrier * modulation_array)
plt.margins(0.1)
Explanation: Interesting, so just multiplying by 2*cos(x) doubles the frequency. Or, to think of it another way, just multiplying by a 90-degree phase-shifted wave of the same frequency (cos(x)... or as it could alternatively be called, sin((pi/2)-x)) doubles the frequency of a wave, but when the frequency doubles, the amplitude halves, and we have to do a scalar doubling in order to restore the original amplitude.
So let's now try to use this to modulate a carrier wave...
End of explanation
samp_rate = 1000  # samples/second
len_in_sec = 1
space_freq = 8  # Hz ("space" = "0")
mark_freq = 16  # Hz ("mark" = "1")
t = np.linspace(0, .5, samp_rate * len_in_sec)
space = np.sin(space_freq * 2 * np.pi * t)
mark = np.sin(mark_freq * 2 * np.pi * t)
modulated_01 = np.append(space, mark)
setup_graph(title='modulation array for "01"', fig_size=(10,4))
plt.plot(modulated_01)
plt.margins(0.1)
Explanation: Excellent, so that worked well.
An easier approach
Another simpler way to do FSK modulation could be to just calculate arrays for each of the different frequencies, and then just concatenate those arrays together. In this case, we wouldn't really be modulating a "carrier" wave so much as building the modulated wave, one bit-worth at a time.
Let's try that approach...
End of explanation
def bfsk_modulate(bit_array, space_freq, mark_freq, baud, sample_rate):
    seconds_per_bit = 1 / baud
    samples_per_bit = sample_rate * seconds_per_bit
    t = np.linspace(0, seconds_per_bit, samples_per_bit)
    space = np.sin(space_freq * 2 * np.pi * t)
    mark = np.sin(mark_freq * 2 * np.pi * t)
    signal = np.array([])
    for bit in bit_array:
        if bit == 0:
            signal = np.append(signal, space)
        elif bit == 1:
            signal = np.append(signal, mark)
    return signal
sig_010110 = bfsk_modulate([0,1,0,1,1,0], 8, 16, 1, 1000)
setup_graph(title='bfsk-modulated "010110"', fig_size=(18,5))
plt.plot(sig_010110)
plt.margins(0.05)
sig_010110 = bfsk_modulate([0,1,0,1,1,0], 400, 800, 5, 3200)
from IPython.display import Audio
Audio(sig_010110, rate=3200)
Explanation: As can be seen, the results are practically equivalent, and this later approach is much less computationally intensive (as it mostly just entails memcpy's).
So let's now build a general FSK algorithm using this approach.
General FSK modulation algorithm
Let's start with BFSK (only 2 tones)...
End of explanation
import scipy
import scipy.io.wavfile
def write_audio_file(filename, filedata, sample_rate):
    scipy.io.wavfile.write(filename, sample_rate, filedata)
write_audio_file('raw_data/bfsk_010110.wav', sig_010110, 3200)
Explanation: Let's see how to write the wave data out to a .wav file
End of explanation
def fsk_modulate(bit_str, bit_freq_map, baud, sample_rate):
    seconds_per_bit = 1 / baud
    samples_per_bit = int(sample_rate * seconds_per_bit)
    t = np.linspace(0, seconds_per_bit, samples_per_bit)
    
    # maps from bit sequence (like "10") to the modulated wave representing that "symbol"
    symbol_map = {bit_seq: np.sin(freq * 2 * np.pi * t) for bit_seq, freq in bit_freq_map.items()}
    signal = np.array([])
    bits_per_symbol = len(list(bit_freq_map.keys())[0])  # Assume all keys are the same length
    for symbol in [bit_str[i:i+bits_per_symbol] for i in range(0, len(bit_str), bits_per_symbol)]:
        symbol_wave = symbol_map[symbol]
        signal = np.append(signal, symbol_wave)
    return signal
bit_freq_map = {
    "00": 5,
    "01": 10,
    "10": 15,
    "11": 20
}
sig_00011011 = fsk_modulate("00011011", bit_freq_map, 1, 5000)
setup_graph(title='fsk-modulated "00011011"', fig_size=(18,5))
plt.plot(sig_00011011)
plt.margins(0.05)
Explanation: General FSK algorithm
End of explanation
bit_freq_map = {
    "00": 500,
    "01": 1000,
    "10": 1500,
    "11": 2000
}
sig_00011011 = fsk_modulate("00011011", bit_freq_map, 8, 5000)
setup_graph(title='fsk-modulated "00011011"', fig_size=(18,5))
plt.plot(sig_00011011)
plt.margins(0.05)
Audio(sig_00011011, rate=5000)
Explanation: Now let's listen to one of audible frequencies (by just multiplying each of the frequencies by 10), and speed it up (by increasing the baud rate).
End of explanation
samp_rate = 1000
len_in_sec = 1
carrier_freq = 20  # Hz
t = np.linspace(0, 1, samp_rate * len_in_sec)
carrier = 1*np.sin(carrier_freq * 2 * np.pi * t)
plt.plot(t, carrier)
zero_crossings = np.where(np.diff(np.sign(carrier)))[0]
zero_crossings
len(zero_crossings)
# Note that in a single wave cycle, there will be 2 zero crossings
frequency_detected = len(zero_crossings) / 2
frequency_detected
Explanation: FSK Demodulation
I want to be able to write a software FSK (Frequency Shift Key) Modem. First, I'll try to figure out the Demodulation step...
FSK Demodulation with Zero-crossing detection
The simplest approach is zero-crossing detection - just see how often the sign changes to determine the frequency in a particular point.
End of explanation
hz_10 = 1*np.sin(10 * 2 * np.pi * t)
hz_20 = 1*np.sin(20 * 2 * np.pi * t)
hz_30 = 1*np.sin(30 * 2 * np.pi * t)
hz_40 = 1*np.sin(40 * 2 * np.pi * t)
[np.dot(carrier, hz_10), np.dot(carrier, hz_20), np.dot(carrier, hz_30), np.dot(carrier, hz_40)]
sum([hz_20[i]*carrier[i] for i in range(len(carrier))])
Explanation: Excellent! So using zero-crossing detection, we accurately determined the frequency to be 20Hz
FSK Demodulation with Multiplying by prospective frequencies
To try to determine the frequency, we can just take the dot product of a number of prospective frequencies, and see which has the largest dot product, which should correspond to the closest approximation of the frequency.
End of explanation
[np.dot(carrier, 1*np.sin(freq * 2 * np.pi * t)) for freq in [19.5, 19.6, 19.7, 19.8, 19.9, 20]]
Explanation: Success. The 20Hz dot product is clearly much larger than the 10, 20, or 40Hz dot products, indicating the frequency is 20Hz. 
It's worth noting, however, that if the frequency is just a bit off, the dot product is still near zero, so there's not a lot of room for error in terms of the frequency with this approach. Only 0.5Hz off has a near-zero dot product.
End of explanation
def dot_between_freqs(f1, f2):
    t = np.linspace(0, 1, samp_rate * len_in_sec)
    f1_samples = 1*np.sin(f1 * 2 * np.pi * t)
    f2_samples = 1*np.sin(f2 * 2 * np.pi * t)
    return np.dot(f1_samples, f2_samples)
center_freq = 20
frequency_diffs = np.linspace(-10, 10, 500)
dots = [dot_between_freqs(center_freq, center_freq+d) for d in frequency_diffs]
setup_graph(title='frequency deviation vs dot product', x_label='frequency deviation (in Hz)', y_label='dot product', fig_size=(14,7))
plt.plot(frequency_diffs, dots)
Explanation: Relationship between frequency difference and dot product?
Out of curiosity, let's graph the relationship between the difference in the frequencies and the dot product...
End of explanation
center_freq = 20000
frequency_diffs = np.linspace(-10, 10, 500)
dots = [dot_between_freqs(center_freq, center_freq+d) for d in frequency_diffs]
setup_graph(title='frequency deviation vs dot product', x_label='frequency deviation (in Hz)', y_label='dot product', fig_size=(14,7))
plt.plot(frequency_diffs, dots)
Explanation: Interesting! So it looks like a reflected damped sine wave. Now, does the variation in the dot product change based on absolute deviation from the center frequency, or is it based on percentage? Let's try a higher frequency with the same deviations...
End of explanation
fft_output = np.fft.rfft(carrier)
_ = plt.plot(fft_output)
len(carrier)
[np.abs(fft_output[10]), np.abs(fft_output[20]), np.abs(fft_output[30])]
Explanation: So it looks like the dot product varies as you move away from the center frequency purely based on absolute difference in Hz, rather than based on percentage of deviation (relative to the frequency).
FSK Demodulation with FFT
The above method is kind of like a manual FFT. Let's see how the normal FFT fares.
End of explanation
fft_output100 = np.fft.rfft(carrier, n=100)
plt.plot(fft_output100)
Explanation: Just as with the manual wave dot products above, you can see that the FFT output corresponding to 20Hz is the largest, indicating the wave is 20Hz. 
Now, it's worth noting that the fft output length defaults to the num_samples/2. But for efficiency, we can constrain the fft to some lower number of output bins...
End of explanation
[np.abs(fft_output100[1]), np.abs(fft_output100[2]), np.abs(fft_output100[3])]
Explanation: And it works, but since the sample size didn't change, we need to divide the frequencies we're looking for by 10 as well (so below, fft_output100[2] corresponds to the 20Hz component.
End of explanation
aprs_msg = "W0HAK>NY5N:>Hello World'"
aprs_msg_bits = ''.join(["{0:b}".format(ord(c)).zfill(8) for c in aprs_msg])
aprs_msg_bits
bit_freq_map = {
    "0": 2200,
    "1": 1200
}
baud = 1200
sample_rate = 44100
aprs_msg_signal = fsk_modulate(aprs_msg_bits, {"0": 2200, "1": 1200}, baud=1200, sample_rate=44100)
len(aprs_msg_signal)
Audio(aprs_msg_signal, rate=44100)
Explanation: Let's do a full demodulation algorithm now
First, let's make a nice testing modulated signal
End of explanation
def fsk_demodulate(raw_signal, bit_freq_map, baud, sample_rate):
    seconds_per_bit = 1 / baud
    samples_per_bit = int(sample_rate * seconds_per_bit)
    t = np.linspace(0, seconds_per_bit, samples_per_bit)
    
    # maps from bit sequence (like "10") to the modulated wave representing that "symbol"
    wave_to_symbol_map = {bit_seq: np.sin(freq * 2 * np.pi * t) for bit_seq, freq in bit_freq_map.items()}
    
    bit_str = ""
    for index in range(0, len(raw_signal), samples_per_bit):
        best_symbol = ""
        highest_dot_abs = 0
        for symbol, symbol_wave in wave_to_symbol_map.items():
            raw_window = raw_signal[index:index+samples_per_bit]
            dot_abs = np.abs(np.dot(symbol_wave[0:len(raw_window)], raw_window))
            if dot_abs > highest_dot_abs:
                best_symbol = symbol
                highest_dot_abs = dot_abs
        bit_str += best_symbol
    return bit_str
demodulated_aprs_msg = fsk_demodulate(aprs_msg_signal, bit_freq_map, baud, sample_rate)
demodulated_aprs_msg
demodulated_aprs_msg == aprs_msg_bits
Explanation: How to do full demodulation algorithm
First, we have to locate the beginning of the frame. We'll postpone that for now, and just assume the beginning of the signal is the beginning of the frame.
Next, we need to split the the signal into the windows to represent each bit.
Finally, we'll determine the bit represented by each window.
End of explanation
<END_TASK>