實現 O(n) 時間的高斯模糊算法

原發於個人我的博客 我的博客網站

本文主要是針對高斯模糊算法進行優化，最後在線性時間內實現高斯模糊效果。固然，該算法並不是本人原創，實現過程當中也借鑑了一些文章和論文，相關連接都在文末貼出

搭配閱讀，我寫了一個簡單的Demo，Demo連接，Demo代碼地址，能夠在Demo中測試各類模糊效果及其耗時

高斯模糊

高斯模糊能夠實現很是平滑的模糊效果，對於每個像素，它會取它自己及其周圍必定範圍內的像素，去計算模糊後的效果，固然，每一個像素在計算過程當中所佔的權重並不相同，距離其自己越近，權重越高——這固然是合理的。權重的計算公式以下：

公式中的σ(sigma)是人爲定值的，值越大，生成的圖像越模糊

固然，並不是全部像素都會互相影響，以當前像素爲圓心，取一個半徑，在這個半徑之內的像素纔會影響該像素的模糊效果。這個半徑與sigma成正比，不過沒有統一的公式，在這個提問下，有回答說NVidia採用的是radius=sigma*3的計算方式，我在這也使用了這個公式

下面是其基本實現，時間複雜度爲O(nrr)，n爲圖片大小

function gaussianBlur(src, dest, width, height, sigma) {
  const radius = Math.round(sigma * 3) // kernel size
  for (let i = 0; i < width; i++) {
    for (let j = 0; j < height; j++) {
      let accumulation = 0
      let weightSum = 0
      for (let dx = -radius; dx <= radius; dx++) {
        for (let dy = -radius; dy <= radius; dy++) {
          const x = Math.min(width - 1, Math.max(0, i + dx))
          const y = Math.min(height - 1, Math.max(0, j + dy))
          // calc weight
          const weight =
            Math.exp(
              -(Math.pow(dx, 2) + Math.pow(dy, 2)) / (2 * Math.pow(sigma, 2))
            ) /
            (Math.PI * 2 * Math.pow(sigma, 2))
          accumulation += src[y * width + x] * weight
          weightSum += weight
        }
      }
      dest[j * width + i] = Math.round(accumulation / weightSum)
    }
  }
}

效果以下，也能夠去 Demo 中嘗試(高斯模糊可能耗時較久，耐心等待吧)

原圖：

高斯模糊後(sigma=5)，5411ms(MBP, 13-inch, 2017, 8G)：

方框模糊

最簡單的方框模糊(box blur)沒有權重，只要在半徑內，權重都是相同的

function simpleBoxBlur(src, dest, width, height, radius) {
  for (let i = 0; i < width; i++) {
    for (let j = 0; j < height; j++) {
      let accumulation = 0
      for (let dx = -radius; dx <= radius; dx++) {
        for (let dy = -radius; dy <= radius; dy++) {
          const x = Math.min(width - 1, Math.max(0, i + dx))
          const y = Math.min(height - 1, Math.max(0, j + dy))
          accumulation += src[y * width + x]
        }
      }
      dest[j * width + i] = Math.round(
        accumulation / Math.pow(2 * radius + 1, 2)
      )
    }
  }
}

radius=15的效果(775ms)：

能夠看出，雖然耗時較高斯模糊短，但其效果是很不平滑的。那有沒有辦法讓方框模糊擁有高斯模糊同樣的品質呢，這篇論文提出了方案

根據論文中所述，能夠經過屢次的方框模糊實現高斯模糊，而這屢次方框模糊的核由如下步驟得出：

• 由下圖公式(n爲模糊次數)計算得出wideal，再計算wl和wu，wl爲第一個小於等於wideal的奇數，wu是第一個大於等於wideal的奇數(很明顯wu = wl+2)(要求奇數的緣由是size=radius*2+1)

• 計算m:

• 前m次用wl做核的大小，其他的n-m次用wu做核的大小

化做代碼：

function genKernelsForGaussian(sigma, n) {
  const wIdeal = Math.sqrt((12 * Math.pow(sigma, 2)) / n + 1)
  let wl = Math.floor(wIdeal)
  if (wl % 2 === 0) {
    wl--
  }
  const wu = wl + 2
  let m =
    (12 * Math.pow(sigma, 2) - n * Math.pow(wl, 2) - 4 * n * wl - 3 * n) /
    (-4 * wl - 4)
  m = Math.round(m)
  const sizes = []
  for (let i = 0; i < n; i++) {
    sizes.push(i < m ? wl : wu)
  }
  return sizes
}

function boxBlur(src, dest, width, height, sigma) {
  const kernels = genKernelsForGaussian(sigma, 3)
  // radius * 2 + 1 = kernel size
  simpleBoxBlur(src, dest, width, height, (kernels[0] - 1) / 2)
  // 注意這裏要顛倒 src 和 dest 的順序
  simpleBoxBlur(dest, src, width, height, (kernels[1] - 1) / 2)
  simpleBoxBlur(src, dest, width, height, (kernels[2] - 1) / 2)
}

效果和高斯模糊基本一致(227ms)：

雖然時間複雜度與高斯模糊相同，但獲得的半徑更小，在sigma=5，n=3的狀況下，radius=[4, 4, 5]，高斯模糊則是15；同時還不用計算複雜的weight，因此從實際結果上看，速度要優於比高斯模糊

水平&垂直模糊

在這篇文章中提到 方框模糊能夠用水平模糊+垂直模糊替代實現。要實現水平/垂直模糊很簡單, 只需考慮水平/垂直線上的像素便可：

// horizontal motion blur
function hMotionBlur(src, dest, width, height, radius) {
  for (let i = 0; i < width; i++) {
    for (let j = 0; j < height; j++) {
      let accumulation = 0
      for (let dx = -radius; dx <= radius; dx++) {
        const x = Math.min(width - 1, Math.max(0, i + dx))
        accumulation += src[j * width + x]
      }
      dest[j * width + i] = Math.round(accumulation / (2 * radius + 1))
    }
  }
}

// vertical motion blur
function vMotionBlur(src, dest, width, height, radius) {
  for (let i = 0; i < width; i++) {
    for (let j = 0; j < height; j++) {
      let accumulation = 0
      for (let dy = -radius; dy <= radius; dy++) {
        const y = Math.min(height - 1, Math.max(0, j + dy))
        accumulation += src[y * width + i]
      }
      dest[j * width + i] = Math.round(accumulation / (2 * radius + 1))
    }
  }
}

應用此優化的模糊算法，時間複雜度能夠優化到O(nr)：

function _mutantBoxBlur(src, dest, width, height, radius) {
  hMotionBlur(dest, src, width, height, radius)
  vMotionBlur(src, dest, width, height, radius)
}

function mutantBoxBlur(src, dest, width, height, sigma) {
  const boxes = genKernelsForGaussian(sigma, 3)
  for (let i = 0; i < src.length; i++) {
    dest[i] = src[i]
  }
  _mutantBoxBlur(src, dest, width, height, (boxes[0] - 1) / 2)
  _mutantBoxBlur(src, dest, width, height, (boxes[1] - 1) / 2)
  _mutantBoxBlur(src, dest, width, height, (boxes[2] - 1) / 2)
}

效果以下(82ms)：

終極優化

以水平模糊爲例(這裏爲了方便，用二維數組表示)

而

則

那很明顯，咱們能夠經過此優化將時間複雜度由O(nr)降到O(n)，最終的代碼:

// horizontal fast motion blur
function hFastMotionBlur(src, dest, width, height, radius) {
  for (let i = 0; i < height; i++) {
    let accumulation = radius * src[i * width]
    for (let j = 0; j <= radius; j++) {
      accumulation += src[i * width + j]
    }
    dest[i * width] = Math.round(accumulation / (2 * radius + 1))
    for (let j = 1; j < width; j++) {
      const left = Math.max(0, j - radius - 1)
      const right = Math.min(width - 1, j + radius)
      accumulation =
        accumulation + (src[i * width + right] - src[i * width + left])
      dest[i * width + j] = Math.round(accumulation / (2 * radius + 1))
    }
  }
}

// vertical fast motion blur
function vFastMotionBlur(src, dest, width, height, radius) {
  for (let i = 0; i < width; i++) {
    let accumulation = radius * src[i]
    for (let j = 0; j <= radius; j++) {
      accumulation += src[j * width + i]
    }
    dest[i] = Math.round(accumulation / (2 * radius + 1))
    for (let j = 1; j < height; j++) {
      const top = Math.max(0, j - radius - 1)
      const bottom = Math.min(height - 1, j + radius)
      accumulation =
        accumulation + src[bottom * width + i] - src[top * width + i]
      dest[j * width + i] = Math.round(accumulation / (2 * radius + 1))
    }
  }
}

function _fastBlur(src, dest, width, height, radius) {
  hFastMotionBlur(dest, src, width, height, radius)
  vFastMotionBlur(src, dest, width, height, radius)
}

function fastBlur(src, dest, width, height, sigma) {
  const boxes = genKernelsForGaussian(sigma, 3)
  for (let i = 0; i < src.length; i++) {
    dest[i] = src[i]
  }
  _fastBlur(src, dest, width, height, (boxes[0] - 1) / 2)
  _fastBlur(src, dest, width, height, (boxes[1] - 1) / 2)
  _fastBlur(src, dest, width, height, (boxes[2] - 1) / 2)
}

效果(20ms)：

總結

前後經過屢次方框模糊、水平+垂直模糊代替方框模糊、對水平/垂直模糊計算過程的優化，逐步將高斯模糊的耗時下降，最後獲得了只與圖片大小相關的O(n)時間複雜度算法

參考

• http://elynxsdk.free.fr/ext-d...
• http://blog.ivank.net/fastest...
• https://www.peterkovesi.com/p...
• https://stackoverflow.com/que...
• https://blog.demofox.org/2015...