shigureni.dev

マーチングスクエア — 場から線を抜く

二値化は値を「上か下か」で塗り分ける。その境目だけを線で抜くと等値線になる。場の値がちょうどあるしきい値になる点を繋ぐと、地形図の等高線と同じ一本の輪郭が出る。塗りでなく境界を描くので、面でなく線で絵を見せる処理になる。

マーチングスクエア。連続な場 f(x, y) の f = iso を満たす点を、格子のセル単位の局所判定だけで繋ぐ。セル4隅が iso の上か下かを4ビットの番号にまとめると 16 通りで、線がそのセルをどう横切るかが番号ごとに決まる。二値化が各画素を「上か下か」で塗り分けるのに対し、こちらは上と下の境目だけを線で抜く。同じしきい値を、面でなく境界として読む。ピクセル単位の3次元版がマーチングキューブで、医療画像の等値面抽出やメタボールの表面生成に使われる。

場は sin と cos を数本足した波。各格子点で値が iso より上か下かだけを見て、上を明るい色、下を暗い色で塗り分ける。二値化と同じ操作で、iso を時間で上下させると明暗の境目が場の上を動く。

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

その境目を線分として引く。セルの4隅 tl tr br bl がそれぞれ iso の上か下かを4ビットの番号にまとめると 16 通りになる。全隅が上（15）か全隅が下（0）のセルは線が通らないので飛ばす。残りは上下が分かれている辺を線が横切る。横切る位置を辺の両端の値の差で線形補間して求め、辺どうしを結ぶ。

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 = 40
  const cellW = w / n
  const cellH = h / n
  let t = 0

const field = (gx, gy) =>
    Math.sin(gx * 0.5 + t) * 0.5 + Math.cos(gy * 0.45 - t) * 0.5 + Math.sin((gx + gy) * 0.3) * 0.4

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

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))
  }

return () => {
    t += 0.02
    const iso = Math.sin(t * 0.7) * 0.7

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

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 から iso までが辺全体の何割かを出し、その割合だけ端から進めた点を返す。場が滑らかなので、隣り合うセルの線分が端で繋がって一本の閉じた輪郭になる。code が 5 と 10 のときだけ対角に二本入るのは、対角の隅どうしが斜めに分かれていて線が二回横切るため。この鞍点は繋ぎ方が二通りあって、ここだけ曖昧さが残る。iso が動くと、輪郭が場の尾根や谷を伝って太ったり痩せたりする。

しきい値を一つでなく等間隔に何本も並べると等高線図になる。場に動くガウスの山をいくつか足して、iso を -0.6 から刻んで同じセル判定を各レベルで回す。一本の輪郭が高さ一つに対応する。

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 = 60
  const cellW = w / n
  const cellH = h / n
  const levels = Array.from({ length: 7 }, (_v, i) => -0.6 + i * 0.4)
  let t = 0

const blobs = Array.from({ length: 4 }, (_v, i) => ({
    x: 12 + i * 12,
    y: 30,
    vx: (Math.random() - 0.5) * 0.3,
    vy: (Math.random() - 0.5) * 0.3,
    a: 0.7 + Math.random() * 0.9,
    s: 0.012 + Math.random() * 0.02,
  }))

const field = (gx, gy) => {
    let v = Math.sin(gx * 0.18 + t) * 0.3 + Math.cos(gy * 0.15 - t) * 0.3
    for (const b of blobs) {
      const d2 = (gx - b.x) ** 2 + (gy - b.y) ** 2
      v += b.a * Math.exp(-d2 * b.s)
    }
    return v
  }

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

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))
  }

return () => {
    t += 0.012
    for (const b of blobs) {
      b.x += b.vx
      b.y += b.vy
      if (b.x < 6 || b.x > n - 6) b.vx = -b.vx
      if (b.y < 6 || b.y > n - 6) b.vy = -b.vy
    }

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

c.fillStyle = colors.paper
    c.fillRect(0, 0, w, h)
    c.lineWidth = 1.5

for (let li = 0; li < levels.length; li++) {
      const iso = levels[li]
      c.strokeStyle = li % 2 === 0 ? colors.tint(30) : colors.tint(120)
      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()
    }
  }
}

ガウスの山が動くと、それを囲む等高線が伸び縮みする。二つの山が近づくと、外側の低いレベルの輪郭が先に繋がって一つの島になる。同じセル判定を iso を変えながら回すだけで、一枚の場が標高ごとの輪郭の束になる。二値化が面で塗り分けるところを、その境界線だけで地形を見せる読み方になる。