メディアネットワーク実験IIA（サウンドプログラミング）

科目名：　メディアネットワーク実験IIA（2021年～）
対象：　メディアネットワークコース3年目
日時：　（11月10日（水） - 11月11日（木）13:00～18:00）
場所：　オンライン（2019年度まではM棟1階計算機室でした．）
レポート提出締切：　（11月25日（木）13:00）
レポート提出先：　slackのダイレクトメッセージで「青木直史」あてに提出すること．
連絡先：　青木 直史（情報エレクトロニクス棟6階6-07）（Tel: 011-706-6532）（E-mail: aoki@ime.ist.hokudai.ac.jp）

目的

　音はマルチメディアコンテンツを構成する重要な要素である．本実験は，Pythonによる具体例を通して，サウンドプログラミングに対する理解を深めることを目的としている．

1．はじめに

　本演習は，Jupyter Notebookを利用し，ブラウザを使ってプログラムを実行しながら進めるものとする．なお，音を確認する場合は各自のイヤフォンまたはヘッドフォンを使うこと．環境のインストールはつぎの手順のとおり．

(1) Anacondaのインストール

https://www.anaconda.com/distribution/
から「Python 3.8」のインストーラをダウンロードして実行．

(2) Jupyter Notebookの起動

Anaconda Promptを開く．

jupyter notebook
と入力して実行．

(3) ブラウザ（Chrome推奨）を利用してプログラミングを行う．

「New」ボタンをクリックして「Python 3」を選択．

音ファイルや画像ファイルを利用するときは「Upload」ボタンをクリックし，必要なファイルをJupyter Notebookにアップロードすること．

2．フィルタ

　フィルタは、特定の帯域の周波数成分だけを選択的に通過させるふるいである．「New」ボタンをクリックし，新しくウィンドウを作成しなさい．つづいて，以下のプログラムを順番にセルに貼りつけ，実行しなさい．高域の周波数成分がカットされることを確認しなさい．

(1)

%matplotlib inline
import numpy as np
from scipy.io import wavfile
from IPython.display import display, Audio

(2)

def LPF(fs, fc, Q):
    fc /= fs
    fc = np.tan(np.pi * fc) / (2.0 * np.pi)
    a = np.zeros(3)
    b = np.zeros(3)
    a[0] = 1.0 + 2.0 * np.pi * fc / Q + 4.0 * np.pi * np.pi * fc * fc
    a[1] = (8.0 * np.pi * np.pi * fc * fc - 2.0) / a[0]
    a[2] = (1.0 - 2.0 * np.pi * fc / Q + 4.0 * np.pi * np.pi * fc * fc) / a[0]
    b[0] = 4.0 * np.pi * np.pi * fc * fc / a[0]
    b[1] = 8.0 * np.pi * np.pi * fc * fc / a[0]
    b[2] = 4.0 * np.pi * np.pi * fc * fc / a[0]
    a[0] = 1.0
    return a, b

(3)

fs = 8000

(4)

length_of_s = int(fs * 1)
s0 = np.zeros(length_of_s)

(5)

for i in range(1, 9):
    for n in range(length_of_s):
        s0[n] += 1.0 * np.sin(2 * np.pi * 440 * i * n / fs)

gain = 0.5 / np.max(np.abs(s0))
s0 *= gain

(6)

for n in range(length_of_s):
    s0[n] = (s0[n] + 1.0) / 2.0 * 65536.0
    if s0[n] > 65535.0:
        s0[n] = 65535.0
    elif s0[n] < 0.0:
        s0[n] = 0.0;
    s0[n] = (s0[n] + 0.5) - 32768

wavfile.write('pulse_train.wav', fs, s0.astype(np.int16))

(7)

Audio('pulse_train.wav')

(8)

s1 = np.zeros(length_of_s)

(9)

fc = 880
Q = 1 / np.sqrt(2)
a, b = LPF(fs, fc, Q)
for n in range(length_of_s):
    for m in range(0, 3):
        if n - m >= 0:
            s1[n] += b[m] * s0[n - m]

    for m in range(1, 3):
        if n - m >= 0:
            s1[n] += -a[m] * s1[n - m]

(10)

master_volume = 0.5
s1 /= np.max(np.abs(s1))
s1 *= master_volume

(11)

for n in range(length_of_s):
    s1[n] = (s1[n] + 1.0) / 2.0 * 65536.0
    if s1[n] > 65535.0:
        s1[n] = 65535.0
    elif s1[n] < 0.0:
        s1[n] = 0.0;
    s1[n] = (s1[n] + 0.5) - 32768

wavfile.write('p1.wav', fs, s1.astype(np.int16))

(12)

Audio('p1.wav')

3．減算合成1（ハイハット）

　原音をフィルタで削ることで音色をつくる音響合成のテクニックを減算合成と呼ぶ。「New」ボタンをクリックし，新しくウィンドウを作成しなさい．つづいて，以下のプログラムを順番にセルに貼りつけ，実行しなさい．白色雑音からハイハットの音をつくることができることを確認しなさい．

(1)

%matplotlib inline
import numpy as np
from scipy.io import wavfile
from IPython.display import display, Audio

(2)

def HPF(fs, fc, Q):
    fc /= fs
    fc = np.tan(np.pi * fc) / (2.0 * np.pi)
    a = np.zeros(3)
    b = np.zeros(3)
    a[0] = 1.0 + 2.0 * np.pi * fc / Q + 4.0 * np.pi * np.pi * fc * fc;
    a[1] = (8.0 * np.pi * np.pi * fc * fc - 2.0) / a[0]
    a[2] = (1.0 - 2.0 * np.pi * fc / Q + 4.0 * np.pi * np.pi * fc * fc) / a[0]
    b[0] = 1.0 / a[0]
    b[1] = -2.0 / a[0]
    b[2] = 1.0 / a[0]
    a[0] = 1.0
    return a, b

(3)

def ADSR(fs, A, D, S, R, gate, duration):
    A = int(fs * A)
    D = int(fs * D)
    R = int(fs * R)
    gate = int(fs * gate)
    duration = int(fs * duration)
    e = np.zeros(duration)
    if A != 0:
        for n in range(A):
            e[n] = 1.0 - np.exp(-5.0 * n / A)

    if D != 0:
        for n in range(A, gate):
            e[n] = S + (1.0 - S) * np.exp(-5.0 * (n - A) / D)

    else:
        for n in range(A, gate):
            e[n] = S

    if R != 0:
        for n in range(gate, duration):
            e[n]= e[gate - 1] * np.exp(-5.0 * (n - gate + 1) / R)

    return e

(4)

def hihat_cymbal_close(fs, velocity, gate):
    duration = 1

    length_of_s = int(fs * duration)
    s0 = np.zeros(length_of_s)

    np.random.seed(0)
    for n in range(length_of_s):
        s0[n] = (np.random.rand() * 2.0) - 1.0

    VCF_A = [0]
    VCF_D = [0]
    VCF_S = [1]
    VCF_R = [0]
    VCF_offset = [10000]
    VCF_depth=[0]

    vcf = ADSR(fs, VCF_A[0], VCF_D[0], VCF_S[0], VCF_R[0], gate, duration)
    for n in range(length_of_s):
        vcf[n] = VCF_offset[0] + vcf[n] * VCF_depth[0]

    s1 = np.zeros(length_of_s)
    Q = 1 / np.sqrt(2)
    for n in range(length_of_s):
        a, b = HPF(fs, vcf[n], Q)
        for m in range(0, 3):
            if n - m >= 0:
                s1[n] += b[m] * s0[n - m]

        for m in range(1, 3):
            if n - m >= 0:
                s1[n] += -a[m] * s1[n - m]

    VCA_A = [0]
    VCA_D = [0.1]
    VCA_S = [0]
    VCA_R = [0.1]
    VCA_offset = [0]
    VCA_depth = [1]

    vca = ADSR(fs, VCA_A[0], VCA_D[0], VCA_S[0], VCA_R[0], gate, duration)
    for n in range(length_of_s):
        vca[n] = VCA_offset[0] + vca[n] * VCA_depth[0]

    for n in range(length_of_s):
        s1[n] *= vca[n]

    gain = velocity / 127 / np.max(np.abs(s1))
    s1 *= gain

    return s1

(5)

def percussion(fs, note_number, velocity, gate):
    if note_number == 36:
        s = bass_drum(fs, velocity, gate)
    elif note_number == 40:
        s = snare_drum(fs, velocity, gate)
    elif note_number == 42:
        s = hihat_cymbal_close(fs, velocity, gate)

    return s

(6)

fs = 44100

(7)

note_number = 42
velocity = 100
gate = 0.1

(8)

s = percussion(fs, note_number, velocity, gate)
length_of_s = len(s)

(9)

for n in range(length_of_s):
    s[n] = (s[n] + 1.0) / 2.0 * 65536.0
    if s[n] > 65535.0:
        s[n] = 65535.0
    elif s[n] < 0.0:
        s[n] = 0.0;
    s[n] = (s[n] + 0.5) - 32768

wavfile.write('hihat.wav', fs, s.astype(np.int16))

(10)

Audio('hihat.wav')

4．減算合成2（バスドラム）

　「New」ボタンをクリックし，新しくウィンドウを作成しなさい．つづいて，以下のプログラムを順番にセルに貼りつけ，実行しなさい．サイン波と白色雑音からバスドラムからの音をつくることができることを確認しなさい．

(1)

%matplotlib inline
import numpy as np
from scipy.io import wavfile
from IPython.display import display, Audio

(2)

def LPF(fs, fc, Q):
    fc /= fs
    fc = np.tan(np.pi * fc) / (2.0 * np.pi)
    a = np.zeros(3)
    b = np.zeros(3)
    a[0] = 1.0 + 2.0 * np.pi * fc / Q + 4.0 * np.pi * np.pi * fc * fc
    a[1] = (8.0 * np.pi * np.pi * fc * fc - 2.0) / a[0]
    a[2] = (1.0 - 2.0 * np.pi * fc / Q + 4.0 * np.pi * np.pi * fc * fc) / a[0]
    b[0] = 4.0 * np.pi * np.pi * fc * fc / a[0]
    b[1] = 8.0 * np.pi * np.pi * fc * fc / a[0]
    b[2] = 4.0 * np.pi * np.pi * fc * fc / a[0]
    a[0] = 1.0
    return a, b

(3)

def ADSR(fs, A, D, S, R, gate, duration):
    A = int(fs * A)
    D = int(fs * D)
    R = int(fs * R)
    gate = int(fs * gate)
    duration = int(fs * duration)
    e = np.zeros(duration)
    if A != 0:
        for n in range(A):
            e[n] = 1.0 - np.exp(-5.0 * n / A)

    if D != 0:
        for n in range(A, gate):
            e[n] = S + (1.0 - S) * np.exp(-5.0 * (n - A) / D)

    else:
        for n in range(A, gate):
            e[n] = S

    if R != 0:
        for n in range(gate, duration):
            e[n]= e[gate - 1] * np.exp(-5.0 * (n - gate + 1) / R)

    return e

(4)

def compressor(fs, x):
    length_of_s = len(x)

    gain = 1.0 / np.max(np.abs(x))
    x *= gain

    threshold = 0.6
    ratio = 1 / 8
    width = 0.2
    gain = 1 / (threshold + (1.0 - threshold) * ratio)

    for n in range(length_of_s):
        if x[n] < 0:
            sign_of_s = -1
        else:
            sign_of_s = 1

        abs_of_s = np.abs(x[n])

        if abs_of_s >= threshold - width / 2 and abs_of_s < threshold + width / 2:
            abs_of_s = abs_of_s + (ratio - 1) * (abs_of_s - threshold + width / 2)*(abs_of_s - threshold + width / 2) / (width * 2)
        elif abs_of_s >= threshold + width / 2:
            abs_of_s = threshold + (abs_of_s - threshold) * ratio

        x[n] = sign_of_s * abs_of_s * gain

    return x

(5)

def bass_drum(fs, velocity, gate):
    duration = 1

    length_of_s = int(fs * duration)
    sa0 = np.zeros(length_of_s)
    sb0 = np.zeros(length_of_s)

    VCO_A = [0]
    VCO_D = [0.1]
    VCO_S = [0]
    VCO_R = [0.1]
    VCO_offset = [55]
    VCO_depth = [50]

    vco = ADSR(fs, VCO_A[0], VCO_D[0], VCO_S[0], VCO_R[0], gate, duration)
    for n in range(length_of_s):
        vco[n] = VCO_offset[0] + vco[n] * VCO_depth[0]

    x = 0
    for n in range(length_of_s):
        sa0[n] = np.sin(2.0 * np.pi * x)
        delta = vco[n] / fs
        x += delta
        if x >= 1:
            x -= 1

    np.random.seed(0)
    for n in range(length_of_s):
        sb0[n] = (np.random.rand() * 2.0) - 1.0

    VCA_A = [0]
    VCA_D = [0.2]
    VCA_S = [0]
    VCA_R = [0.2]
    VCA_offset = [0]
    VCA_depth = [1]

    vca = ADSR(fs, VCA_A[0], VCA_D[0], VCA_S[0], VCA_R[0], gate, duration)
    for n in range(length_of_s):
        vca[n] = VCA_offset[0] + vca[n] * VCA_depth[0]

    for n in range(length_of_s):
        sb0[n] *= vca[n]

    s0 = np.zeros(length_of_s)
    for n in range(length_of_s):
        s0[n] = sa0[n] * 0.6 + sb0[n] * 0.4

    VCF_A = [0]
    VCF_D = [0.04]
    VCF_S = [0]
    VCF_R = [0.04]
    VCF_offset = [200]
    VCF_depth = [4000]

    vcf = ADSR(fs, VCF_A[0], VCF_D[0], VCF_S[0], VCF_R[0], gate, duration)
    for n in range(length_of_s):
        vcf[n] = VCF_offset[0] + vcf[n] * VCF_depth[0]

    s1 = np.zeros(length_of_s)
    Q = 1 / np.sqrt(2)
    for n in range(length_of_s):
        a, b = LPF(fs, vcf[n], Q)
        for m in range(0, 3):
            if n - m >= 0:
                s1[n] += b[m] * s0[n - m]

        for m in range(1, 3):
            if n - m >= 0:
                s1[n] += -a[m] * s1[n - m]

    VCA_A = [0]
    VCA_D = [0.3]
    VCA_S = [0]
    VCA_R = [0.3]
    VCA_offset = [0]
    VCA_depth = [1]

    vca = ADSR(fs, VCA_A[0], VCA_D[0], VCA_S[0], VCA_R[0], gate, duration)
    for n in range(length_of_s):
        vca[n] = VCA_offset[0] + vca[n] * VCA_depth[0]

    for n in range(length_of_s):
        s1[n] *= vca[n]

    s1 = compressor(fs, s1)

    gain = velocity / 127 / np.max(np.abs(s1))
    s1 *= gain

    return s1

(6)

def percussion(fs, note_number, velocity, gate):
    if note_number == 36:
        s = bass_drum(fs, velocity, gate)
    elif note_number == 40:
        s = snare_drum(fs, velocity, gate)
    elif note_number == 42:
        s = hihat_cymbal_close(fs, velocity, gate)

    return s

(7)

fs = 44100

(8)

note_number = 36
velocity = 100
gate = 0.1

(9)

s = percussion(fs, note_number, velocity, gate)
length_of_s = len(s)

(10)

for n in range(length_of_s):
    s[n] = (s[n] + 1.0) / 2.0 * 65536.0
    if s[n] > 65535.0:
        s[n] = 65535.0
    elif s[n] < 0.0:
        s[n] = 0.0;
    s[n] = (s[n] + 0.5) - 32768

wavfile.write('bass_drum.wav', fs, s.astype(np.int16))

(11)

Audio('bass_drum.wav')

5．減算合成3（スネアドラム）

　「New」ボタンをクリックし，新しくウィンドウを作成しなさい．つづいて，以下のプログラムを順番にセルに貼りつけ，実行しなさい．矩形波と白色雑音からスネアドラムからの音をつくることができることを確認しなさい．

(1)

%matplotlib inline
import numpy as np
from scipy.io import wavfile
from IPython.display import display, Audio

(2)

def LPF(fs, fc, Q):
    fc /= fs
    fc = np.tan(np.pi * fc) / (2.0 * np.pi)
    a = np.zeros(3)
    b = np.zeros(3)
    a[0] = 1.0 + 2.0 * np.pi * fc / Q + 4.0 * np.pi * np.pi * fc * fc
    a[1] = (8.0 * np.pi * np.pi * fc * fc - 2.0) / a[0]
    a[2] = (1.0 - 2.0 * np.pi * fc / Q + 4.0 * np.pi * np.pi * fc * fc) / a[0]
    b[0] = 4.0 * np.pi * np.pi * fc * fc / a[0]
    b[1] = 8.0 * np.pi * np.pi * fc * fc / a[0]
    b[2] = 4.0 * np.pi * np.pi * fc * fc / a[0]
    a[0] = 1.0
    return a, b

(3)

def ADSR(fs, A, D, S, R, gate, duration):
    A = int(fs * A)
    D = int(fs * D)
    R = int(fs * R)
    gate = int(fs * gate)
    duration = int(fs * duration)
    e = np.zeros(duration)
    if A != 0:
        for n in range(A):
            e[n] = 1.0 - np.exp(-5.0 * n / A)

    if D != 0:
        for n in range(A, gate):
            e[n] = S + (1.0 - S) * np.exp(-5.0 * (n - A) / D)

    else:
        for n in range(A, gate):
            e[n] = S

    if R != 0:
        for n in range(gate, duration):
            e[n]= e[gate - 1] * np.exp(-5.0 * (n - gate + 1) / R)

    return e

(4)

def compressor(fs, x):
    length_of_s = len(x)

    gain = 1.0 / np.max(np.abs(x))
    x *= gain

    threshold = 0.6
    ratio = 1 / 8
    width = 0.2
    gain = 1 / (threshold + (1.0 - threshold) * ratio)

    for n in range(length_of_s):
        if x[n] < 0:
            sign_of_s = -1
        else:
            sign_of_s = 1

        abs_of_s = np.abs(x[n])

        if abs_of_s >= threshold - width / 2 and abs_of_s < threshold + width / 2:
            abs_of_s = abs_of_s + (ratio - 1) * (abs_of_s - threshold + width / 2)*(abs_of_s - threshold + width / 2) / (width * 2)
        elif abs_of_s >= threshold + width / 2:
            abs_of_s = threshold + (abs_of_s - threshold) * ratio

        x[n] = sign_of_s * abs_of_s * gain

    return x

(5)

def snare_drum(fs, velocity, gate):
    duration = 1

    length_of_s = int(fs * duration)
    sa0 = np.zeros(length_of_s)
    sb0 = np.zeros(length_of_s)

    VCO_A = [0]
    VCO_D = [0]
    VCO_S = [1]
    VCO_R = [0]
    VCO_offset = [150]
    VCO_depth = [0]

    vco = ADSR(fs, VCO_A[0], VCO_D[0], VCO_S[0], VCO_R[0], duration, duration)
    for n in range(length_of_s):
        vco[n] = VCO_offset[0] + vco[n] * VCO_depth[0]

    x = 0
    for n in range(length_of_s):
        if x < 0.5:
            sa0[n] = 1
        else:
            sa0[n] = -1

        delta = vco[n] / fs

        if 1 - delta <= x and x < 1:
            t = (x - 1) / delta
            d = t * t + 2 * t + 1
            sa0[n] += d
        elif 0 <= x and x < delta:
            t = x / delta
            d = -t * t + 2 * t - 1
            sa0[n] += d

        if 0.5 - delta <= x and x < 0.5:
            t = (x - 0.5) / delta
            d = t * t + 2 * t + 1
            sa0[n] -= d
        elif 0.5 <= x and x < 0.5 + delta:
            t = (x - 0.5) / delta
            d = -t * t + 2 * t - 1
            sa0[n] -= d

        x += delta
        if x >= 1:
            x -= 1

    np.random.seed(0)
    for n in range(length_of_s):
        sb0[n] = (np.random.rand() * 2.0) - 1.0

    VCA_A = [0]
    VCA_D = [0]
    VCA_S = [1]
    VCA_R = [0]
    VCA_offset = [0]
    VCA_depth = [1]

    vca = ADSR(fs, VCA_A[0], VCA_D[0], VCA_S[0], VCA_R[0], duration, duration)
    for n in range(length_of_s):
        vca[n] = VCA_offset[0] + vca[n] * VCA_depth[0]

    for n in range(length_of_s):
        sb0[n] *= vca[n]

    s0 = np.zeros(length_of_s)
    for n in range(length_of_s):
        s0[n] = sa0[n] * 0.3 + sb0[n] * 0.7

    VCF_A = [0]
    VCF_D = [0.1]
    VCF_S = [0]
    VCF_R = [0.1]
    VCF_offset = [8000]
    VCF_depth = [-7800]

    vcf = ADSR(fs, VCF_A[0], VCF_D[0], VCF_S[0], VCF_R[0], gate, duration)
    for n in range(length_of_s):
        vcf[n] = VCF_offset[0] + vcf[n] * VCF_depth[0]

    s1 = np.zeros(length_of_s)
    Q = 1 / np.sqrt(2)
    for n in range(length_of_s):
        a, b = LPF(fs, vcf[n], Q)
        for m in range(0, 3):
            if n - m >= 0:
                s1[n] += b[m] * s0[n - m]

        for m in range(1, 3):
            if n - m >= 0:
                s1[n] += -a[m] * s1[n - m]

    VCA_A = [0]
    VCA_D = [0.2]
    VCA_S = [0]
    VCA_R = [0.2]
    VCA_offset = [0]
    VCA_depth = [1]

    vca = ADSR(fs, VCA_A[0], VCA_D[0], VCA_S[0], VCA_R[0], gate, duration)
    for n in range(length_of_s):
        vca[n] = VCA_offset[0] + vca[n] * VCA_depth[0]

    for n in range(length_of_s):
        s1[n] *= vca[n]

    s1 = compressor(fs, s1)

    gain = velocity / 127 / np.max(np.abs(s1))
    s1 *= gain

    return s1

(6)

def percussion(fs, note_number, velocity, gate):
    if note_number == 36:
        s = bass_drum(fs, velocity, gate)
    elif note_number == 40:
        s = snare_drum(fs, velocity, gate)
    elif note_number == 42:
        s = hihat_cymbal_close(fs, velocity, gate)

    return s

(7)

fs = 44100

(8)

note_number = 40
velocity = 100
gate = 0.1

(9)

s = percussion(fs, note_number, velocity, gate)
length_of_s = len(s)

(10)

for n in range(length_of_s):
    s[n] = (s[n] + 1.0) / 2.0 * 65536.0
    if s[n] > 65535.0:
        s[n] = 65535.0
    elif s[n] < 0.0:
        s[n] = 0.0;
    s[n] = (s[n] + 0.5) - 32768

wavfile.write('snare_drum.wav', fs, s.astype(np.int16))

(11)

Audio('snare_drum.wav')

6．自動演奏

　「New」ボタンをクリックし，新しくウィンドウを作成しなさい．つづいて，以下のプログラムを順番にセルに貼りつけ，実行しなさい．ドラムセットの演奏が聞こえてくることを確かめてください．

(1)

%matplotlib inline
import numpy as np
from scipy.io import wavfile
from IPython.display import display, Audio

(2)

def LPF(fs, fc, Q):
    fc /= fs
    fc = np.tan(np.pi * fc) / (2.0 * np.pi)
    a = np.zeros(3)
    b = np.zeros(3)
    a[0] = 1.0 + 2.0 * np.pi * fc / Q + 4.0 * np.pi * np.pi * fc * fc
    a[1] = (8.0 * np.pi * np.pi * fc * fc - 2.0) / a[0]
    a[2] = (1.0 - 2.0 * np.pi * fc / Q + 4.0 * np.pi * np.pi * fc * fc) / a[0]
    b[0] = 4.0 * np.pi * np.pi * fc * fc / a[0]
    b[1] = 8.0 * np.pi * np.pi * fc * fc / a[0]
    b[2] = 4.0 * np.pi * np.pi * fc * fc / a[0]
    a[0] = 1.0
    return a, b

(3)

def HPF(fs, fc, Q):
    fc /= fs
    fc = np.tan(np.pi * fc) / (2.0 * np.pi)
    a = np.zeros(3)
    b = np.zeros(3)
    a[0] = 1.0 + 2.0 * np.pi * fc / Q + 4.0 * np.pi * np.pi * fc * fc;
    a[1] = (8.0 * np.pi * np.pi * fc * fc - 2.0) / a[0]
    a[2] = (1.0 - 2.0 * np.pi * fc / Q + 4.0 * np.pi * np.pi * fc * fc) / a[0]
    b[0] = 1.0 / a[0]
    b[1] = -2.0 / a[0]
    b[2] = 1.0 / a[0]
    a[0] = 1.0
    return a, b

(4)

def ADSR(fs, A, D, S, R, gate, duration):
    A = int(fs * A)
    D = int(fs * D)
    R = int(fs * R)
    gate = int(fs * gate)
    duration = int(fs * duration)
    e = np.zeros(duration)
    if A != 0:
        for n in range(A):
            e[n] = 1.0 - np.exp(-5.0 * n / A)

    if D != 0:
        for n in range(A, gate):
            e[n] = S + (1.0 - S) * np.exp(-5.0 * (n - A) / D)

    else:
        for n in range(A, gate):
            e[n] = S

    if R != 0:
        for n in range(gate, duration):
            e[n]= e[gate - 1] * np.exp(-5.0 * (n - gate + 1) / R)

    return e

(5)

def compressor(fs, x):
    length_of_s = len(x)

    gain = 1.0 / np.max(np.abs(x))
    x *= gain

    threshold = 0.6
    ratio = 1 / 8
    width = 0.2
    gain = 1 / (threshold + (1.0 - threshold) * ratio)

    for n in range(length_of_s):
        if x[n] < 0:
            sign_of_s = -1
        else:
            sign_of_s = 1

        abs_of_s = np.abs(x[n])

        if abs_of_s >= threshold - width / 2 and abs_of_s < threshold + width / 2:
            abs_of_s = abs_of_s + (ratio - 1) * (abs_of_s - threshold + width / 2)*(abs_of_s - threshold + width / 2) / (width * 2)
        elif abs_of_s >= threshold + width / 2:
            abs_of_s = threshold + (abs_of_s - threshold) * ratio

        x[n] = sign_of_s * abs_of_s * gain

    return x

(6)

def hihat_cymbal_close(fs, velocity, gate):
    duration = 1

    length_of_s = int(fs * duration)
    s0 = np.zeros(length_of_s)

    np.random.seed(0)
    for n in range(length_of_s):
        s0[n] = (np.random.rand() * 2.0) - 1.0

    VCF_A = [0]
    VCF_D = [0]
    VCF_S = [1]
    VCF_R = [0]
    VCF_offset = [10000]
    VCF_depth=[0]

    vcf = ADSR(fs, VCF_A[0], VCF_D[0], VCF_S[0], VCF_R[0], gate, duration)
    for n in range(length_of_s):
        vcf[n] = VCF_offset[0] + vcf[n] * VCF_depth[0]

    s1 = np.zeros(length_of_s)
    Q = 1 / np.sqrt(2)
    for n in range(length_of_s):
        a, b = HPF(fs, vcf[n], Q)
        for m in range(0, 3):
            if n - m >= 0:
                s1[n] += b[m] * s0[n - m]

        for m in range(1, 3):
            if n - m >= 0:
                s1[n] += -a[m] * s1[n - m]

    VCA_A = [0]
    VCA_D = [0.1]
    VCA_S = [0]
    VCA_R = [0.1]
    VCA_offset = [0]
    VCA_depth = [1]

    vca = ADSR(fs, VCA_A[0], VCA_D[0], VCA_S[0], VCA_R[0], gate, duration)
    for n in range(length_of_s):
        vca[n] = VCA_offset[0] + vca[n] * VCA_depth[0]

    for n in range(length_of_s):
        s1[n] *= vca[n]

    gain = velocity / 127 / np.max(np.abs(s1))
    s1 *= gain

    return s1

(7)

def bass_drum(fs, velocity, gate):
    duration = 1

    length_of_s = int(fs * duration)
    sa0 = np.zeros(length_of_s)
    sb0 = np.zeros(length_of_s)

    VCO_A = [0]
    VCO_D = [0.1]
    VCO_S = [0]
    VCO_R = [0.1]
    VCO_offset = [55]
    VCO_depth = [50]

    vco = ADSR(fs, VCO_A[0], VCO_D[0], VCO_S[0], VCO_R[0], gate, duration)
    for n in range(length_of_s):
        vco[n] = VCO_offset[0] + vco[n] * VCO_depth[0]

    x = 0
    for n in range(length_of_s):
        sa0[n] = np.sin(2.0 * np.pi * x)
        delta = vco[n] / fs
        x += delta
        if x >= 1:
            x -= 1

    np.random.seed(0)
    for n in range(length_of_s):
        sb0[n] = (np.random.rand() * 2.0) - 1.0

    VCA_A = [0]
    VCA_D = [0.2]
    VCA_S = [0]
    VCA_R = [0.2]
    VCA_offset = [0]
    VCA_depth = [1]

    vca = ADSR(fs, VCA_A[0], VCA_D[0], VCA_S[0], VCA_R[0], gate, duration)
    for n in range(length_of_s):
        vca[n] = VCA_offset[0] + vca[n] * VCA_depth[0]

    for n in range(length_of_s):
        sb0[n] *= vca[n]

    s0 = np.zeros(length_of_s)
    for n in range(length_of_s):
        s0[n] = sa0[n] * 0.6 + sb0[n] * 0.4

    VCF_A = [0]
    VCF_D = [0.04]
    VCF_S = [0]
    VCF_R = [0.04]
    VCF_offset = [200]
    VCF_depth = [4000]

    vcf = ADSR(fs, VCF_A[0], VCF_D[0], VCF_S[0], VCF_R[0], gate, duration)
    for n in range(length_of_s):
        vcf[n] = VCF_offset[0] + vcf[n] * VCF_depth[0]

    s1 = np.zeros(length_of_s)
    Q = 1 / np.sqrt(2)
    for n in range(length_of_s):
        a, b = LPF(fs, vcf[n], Q)
        for m in range(0, 3):
            if n - m >= 0:
                s1[n] += b[m] * s0[n - m]

        for m in range(1, 3):
            if n - m >= 0:
                s1[n] += -a[m] * s1[n - m]

    VCA_A = [0]
    VCA_D = [0.3]
    VCA_S = [0]
    VCA_R = [0.3]
    VCA_offset = [0]
    VCA_depth = [1]

    vca = ADSR(fs, VCA_A[0], VCA_D[0], VCA_S[0], VCA_R[0], gate, duration)
    for n in range(length_of_s):
        vca[n] = VCA_offset[0] + vca[n] * VCA_depth[0]

    for n in range(length_of_s):
        s1[n] *= vca[n]

    s1 = compressor(fs, s1)

    gain = velocity / 127 / np.max(np.abs(s1))
    s1 *= gain

    return s1

(8)

def snare_drum(fs, velocity, gate):
    duration = 1

    length_of_s = int(fs * duration)
    sa0 = np.zeros(length_of_s)
    sb0 = np.zeros(length_of_s)

    VCO_A = [0]
    VCO_D = [0]
    VCO_S = [1]
    VCO_R = [0]
    VCO_offset = [150]
    VCO_depth = [0]

    vco = ADSR(fs, VCO_A[0], VCO_D[0], VCO_S[0], VCO_R[0], duration, duration)
    for n in range(length_of_s):
        vco[n] = VCO_offset[0] + vco[n] * VCO_depth[0]

    x = 0
    for n in range(length_of_s):
        if x < 0.5:
            sa0[n] = 1
        else:
            sa0[n] = -1

        delta = vco[n] / fs

        if 1 - delta <= x and x < 1:
            t = (x - 1) / delta
            d = t * t + 2 * t + 1
            sa0[n] += d
        elif 0 <= x and x < delta:
            t = x / delta
            d = -t * t + 2 * t - 1
            sa0[n] += d

        if 0.5 - delta <= x and x < 0.5:
            t = (x - 0.5) / delta
            d = t * t + 2 * t + 1
            sa0[n] -= d
        elif 0.5 <= x and x < 0.5 + delta:
            t = (x - 0.5) / delta
            d = -t * t + 2 * t - 1
            sa0[n] -= d

        x += delta
        if x >= 1:
            x -= 1

    np.random.seed(0)
    for n in range(length_of_s):
        sb0[n] = (np.random.rand() * 2.0) - 1.0

    VCA_A = [0]
    VCA_D = [0]
    VCA_S = [1]
    VCA_R = [0]
    VCA_offset = [0]
    VCA_depth = [1]

    vca = ADSR(fs, VCA_A[0], VCA_D[0], VCA_S[0], VCA_R[0], duration, duration)
    for n in range(length_of_s):
        vca[n] = VCA_offset[0] + vca[n] * VCA_depth[0]

    for n in range(length_of_s):
        sb0[n] *= vca[n]

    s0 = np.zeros(length_of_s)
    for n in range(length_of_s):
        s0[n] = sa0[n] * 0.3 + sb0[n] * 0.7

    VCF_A = [0]
    VCF_D = [0.1]
    VCF_S = [0]
    VCF_R = [0.1]
    VCF_offset = [8000]
    VCF_depth = [-7800]

    vcf = ADSR(fs, VCF_A[0], VCF_D[0], VCF_S[0], VCF_R[0], gate, duration)
    for n in range(length_of_s):
        vcf[n] = VCF_offset[0] + vcf[n] * VCF_depth[0]

    s1 = np.zeros(length_of_s)
    Q = 1 / np.sqrt(2)
    for n in range(length_of_s):
        a, b = LPF(fs, vcf[n], Q)
        for m in range(0, 3):
            if n - m >= 0:
                s1[n] += b[m] * s0[n - m]

        for m in range(1, 3):
            if n - m >= 0:
                s1[n] += -a[m] * s1[n - m]

    VCA_A = [0]
    VCA_D = [0.2]
    VCA_S = [0]
    VCA_R = [0.2]
    VCA_offset = [0]
    VCA_depth = [1]

    vca = ADSR(fs, VCA_A[0], VCA_D[0], VCA_S[0], VCA_R[0], gate, duration)
    for n in range(length_of_s):
        vca[n] = VCA_offset[0] + vca[n] * VCA_depth[0]

    for n in range(length_of_s):
        s1[n] *= vca[n]

    s1 = compressor(fs, s1)

    gain = velocity / 127 / np.max(np.abs(s1))
    s1 *= gain

    return s1

(9)

def percussion(fs, note_number, velocity, gate):
    if note_number == 36:
        s = bass_drum(fs, velocity, gate)
    elif note_number == 40:
        s = snare_drum(fs, velocity, gate)
    elif note_number == 42:
        s = hihat_cymbal_close(fs, velocity, gate)

    return s

(10)

score = np.array([[1, 1920, 42, 100, 96],
                  [1, 2160, 42, 100, 96],
                  [1, 2400, 42, 100, 96],
                  [1, 2640, 42, 100, 96],
                  [1, 2880, 42, 100, 96],
                  [1, 3120, 42, 100, 96],
                  [1, 3360, 42, 100, 96],
                  [1, 3600, 42, 100, 96],
                  [1, 3840, 42, 100, 96],
                  [1, 4080, 42, 100, 96],
                  [1, 4320, 42, 100, 96],
                  [1, 4560, 42, 100, 96],
                  [1, 4800, 42, 100, 96],
                  [1, 5040, 42, 100, 96],
                  [1, 5280, 42, 100, 96],
                  [1, 5520, 42, 100, 96],
                  [2, 1920, 36, 100, 96],
                  [2, 2880, 36, 100, 96],
                  [2, 3840, 36, 100, 96],
                  [2, 4800, 36, 100, 96],
                  [3, 2400, 40, 100, 96],
                  [3, 3360, 40, 100, 96],
                  [3, 4320, 40, 100, 96],
                  [3, 5280, 40, 100, 96]])

division = 480
tempo = 120
number_of_track = 3
end_of_track = 6
number_of_note = score.shape[0]

(11)

fs = 44100
length_of_s = int(fs * (end_of_track + 2))
track = np.zeros((length_of_s, number_of_track))
s = np.zeros(length_of_s)

(12)

for i in range(number_of_note):
    j = int(score[i, 0] - 1)
    onset = (score[i, 1] / division) * (60 / tempo)
    note_number = score[i, 2]
    velocity = score[i, 3]
    gate = (score[i, 4] / division) * (60 / tempo)
    x = percussion(fs, note_number, velocity, gate)
    offset = int(fs * onset)
    length_of_x = len(x)
    for n in range(length_of_x):
        track[offset + n, j] += x[n]

(13)

for j in range(number_of_track):
    for n in range(length_of_s):
        s[n] += track[n, j]

(14)

master_volume = 0.5
s /= np.max(np.abs(s))
s *= master_volume

(15)

for n in range(length_of_s):
    s[n] = (s[n] + 1.0) / 2.0 * 65536.0
    if s[n] > 65535.0:
        s[n] = 65535.0
    elif s[n] < 0.0:
        s[n] = 0.0;
    s[n] = (s[n] + 0.5) - 32768

wavfile.write('p5.wav', fs, s.astype(np.int16))

(16)

Audio('p5.wav')

7．レポートについて

　下記の課題について，レポートを作成しなさい．

（課題1） ウェブサイトなどから楽器音の音データを入手し，波形，周波数特性，スペクトログラムを表示しなさい．これらをお手本とすることで，楽器音を合成するプログラムをつくりなさい．どのような手順で音をつくったのか，自分のプログラムについて説明してください．

（課題2） 自分でつくった楽器音を使ってカノン以外の音楽データを作成し、そのWAVEファイルをメールに添付して提出してください．

8．創成課題

1．自分でつくった楽器音を使ってカノン以外の音楽データを作成しなさい．少なくともふたつ以上の楽器音をつくり，音楽データを作成すること．

2．創成課題の発表会に向けて，プレゼンテーションのための資料を作成しなさい．かならず，音楽データのデモをまじえること．どのような手順で音をつくったのか，具体的な説明があるものを高く評価します．

Last Modified: April 1 12:00 JST 2021 by Naofumi Aoki
E-mail: aoki@ime.ist.hokudai.ac.jp