動かして学ぶバイオメカニクス #15 〜3Dでスティックピクチャを動かす〜
今回は数理的な話題はお休みとし，アニメーションを動かすようにインタラクティブにデータを表示する方法を示す．ロジックを検証する方法，考えを可視化する方法としてこのようなツールも重要である．アイデアやPythonの知識次第で，いろいろな表示方法が可能だろう．ここで紹介する方法はひとつの方法に過ぎないが，読者の便利なツールとして発展させるきっかけになってくれればと思う．
はじめに
で示してきたように，角度やトルクのグラフを見ていても，実際の現象を想像するのは困難で，得られた計測データを実際の動きと比較して確認できなくては，はじまらない．また，実際にこのシリーズで過去に提示したコードにも誤りがありご迷惑をおかけしたが，理論からもコードからも，多角的に検証するすべが必要である．推定した重心位置が正しいかどうかなどは，３次元空間でその位置を描かないとよくわからないだろう．
そこで，今回はステックピクチャとグラフをアニメーションで確認する方法を提示する．ここで示した方法は，特別なものではなく，ありふれた技術でこんな事もできるという程度の参考例程度の内容だが，読者の方でよりよい可視化を構築していただきたい．
なお，Pythonについては，最後に取り上げる書籍や
などを参照するとよいだろう．ここではPythonの解説はあまり深入りしない．
確認した動作環境とコードについて
いつものように，以下のGoogle ColaboratoryのリンクからブラウザでPythonコードを実行できるようにした．なお，ログインするためには，Googleアカウントが必要となるので注意をされたい．
また，使用する飛び降り動作のモーションキャプチャのサンプルデータ（sample_optitrack_221027.csv）とフォースプレートのサンプルデータ（sample_fp_221027.csv），投球動作のモーションキャプチャのサンプルデータ（sample_optitrack_sec15.csv）は，このリンクからダウンロードして，Google Driveをマウントし，マウントしたドライブにそれらをコピーしてご使用いただきたい．飛び降り動作のデータは，第６〜12章で使用したデータと同じであるので，すでにダウンロード済みの方は，サイドダウンロードする必要はない．投球動作のデータは今回新しく使用するデータで必ずダウンロードする必要がある．
ここでのコードは，JupyterLab（デスクトップアプリ），Google Colaboratoryで検証したが，とくにグラフィックスの表示結果やコードは，その二つでも若干動作が異なる．動作させるPython環境に依存し，少々修正が必要となるかも知れない．とりあえずJuputerLabのアプリやGoogle Colabratoryで確認したが，環境によってグラフィックスのコードの調整をお願いしたい．ここで示したコードは，JupyterLabで動作するが，Colabでは若干修正している．
なお，Google ColaboratoryやGoogle Driveの使用方法は，
を参照されたい．
スティックピクチャを動かす
使用する関数は第12章
と共通なので，ここでは示さないが，新しく定義した関数を以下に示す．
パッケージの読み込み
from ipywidgets import interact, IntSlider, FloatSliderこのパッケージはスライダーバーを使うアニメーションで使用する．動作環境によっては異なるパッケージと関数が必要となるかも知れない．これは，Jupyter系の環境で動作するようだ．
Animation 1（飛び降りバランス動作）
準備・データの格納
今回の計算に必要なBodyLinkクラスの配列b_linkの各インスタンス変数にデータを格納する．
親子関係（child_val），兄弟姉妹関係（sister_val），遠位端（distal_val）
[b_link[i[0]].set_child_sister(b_link[i[1]], b_link[i[2]]) for i in child_sister_list];
[b_link[i[0]].set_distal_end(b_link[i[1]]) for i in distal_end_list];linkの原点の位置ベクトル（p）
手，前腕の姿勢回転行列を計算するためだけに必要なマーカの位置ベクトルもインスタンス変数pに格納
# 関節位置をBodyLinkのpに格納
set_link_p_list(b_link, range(1,21), file_path_op_15, .00005, 360.)
# 姿勢回転行列を計算するために必要なマーカ位置ををBodyLinkのl_id: 21~24のpに格納
set_rot_p_list(b_link, range(21,27), file_path_op_15, .00005, 360.)重心の位置ベクトル（cg）
各リンクの重心の加速度ベクトルを格納する前に，各リンクの重心の位置ベクトルcgを計算
# 各linkの重心を，self.cgに格納
for n in range(1,16):
b_link[n].cg = link_cg(b_link, n, file_path_op_15, .00005, 360.)解析例１
まず以下のコマンドで「3Dグラフの領域」を作る．
# ３Dグラフ領域の作成
fig = plt.figure(figsize=(10, 10)) # figsizeで大きさを指定
ax = fig.add_subplot(projection='3d') # 3Dの図を指定
そこに，３つの座標[X1,Y1, Z1]，座標[X2, Y2, Z2]，座標[X3, Y3, Z3] を線で結ぶようなグラフは
ax.plot([X1, X2, X3], [Y1, Y2, Y3], [Z1, Z2, Z3])
で描けるので，時刻 t での，右脚の大腿，下腿，足部の各関節（右股関節:10，右膝関節：11，右足関節：12）の点 b_link[].pの各座標を線で結ぶようなグラフは
ax.plot([b_link[10].p[t][0], b_link[11].p[t][0], b_link[12].p[t][0]], [b_link[10].p[t][1], b_link[11].p[t][1], b_link[12].p[t][1]])
と書ける．そこで
# ３Dグラフ領域の作成
fig = plt.figure(figsize=(10, 10))
ax = fig.add_subplot(projection='3d')
# データの準備
time = 100
x = [b_link[10].p[time][0], b_link[11].p[time][0], b_link[12].p[time][0]]
y = [b_link[10].p[time][1], b_link[11].p[time][1], b_link[12].p[time][1]]
z = [b_link[10].p[time][2], b_link[11].p[time][2], b_link[12].p[time][2]]
print(x)
print(y)
print(z)
# グラフの描画
ax.plot(x, y, z)
ax.set_aspect('equal')
とるすることで，時刻（フレーム数） time = 100 の，股関節b_link[10]と，膝関節b_link[11]と，足関節b_link[12]を結ぶ直線が描ける．
スティックピクチャ
ここで，b_linkから，ステックピクチャを描く座標を作る関数を作る．
まず，スティックピクチャを描く関節点の順番を配列に収める．
なお，ここではスティックピクチャは以下の上肢，頭部から体幹まで，下肢の３つの線で描く．
upper_arm_list = [17, 6, 5, 4, 3, 7, 8, 9, 18] # 腕の関節位置
body_list = [16, 3, 2, 1] # 頭から腰までの関節位置
lower_body_list = [19, 12,11, 10, 13, 14, 15,20] # 下肢の関節位置以下は，for文の内包表記（補足１参照．内包表記については，note.nkmk.meなども参照）を使用し，時刻（フレーム数）timeのスティックピクチャ用データを出力する関数である．
# 22.12.29修正
# 23.1.2修正
# 上記で定義した関節位置のリストを利用し
# xyz３成分を出力
def arm_list_xyz(body_link, time, joint_list):
return np.array([body_link[num].p[time] for num in joint_list]).T以下は，この関数arm_list_xyz()を使用し，時刻（フレーム数） 100時点の上肢の各位置座標をx, y, zの順番で配列に収める例．
arm_list_xyz(b_link, 100, upper_arm_list)
>>>
array([[-0.73957689, -0.69166199, -0.59685296, -0.43669012, -0.26924387,
-0.09797377, 0.066447 , 0.15651652, 0.21985834],
[ 0.82186204, 0.89280468, 0.96983795, 1.03769905, 1.05740791,
1.06835112, 1.03603764, 0.99535049, 0.95483711],
[ 1.33222412, 1.39054271, 1.63962887, 1.87142728, 1.91858204,
1.87472344, 1.63800472, 1.38685241, 1.33298292]])以下のstick()関数は，さらにこれを利用し，スティックピクチャを描く関数で
def stick(ax, body_link, time, lw=1, xrange=(-.2, 1.75), yrange=(-.5, 1.), zrange=(0, 2.1)):
ax.plot(*arm_list_xyz(body_link, time, upper_arm_list), linewidth=lw, c='C0')
ax.plot(*arm_list_xyz(body_link, time, body_list), linewidth=lw, c='C1')
ax.plot(*arm_list_xyz(body_link, time, lower_body_list), linewidth=lw, c='C2')
ax.set_xlabel('x')
ax.set_ylabel('y')
ax.set_zlabel('z')
ax.set_xlim(*xrange)
ax.set_ylim(*yrange)
ax.set_zlim(*zrange)
ax.set_aspect('equal') # Colabでは動作しないこれを利用し，
fig = plt.figure(figsize=(10, 10))
ax = fig.add_subplot(projection='3d')
frame = 100 # データのフレーム数（時間）
stick(ax, b_link, frame, 2, (-1.25, .75), (0., 2.), (0, 2))
のように，時刻（frame）100のときのスティックピクチャを描いた．これを，スライダーバーで動かすアニメーションにする．
def stick_plot(body_link):
#upper_arm_list = [17, 6, 5, 4, 3, 7, 8, 9, 18] # 腕の関節位置
#body_list = [16, 3, 2, 1] # 頭から腰までの関節位置
#lower_body_list = [19, 12,11, 10, 13, 14, 15,20] # 下肢の関節位置
step_frame = 10 # アニメーションの際の描画でスキップするフレーム数
q = len(body_link[10].p) // step_frame # フレーム数をstep_sizeで割った，整数の商
max_frame = q * step_frame # 描画する最大フレーム数
step_time = .2
sf = 360.
q = (len(body_link[10].p) / sf) // step_time # フレーム数をstep_sizeで割った，整数の商
max_time = q * step_time
@interact(time = FloatSlider(min=0, max=max_time, step=step_time, value=0, continuous_update=True), # 描画する時間 [s]
elev = IntSlider(min=-10, max=90, step=5, value=0, continuous_update=True), # ViewPoint（視点）の仰角
azim = IntSlider(min=-100, max=180, step=5, value=-100, continuous_update=True)) # ViewPointの方位角
def plot_3d(elev, azim, time):
fig = plt.figure(figsize=(10, 10))
ax = fig.add_subplot(projection='3d')
ax.view_init(elev=elev, azim=azim) # 視点をここで設定．
time2 = int(time * sf)
stick(ax, b_link, time2, 2, (-1.25, .75), (0., 2.), (0, 2.2))ここで，スライダーバーは@interactで３つの変数（time, elev, azim）定義する．timeは時間軸を変え，アニメーションを動かすためにtimeのバーを動かす必要がある．自動再生は異なるanimation.FuncAnimationなどを利用するとよいだろう（最後の「アニメーション」で述べる）．elevとazimは３Dビューの視点を調整するバーである．関数の中で，初期値や終端値などを設定している．
その変数の定義後に，実際に動かす関数plot_3d(elev, azim, time)を，この関数内で記述することで，スライダーバーで動かすことができる．
Interactの詳細は，
などをご覧になっていただきたい．
そこで，これを使って以下のように実行すると
stick_plot(b_link)
のように左上のスライダーバーで，時間（アニメーション）と視点を動かすことができる．動作はGoogle Colabでご確認いただきたい．ただし，Colab上では反応が遅いことお許し頂きたい．
Animation ２（投球動作）
ここでは，いつもと異なる投球動作のデータを利用し，グラフと連動する例を示す．130 km/hぐらいの球速の投球データを活用し，過去の章を参考に，ここで示した以外の計算を試してみると良いだろう．
準備・データの格納
ここで，別途，BodyLinkの新しい配列b_link2を作り，Animation 1と同じようにここでの計算に必要なデータを格納する．
b_link2 = [0] * 27 # 配列の初期化
b_link2[0] = BodyLink(0, '', '', .0, .0) # 0のときストップ
b_link2[1] = BodyLink(1, 'Hip', 'Bone', .187, (1.-.609)) # 腰（下体幹）
b_link2[2] = BodyLink(2, 'Chest', 'Bone', .302, (1.-.428)) # 胸（上体幹）
b_link2[3] = BodyLink(3, 'Neck', 'Bone', .069, (1.-.821)) # 頭部
b_link2[4] = BodyLink(4, 'RUArm', 'Bone', .027, .529) # 右上腕
b_link2[5] = BodyLink(5, 'RFArm', 'Bone', .016, .415) # 右前腕
b_link2[6] = BodyLink(6, 'RHand', 'Bone', .006, .891) # 右手
b_link2[7] = BodyLink(7, 'LUArm', 'Bone', .027, .529) # 左上腕
#b_link2[8] = BodyLink(8, 'LFArm', 'Bone', .016, .415) # 左前腕
b_link2[8] = BodyLink(8, 'LELB', 'Bone Marker', .016, .415) # 左前腕
b_link2[9] = BodyLink(9, 'LHand', 'Bone', .006, .891) # 左手
b_link2[10] = BodyLink(10, 'RThigh', 'Bone', .110, .475) # 右大腿
b_link2[11] = BodyLink(11, 'RShin', 'Bone', .051, .406) # 右下腿
b_link2[12] = BodyLink(12, 'RFoot', 'Bone', .011, (1.-.595)) # 右足
b_link2[13] = BodyLink(13, 'LThigh', 'Bone', .110, .475) # 左大腿
b_link2[14] = BodyLink(14, 'LShin', 'Bone', .051, .406) # 左下腿
b_link2[15] = BodyLink(15, 'LFoot', 'Bone', .011, (1.-.595)) # 左足
# 各リンクの遠位端も定義（遠位端の位置ベクトルを計算するためだけに必要な変数
b_link2[16] = BodyLink(16, 'Head', 'Bone', .0, .0) # 3 Neckの遠位端
b_link2[17] = BodyLink(17, 'RFIN', 'Bone Marker',.0, .0) # 6 RHandの遠位端
b_link2[18] = BodyLink(18, 'LFIN', 'Bone Marker', .0, .0) # 9 LHandの遠位端
b_link2[19] = BodyLink(19, 'RToe', 'Bone', .0, .0) # 12 RFootの遠位端
b_link2[20] = BodyLink(20, 'LToe', 'Bone', .0, .0) # 15 LFootの遠位端
# 手，前腕の姿勢回転行列を計算するためだけに必要な変数
b_link2[21] = BodyLink(21, 'RWRA', 'Bone Marker', .0, .0) # 右手首橈骨側
b_link2[22] = BodyLink(22, 'RWRB', 'Bone Marker', .0, .0) # 右手首尺骨側
b_link2[23] = BodyLink(23, 'LWRA', 'Bone Marker', .0, .0) # 左手首橈骨側
b_link2[24] = BodyLink(24, 'LWRB', 'Bone Marker', .0, .0) # 左手首尺骨側
b_link2[25] = BodyLink(25, 'RFHD', 'Bone Marker', .0, .0) # 右前頭
b_link2[26] = BodyLink(26, 'LFHD', 'Bone Marker', .0, .0) # 左前頭新しいBodyLinkクラスの配列b_link2の各インスタンス変数にデータを格納する．
親子関係（child_val），兄弟姉妹関係（sister_val），遠位端（distal_val）
[b_link2[i[0]].set_child_sister(b_link2[i[1]], b_link2[i[2]]) for i in child_sister_list];
[b_link2[i[0]].set_distal_end(b_link2[i[1]]) for i in distal_end_list];linkの原点の位置ベクトル（p）
手，前腕の姿勢回転行列を計算するためだけに必要なマーカの位置ベクトルもインスタンス変数pに格納
# 関節位置をBodyLinkのpに格納
set_link_p_list(b_link2, range(1,21), file_path_op_15_2, .00005, 360.)
# 姿勢回転行列を計算するために必要なマーカ位置ををBodyLinkのl_id: 21~24のpに格納
set_rot_p_list(b_link2, range(21,27), file_path_op_15_2, .00005, 360.)重心の位置ベクトル（cg）
各リンクの重心の加速度ベクトルを格納する前に，各リンクの重心の位置ベクトルcgを計算
# 各linkの重心を，self.cgに格納
for n in range(1,16):
b_link2[n].cg = link_cg(b_link2, n, file_path_op_15_2, .00005, 360.)解析例
stick()関数で，ある時刻のスティックピクチャを描くと
fig = plt.figure(figsize=(10, 10))
ax = fig.add_subplot(projection='3d')
stick(ax, b_link2, 100)
これを利用して，時間経過するスティックピクチャを重ね書きする関数
def stick_plot3(body_link):
upper_arm_list = [17, 6, 5, 4, 3, 7, 8, 9, 18] # 腕の関節位置
body_list = [16, 3, 2, 1] # 頭から腰までの関節位置
lower_body_list = [19, 12,11, 10, 13, 14, 15,20] # 下肢の関節位置
x_range=(-.2, 1.75) # x軸描画範囲
y_range=(-.5, 1.) # y軸描画範囲
z_range=(0, 2.1) # z軸描画範囲
step_frame = 10 # アニメーションの際の描画でスキップするフレーム数
q = len(body_link[10].p) // step_frame # フレーム数をstep_sizeで割った，整数の商
max_frame = q * step_frame # 描画する最大フレーム数
# 時間（s）に変換
step_time = .005 # 5 ms毎に描画
sf = 360. # サンプリングレート
q = (len(body_link[10].p) / sf) // step_time # 時間をstep_sizeで割った，整数の商
max_time = q * step_time
fig = plt.figure(figsize=(10, 10))
ax = fig.add_subplot(projection='3d')
for time2 in range(0, max_frame, step_frame):
stick(ax, body_link, time2, 1, x_range, y_range, z_range)を作って，描くと
stick_plot3(b_link2)
さらに，これに，@interactを利用してアニメーションする関数stick_plot2()を
def stick_plot2(body_link):
upper_arm_list = [17, 6, 5, 4, 3, 7, 8, 9, 18] # 腕の関節位置
body_list = [16, 3, 2, 1] # 頭から腰までの関節位置
lower_body_list = [19, 12,11, 10, 13, 14, 15,20] # 下肢の関節位置
x_range=(-.2, 1.75) # x軸描画範囲
y_range=(-.5, 1.) # y軸描画範囲
z_range=(0, 2.1) # z軸描画範囲
step_frame = 10 # アニメーションの際の描画でスキップするフレーム数
q = len(body_link[10].p) // step_frame # フレーム数をstep_sizeで割った，整数の商
max_frame = q * step_frame # 描画する最大フレーム数
# 時間（s）に変換
step_time = .005 # 5 ms毎に描画
sf = 360. # サンプリングレート
q = (len(body_link[10].p) / sf) // step_time # 時間をstep_sizeで割った，整数の商
max_time = q * step_time
@interact(time = FloatSlider(min=0, max=max_time, step=step_time, value=0, continuous_update=True), # 描画する時間 [s]
elev = IntSlider(min=-10, max=90, step=5, value=0, continuous_update=True), # ViewPoint（視点）の仰角
azim = IntSlider(min=-100, max=180, step=5, value=-100, continuous_update=True)) # ViewPointの方位角
def plot_3d(elev, azim, time):
fig = plt.figure(figsize=(10, 10))
ax = fig.add_subplot(projection='3d')
ax.view_init(elev=elev, azim=azim)
time2 = int(time * sf)
[stick(ax, body_link, time2, .1, x_range, y_range, z_range) for time2 in range(0, max_frame, step_frame)] # 細い線で重ね描画
stick(ax, body_link, time2, 2, x_range, y_range, z_range) # time2の時間だけを太線でアニメーションのように作り，
stick_plot2(b_link2)で実行すると，
を実行できる（動作はGoogle Colabで確認されたい．ただし，Colabだと，ゆっくりとしか動かないだろう）．
次に，推定した床反力のグラフといっしよにアニメーションを行う．以下のように床反力のグラフを描く事ができるが，
t_list = np.arange(len(b_link2[1].p))/360.
plt.plot(t_list, smoothing_spline(b_link2[1].Resultant_Cg(86., True), .00005, 360., 2)-[0., 0., -9.8])
plt.xlabel('Time [s]')
plt.ylabel('Force [N]')
plt.legend(['x', 'y', 'z'])
これも，一緒に以下の関数で動かしてみる．
def stick_graph_plot2(body_link):
upper_arm_list = [17, 6, 5, 4, 3, 7, 8, 9, 18] # 腕の関節位置
body_list = [16, 3, 2, 1] # 頭から腰までの関節位置
lower_body_list = [19, 12,11, 10, 13, 14, 15,20] # 下肢の関節位置
x_range=(-.2, 1.75) # x軸描画範囲
y_range=(-.5, 1.) # y軸描画範囲
z_range=(0, 2.1) # z軸描画範囲
step_frame = 10 # アニメーションの際の描画でスキップするフレーム数
q = len(body_link[10].p) // step_frame # フレーム数をstep_sizeで割った，整数の商
max_frame = q * step_frame # 描画する最大フレーム数
# 時間（s）に変換
step_time = .0025 # 2.5 ms毎に描画
sf = 360. # サンプリングレート
q = (len(body_link[10].p) / sf) // step_time # 時間をstep_sizeで割った，整数の商
max_time = q * step_time
cg = b_link2[1].Resultant_Cg(86., True) # 重心位置
cg_acc = smoothing_spline(b_link2[1].Resultant_Cg(86., True), .00005, 360., 2)-[0., 0., -9.8] # 重心位置の加速度＋重力加速度
t_list = np.arange(len(b_link2[1].p))/360.
m = 86. #体重
@interact(time = FloatSlider(min=0, max=max_time, step=step_time, value=0, continuous_update=True), # 描画する時間 [s]
elev = IntSlider(min=-10, max=90, step=5, value=0, continuous_update=True), # ViewPoint（視点）の仰角
azim = IntSlider(min=-180, max=180, step=5, value=-90, continuous_update=True)) # ViewPointの方位角
def plot_3d(elev, azim, time):
fig = plt.figure(figsize=(10, 15))
ax = fig.add_subplot(projection='3d')
ax.view_init(elev=elev, azim=azim)
time2 = int(time * sf)
[stick(ax, body_link, time2, .1, x_range, y_range, z_range) for time2 in range(0, max_frame, step_frame)] # 細い線で重ね描画
stick(ax, body_link, time2, 2, x_range, y_range, z_range) # time2の時間だけを太線でアニメーション
ax.scatter(*cg[time2], s=50, c='k') # 重心
scale = .05
ax.plot([cg[time2, 0], cg[time2, 0]+scale*cg_acc[time2, 0]],
[cg[time2, 1], cg[time2, 1]+scale*cg_acc[time2, 1]],
[cg[time2, 2], cg[time2, 2]+scale*cg_acc[time2, 2]], c='r', linewidth=2) # 重心に作用する力（＝床反力の合力）のベクトル
# 床反力のグラフ
ax2 = fig.add_subplot(5,5,(2,4))
ax2.plot(t_list, m * cg_acc)
ax2.axvline(x=time2/sf, linestyle='--', c='k')
plt.xlabel('Time [s]')
plt.ylabel('Force [N]')
plt.legend(['x', 'y', 'z'])これを使用すると
stick_graph_plot2(b_link2)
のように，グラフの時間軸も移動させ，連動させることもできる（動作の確認はGoogle Colabや以下の動画で）．
アニメーション
以下の例は，スライダーバーで動かすのではなくアニメーション化した．まず，さらに必要なパッケージを呼び込み，
from matplotlib.animation import FuncAnimation # アニメーション用
from IPython.display import HTML
import pathlib # ファイルパス指定用ここで，アニメーションで使用する重心と，床反力（重心の加速度から計算）を計算しておき，
# ここで，アニメーションで使用する重心と，床反力（重心の加速度から計算）を計算しておく
cg = b_link2[1].Resultant_Cg(86., True) # 重心位置
cg_acc = smoothing_spline(b_link2[1].Resultant_Cg(86., True), .00005, 360., 2)-[0., 0., -9.8] # 重心位置の加速度＋重力加速度
t_list = np.arange(len(b_link2[1].p))/360.以下で，実行する．
ここで，MatplotlibのFuncAnimation()を利用しているが，利用方法は検索して調べていただきたい．
# アニメーションの初期設定
upper_arm_list = [17, 6, 5, 4, 3, 7, 8, 9, 18] # 腕の関節位置
body_list = [16, 3, 2, 1] # 頭から腰までの関節位置
lower_body_list = [19, 12,11, 10, 13, 14, 15,20] # 下肢の関節位置
x_range = (-.2, 1.75) # x軸描画範囲
y_range = (-.5, 1.) # y軸描画範囲
z_range = (0., 2.1) # z軸描画範囲
p_elev = 20
p_azim = -110
f_len = len(b_link2[10].p)
step_frame = 10 # アニメーションの際の描画でスキップするフレーム数
q = f_len // step_frame # フレーム数をstep_sizeで割った，整数の商
max_frame = q * step_frame # 描画する最大フレーム数
sf = 360. # サンプリングレート
m = 86. #体重
# ３Dグラフ領域の作成
fig = plt.figure(figsize=(10, 12))
ax = fig.add_subplot(projection='3d')
ax.view_init(elev=p_elev, azim=p_azim)
ax2 = fig.add_subplot(5,5,(2,4))
# アニメーション更新用関数の作成
def plot_3d(frame):
ax.cla()
ax2.cla()
[stick(ax, b_link2, time, .1, x_range, y_range, z_range) for time in range(0, max_frame, 10)] # 細い線で重ね描画
frame_new = step_frame * frame
stick(ax, b_link2, frame_new, 2, x_range, y_range, z_range) # time2の時間だけを太線でアニメーション
ax.scatter(*cg[frame_new], s=50, c='k') # 重心
scale = .05
ax.plot([cg[frame_new, 0], cg[frame_new, 0] + scale * cg_acc[frame_new, 0]],
[cg[frame_new, 1], cg[frame_new, 1] + scale * cg_acc[frame_new, 1]],
[cg[frame_new, 2], cg[frame_new, 2] + scale * cg_acc[frame_new, 2]], c='r', linewidth=2) # 重心に作用する力（＝床反力の合力）のベクトル
# 床反力のグラフ
ax2.plot(t_list, m * cg_acc)
ax2.axvline(x=frame_new / sf, linestyle='--', c='k')
ax2.set_xlabel('Time [s]')
ax2.set_ylabel('Force [N]')
ax2.legend(['x', 'y', 'z'])
# アニメーションのオブジェクトを作成
anim = FuncAnimation(fig, plot_3d, frames=q, interval=50)
# Jupyter内でアニメーションの表示
plt.close()
HTML(anim.to_jshtml())
# Gif/MP4に保存する場合
#path_dir = pathlib.Path('/Users/.../sample_data/')
#path_gif = path_dir.joinpath('animation.gif')
#path_mp4 = path_dir.joinpath('animation.mp4')
#anim.save(path_mp4)
#plt.close()ここで，
anim = FuncAnimation(fig, plot_3d, frames=q, interval=50)は，アニメーションを実行する関数だが，最初の引数figは３Dのグラフを各領域を，第２引数plot_3dは更新して描いていくグラフの関数名を与えている．別途定義しているその関数plot_3d(frame)の引数frameは，更新する際の時間（フレーム）の変数である．この引数frameに応じて更新関数plot_3d()を書き換えていく．FuncAnimation()の第３引数は繰り返し，一つずつ増加するframeの最大値で，更新関数plot_3d(frame)を呼び出す回数を指定し，その値で繰り返しを終了する．第４引数intervalは更新する時間（ms で指定．デフォルトは 200 ms）を与える．
出力は，前述のスライダーバーと似ているので，ここでは示さない．なお，実行してから，アニメーションが出力されまでは，かなり時間を要する．gifへの出力は最後のコメントアウトした部分（#）を参照していただきたい．mp4への出力もできる（拡張子を.mp4とするだけ）．
ここでは，床反力のグラフ（破線が時刻を示している），重心位置（黒），推定した床反力（赤線）をスライダーバーで動かす例を示したが，ご自身でトルクや必要な情報を計算したり，色を変えたり，線の太さを変えたり，自由に変更していただければと思う．アニメーションの方法もPython環境に依存する．上記コードはJupyterLabやGoogle Colabでの動作を確認しているが，他の方法もいろいろとあるので，検索して試してみていただきたい．
Pythonが動作すれば，できることにも限りはあるが，高価なソフトなどを購入せずとも逆動力学計算や可視化も自由に行える．平滑化スプラインではなくバターワースフィルタを使いたい場合は，そこだけ変えればよいだろう．過去に示した章を参考に，コードとサンプルデータを活用してみてはいかがであろうか．行列の扱いもMatlabと似ていて，乗り換えるのもそれほど敷居は高くないだろう．
他人のプログラムを理解するのは，なかなかつらいので，自分ならこう動かすと思い切って書き換え，設計者の立場で作り変えてしまうのが良い．このマガジンでは筆者の力不足でコードの説明が足りないが，以下の補足に少しばかり，for文の内包表記と，ここで使用しているクラスの復習を行った．
過去にも紹介したが，また多くの書籍を眺めたわけではないが，筆者のPyhtonの教科書のおすすめは
である．Pythonを動かすために，他に環境構築も必要で，この１冊ですべて実現できるわけではないが，クラスやここで紹介している内包表記なども記述され，平易で丁寧な説明である．
補足
補足１：内包表記
Pythonの繰り返しのfor文はたとえば，
s_list = [17, 6, 5] # サンプルの配列
for num in s_list:
print(num)
>>>
17
6
5と書くことで，タスクを繰り返し計算することができる．一般的な表現方法だ．これをさらに内包表記を利用すると
s_list = [17, 6, 5] # サンプルの配列
[print(num) for num in s_list]と書き換えることができる．つまり，
for 任意の変数名 in 繰り返す（代入する）オブジェクト:
評価する式と書くところを
[評価する式 for 任意の変数名 in 繰り返す（代入する）オブジェクト]と１行で書くことができる．記述は簡潔になる．そこで，今回の例題で使用している配列b_linkを用いたfor文の内包表記の具体例を示す．
なお，b_linkはクラスBodyLinkの変数だが，たとえばb_link[10].pには右大腿（番号10）の関節の位置ベクトル（p）のxyzの時系列データが配列で収められている．その100番目の時刻の座標を抽出する際には，b_link[10].p[100]と記述すれば良い．
ここで使用している多関節構造を表現しているクラスBodyLinkを理解ない限り，このことは曖昧になるので，まだ呼んでいない方は，第2章
などを復習していただき，クラスを利用してどのように定義を与えているか理解していただきたい．もちろん無理にクラスを使用しなくても良いが，力学計算のコードは長くなってしまうかも知れない．あとはnumpyの配列の取り扱いになれることだ．
そこで，
upper_arm_list = [17, 6, 5, 4, 3, 7, 8, 9, 18] # 上肢の関節位置
for num in upper_arm_list:
print(b_link[num].p[100])
>>>
[-0.73957689 0.82186204 1.33222412]
[-0.69166199 0.89280468 1.39054271]
[-0.59685296 0.96983795 1.63962887]
[-0.43669012 1.03769905 1.87142728]
[-0.26924387 1.05740791 1.91858204]
[-0.09797377 1.06835112 1.87472344]
[0.066447 1.03603764 1.63800472]
[0.15651652 0.99535049 1.38685241]
[0.21985834 0.95483711 1.33298292]のように，b_link[num].p[100]，すなわち配列b_linkの関節の位置ベクトルpの時刻100の座標をfor文で出力するところを以下で内包表記に書き換える．なお，ここで，numに，上肢の関節位置upper_arm_list = [17, 6, 5, 4, 3, 7, 8, 9, 18]の数字をfor文で繰り返し代入している．そこで，
[b_link[num].p[100] for num in upper_arm_list]と書けば良い．関数arm_list_xyz()はこれを利用している．
【著作権・転載・免責について】
権利の帰属
本ホームページで提示しているソフトウェアならびにプログラムリストは，スポーツセンシング社の著作物であり，スポーツセンシング社に知的所有権がありますが，自由にご利用いただいて構いません．
本ページに掲載されている記事，ソフトウェア，プログラムなどに関する著作権および工業所有権については，株式会社スポーツセンシングに帰属するものです．非営利目的で行う研究用途に限り，無償での使用を許可します．
転載
本ページの内容の転載については非営利目的に限り，本ページの引用であることを明記したうえで，自由に行えるものとします．
免責
本ページで掲載されている内容は，特定の条件下についての内容である場合があります． ソフトウェアやプログラム等，本ページの内容を参照して研究などを行う場合には，その点を十分に踏まえた上で，自己責任でご利用ください．また，本ページの掲載内容によって生じた一切の損害については，株式会社スポーツセンシングおよび著者はその責を負わないものとします．
【プログラムの内容について】
プログラムや内容に対する質問に対しては，回答できないことのほうが多くなると思いますが，コメントには目は通します．回答は必要最低限にとどめますので，返信はあまり期待しないでいただけると幸いです，
「動かして学ぶ」という大それたタイトルをつけたものの，また，きれいなプログラムに対するこだわりはあるものの，実際のプログラミングのスキルは決して高くありません．最下部の方のコメント欄によるプログラムの間違いのご指摘は歓迎します．できるだけ反映します．
【解析・受託開発について】
スポーツセンシングでは，豊富な知見を持つ，研究者や各種エンジニアが研究・開発のお手伝いをしております．研究・開発でお困りの方は，ぜひスポーツセンシングにご相談ください．
【例】
・データ解析の代行
・受託開発
（ハードウェア、組込みソフトウェア、PC/モバイルアプリ）
・測定システム構築に関するコンサルティング など
その他，幅広い分野をカバーしておりますので，まずはお気軽にお問い合わせください．
【データの計測について】
スポーツセンシング社のスタジオで，フォースプレートやモーションキャプチャを利用した計測も行えます．出力されるデータと，ここで示したプログラム（入力データの取り込み関数を少々改変する必要があるが）で，同様な解析を行えますので，まずはお気軽にお問い合わせください．
株式会社スポーツセンシング
【ホームページ】sports-sensing.com
【Facebook】sports.sensing
【Twitter】Sports_Sensing
【メール】support@sports-sensing.com