shigureni.dev

等値面 — マーチングスクエア

レイマーチが光線で面を探すのに対し、場の値から直接、面の輪郭線を抜く手がある。各セルの 4 隅で値がしきい値より上か下かを判定すると、4 隅 × 2 値で 16 通りのパターンになる。そのパターンごとに、どの辺とどの辺を線で結ぶかが決まっている。値が分かれている辺の上で交点を線形補間して求め、辺どうしを結ぶ。

マーチングスクエア。連続な場 f(x, y) の f = iso を満たす曲線を、格子の各セルの局所判定だけで抜く手法。セル4隅の上下だけで線がそのセルを横切るかと通る辺が決まるので、曲線追跡なしに並列に処理できる。三次元の等値「面」を抜くマーチングキューブの二次元版で、医療画像の表面抽出やメタボールの輪郭生成に使われる。

場は metaballs の影響の合計にする。各球が距離 d の点に r² / d² の影響を投げ、その合計が (x, y) の値になる。球がいくつか漂っていて、各セルでその合計を明るさに塗る。しきい値 1.0 を境に、超えた側を明るく、下を暗く塗り分ける。

const influence = (px, py, blobs) => {
  let sum = 0
  for (const b of blobs) {
    const dx = px - b.x
    const dy = py - b.y
    sum += (b.r * b.r) / (dx * dx + dy * dy + 1)
  }
  return sum
}

export default (c) => {
  const w = c.canvas.width
  const h = c.canvas.height
  const n = 64
  const cellW = w / n
  const cellH = h / n
  let t = 0

const blobs = [
    { x: w * 0.3, y: h * 0.4, r: 70, vx: 2.1, vy: 1.3 },
    { x: w * 0.7, y: h * 0.5, r: 70, vx: -1.7, vy: 2.0 },
    { x: w * 0.5, y: h * 0.7, r: 56, vx: 1.5, vy: -2.2 },
  ]

return () => {
    for (const b of blobs) {
      b.x += b.vx
      b.y += b.vy
      if (b.x < b.r || b.x > w - b.r) b.vx = -b.vx
      if (b.y < b.r || b.y > h - b.r) b.vy = -b.vy
      b.x = Math.max(b.r, Math.min(w - b.r, b.x))
      b.y = Math.max(b.r, Math.min(h - b.r, b.y))
    }

for (let gy = 0; gy < n; gy++) {
      for (let gx = 0; gx < n; gx++) {
        const v = influence((gx + 0.5) * cellW, (gy + 0.5) * cellH, blobs)
        c.fillStyle = v > 1 ? colors.tint(210) : colors.tint(30)
        c.fillRect(gx * cellW, gy * cellH, cellW + 1, cellH + 1)
      }
    }
  }
}

明るい側と暗い側の境目が、抜きたい等値線の通り道になる。格子で塗ると境目は階段状で、まだ線にはなっていない。

その境目を線分として引く。格子点ごとに影響の合計を一度だけ計算しておき、各セルの4隅 tl tr br bl がそれぞれしきい値の上か下かを4ビットの番号 code にまとめる。全隅が上（15）か全隅が下（0）のセルは線が通らないので飛ばす。残りは上下が分かれている辺を線が横切る。その横切る位置を辺の両端の値の差で線形補間して求め、辺どうしを結ぶ。

const influence = (px, py, blobs) => {
  let sum = 0
  for (const b of blobs) {
    const dx = px - b.x
    const dy = py - b.y
    sum += (b.r * b.r) / (dx * dx + dy * dy + 1)
  }
  return sum
}

const cross = (pa, pb, va, vb, iso) => {
  const d = vb - va
  const f = Math.abs(d) < 1e-6 ? 0.5 : (iso - va) / d
  return pa + (pb - pa) * Math.min(1, Math.max(0, f))
}

const seg = (c, ax, ay, bx, by) => {
  c.moveTo(ax, ay)
  c.lineTo(bx, by)
}

export default (c) => {
  const w = c.canvas.width
  const h = c.canvas.height
  const n = 56
  const cellW = w / n
  const cellH = h / n
  const iso = 1

const blobs = [
    { x: w * 0.3, y: h * 0.4, r: 70, vx: 2.1, vy: 1.3 },
    { x: w * 0.7, y: h * 0.5, r: 70, vx: -1.7, vy: 2.0 },
    { x: w * 0.5, y: h * 0.7, r: 56, vx: 1.5, vy: -2.2 },
  ]

const grid = new Float32Array((n + 1) * (n + 1))

for (let gy = 0; gy <= n; gy++) {
      for (let gx = 0; gx <= n; gx++) {
        grid[gy * (n + 1) + gx] = influence(gx * cellW, gy * cellH, blobs)
      }
    }

c.fillStyle = colors.paper
    c.fillRect(0, 0, w, h)
    c.strokeStyle = colors.tint(40)
    c.lineWidth = 2
    c.beginPath()

for (let gy = 0; gy < n; gy++) {
      for (let gx = 0; gx < n; gx++) {
        const tl = grid[gy * (n + 1) + gx]
        const tr = grid[gy * (n + 1) + gx + 1]
        const br = grid[(gy + 1) * (n + 1) + gx + 1]
        const bl = grid[(gy + 1) * (n + 1) + gx]
        let code = 0
        if (tl > iso) code |= 8
        if (tr > iso) code |= 4
        if (br > iso) code |= 2
        if (bl > iso) code |= 1
        if (code === 0 || code === 15) continue

const x = gx * cellW
        const y = gy * cellH
        const top = cross(x, x + cellW, tl, tr, iso)
        const right = cross(y, y + cellH, tr, br, iso)
        const bot = cross(x, x + cellW, bl, br, iso)
        const left = cross(y, y + cellH, tl, bl, iso)

if (code === 1 || code === 14) seg(c, x, left, bot, y + cellH)
        else if (code === 2 || code === 13) seg(c, bot, y + cellH, x + cellW, right)
        else if (code === 3 || code === 12) seg(c, x, left, x + cellW, right)
        else if (code === 4 || code === 11) seg(c, top, y, x + cellW, right)
        else if (code === 6 || code === 9) seg(c, top, y, bot, y + cellH)
        else if (code === 7 || code === 8) seg(c, x, left, top, y)
        else if (code === 5) {
          seg(c, x, left, top, y)
          seg(c, bot, y + cellH, x + cellW, right)
        } else if (code === 10) {
          seg(c, top, y, x + cellW, right)
          seg(c, x, left, bot, y + cellH)
        }
      }
    }

c.stroke()
  }
}

球を囲む滑らかな閉じた輪郭が出て、球が近づくと二つの輪が一本に繋がる。cross が辺を横切る位置で、(iso - va) / d は辺の片端 va からしきい値までが辺全体の何割かを返し、その割合だけ端から進めた点が交点になる。場が滑らかなので隣り合うセルの線分が端で繋がり、一本の連続した輪郭になる。code が 5 と 10 のときだけ対角に二本入るのは、上の隅と下の隅が斜めに分かれていて線が二回横切るため。

塗りが面で見せるところを、マーチングスクエアは同じしきい値の輪郭だけで抜く。r² / d² を足し合わせた影響の合計という場に、しきい値で線を引いている。smin で繋いだ円が距離で融合を作るのと別ルートで、影響を足し合わせる metaballs も同じ等値面の読み方に届く。