Pythonで〇〇が立った！〜scipyで慣性モーメント最小化問題を解く〜

7995 ワード
Python scipy numpy Python テキストリンク

teratailの質問回答で、図形が立つとは何ぞや、という質問に答えてみましたので、一般化して記録しておきます。

皆さんおなじみの左のような絵、走っている状況であるのと漫画風のデフォルメが効いているので、自然に感じます。しかし、真ん中のように2値化してみると、漫画効果が薄れて、なんだかバランスが悪く、倒れてしまいそうに感じます。それを、右のように「立たせる」ための角度を、理論的に求めて、立たせてみよう、という試みです。

アイデアは次のようなものです。

  • まずは図形の重心が中心になるように並行移動する
  • そのあと、中心を軸に回転させて、最も「バランスよく立った」状態にする
  • 「バランスよく立った」を、中心縦軸を固定して「まわした」時の慣性モーメントが最小になる、と定義する

慣性モーメントを最小化することで、左右をなるべく縮めた位置が求められます。結果的に図形が「立つ」ことになります。フィギュアスケートの選手が、腕をたたむことで、安定して高速回転できるのと、原理は同じです。

重心を求めること、並行移動、回転、目的関数の最小化にはPythonの数値計算ライブラリであるscipyの関数を使いました。scipyの勉強にもなりますね。目的関数である慣性モーメントを算出する部分のみnumpyで手作りしました。なお、上の図はレイアウトの関係上縦長の枠にしていますが、実際に計算した図形は、回転してはみ出さないように、ほぼ正方形を枠としています。

以下がコードです。-14.12度ほどの回転をさせることで、気持ちよく「立たせる」ことができました。

import numpy as np
from scipy import ndimage, optimize
import matplotlib.pyplot as plt
import cv2

# 画像読み込みと2値化、白黒反転
img = cv2.imread('sample.jpg', 0)
_, img = cv2.threshold(cv2.bitwise_not(img), 1, 1, cv2.THRESH_BINARY)

# 平行移動して重心を中心にする
center = ndimage.center_of_mass(img)
img = ndimage.shift(img, np.array(img.shape)/2-np.array(center))

# 2値化した図形を表示
plt.imshow(img)
plt.show()

# 図形を立てた時の、中心軸に対する慣性モーメント
# mass:縦軸ピクセル数、radius:中心軸からの左右のピクセル距離
# 慣性モーメント = Σ mass * radius^2
def inertia(img):
    mass = img.sum(axis=0)
    radius = np.abs(np.arange(-len(mass)//2, len(mass)//2+1))
    if len(mass) != len(radius):
        radius = radius[radius != 0] - 0.5  # 中心位置の補正
    return (mass * radius * radius).sum()

# 回転させて慣性モーメントを測定する=最小化関数
def rotated_inertia(degree, img):
    return inertia(ndimage.rotate(img, degree, reshape=False))

# 最小化関数の最小値を求めることで、-90〜90度の回転時の最小慣性モーメントを得る
res = optimize.minimize_scalar(
    rotated_inertia, bounds=[-90,90], args=(img), method='Bounded')

# 最小値の時の回転角度を求める
print(res.x)

# 最小値の時の回転を実際に行う
img = ndimage.rotate(img, res.x, reshape=False)

# 「立った」図形を表示
plt.imshow(img)
plt.show()