<jupyter_start><jupyter_text>Lab 2. Sounds and HarmonicsIn this lab we will produce some sounds with their harmonics, in order to study their acoustic effects<jupyter_code>import numpy as np
import matplotlib.pyplot as plt
import scipy.signal as ss
from math import pi, sin
from matplotlib.ticker import LogFormatter
import simpleaudio as sa
import ipywidgets as widgets
fs = 44100        # 44100 samples per second
samples = 2**16   # for the FFT<jupyter_output><empty_output><jupyter_text>Let's define a function to map the harmonic index (1 = fundamental, 2 = II harmonic, etc.) to the relative amplitude, according to different models.For the _harp_ model, we use the coefficients derived in the first lesson:$$C_n \propto \frac{1}{n^2} \sin(n \pi k), \quad k \in \ ]0, 0.5]$$Where $k$ represents where the string is pinched, from the edge (0) to the middle (0.5) of the string.<jupyter_code>def htoamp(h, model, k=0.5):
    if model == 'flat':
        return 1
    if model == 'harp':
        return 1 / h ** 2 * sin(h * pi * k)
    if model == 'even':
        return 1 / h ** (2 if h % 2 > 0 else 1.1)   # even harmonics decade faster
    if model == 'odd':
        return 1 / h ** (2 if h % 2 == 0 else 1.1)    # odd harmonics decade faster
    if model == 'verotta':
        return 1 / h
    if model == 'sethares':
        return 0.8**h
    return 0
htoamp(1, 'harp', 3)
class Chord:
    def __init__(self, f, duration, harmonics, onlyharm, hmodel, k_harp):
        self.f = f
        self.duration = duration
        self.harmonics = harmonics
        self.onlyharm = onlyharm
        self.hmodel = hmodel
        self.k_harp = k_harp
    
    def note(self):
        '''Generate a np array for a given frequency. If onlyharm = True, generates only the harmonics'''
        # generate array with duration*sample_rate steps, ranging between 0 and duration
        t = np.linspace(0, self.duration, int(self.duration * fs), False)

        if self.onlyharm:
            # generate only the requested harmonic
            w = np.sin(2 * np.pi * self.f * self.harmonics * t) * htoamp(self.harmonics, self.hmodel, self.k_harp) * \
                np.exp(- self.harmonics * t * 3/self.duration)
            return w

        # generate the fundamental wave
        w = np.sin(2 * np.pi * self.f * t) * htoamp(1, self.hmodel, self.k_harp) * np.exp(- t * 3/self.duration)
        # add the harmonics with given weights
        if self.harmonics > 1:
            for h in range(self.harmonics-1):   # h = 0..harmonics-2
                w += np.sin(2 * np.pi * self.f * (h+2) * t) * htoamp(h+2, self.hmodel, self.k_harp) * \
                     np.exp(- (h+2) * t * 3/self.duration)
        return w

    def play(self, audio):
        '''Play a given signal to the audio card'''
        # ensure that highest value is in 16-bit range
        playable = audio * (2**15 - 1) / np.max(np.abs(audio))
        if self.onlyharm:
            # scale down so to hear the effect of the scaled single harmonic
            playable *= htoamp(self.harmonics, self.hmodel, self.k_harp)
        # stop any ongoing play
        #sa.stop_all()
        # convert to 16-bit data and play
        return sa.play_buffer(playable.astype(np.int16), 1, 2, fs)

    def plottimefreq(self, s):
        # plot the time series
        plt.subplots(figsize=(15, 5)) 
        ax = plt.subplot(1, 2, 1)
        plt.plot(s)
        plt.xlim(0)
        ax.set_xlabel('t [ms]')
        maxx = int(self.duration*10 + 1)*100
        ax.set_xticks(np.arange(0, maxx*fs/1000, maxx*fs/10000, dtype=int))
        ax.set_xticklabels(np.arange(0, maxx, maxx/10, dtype=int))
        plt.grid(which='major')
        plt.title("Wave packet")

        # also compute and plot power spectrum
        ax = plt.subplot(1, 2, 2)
        s = np.pad(s, (0, samples-s.size), mode='constant')
        W = np.abs(np.fft.fft(s) ** 2)
        f = np.fft.fftfreq(s.size, 1/fs)
        plt.plot(f, W)
        plt.xlim(20, 15000)
        #formatter = LogFormatter(labelOnlyBase=False, minor_thresholds=(1, 0.1))
        #ax.get_xaxis().set_minor_formatter(formatter)
        ax.set_xlabel('f [Hz]')
        plt.ylim(1E-2)
        plt.xscale('log')
        plt.yscale('log')
        plt.grid(which='both')
        plt.title("Power spectrum (log/log)")
        plt.show()

@widgets.interact(
    frequency=widgets.IntSlider(min=100, max=500, value=220, continuous_update=False),
    duration=widgets.FloatSlider(min=0.01, max=1.48, value=0.7, step=0.01, continuous_update=False),
    harmonics=widgets.IntSlider(min=1, max=30, value=0, continuous_update=False),
    only_harm=widgets.Checkbox(description='only harmonic',value=True),
    h_model=widgets.RadioButtons(options=['flat', 'harp', 'even', 'odd', 'verotta']),
    k_harp=widgets.FloatSlider(min=0.01, max=0.5, value=0.5, step=0.01, continuous_update=False),
    window=widgets.RadioButtons(options=['hanning', 'rect', 'gaussian'])
    )
def interactiveplay(frequency, duration, harmonics, only_harm, h_model, k_harp, window):
    # generate signal for the given chord
    c = Chord(frequency, duration, harmonics, only_harm, h_model, k_harp)
    s = c.note()
    # use a window to smooth begin and end
    if window == 'hanning':
        s *= np.hanning(s.size)
    elif window == 'gaussian':
        s *= ss.gaussian(s.size, 5000*duration)
    c.plottimefreq(s)
    c.play(s)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.12946180876641392, 'fraction_numerical': 0.02903643425189569, 'mean_word_length': 2.6270959087860497, 'pattern_counts': {'":': 0, '<': 9, '<?xml version=': 0, '>': 11, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '2', 'dir': 'jupyter-structured-clean-dedup', 'id': '664667', 'license': 'MIT', 'n_tokens_mistral': 1946, 'n_tokens_neox': 1835, 'n_words': 515, 'path': 'PhysicsAndMusic/Lab2 - Harmonics.ipynb', 'repo_name': 'glpatcern/inotebooks', 'size': '5407'}