GPU 加速下的图像视觉
越來越多的移動計算設備都開始攜帶照相機鏡頭，這對于攝影界來說是一個好事情，不僅如此攜帶鏡頭也為這些設備提供了更多的可能性。除了最基本的拍攝功能，結合合適的軟件這些更為強大的硬件設備可以像人腦一樣理解它看到了什么。
僅僅具備一點點的理解能力就可以催生一些非常強大的應用，比如說條形碼識別，文檔識別和成像，手寫文字的轉化，實時圖像防抖，增強現實等。隨著處理能力變得更加強大，鏡頭保真程度更高，算法效率更好，機器視覺 (machine vision) 這個技術將會解決更加重大的問題。
有些人認為機器視覺是個非常復雜的領域，是程序員們的日常工作中絕不會遇到的。我認為這種觀點是不正確的。我發起了一個開源項目?GPUImage，其實在很大程度上是因為我想探索一下高性能的機器視覺是怎么樣的，并且讓這種技術更易于使用。
GPU 是一種理想的處理圖片和視頻的設備，因為它是專門為并行處理大量數據而生的，圖片和視頻中的每一幀都包含大量的像素數據。在某些情況下 GPU 處理圖片的速度可以是 CPU 的成千上百倍。
在我開發 GPUImage 的過程中我學到了一件事情，那就是即使是圖片處理這樣看上去很復雜的工作依然可以分解為一個個更小更簡單的部分。這篇文章里我想將一些機器視覺中常見的過程分解開來，并且展示如何在現代的 GPU 設備上讓這些過程運行地更快。
以下的每一步在 GPUImage 中都有完整的實現，你可以下載包含了 OS X 和 iOS 版本的示例工程 FilterShowcase，在其中體驗一下各個功能。此外，這些功能都有基于 CPU (有些使用了 GPU 加速) 的實現，這些實現是基于 OpenCV 庫的，在另一片文章中 Engin Kurutepe 詳細地講解了這個庫。
索貝爾 (Sobel) 邊界探測
我將要描述的第一種操作事實上在濾鏡方面的應用比機器視覺方面更多，但是從這個操作講起是比較合適的。索貝爾邊界探測用于探測一張圖片中邊界的出現位置，邊界是指由明轉暗的突然變化或者反過來由暗轉明的區域[1]。在被處理的圖片中一個像素的亮度反映了這個像素周圍邊界的強度。
下面是一個例子，我們來看看同一張圖片在進行索貝爾邊界探測之前和之后：
?
正如我上面提到的，這項技術通常用來實現一些視覺效果。如果在上面的圖片中將顏色進行反轉，最明顯的邊界用黑色代表而不是白色，那么我們就得到了一張類似鉛筆素描效果的圖片。
那么這些邊界是如何被探測出來的？第一步這張彩色圖片需要減薄成一張亮度 (灰階) 圖。Janie Clayton 在她的文章中解釋了這一步是如何在一個片斷著色器?(fragment shader) 中完成的。簡單地說這個過程就是將每個像素的紅綠藍部分加權合為一個代表這個像素亮度的值。
有的視頻設備和相機提供的是 YUV 格式的圖片，而不是 RGB 格式。YUV 這種色彩格式已經將亮度信息 (Y) 和色度信息 (UV) 分離，所以如果原圖片是這種格式，顏色轉換這個步驟就可以省略，直接用其中亮度的部分就可以。
圖片一旦減薄到僅剩亮度信息，一個像素周圍的邊界強度就可以由它周圍 3*3 個臨近像素計算而得。在一組像素上進行圖片處理的計算過程涉及到一個叫做卷積矩陣 (參考：convolution matrix) 的東西。卷積矩陣是一個由權重數據組成的矩陣，中心像素周圍像素的亮度乘以這些權重然后再相加就能得到中心像素的轉化后數值。
圖片上的每一個像素都要與這個矩陣計算出一個數值。在處理的過程中像素的處理順序是無關緊要的，所以這種計算很容易并行運行。因此，這個計算過程可以通過一個片斷著色器的方式運行在可編程的 GPU 上，來極大地提高處理效率。正如在?Janie 的文章中所提到的，片斷著色器是一些 C 語言風格的程序，運行在 GPU 上可以進行一些非常快速的圖片處理。
下面這個是索貝爾算子的水平處理矩陣：
| ?1 | 0 | +1 |
| ?2 | 0 | +2 |
| ?1 | 0 | +1 |
為了進行某一個像素的計算，每一個臨近像素的亮度信息都要讀取出來。如果要處理的圖片已經被轉化為灰階圖，亮度可以從紅綠藍任意顏色通道中抽樣。臨近像素的亮度乘以矩陣中對應的權重，然后加到最終值里去。
在一個方向上尋找邊界的過程是這樣的：轉化之后對比一個像素左右兩邊像素的亮度差。如果當前這個像素左右兩邊的像素亮度相同也就是說在圖片上是一個柔和的過度，它們的亮度值和正負權重會相互抵消，于是這個區域不會被判定為邊界。如果左邊像素和右邊像素的亮度差別很大也就是說是一個邊界，用其中一個亮度減去另一個，這種差異越大這個邊界就越強 (越明顯)。
索貝爾過程有兩個步驟，首先是水平矩陣進行，同時一個垂直矩陣也會進行，這個垂直矩陣中的權重如下
| ?1 | ?2 | ?1 |
| 0 | 0 | 0 |
| +1 | +2 | +1 |
兩個方向轉化后的加權和會被記錄下來，它們的平方和的平方根也會被計算出來。之所以要進行平方是因為計算出來的值可能是正值也可能是負值，但是我們需要的是值的量級而不關心它們的正負。有一個好用內建的 GLSL 函數能夠幫助我們快速完成這個過程。
最終計算出來的這個值會用來作為輸出圖片中像素的亮度。因為索貝爾算子會突出顯示兩邊像素亮度的不同的地方，所以圖片中由明轉暗或者相反的突然轉變會成為結果中明亮的像素。
索貝爾邊界探測有一些相似的變體，例如普里維特 (Prewitt) 邊界探測[2]。普里維特邊界探測會在橫向豎向矩陣中使用不同的權重，但是它們運作的基本過程是一樣的。
作為索貝爾邊界探測如何用代碼實現的一個例子，下面是用 OpenGL ES 進行索貝爾邊界探測的片斷著色器：
precision mediump float;varying vec2 textureCoordinate; varying vec2 leftTextureCoordinate; varying vec2 rightTextureCoordinate;varying vec2 topTextureCoordinate; varying vec2 topLeftTextureCoordinate; varying vec2 topRightTextureCoordinate;varying vec2 bottomTextureCoordinate; varying vec2 bottomLeftTextureCoordinate; varying vec2 bottomRightTextureCoordinate;uniform sampler2D inputImageTexture;void main() {float bottomLeftIntensity = texture2D(inputImageTexture, bottomLeftTextureCoordinate).r;float topRightIntensity = texture2D(inputImageTexture, topRightTextureCoordinate).r;float topLeftIntensity = texture2D(inputImageTexture, topLeftTextureCoordinate).r;float bottomRightIntensity = texture2D(inputImageTexture, bottomRightTextureCoordinate).r;float leftIntensity = texture2D(inputImageTexture, leftTextureCoordinate).r;float rightIntensity = texture2D(inputImageTexture, rightTextureCoordinate).r;float bottomIntensity = texture2D(inputImageTexture, bottomTextureCoordinate).r;float topIntensity = texture2D(inputImageTexture, topTextureCoordinate).r;float h = -bottomLeftIntensity - 2.0 * leftIntensity - topLeftIntensity + bottomRightIntensity + 2.0 * rightIntensity + topRightIntensity;float v = -topLeftIntensity - 2.0 * topIntensity - topRightIntensity + bottomLeftIntensity + 2.0 * bottomIntensity + bottomRightIntensity;float mag = length(vec2(h, v));gl_FragColor = vec4(vec3(mag), 1.0); }上面這段著色器中中心像素周圍的像素都有用戶定義的名稱，是由一個自定義的頂點著色器提供的，這么做可以優化減少對移動設備環境的依賴。從 3*3 網格中抽樣出這些命名了的像素，然后用自定義的代碼來進行橫向和縱向索貝爾探測。為簡化計算權重為 0 的部分會被忽略。GLSL 函數?length()?計算出水平和垂直矩陣轉化后值的平方和的平方根。然后這個代表量級的值會被拷貝進輸出像素的紅綠藍通道中，這樣就可以用來代表邊界的明顯程度。
坎尼 (Canny) 邊界探測
索貝爾邊界探測可以給你一張圖片邊界強度的直觀印象，但是并不能明確地說明某一個像素是否是一個邊界。如果要判斷一個像素是否是一個邊界，你要設定一個類似閾值的東西，亮度高于這個閾值的像素會被判定為邊界的一部分。然而這樣并不是最理想的，因為這樣的做法判定出的邊界可能會有好幾個像素寬，并且不同的圖片適合的閾值不同。
這里你更需要一種叫做坎尼邊界探測[3]的邊界探測方法?？材徇吔缣綔y可以在一張圖片中探測出連貫的只有一像素寬的邊界：
坎尼邊界探測包含了幾個步驟。和索貝爾邊界探測以及其他我們接下來將要討論的方法一樣，在進行邊界探測之前首先圖片需要轉化成亮度圖。一旦轉化為灰階亮度圖緊接著進行一點點的高斯模糊，這么做是為了降低傳感器噪音對邊界探測的影響。
一旦圖片已經準備好了，邊界探測就可以開始進行。這里的 GPU 加速過程原本是在 Ensor 和 Hall 的文章 "GPU-based Image Analysis on Mobile Devices" [4]中所描述的。
首先，一個給定像素的邊界強度和邊界梯度要確定下來。邊界梯度是指亮度發生變化最大的方向，也是邊界延伸方向的垂直方向。
為了尋找邊界梯度，我們要用到上一章中的索貝爾矩陣。索貝爾轉化得到的橫豎值加合后就是邊界梯度的強度，這個值會編碼進輸出像素的紅色通道。然后橫向豎向索貝爾結果值會與八個方向 (對應一個中心像素周圍的八個像素) 中的一個結合起來，一個方向上 X 部分值會作為輸出像素的綠色通道值，Y 部分則作為藍色通道值。
這個方法使用的著色器和索貝爾邊界探測使用的類似，只是最后一個計算步驟用下面這段代碼：
vec2 gradientDirection;gradientDirection.x = -bottomLeftIntensity - 2.0 * leftIntensity - topLeftIntensity + bottomRightIntensity + 2.0 * rightIntensity + topRightIntensity;gradientDirection.y = -topLeftIntensity - 2.0 * topIntensity - topRightIntensity + bottomLeftIntensity + 2.0 * bottomIntensity + bottomRightIntensity;float gradientMagnitude = length(gradientDirection);vec2 normalizedDirection = normalize(gradientDirection);normalizedDirection = sign(normalizedDirection) * floor(abs(normalizedDirection) + 0.617316); // Offset by 1-sin(pi/8) to set to 0 if near axis, 1 if awaynormalizedDirection = (normalizedDirection + 1.0) * 0.5; // Place -1.0 - 1.0 within 0 - 1.0gl_FragColor = vec4(gradientMagnitude, normalizedDirection.x, normalizedDirection.y, 1.0);為確保坎尼邊界一像素寬，只有邊界中強度最高的部分會被保留下來。因此，我們需要在每一個切面邊界梯度的寬度之內尋找最大值。
這就是我們在上一步中算出的梯度方向起作用的地方。對每一個像素，我們根據梯度值向前和向后取出最近的相鄰像素，然后對比他們的梯度強度 (邊界明顯程度)。如果當前像素的梯度強度高于梯度方向前后的像素我們就保留當前像素。如果當前像素的梯度強度低于任何一個臨近像素，我們就不再考慮這個像素并且將他變為黑色。
執行這個步驟的著色器如下：
precision mediump float;varying highp vec2 textureCoordinate;uniform sampler2D inputImageTexture; uniform highp float texelWidth; uniform highp float texelHeight; uniform mediump float upperThreshold; uniform mediump float lowerThreshold;void main() {vec3 currentGradientAndDirection = texture2D(inputImageTexture, textureCoordinate).rgb;vec2 gradientDirection = ((currentGradientAndDirection.gb * 2.0) - 1.0) * vec2(texelWidth, texelHeight);float firstSampledGradientMagnitude = texture2D(inputImageTexture, textureCoordinate + gradientDirection).r;float secondSampledGradientMagnitude = texture2D(inputImageTexture, textureCoordinate - gradientDirection).r;float multiplier = step(firstSampledGradientMagnitude, currentGradientAndDirection.r);multiplier = multiplier * step(secondSampledGradientMagnitude, currentGradientAndDirection.r);float thresholdCompliance = smoothstep(lowerThreshold, upperThreshold, currentGradientAndDirection.r);multiplier = multiplier * thresholdCompliance;gl_FragColor = vec4(multiplier, multiplier, multiplier, 1.0); }其中?texelWidth?和?texelHeight?是要處理的圖片中臨近像素之間的距離，lowerThreshold?和?upperThreshold?分別設定了我們預期的邊界強度上下限。
在坎尼邊界探測的最后一步，邊界上出現像素間間隔的地方要被填充，出現間隔是因為有一些點不在閾值范圍之內或者是因為非最大值轉化沒有起作用。這一步會完善邊界使邊界連續起來。
在最后一步中需要考慮一個中心像素周圍的所有像素。如果這個中心像素是最大值，上一步中非最大值轉化就不會影響它，它依然是白色。如果它不是最大值，就會變成黑色。對于中間的灰色像素，會考察它周圍像素的信息。凡是與超過一個白色像素挨著的都會變為白色，相反就會變成黑色。這樣就可以將邊界分離的部分接合起來。
正如你所看到的，坎尼邊界探測會比索貝爾邊界探測更復雜一些，但是它會探測出一條物體邊界的干凈線條。這是線條探測、輪廓探測或者其他圖片分析很好的起點。同時也可以被用來生成一些有趣的美學效果。
哈里斯 (Harris) 邊角探測
雖然利用上一章中的邊界探測技術我們可以獲取關于圖片邊界的信息，我們會得到一張可以直觀觀察到邊界所在位置的圖片，但是并沒有更高層面有關圖片中所展示內容的信息。為了得到這些信息，我們需要一個可以處理場景中的像素然后返回場景中所展示內容的描述性信息的算法。
進行物體探測和匹配時一個常見的出發點是特征探測。特征是指一個場景中具有特殊意義的點，這些點可以唯一的區分出一些結構或者物體。由于邊角的出現往往意味著亮度或者顏色的突然變化，所以邊角常常會作為特征的一種。
在 Harris 和 Stephens 的文章 "A Combined Corner and Edge Detector."[5] 中他們提出一個邊角探測的方法。這個命名為哈里斯邊角探測的方法采用了一個多步驟的方法來探測場景中的邊角。
像我們已經討論過的其他方法一樣，圖片首先需要減薄到只剩亮度信息。通過索貝爾矩陣，普里維特矩陣或者其他相關的矩陣計算出一個像素 X 和 Y 方向上的梯度值，計算出的值并不會合并為邊界的量級。而是將 X 梯度傳入紅色部分，Y 梯度傳入綠色部分，X 與 Y 梯度的乘積傳入藍色部分。
然后對上述計算結果進行一個高斯模糊。從模糊后的照片中取出紅綠藍部分中的編碼過的值，并將值帶入一個計算像素是邊角點可能性的公式：
R = Ix2?× Iy2?? Ixy?× Ixy?? k × (Ix2?+ Iy2)2
其中 Ix?是 X 方向梯度值 (模糊后圖片中紅色部分)，Iy?是 Y 梯度值 (綠色部分)，Ixy?是 XY 值的乘積 (藍色部分)，k 是一個靈敏性常數，R 是計算出來的這個像素是邊角的確定程度。Shi，Tomasi[6] 和 Noble[7] 提出過這種計算的另一種實現方法但是結果其實是十分接近的。
在公式中你可以會覺得頭兩項會抵消掉。但這就是前面高斯模糊那一步起作用的地方。通過在一些像素上分別模糊 X、Y 和 XY 的乘積，在邊角附近就會出現可以被探測到的差異。
我們從 Stack Exchange 信號處理分站中的一個問題中取來一張測試圖片：
經過前面的計算過程得到的結果如下圖：
為了找出邊角準確的位置，我們需要選出極點 (一個區域內亮度最高的地方)。這里需要使用一個非最大值轉化。和我們在坎尼邊界探測中所做的一樣，我們要考察一個中心像素周圍的臨近像素 (從一個像素半徑開始，半徑可以擴大)，只有當中心像素的亮度高于它所有臨近像素時才保留他，否則就將這個像素變為黑色。這樣一來最后留下的就應該是一片區域中亮度最高的像素，也就是最可能是邊角的地方。
通過這個過程，我們現在可以從圖片中看到任意不是黑色的像素都是一個邊角所在的位置：
目前我是使用 CPU 來進行點的提取，這可能會是邊角探測的一個瓶頸，不過在 GPU 上使用柱狀圖金字塔[8]可能會加速這個過程。
哈里斯邊角探測只是在場景中尋找邊角的方法之一。"Machine learning for high-speed corner detection,"[9] 中 Edward Rosten 的 FAST 邊角探測方法是另一個性能更好的方法，甚至可能超越基于 GPU 的哈里斯探測。
霍夫 (Hough) 變換線段探測
筆直的線段是另一種我們會在一個場景需要探測的常見的特征。尋找筆直的線段可以幫助應用進行文檔掃描和條形碼讀取。然而，傳統的線段探測方法并不適合在 GPU 上實現，特別是在移動設備的 GPU 上。
許多線段探測過程都基于霍夫變換，這是一項將真實世界笛卡爾直角坐標空間中的點轉化到另一個坐標空間中去的技術。轉化之后在另一個坐標空間中進行計算，計算的結果又轉化回正?？臻g代表線段的位置或者其他特征信息。不幸的是，許多已經提出的計算方法都不適合在 GPU 上運行，因為它們在特性上就不太可能充分地并行執行，并且都需要大量的數學計算，比如在每個像素上進行三角函數計算。
2011年，Dubská 等人 [10] [11] 提出了一種更簡單并更有效的坐標空間轉換方法和分析方法，這種方法更合適在 GPU 上運行。他們的方法依賴與一個叫做平行坐標空間的概念，聽上去很抽象但是我會展示出它其實很容易理解。
我們首先選擇一條線段和線段上的三個點：
要將這條線段轉化到平行坐標空間去，我們需要畫出三個平行的垂直軸。在中間的軸上，我們選取三個點在 X 軸上的值，也就是 1，2，3 處畫一個點。在左邊的軸上，我們選取三個點在 Y 軸上的值，在 3，5，7 處畫一個點。在右邊的軸上我們做同樣的事情，但是取 Y 軸的負值。
接下來我們將代表 Y 軸值的點和它對應的 X 軸值連接起來。連接后的效果像下圖：
你會注意到在右邊的三條線會相交于一點。這個點的坐標值代表了在真實空間中線段的斜率和截距。如果我們用一個向下斜的線段，那么相交會發生在圖的左半邊。
如果我們取交點到中間軸的距離作為 u (在這個例子中是 2)，取豎直方向到 0 的距離作為 v (這里是 1/3)，將軸之間的距離作為 d (這個例子中我使用的距離是 6)，我們可以用這樣的公式計算斜率和截距
斜率 = ?1 + d/u
截距 = d × v/u
斜率是 2，截距是 1，和上面我們所畫的線段一致。
這種簡單有序的線段繪畫非常適合 GPU 進行，所以這種方法是一種利用 GPU 進行線段探測理想的方式。
探測線段的第一步是尋找可能代表一個線段的點。我們尋找的是位于邊界位置的點，并且我們希望將需要分析的點的數量控制在最少，所以之前談論的坎尼邊界探測是一個非常好的起點。
進行邊界探測之后，邊界點被用來在平行坐標空間進行畫線。每一個邊界點會畫兩條線，一條在中間軸和左邊軸之間，另一條在中間軸和右邊軸之間。我們使用一種混合添加的方式使線段的交點變得更亮。在一片區域內最亮的點代表了線段。
舉例來說，我們可以從這張測試圖片開始：
下面是我們在平行坐標空間中得到的 (我已經將負值對稱過來使圖片高度減半)
圖中的亮點就是我們探測到線段的地方。進行一個非最大值轉化來找到區域最值并將其他地方變為黑色。然后，點被轉化回線段的斜率和截距，得到下面的結果：
我必須指出在 GPUImage 中這個非最大值轉換過程是一個薄弱的環節。它可能會導致錯誤的探測出線段，或者在有噪點的地方將一條線段探測為多條線段。
正如之前所提到的，線段探測有許多有趣的應用。其中一種就是條形碼識別。有關平行坐標空間轉換有趣的一點是，在真實空間中平行的線段轉換到平行坐標空間中后是垂直對齊的一排點。不論平行線段是怎樣的都一樣。這就意味著你可以通過一排有特定順序間距的點來探測出條形碼無論條形碼是怎樣擺放的。這對于有視力障礙的手機用戶進行條形碼掃描是有巨大幫助的，畢竟他們無法看到盒子也很難將條形碼對齊。
對我而言，這種線段探測過程中的幾何學優雅是令我感到十分著迷的，我希望將它介紹給更多開發者。
小結
這些就是在過去幾年中發展出來的機器視覺方法中的幾個，它們僅僅是適合在 GPU 上工作的方法中的一部分。我個人認為在這個領域還有著令人激動的開創性工作要去做，這將會誕生可以提高許多人生活質量的應用。希望這篇文章至少為你提供了一個機器視覺領域簡要的總體介紹，并且展示了這個領域并不像許多開發者想象的那樣無法進入。
參考文獻
話題 #21 下的更多文章
原文?GPU-Accelerated Machine Vision
總結
以上是生活随笔為你收集整理的GPU 加速下的图像视觉的全部內容，希望文章能夠幫你解決所遇到的問題。
- 上一篇： Caffe使用step by step：
- 下一篇： 基于 OpenCV 的人脸识别