Python 初級者の学生さんに課した地獄のレッスン記録

**【登場人物】
M:　(＝メンター＝京都のわし)　東京から見えてるかな？
Q:　(＝東京の学生さん＝Python初級講座みたいなやつの修了者)**

```python
print("見えた")
```
    見えた
    
M: Pythonはレッスン動画とか参考書とかいっぱいあるけど、ちょっとつまづくと先に進めない（あるある）ので、それらを利用するのはいいけど、この場所を併用するといいかと。 

**ナレーション：
これは、Google colaboratory の共有機能を使って実施した Python プログラミングの地獄のマンツーマンレッスンの記録である。目的は基本的にはプログラミングの習得なんだが、ネタとして音楽の音階や和音に関するディープなウンチクが展開される。**

M:　例題：　y = sin(x) （0 ≦ x ≦ 4π）のグラフは、描けますか？

本質的じゃない課題なので、ネタばらし、しちゃいます。（一つの方法として）


```python
import numpy as np
# numpy は math よりちょっと賢い数値計算ライブラリで、行列演算とかもできる。
# numpy は長い名前なので、以後、使いやすいように np とする。
import matplotlib.pyplot as plt
# colab には、ita が最初から入ってないので、matplotlib というグラフィックライブラリを使用

x = np.arange(0, 4 * np.pi, 0.1)  # 0 から 4π　まで　0.1 刻みのデータ（配列みたいなもの）を得る。
y = np.sin(x)
plt.plot(x, y)
```
![output_7_1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/216556/54ab3e5b-89bd-8235-07bf-bb69ae2ec0e6.png)

```python
# Q:
print(x)
```

    [ 0.   0.1  0.2  0.3  0.4  0.5  0.6  0.7  0.8  0.9  1.   1.1  1.2  1.3
      1.4  1.5  1.6  1.7  1.8  1.9  2.   2.1  2.2  2.3  2.4  2.5  2.6  2.7
      2.8  2.9  3.   3.1  3.2  3.3  3.4  3.5  3.6  3.7  3.8  3.9  4.   4.1
      4.2  4.3  4.4  4.5  4.6  4.7  4.8  4.9  5.   5.1  5.2  5.3  5.4  5.5
      5.6  5.7  5.8  5.9  6.   6.1  6.2  6.3  6.4  6.5  6.6  6.7  6.8  6.9
      7.   7.1  7.2  7.3  7.4  7.5  7.6  7.7  7.8  7.9  8.   8.1  8.2  8.3
      8.4  8.5  8.6  8.7  8.8  8.9  9.   9.1  9.2  9.3  9.4  9.5  9.6  9.7
      9.8  9.9 10.  10.1 10.2 10.3 10.4 10.5 10.6 10.7 10.8 10.9 11.  11.1
     11.2 11.3 11.4 11.5 11.6 11.7 11.8 11.9 12.  12.1 12.2 12.3 12.4 12.5]
<font color="Blue">
Q: ↑このグラフでｘはとびとびなのに、グラフが滑らかにつながるのはどうして？
</font>

M: 【解答0001】=> Good!
matplotlib.pyplot.plot() は、じつは滑らかなグラフではなくて、
座標と座標を順次、線分でつなぐ、つまり「折れ線グラフ」を描くメソッドようじゃな。
x のキザミ（tick）を大きく(1とかに) してみるとこんなかんじ。
↓

```python
 x = np.arange(0, 4 * np.pi, 1)
 y = np.sin(x)
 print(x)
 plt.plot(x, y)
```
![kakukaku.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/216556/e11637b3-622c-9c99-474b-240f9d7f8e3a.png)
<font color="Blue">
Q:　なるほど。
</font>

M:【課題0001】上記の sin 関数のグラフを
修正して、周期を半分にしてみましょう。
（定義域は変えない）

```python
#　Q:【解答0001】
x = np.arange(0, 4 * np.pi, 0.1)
y = np.sin(2 * x)
plt.plot(x, y)
```
![output_11_1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/216556/81225e63-2fc3-d6df-e358-ccf7c07d3125.png)

M:【課題0002】課題0001の元のグラフに、周期・振幅ともに半分にしたデータを加算してしみましょう。 

```python
#　Q:【解答0002】(データの加算とは？重ね合わせということかな)
x = np.arange(0, 4 * np.pi, 0.1) 
y1 = np.sin(x)
y2 = 0.5 * (np.sin(2 * x))

y = y1 + y2
plt.plot(x, y)
```

![output_15_1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/216556/8c307634-b96b-4823-02f2-2c13f290a6a3.png)


M:【課題0003】課題0001の元のグラフに、
周期・振幅ともに 1/i にした複数の
データ系列を加算してしみましょう。 
ただし、i は 2 から 9 までの整数とします。 
それができたら、i を 99 までにしてみましょう。

```python
#　M:【解答0003】①
# いきなりヒント
import numpy as np
import matplotlib.pyplot as plt 
x = np.arange(0, 4 * np.pi, 0.1) 
y = np.zeros(len(x)) 
# y は x と同じ長さのオールゼロの系列
# len は length の略
print(y) 
# ここで y に 9 個の系列を
# 順番に足していけな良いのじゃな。
#（いくつかの方法あり）

for i in range(1,10):
  y_temp = (np.sin(i * x)) / i
  y += y_temp


plt.plot(x, y)
```

    [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.
     0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.
     0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.
     0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.
     0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.
     0. 0. 0. 0. 0. 0.]
    ![output_18_1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/216556/a4797ca9-5a0f-e2c9-9bdb-2487cfec49b9.png)

![output_17_2.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/216556/4aa965d1-c4a9-0401-c9ce-71af506ab79f.png)

```python
#　Q:【解答0003】②
y2 = np.zeros(len(x))　# ＝後で追加した行
for i in range(1,100):
  y_temp = (np.sin(i * x)) / i
  y2 += y_temp

plt.plot(x, y2)
```

![output_18_1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/216556/9be72327-f367-5e55-e288-e279bcf06415.png)

M:　ふむ、素晴らしい！
あ、ちょっと待った！
②は惜しい！
ループの前にyをゼロで初期化しないと、①の結果を引き継いでしまうぞ！
だって、振幅大きすぎ！
ところで、y_temp は使わなくても出来るんではないかな？
というわけで②だけ修正してみよう。

で、それが出来たら、
【課題0004】課題0003②で、iを奇数だけでやってみよう！
<font color="Blue">
Q:　ほんまやね、②訂正しました！
確かにy_tempを使わなくてもできますね！
フーリエ変換の話につながるのかな。物理か数学でやったことある！
</font>

```python
#　Q:【解答0004】
y_0004 = np.zeros(len(x))
for i in range(1,100, 2):
  y_0004 += (np.sin(i * x)) / i

plt.plot(x, y_0004)
```
![output_21_1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/216556/d289c817-f0b6-f3db-b677-9d101c7d733b.png)

M:　物理、数学、そして音楽につながるのじゃ。  
ルート音をドとすると、２倍音はオクターブ上のド、３倍音はその上のソ、４倍音は２オクターブ上のド、５倍音はその上のミ、つまりCメジャーコードの構成音となる。（ただしド以外は近似値）  
また、ギターで１２フレットの位置でハーモニクス奏法すると２倍音が出て、大体７フレットの位置だと３倍音、５フレットの位置だと４倍音が出る。（なんとかいう芸人がやってたのは７フレットだったか。）  
実はギターの音は鋸歯状波（課題0003）に近いので、普通に開放弦を弾くと整数倍の倍音を含んで鳴っている。  
鋸歯状波は、ギター以外に、バイオリンやトランペットなど多くの楽器の倍音構成に近い。 それに対して、矩形波（課題0004）、つまり奇数倍音を多く持つ波形は、リコーダーやフルートやオーボエなど、木管楽器系の音色となるのじゃ。 
奇数倍音を持つ波形は矩形波だけじゃなくて三角波などもある。 
https://ja.wikipedia.org/wiki/%E4%B8%89%E8%A7%92%E6%B3%A2_(%E6%B3%A2%E5%BD%A2)  
一般に奇数倍音だけを持つ波形は周期の半分が左右対称になるのじゃな。 

M:【課題0005】ギターのナットの x 座標を 0 、ブリッジの x 座標を 1 とすると、 各フレットの x 座標は、 
 
*  12 フレット： 1/2 , 
*  24 フレット： 1 - ((1/2)の2乗) = 3/4 (24フレットは無いけど) 
*   1 フレット： 1 - ((1/2)の12乗根)  
*   2 フレット： 1 - (((1/2)の12乗根)の2乗)  
 
などとなります。（・・ホンマかいな？） 
 
次のプログラムは、ギターのナットから12フレットまでの位置をプロットしようという意図(y 座標は 0 に固定)だったのですが、 
a と x の求め方が間違っていて、フレットが等間隔になってしまっています。  
 
①本物らしいフレットの位置に直せますか？
 
②７フレット付近でやってたハーモニクス奏法は、本当は微妙にナット側かブリッジ側か、どちらにずらすべきでしょうか？


```python
# M:
import numpy as np
import matplotlib.pyplot as plt
 
a = 1 / 12
n = 13
x = [1 - a * i for i in range(n)] 
y = [0] * n 
print(x)
plt.plot(x, y, "o")
```

    [1.0, 0.9166666666666666, 0.8333333333333334, 0.75, 0.6666666666666667, 0.5833333333333334, 0.5, 0.41666666666666674, 0.33333333333333337, 0.25, 0.16666666666666674, 0.08333333333333337, 0.0]

![output_24_2.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/216556/0d8a1b43-fe62-3b74-b3d9-9743922a143e.png)

```python
# Q:
import numpy as np
import matplotlib.pyplot as plt
 
a = (1/2)**(1 / 12)
n = 13
x = [1 - a ** i for i in range(n)] 
y = [0] * n 
print(x)
plt.plot(x, y, "o")
```

    [0.0, 0.05612568731830647, 0.10910128185966061, 0.1591035847462854, 0.2062994740159002, 0.25084646156165913, 0.2928932188134523, 0.3325800729149827, 0.3700394750525633, 0.40539644249863926, 0.43876897584531327, 0.4702684528203521, 0.4999999999999998]

![output_25_2.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/216556/7ea79368-19c5-a841-4895-4bab3b635f74.png)
<font color="Blue">
s: 3倍音がでるのは, x = 0.333...
対して7フレットの位置は x = 0.3325...
よって、微妙にブリッジ側にずらすべき
</font>
M: たいへんよくできました！ 
さて、なぜ７フレットが弦のちょうど 1/3 の長さではないのか？ 

答えは、「そうなる理由がないから。」 

7 フレットの位置 1-(1/2)**(7/12) は、おそらく無理数になる（つまり　1/3 とかにはならない）んじゃないかな？ 

たぶん1/3 になる証明はできないでしょ。 

歴史的には、バッハさんとかの時代？に楽曲を自由に転調するために、音階の周波数を等比数列として設計する必要があった。 

つまり隣り合う周波数の比は一定になるようにした。
　ド ： ド# = ド# : レ = レ : レ# = ・・・
で、たまたま１オクターブを１２個の等比数列に分ける（これを平均律という）と、ところどころ、気持ちいい音程があることに気づいた。 

例えばドミソは一緒に鳴らすと何となく気持ちいいい（調和している）
よく調べると、ドミソの周波数の比が「だいたい」4:5:6 になっていた、というわけじゃ。 

それに対して、無理やりドミソを 4:5:6 とかキッチリした比率（つまり有理数かな）とするような音階の設計もあり、これを「純正調」という。(ルーツはピタゴラスさんとか。)   ただ、純正調では転調する（あるいはカポとかつける）とチューニングし直す必要がある、のでほとんど使われない。

M:【課題0006】周波数, 振幅, 秒数, データの最小間隔, それぞれを与えて任意のサイン波のデータを生成する関数 sin_wave を実装してみましょう。  
（ここでいう関数は、x から y を求める数学の関数ではないよ。）

```python
# M:
import numpy as np
import matplotlib.pyplot as plt

def sin_wave(frequency, amplitude, duration, tick):
    # sine_wave: 任意のサイン波のデータを生成する関数
    #   frequency: 周波数(Hz),  
    #   amplitude: 振幅(V, すなわちボルト), 
    #   duration : データの秒数, 
    #   tick     : データの最小間隔(秒), 

    # ここに実装する。
    # s:
    x = np.arange(0, duration, tick)
    y = amplitude * (np.sin(2 * np.pi * x * frequency))

    return x, y

x, y = sin_wave(frequency=10, amplitude=3, duration=2, tick=0.01)
plt.plot(x, y)

```

![output_29_1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/216556/a10c2084-fc3c-7522-0d2e-ee9ee0abe67f.png)

M: 課題0006=大体Good! 

 一般に、結果が定数となる計算はループの外に出したほうがよい。
たとえば 2 * np.pi * x[i] * frequency を i を回す for ループの中で使いたい場合、


```python
npf = 2 * np.pi *  frequency
for i in range(0,duration,tick):
    y[i] = amplitude * (np.sin(npf * x[i]))
```

などとする方が高速になる。 

ただし、今回のように for を使わず numpy のベクトル演算を使うような場合、numpy の実装がかしこければ速度はほとんど変わらないかもしれない。

M:【課題0007】sin_wave関数を使って、Aのメジャーコードの基本形の３つの構成音とそれらを重ね合わせた合計４本の系列をそれぞれ別の色のグラフ（plot の colorパラメータをそれぞれ "r","g","b","k" とする）にして、ひとつの座標系に重ねて描きこんでみましょう。 

ただし、時間は0.02秒、データの最小間隔は１万分の１秒、Aの基音の周波数は220Hz、構成音の振幅は全て１Vとします。

<font color="Blue">Q: 勘でやった</font>

```python
# Q:
x, y1 = sin_wave(frequency=220, amplitude=1, duration=0.02, tick=0.0001)
plt.plot(x, y1, color="r")
x, y2 = sin_wave(frequency=220*2**(1/3), amplitude=1, duration=0.02, tick=0.0001)
plt.plot(x, y2, color="g")
x, y3 = sin_wave(frequency=220*2**(7/12), amplitude=1, duration=0.02, tick=0.0001)
plt.plot(x, y3, color="b")

plt.plot(x, y1+y2+y3, color="k")
```
![output_34_1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/216556/e5db2943-a08f-64e9-e3cd-c5eaa08a47c1.png)

M: なるほど、アリですね。・・平均律か。（そりゃそうかw）  
(2\*\*(1/12))\*\*4 が 2\*\*(1/3) になるとは気が付かんかった・・（文系爆）  
次の解答例はメモリーを節約する必要がある場合。  
参考までにさらっと。

```python
# M: 解答例
plt.figure(figsize=(40, 10), dpi=72)
root = 2 ** (1 / 12)
for i, s in enumerate([[0, "r"], [4, "g"], [7, "b"]]):
    x, y = sin_wave(220 * (root ** s[0]), 1, 0.05, 0.0001)
    plt.plot(x, y, color=s[1], lw=2)
    y_sum = y if i == 0 else y_sum + y
plt.plot(x, y_sum, color="k", lw=5)
plt.plot([0, 0.05], [0, 0], color="k")
```
![output_36_1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/216556/9a42d427-204a-4208-2a2a-920dd73c8650.png)

M: 【課題0008】さて、課題0007 は平均律だったけど、純正調だと、どうなるでしょう？  
せっかくなので、上の解答例をベースでやってみては？  
あと、x 軸 ( y=0 の線) も描いてみて。

```python
# Q:
plt.figure(figsize=(40, 10), dpi=72)
for i, s in enumerate([[1, "r"], [1.25, "g"], [1.5, "b"]]):
    x, y = sin_wave(220 * (s[0]), 1, 0.05, 0.0001)
    plt.plot(x, y, color=s[1], lw=2)
    y_sum = y if i == 0 else y_sum + y
plt.plot(x, y_sum, color="k", lw=5)
#y = [0] * len(x)
#plt.plot(x, y, color="k")
plt.plot([0, 0.05], [0, 0], color="k")
```
![output_38_1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/216556/7c6ea91a-73d7-e0f4-0ee8-73e3bebd8b71.png)

M: 上の２つのグラフの表示サイズを変えてみた。 
（PCで見るほうがわかりやすいかも） 

下の純正調だと各系列がちょうど y=0 でクロスしてる点が２つあるのがわかるね。 

あと、x軸だけど、１本の線を引くだけだから、メモリーを節約するために２つの座標だけにしたよ。 
（そもそも x 軸を引くメソッドがあればいいんだけど、未確認）


で、綺麗にハモるということは、ゼロクロスすることが重要ではなくて「合成した波形が短い周期で同じパターンを繰り返す」ことだと思う。 

時間を0.05秒に伸ばしてみたら、後者の純正調のほうは全く同じパターンを繰り返してるのに対して、前者の平均律では。次第にパターンが変わっていくように見える。 

（周期の最小公倍数の時間までいくと、同じ波形を繰り返すんだね。）

さて、お待ちかねのオクターブの話題。 なぜ、１オクターブ間隔ありきで、その間が周波数の等比数列となる音階を設計することになったのか？ 逆に言うと、たとえば、２オクターブを単位としてその中を７個の音階に分けたりしてははいけなかったのか？
さらには、たとえば１Hzから1KHzまでの間を周波数の「等差」数列で音階を作ってはいけなかったのか？ それは音色と音程の関係にある。 鋸歯状波でみたように、サイン波以外の楽音（音程のある音）にはふつう、基音の整数倍の周波数の「倍音」が含まれ、これが「音色」を形成する。 この倍音の最低の高さが普通、２倍の周波数、すなわち１オクターブ上である。 で、課題0002 では、基音と１オクターブ上のサイン波を合成して、「音色」を調整した。この高いほうの音の比率をさらに上げていくと、だんだん１オクターブ上の音に近づく。 つまりドからレなどのように、だんだんと周波数を上げるのではなく、倍音の比率を変えることにより、いつの間にか１オクターブ上の音になっていくのである。 すなわちオクターブ違いの音というのは、音色が違うだけの（元の音程が同じ）音とも言える。たとえばカラオケで女声と男声がオクターブ違いの同じメロディーをデュエットするのは、音色が違うが同じ音程で歌っている、ともいえる。 このようにオクターブ違いの音は「もともと同じ音程」なので、１オクターブという概念を音階の設計から外すわけにはいかないのである。・・・って説明になってるかな？ （さて矩形波の奇数倍音は、どう説明したらいいのか？　まぁ「ド」に対する「ソ」も外せない音程ということで・・）

**いったん修了**