３Dプリンタで暇つぶし

写真1枚から歯形の「乱れ」を見つける ― r(θ)の数値微分でPETG/POMギヤを読む

この記事は、これまでの歯形解析（gear_analyzerの開発、A面/B面の印刷方向による差、偶数歯・奇数歯のフロー補正など）の続きにあたる。初めて読む方向けに、必要な前提は都度カッコ書きで補っている。

きっかけ

将棋AIがなぜ人間に勝てるのかという雑談から始まった。評価関数の精緻化、力任せの探索、人間の固定観念からの解放、終盤の精度。この4つの柱を並べてみると、そのままPETGギヤの寸法最適化に応用できそうに見えた。

だが実際に検討すると、すぐに壁にぶつかった。将棋AIは「自己対戦」というほぼ無料で正確なシミュレータを持っている。一手指せば即座に正しい結果が返る。一方ギヤの最適化では、シミュレーション1回＝実際に印刷してベンチ測定するという、時間もコストもかかる工程になる。FDM特有の層間強度や、プリンタ個体差といった物理は数値モデル化が難しく、AIが「頭の中で」大量に試行することができない。将棋AIの強さを支えている「安価で正確なシミュレータ」は、ギヤ開発には存在しない。

代わりに長年頼ってきたのは、理想インボリュート曲線と歯型写真を目視で見比べる方法だった。「綺麗な戦闘機は性能が良い」という経験則と同じ構造で、機能的に優れた形状は結果として視覚的にも滑らかに見える。ではこの目視評価は、どこまで数値で裏付けられるのか。今回はそれを試した記録である。

なぜ数値微分なのか

「曲線を描くのに、馬と鹿をまちがえたら、馬鹿になる」という地口が、そのまま今回のテーマになった。実際、歯形解析ツール（以下gear_analyzer。USB顕微鏡で撮った歯車写真から圧力角や歯先径を測定するために自作したツール）の開発でも、輪郭のフィッティング手法をあるバージョンから次のバージョンへ強化した際、狙いとは逆に圧力角の算出値がずれるという現象が起きたことがある。曲線当てはめを精緻化したつもりが、歯形が実際に持っている別の性質（根元付近の輪郭の交差パターン）を無視した結果、数値だけが空回りしていた。

目視でこの種の「馬鹿」を避けられているとしたら、それは長年の経験が培った勘による部分が大きい。ならばその勘が見ている「何か」を、もっと直接的な数値に置き換えられないか。今回試したのが、歯形の輪郭を極座標 r(θ) として捉え、それを数値微分するという方法である。

手法：写真から二階微分まで

USB顕微鏡で撮影した歯車の写真1枚から、次の手順で解析した。実際に使ったコードの骨格は以下の通り（OpenCV + scipy）。

import cv2
import numpy as np
from scipy.signal import savgol_filter

# 1. 二値化して輪郭抽出
gray = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)
_, mask = cv2.threshold(gray, 40, 255, cv2.THRESH_BINARY)  # 閾値40
mask = cv2.morphologyEx(mask, cv2.MORPH_OPEN, np.ones((3, 3), np.uint8))
mask = cv2.morphologyEx(mask, cv2.MORPH_CLOSE, np.ones((3, 3), np.uint8))
contours, hierarchy = cv2.findContours(mask, cv2.RETR_CCOMP, cv2.CHAIN_APPROX_NONE)

# 2. シャフト穴の重心を中心に、外形輪郭を極座標へ
#    (outer=最大面積の輪郭、hole=その子輪郭のうち最大のもの)
cx, cy = hole[:, 0].mean(), hole[:, 1].mean()
dx, dy = outer[:, 0] - cx, -(outer[:, 1] - cy)  # 画像はy下向きなので反転
theta, r = np.arctan2(dy, dx), np.hypot(dx, dy)

# 3. 等間隔角度グリッドへ再サンプリングし、Savitzky-Golayで平滑化
theta_grid = np.linspace(-np.pi, np.pi, 2000, endpoint=False)
r_grid = np.interp(theta_grid, theta_sorted, r_sorted)  # 周期境界は前後に複製して処理
r_smooth = savgol_filter(r_grid, window_length=31, polyorder=3)

# 4. 一階微分・二階微分（数値微分）
dth = theta_grid[1] - theta_grid[0]
dr = np.gradient(r_smooth, dth)
d2r = np.gradient(dr, dth)

r(θ)そのものは理想的な山型（歯元→歯先→歯元）を描く。だが微分を重ねるほどノイズが増幅されるのは教科書通りだった。特に二階微分は、傾きが緩やかな区間では曲率 κ の近似値として振る舞う（κ ≈ d²r/dθ² / (1+(dr/dθ)²)^1.5 の分子が支配的になるため）。理想的なインボリュート曲線であればd²r/dθ²はなだらかに変化するはずで、鋭いスパイクは「本来なめらかであるべき場所になめらかでない何かがある」ことの直接的な証拠になる。

ピッチ円付近に現れた異常

解析したのは、8T・モジュール0.5相当のPETGギヤ1個体（2026年9月15日撮影）。写真1枚から1歯を切り出し、B面（以前G-code解析で確認した、印刷時のノズルが歯先を過ぎたあと通る側）に絞って見ると、ピッチ円から歯先までの区間（この個体では概ね-145.6°〜-135.7°）のうち、ピッチ円のすぐ外側、-143.7°付近だけに二階微分の突出したスパイクが現れた。

r(θ)（上段）を見る限り目立った異常はないが、一階微分（中段）では滑らかな減少カーブに一度だけ明確な凹みがあり、二階微分（下段）ではその位置に約+25500から-18400まで跳ぶスパイクが立っている。

最初はこれを印刷条件のせいだと疑った。

まずフロー補正の切替点を疑った（歯車では偶数歯・奇数歯でフロー倍率を変える補正を使うことがあるため）。しかし今回のギヤは圧力角のみを補正しており、フロー補正自体を使っていないと確認できたため除外。
次にZシームの位置を疑った。スライサーのシーム（層の継ぎ目）が同じ場所に集中すれば段差になり得る。しかしシームは歯の根元に設定しており、ピッチ円付近に出る理由にならないため除外。
印刷条件を一つずつ消していってなお異常が残ったため、インボリュート曲線と歯元のフィレット曲線が接続する、設計上避けられない曲率不連続点という仮説が最後に残った。

これを後押ししたのが対称性だった。歯先付近は微分値の暴れが小さいのに対し、歯先から見て±10〜13°付近ではA面・B面の両方に同時に暴れの幅が広がる区間が現れた。もし原因が印刷や測定のノイズなら、左右のフランクでたまたま同じ角度・同じ距離のところで揃って荒れる理由がない。左右対称に現れたという事実が、これが偶然ではなく歯形に内在する特徴だと考える最大の根拠になった。

PETGとPOMを重ねる

同じ手法を市販POMギヤの写真（2026年3月19日撮影、同じく8T相当）にも適用し、歯先を基準に位置を揃えて重ねた。

r(θ)自体はほぼ同じ形の滑らかな曲線を描くが、二階微分で見ると異常が出る位置がずれている。POMの異常はピッチ円のすぐ外側（-10°〜-9.5°付近）、PETGの異常はB面の中間あたり（-8°〜-7.5°付近）で発生しており、振幅もPOM（-64000〜+21000）の方がPETG（-18000〜+26000）より大きい。同じ場所で起きているわけではなく、素材や印刷方式による歯面形成メカニズムの違いを反映している可能性がある。

ただしこれは1枚の写真・1本の歯だけの比較であり、撮影倍率が同一である保証もない。この時点で「性能差」を語るのは早計で、曲率解析の指標と実際の消費電流・騒音との相関はまだ一度も検証していない。ここは今後の宿題として残る。

平滑化のせいではないか、を疑う

ここで一つ気になる点が出てきた。二階微分を全区間で見ると、実はピッチ円付近に限らず、どこもそれなりに凸凹している。だとすると「ピッチ円付近だけが特別」という見立て自体が、Savitzky-Golayという特定の平滑化手法の癖で生まれた見かけの現象ではないか、という疑いが生じる。これを二つの別の手法で検証した。

円弧近似（ArcWelder風）。3Dプリンタのスライサーが生成するG-codeは標準では直線移動（G1）の集合でしかなく、曲線も細かい多角形として近似されている。これを後から円弧（G2/G3）に変換するArc Welderという仕組みがあるが、それを模して、連続する点に貪欲に円を当てはめ、許容誤差内で区間をつないでいく処理を実装した。

def fit_circle(xs, ys):
    # 代数的円フィット（Kasa法）
    A = np.c_[2*xs, 2*ys, np.ones(len(xs))]
    b = xs**2 + ys**2
    (a, bb, c), *_ = np.linalg.lstsq(A, b, rcond=None)
    return a, bb, np.sqrt(c + a**2 + bb**2)

結果、許容誤差を0.15px〜2.0pxまで変えても、二階微分の荒れはほとんど改善しなかった（生データとほぼ同じ-55万〜+34万のオーダー）。原因は、隣り合う円弧同士の接線方向を一致させる制約を入れていないため、円弧の継ぎ目で傾きが不連続に跳ぶこと。実際のArc Welderも同じ構造で、G1の集まりを円弧に変換するだけであり、接続点での接線連続性までは保証しない。

B-スプライン（5次）。有限差分を使わず、スプラインの導関数を直接計算できる。

from scipy.interpolate import UnivariateSpline
spl = UnivariateSpline(theta_tooth, r_tooth_raw, k=5, s=s_val)  # s_val = 係数 × 点数 × 残差分散
r_bs, dr_bs, d2r_bs = spl(theta_tooth), spl.derivative(1)(theta_tooth), spl.derivative(2)(theta_tooth)

平滑化パラメータsを「係数×点数×残差分散」の形で3段階（0.3倍・1.0倍・3.0倍）試したところ、0.3倍（節点数52、ほぼ生データを補間）では暴れが1階微分で3200、2階微分で120万を超える爆発が起きたが、1.0倍・3.0倍（節点数23前後）ではSG平滑化とほぼ重なる結果になり、±10〜13°付近の特徴も同じ場所に同じ形で現れた。

有限差分（SG＋np.gradient）と解析的微分（B-スプライン）という、全く異なる数学的基盤で同じ特徴が再現されたことは、この特徴が特定の微分手法のアーティファクトではないという、追加の裏付けになった。

ただし一つ落とし穴があった。B面だけ（56点）に範囲を絞って、全歯（426点）で使ったのと同じ倍率のsをそのまま当てはめたところ、節点数がわずか2つまで潰れ、ノイズと一緒に肝心のスパイクまで消えてしまった。

「滑らかになった」ことと「正しく捉えられている」ことは別問題で、平滑化パラメータは対象データの点数・密度ごとに調整し直す必要がある。円弧近似・B-スプラインという二つの寄り道は、結果としてこの一点を繰り返し確認する検証になった。

gear_analyzerへの実装

以上の手順を、既存のgear_analyzerに機能として組み込んだ。このツールは既に歯先方向を基準に正規化したmm座標（tooth_contours、各歯についてx,y配列を保持）を持っていたため、新たな画像処理は不要で、既存の解析結果をそのまま数値微分に回すだけで済んだ。サイドバーに「微分グラフ」ボタンを追加し、歯番号とA面／B面／両面を切り替えながらr(θ)・dr/dθ・d²r/dθ²の3段グラフを確認できるようにしている。

余談：直径が縮んだ話

本題とは別に、圧力角のみを補正し直径全体は無補正のまま印刷した8Tギヤ（設計直径5.1mm）が、実測4.87mmだったという事例も同じ日の作業で出てきたので記録しておく。無補正である以上、これは失敗ではなく「このプリンタ・この材料・この形状における素の収縮率（約4.5%）」を示す基準データと言える。フローを変えていなくても、スライサーのライン幅設定と実効ビード幅のずれだけで直径が系統的に縮むことがある点は、見落としやすい要因として記録しておく価値がある。単純な円柱を同条件で印刷して比較すれば、収縮率由来か、小フィーチャー特有の誤差かを切り分けられるはずで、これは次回試したい。

今回使ったパラメータ（再現用メモ）

項目	値
対象PETGギヤ写真	`WIN_20260915_20_44_15_Pro.jpg`
対象POMギヤ写真	`WIN_20260319_04_35_17_Pro.jpg`
二値化閾値	40
モルフォロジー処理	OPEN→CLOSE、カーネル3×3
角度再サンプリング数	2000点（1歯あたり426点前後）
Savitzky-Golay	窓幅31、次数3
ピッチ円半径の概算	歯元からの比率 = 歯元丈:歯末丈 = 1.25:1（標準比率、実モジュールでは未較正）
円弧近似の許容誤差	0.15px／0.5px／1.0px／2.0px
B-スプライン	次数k=5、s = 係数(0.3/1.0/3.0) × 点数 × 残差分散

この表と本文中のコード片があれば、同じ手順は再現できるはずである。次に試すべきことは、①異なる個体・複数歯での再現性確認、②曲率指標と消費電流・騒音の相関検証、③単純な円柱での収縮率の切り分け、の3つ。

忘備録なんで、qiitaにグラフを張り付けてあるから、qiitaも読むこと。