乳揺れは単振動しているのか?

高校物理で習ったことですが、ばねに付いたおもりは単振動する。
数式で書くと正弦波(sinカーブ)に従って、位置が変わっていく。

f:id:tottokotokoroten:20210813130140p:plain

振動、つまり揺れるといえば……みんな大好き、おっぱいですね。
乳揺れも、単振動しているのだろうか??

この疑問に答えるため、かんたんなプログラミングで乳揺れの数式化を試みた。
今回はPythonopenCVというライブラリを主に用いた。なぜPythonか?そりゃあ、おっぱいそんだからに決まっているでしょう。

流れとしては、乳揺れの動画から、おっぱいの位置を検出、その座標変化を求め、これを正弦波にフィッティングできるかを検証する、という戦略。

今回は便宜上、おっぱいの位置は、乳首の位置として定義した(乳首は色の濃さをもとに検出できるので)
検証用動画には「boob physics nipple」でGoogle GIF画像検索して出てきた以下のGIF画像を用いた。
この画像が解析に適していそうな点としては、乳首の色が比較的濃いこと、おっぱい以外にあまり動きがないこと、が挙げられる。

f:id:tottokotokoroten:20210812224638g:plain

import cv2
from matplotlib import pyplot as plt
import numpy as np
import csv

gif = cv2.VideoCapture('boob.gif') # gif画像の読み込み
nflames = gif.get(cv2.CAP_PROP_FRAME_COUNT)

# 1フレームごとに画像取得
images = []
i = 0
while True:
    is_success, img = gif.read()
    if not is_success:
        break
    images.append(img)
    i += 1

さて、GIF画像は3色のチャネルから成っているので、その中で最も乳首のコントラストが良いチャネルを選ぶ。

f:id:tottokotokoroten:20210812223120p:plain

一番右が良さそうである。

次に、この画像をもとに、乳首の検出が可能か試みた。

左から、①元画像、②乳首付近の切り抜き、③それを平滑化処理したもの、④それを閾値で切り二値化したもの

f:id:tottokotokoroten:20210812223620p:plain

④でなんとなく乳首付近を切り抜くことはできているが、バックグラウンドノイズ(乳の陰影)が邪魔をして、右下の影の部分も切り出されてしまっている。

そこで次に、時間方向に動画を平均化した一枚画像を作り、各フレームについて平均画像との差分をとることで、乳首を際立たせることにした。

左から、①元画像、②平均化画像、③それらの差分

f:id:tottokotokoroten:20210812224037p:plain

③を見ると、おお! なんとなく乳首が際立っていそうだ!

以下は乳首付近の拡大図。左から、①元画像、②平均化画像、③差分&平滑化したもの、④閾値で切ったもの

f:id:tottokotokoroten:20210812224146p:plain

きれいに乳首検出できました。やったね!

そこで、この乳首の重心を乳首の位置として定義し、動画の各フレームについて乳首の位置を検出しました。

# 平均化画像を取得
imgarray = np.array(images)
imgmean = imgarray[:,:,:,2].mean(axis=0) - imgarray[:,:,:,2].mean(axis=None)
imgmean0 = imgmean[60:140,120:180]

# main
j = 0
x = []
y = []
xx = []
yy = []
threshold = 80 #二値化の閾値
while True:
    img0 = imgarray[j,60:140,120:180,2] #ROIを切り取り
    img1 = img0-imgmean0 #平均画像との差分
    img2 = cv2.blur(img1, (3,3)) #3x3カーネルで平滑化
    ret, img3 = cv2.threshold(img2,threshold,255,cv2.THRESH_BINARY) #閾値で二値化
    img4 = img3.astype(np.uint8) #uint8に変換
    nlabels, labels, stats, centroids = cv2.connectedComponentsWithStats(img4) #オブジェクト検出

    if j==nflames-1: #最後のフレームで処理を停止
        break
    y.append(centroids[0][0]) #ROI中の相対座標
    x.append(centroids[0][1])
    yy.append(centroids[0][0]+120) #画像全体中での座標
    xx.append(centroids[0][1]+60)
    j += 1

こうして検出した乳首の軌跡を元画像に重ねると、以下になります。

f:id:tottokotokoroten:20210813113959p:plain

# 元画像と軌跡を並べて表示
fig = plt.figure()
ax4 = fig.add_subplot(1,2,1)
ax4.imshow(images[8][:,:,2])
ax4.plot(yy[2:12],xx[2:12], color="red")
ax4.set_xlim([0, 299])
ax4.set_ylim([200, 0])

これを見て思ったこと……乳首あんまり動いてなくね??

そう、ここで改めて画像を見返すと、まさにその通り……素材選定を誤ったことに気づきました。
このGIF動画は全部で31フレームと短いうえ、乳首が指でつままれている時間も長いので、実際に乳揺れしている瞬間はほんの一瞬で、たった11フレームしかない。上記の軌跡はその11フレームだけなので、このような微々たる動きになってしまうのだ。

さて、とにもかくにもここで、本来の目的である「乳揺れは単振動しているのか?」を検証するため、乳首のy座標の時間変化を出してみた。

f:id:tottokotokoroten:20210813114842p:plain

# y座標の時間変化をプロット
fig = plt.figure(figsize=(3,4))
ax5 = fig.add_subplot(ylabel='y',xlabel='flame')
ax5.plot(xx[2:12])
ax5.set_ylim([200, 0])

うーん、せいぜい1周期しかない。

かんたんにエクセルで近似してみました。

f:id:tottokotokoroten:20210813123532p:plain

振幅は減衰していくはずなので、もっとパラメータをいじる必要がある。

ちなみに、多項式近似してみると…

3次

f:id:tottokotokoroten:20210813123651p:plain

4次

f:id:tottokotokoroten:20210813123703p:plain

おっ、これは…

結論:この乳揺れは4次元
 

次はもっと長い時間揺れてる動画でやります!!

あと、あんまり乳首がモロに写ってる画像って公開しづらいので、できれば乳首に頼らず水着グラビアとかでやりたいなぁ

30歳男の夏休みの自由研究でした