ベジエ曲面のスライス断面抽出

任意の平面でベジエ曲面を切ったときの断面（等高線）を、ベジエクリッピング法の考え方を使って求める。
このページでは bezier_surface_slice.html が内部で行っている4段階の処理を、図とコードで順に説明する。

00概要・全体パイプライン

入力は4×4個の制御点からなる1枚のベジエパッチと、1枚の切断平面。出力はその平面上に落ちる断面を近似する線分（ポリライン）の集合。この2つの間を、次の4段階が繋いでいる。

符号が混在する子パッチに再帰 片側 ベジエパッチ (4×4) 凸包判定 f = n·p − d 棄却（枝刈り） 4分割 de Casteljau 末端パッチ depth = maxDepth マーチング スクエア 断面 線分列
図0. 全体パイプライン。片側判定で棄却できないパッチだけを4分割して再帰し、既定の深度に達した末端パッチだけをマーチングスクエアで実際に線分化する。

肝心なのは、4分割の再帰と棄却判定を行っているのが1〜2、実際に交線を取り出しているのが3という役割分担。1〜2はどこに交差の可能性があるかを絞り込むだけで、実座標は生成しない。3で初めて座標を持った線分が生まれる。

添付ファイルとの関係
元にした bezier_clip4.html の clipBezier() は、視錐台6平面に対する凸包判定＋4分割再帰で「画面内に入るパッチ」を絞り込むものだった。本ツールはこの枝刈り＋再帰分割という骨格をそのまま流用し、判定対象を「6平面のAND」から「1枚の平面」に置き換え、さらに末端パッチでの交線抽出（マーチングスクエア）を新たに追加している。

01距離関数と凸包性質

交差判定の入口。切断平面からの符号付き距離を求め、パッチの制御点を調べるだけで「交差の可能性があるか」を安価に判定する。

切断平面を単位法線 n とオフセット d で表すと、曲面上の点 p の平面からの符号付き距離は次の一次式になる。

f(p) = n·p − d　　// f>0: 法線側　f<0: 反対側　f=0: 平面上

ここで使うのがベジエ曲面の凸包性質：曲面上のどんな点 S(u,v) も、必ずその16個の制御点が作る凸包の内側に収まる、という性質。

凸包性質（Convex Hull Property）
ベジエ基底関数（バーンスタイン基底）は常に非負かつ総和が1になるため、曲面上の任意の点は制御点の凸結合（重み付き平均）として表される。したがって曲面は制御点群の凸包を決してはみ出さない。

この性質から、次の判定がそのまま成り立つ。

  • 全制御点が f>0（または全部 f<0）→ 凸包全体が平面の片側に収まる → 曲面はこの範囲では絶対に平面と交差しない
  • f の符号が混在 → 凸包が平面をまたぐ → 曲面が交差する可能性がある（確定ではない）ので次の段階へ
f = 0（切断面） 全制御点が片側 → 棄却 f = 0（切断面） 符号が混在 → 要再分割
図1. わかりやすさのため2次元の曲線・4制御点で図示しているが、実装は4×4=16個の3D制御点全部に対して同じ判定を行う。オレンジの点は反対側にある制御点。
制御点・曲線（同じ側） 反対側の制御点／切断面

実装では、パッチの16点を1つの配列に平坦化し、全点の f の符号を1回のループで走査するだけ。

// P: 4x4 の制御点配列, plane: {n, d}
const pts = P[0].concat(P[1], P[2], P[3]);
let allPos = true, allNeg = true;
for (let i = 0; i < pts.length; i++) {
  const f = vdot(plane.n, pts[i]) - plane.d;
  if (f >  1e-9) allNeg = false;
  if (f < -1e-9) allPos = false;
}
if (allPos || allNeg) return; // ★ 凸包が片側に完全収容 → 交差なし（枝刈り）

02再帰的な4分割

符号が混在した「わからないパッチ」を de Casteljau 分割で4つに割り、それぞれへ同じ判定を繰り返す。棄却できた領域は二度と調べない。

符号が混在した場合、パッチを u=0.5・v=0.5 で4分割する（添付ファイルの split4() をそのまま使用）。4つの子パッチはそれぞれ元のUV範囲の 1/4 を担当し、深度 depth+1 として再び①の判定にかけられる。

split4(P) → [左上, 左下, 右上, 右下] 　各子パッチ: u ∈ [0,½] or [½,1]　　v ∈ [0,½] or [½,1]

これを既定の深度 maxDepth に達するまで繰り返す。片側判定で棄却されたパッチは以後まったく処理されないので、実際に細かく分割されるのは切断面の近くだけになる（下図でオレンジ色の部分だけがマーチングスクエアまで進む）。

depth 1 で棄却 depth 1: 混在 → 再分割 depth 2 = maxDepth: 末端パッチ（マーチングスクエアへ）
図2. 概念図。灰色（pruned）は片側判定で棄却された領域で以後分割されない。青の切断曲線が通る領域だけがオレンジ（末端パッチ）まで細分化され、③のマーチングスクエアに渡される。
棄却（片側） 末端パッチ（depth=maxDepth） 真の交線（説明用） 
function slicePatchRecursive(P, plane, depth, maxDepth, leafGridN, segsOut) {
  // …① 凸包の片側判定（前節）…
  if (depth >= maxDepth) {
    // ③ 末端パッチ: マーチングスクエアへ（次節）
    return;
  }
  const [tl, bl, tr, br] = split4(P);   // de Casteljau で4分割
  slicePatchRecursive(tl, plane, depth+1, maxDepth, leafGridN, segsOut);
  slicePatchRecursive(bl, plane, depth+1, maxDepth, leafGridN, segsOut);
  slicePatchRecursive(tr, plane, depth+1, maxDepth, leafGridN, segsOut);
  slicePatchRecursive(br, plane, depth+1, maxDepth, leafGridN, segsOut);
}
添付ファイルとの対応
この再帰の骨格は bezier_clip4.html の clipBezier() と完全に同じ形（凸包判定 → IN/OUT/PARTの3値分岐 → PARTだけsplit4()で4分割再帰）。違いは判定対象を「視錐台6平面」から「1枚の切断平面」に変えた点と、末端到達時の処理（元のファイルはそのままパッチを採用するだけだったが、本ツールはここでさらにマーチングスクエアに進む）だけ。

03マーチングスクエア法で線分を取り出す

①②は「どこに交差があるか」を絞り込むだけで座標を作らない。実際の交線（座標つき線分）を作るのはここだけ。

再帰が maxDepth に達した末端パッチは、もう十分小さいとみなし、その内部を leafGridN × leafGridN 個の格子に分けて実座標を評価する。

grid[i][j] = evalPatch(P, i/K, j/K) f[i][j] = n · grid[i][j] − d　　　（K = leafGridN）

隣接する4点（1セル）ごとに、4隅の符号パターンから交線がどの辺を横切るかを判定する——これがマーチングスクエア法。パターンは4隅の +/− の組み合わせで16通りあるが、本質的には次の3タイプに集約できる。

① 角が1つだけ違う ② 隣同士が同符号 ③ 対角が同符号（鞍点・要注意）
図3. 4隅の符号パターン（16通り）は本質的にこの3種に集約される。緑の点が線形補間で求めた辺上の交点、太線が抽出される線分。③の鞍点だけは2通りの繋ぎ方があり得るため注意が必要（後述）。
f>0 の角 f<0 の角 辺上の交点／抽出される線分

辺の両端で符号が反転していれば、その辺上のどこかに f=0 の点がある。位置は3D座標を直接線形補間して求める（座標の線形補間であって、距離関数の再計算はしない）。

function edgeInterp(pA, fA, pB, fB) {
  const t = fA / (fA - fB);     // f=0 になる位置（0..1）
  return vlerp(pA, pB, t);   // 3D座標を直接補間
}

1セルにつき出てくる線分は0〜2本。0本(全部同符号)/1本(①②)はそのまま処理できるが、③の鞍点（対角が同符号）だけは「どちらの対角同士を繋ぐか」が2通りあって決められない——古典的なマーチングスクエアのあいまいさである。ここではセル中心を追加サンプリングし、中心の符号がどちらの対角グループと同じかで自然な繋ぎ方を選ぶ。

const center = (f00 + f10 + f11 + f01) / 4;
case 5: // (00,11)が+ / (10,01)が− の鞍点
  if (center > 0) { /* 00-11 側で繋ぐ2本 */ }
  else          { /* 10-01 側で繋ぐ2本 */ }
  break;

末端パッチ内の全セルについてこれを行い、得られた線分をすべて集めたものが、そのスライス面での断面全体になる。

04精度・限界・パラメータ

出力の正体、精度を決める2つのパラメータ、そして「本物の」ベジエクリッピング法との違いをまとめる。

出力は線分（ポリライン）であって曲線ではない

①②で絞り込み、③で座標を作る——この設計の帰結として、最終的に得られるのは [[p1,p2], [p3,p4], …] という3D線分の配列である。真の交線をスプラインなどで曲線としてフィッティングする後処理は行っていない。見た目が滑らかになるのは、線分が十分細かく敷き詰められて折れ目が視認できなくなるからにすぎない。

粗い（depth・leafGridN 小）— 折れ目が見える 細かい（depth・leafGridN 大）— 曲線に見える
図4. どちらも中身は直線分の集合。点線が仮想的な「真の交線」、太線が実際に出力される折れ線。右のように十分細かければ目には曲線に見えるが、データとしては常にポリライン。

2つのパラメータ

  • maxDepth（再帰分割深度）— ②の4分割を何段まで許すか。大きいほど末端パッチが細かく切断面に密着するが、パッチ数は最悪 4^maxDepth 倍に増えうる（実際は片側判定の枝刈りでその大半は発生しない）。
  • leafGridN（末端グリッド解像度）— 末端パッチ内をさらに何×何のマーチングスクエア格子に分けるか。大きいほど1本の末端パッチから出る線分が増え、境界の曲率をより滑らかに追える。

両者はどちらも「線分をどれだけ細かく刻むか」を制御する点で似ているが、maxDepthは境界へ収縮していく段数、leafGridNは収縮しきった末端パッチの中の解像度という役割分担になっている。

古典的な Bezier Clipping 法との違い

誠実な注記
Sederberg & Nishita が提案した本来の Bezier Clipping は、距離関数 f(u,v) 自体を1次元のベジエ曲線（制御多角形）とみなし、その凸包を u 軸（または v 軸）に投影して UV区間そのものを直接収縮させる（fat line 法）。1回の反復で区間が数分の一に縮むため収束が速い。

本実装（および元にした bezier_clip4.html）は、それよりシンプルな「片側判定による枝刈り＋一様な4分割再帰（quadtree）」という変種である。実装は単純だが、区間を積極的に収縮させるわけではないため、同じ精度を得るのに古典的な手法より多くの分割段数を要する場合がある。用途（インタラクティブな可視化）には十分な速度が出ているが、この違いは把握しておく価値がある。

退化ケース：接線状の切断面

対称な形状（トーラスなど）で切断面がちょうど制御パッチの境界と一致すると、全制御点の f が厳密に0になり「両側判定」が同時に真になって交差なしと誤判定されることがある。これを避けるため、スライス面群の中心位置に span×0.0007 程度のごく小さいオフセットを加え、厳密な対称一致を避けている。