卷積神經網路中卷積的OpenCL實現
卷積神經網路中卷積的OpenCL實現
==============================================================
目錄
==============================================================
本文主要介紹深度學習中卷積，以及卷積核的OpenCL實現。
另外：新版對齊方式都沒了，莫名其妙。
1、卷積
在CNN（Convolutional Neural Network, 卷積神經網路）中，卷積是個必不可少的部件，該網路也因卷積得名。CNN的卷積層中，有許多濾波器，這些濾波器也稱為卷積核，用於卷積運算，提取資料特徵，這些濾波器的係數由訓練得到。
在數字訊號處理中，訊號x(n)經過系統，相當於訊號與系統的衝激響應函式h(n)進行卷積，得到輸出y(n)，其定義如下式。
上式中的符號“*”是卷積的意思。
可以看出，訊號x(n)經過系統，與衝擊響應函式序列經過縱軸翻折過的序列對應相乘。
而在卷積神經網路中卷積操作的原理即卷積核與對應的點對應相乘，核心即是乘累加，
卷積的方式比較多，假定濾波器的滑動步長stride（用s表示），以及濾波器的填充方式（填充長度用p表示），原始特徵圖W*H（寬W，高H），濾波器的尺寸用K*K表示，則輸出特徵圖的尺寸為
上式中括號為下取整。
這個關係暫時沒想通也沒關係，可以搜尋卷積神經網路卷積瞭解一下，卷積的核心即乘累加，如下為卷積圖解。
圖1. 卷積示例
上圖中左側為特徵圖與卷積核卷積，右側為輸出，其中卷積的步長stride為1，濾波器滑動方向先進行行處理，再進行列處理。
如右側的最左上角的數 3= 2*1 + 3*(-1) + 4 + 6*0；
隨後濾波器向右滑動1個單位，此時濾波器對應的資料為[3,7,6,3]，所以3*1 + 7*(-1) + 6*1 +3*0 = 2，其他同理。
行處理完成，濾波器向下移動一個單位，接著處理下一行，此時與之對應的數為[4,6,12,0]，所以4*1 + 6*(-1) + 12*1 + 0*0 = 10。其他同理。
圖2. 特徵圖填充，卷積示例
圖2中的卷積與圖1不同在於，是否對原始特徵圖進行填充，填充可以保留特徵圖的邊緣資訊，不至於在卷積過程中丟失邊緣資訊。以上是其中一種填充方式，在特徵圖的周圍環繞填充0，本例中忽略頂部及左側的0值，保留右側以及底層的邊緣資訊，使輸出特徵圖大小與原始的特徵圖大小相等。
卷積以及更多填充方式的詳細介紹可以參看如下連結：
Convolution arithmetic tutorial：
A guide to convolution arithmetic for deep learning：
2、卷積與影象處理
瞭解了卷積之後，如何在OpenCL中實現它呢？首先，一個單程序程式，使用一個卷積核對一張特徵圖進行處理，與多個並行程序，使用共享的卷積核對一張特徵圖進行處理相對比，顯然後者計算速度快些。這就是加速，正是時代潮流上的加速。
假設有一個已經訓練好的CNN網路用於識別著名的Lena，原始圖片如下，寬高分別為W和H：
圖 3. 名模Lena
卷積核如下，寬高均為K:
圖 4. K*K卷積核
現對圖片進行周圍零填充，保持經過卷積後的輸出影象與原始影象大小一致，假設保留右側和底側邊緣資訊，只對這兩側進行填充。（為什麼保留右側及底側，其實如果你保留頂部和左側邊緣，又想輸出影象大小一致，可能右側和底側的邊緣資訊也保留不了，在這我單純只是為了好運算）。
圖 5. 填充後的圖大小(W+K-1)*(W+K-1)
如果使用多個程序處理卷積操作，每個程序輸入對應的濾波器係數和影象資料，對於每個程序，濾波器係數一致，可以共享；資料來自影象的不同位置，資料塊大小為K*K。在深度學習中，訓練完成的網路抽象了資料特徵，得到不同特徵的卷積核（濾波器），所以為了直觀，本文將使用Sobel運算元作為卷積核，K=3，可以用於影象邊緣檢測，生成邊緣特徵圖。
Sobel運算元的原理參考：
圖 6. 水平梯度和垂直梯度卷積因子
3、卷積的OpenCL實現
Sobel運算元與影象卷積的Matlab程式，我已放入Github中，如下：
該程式碼與OpenCL中需要實現的卷積邏輯一致：
（1）影象資料與卷積核在儲存器中是連續的，一維陣列；
（2）輸出特徵圖也是連續的一維陣列。
圖 7. Matlab輸出水平梯度特徵圖
圖 8. Matlab輸出垂直梯度特徵圖
圖 9. Matlab輸出影象卷積圖
在OpenCL中將使用W*H個工作項來完成卷積操作，工作項ID為（x, y），x和y取值分別為[0, W-1]和[0, H-1]。
濾波器的座標值為（kj, ki），kj和ki取值範圍均為[0, 1, 2]。
座標為（kj, ki）卷積核對應的資料座標示意圖如下：
圖 10. 濾波器座標與影象資料座標
所以，濾波器座標(kj, ki)對應的儲存索引值為[kj + K*ki]；
影象資料對應的儲存索引值為[( y + ki ) * Wn + x + kj]，其中Wn=W+K-1；因為Matlab的首個數據索引為1，但在C語言中首個索引值為0，所以該處與上文Matlab程式碼中索引值略微不同。
以下為OpenCL的核函式Conv2D.cl:
/**********************************************
function: 2D Convolution of CNN
date : 2018/07/21
**********************************************/
__kernel void Conv2D( __global int * image_in, //image input
__global int * filter_in, //filter input
int K, //filter kernel size
__global int * image_out) //feature map output
{
int W; //work group global size
int Wn; //padded image width
int x; //global id x
int y; //global id y
int ki, kj; //filter coordinate,(kj, ki)
int sum = 0; //multiply and sum of filter and data
W = get_global_size(0);
x = get_global_id(0);
y = get_global_id(1);
Wn = W + (K - 1);
for(ki=0; ki<K; ki++)
for(kj=0; kj<K; kj++)
{
sum = sum + filter_in[ki*K + kj] * image_in[Wn*(y+ki) + x + kj];
}
image_out[y*W + x] = sum;
}OpenCL主函式main.cpp:
圖 11. 隨機數Sobel水平梯度卷積的結果顯示
將輸出影象匯入到Matlab中顯示即可作為驗證，本程式碼中沒有讀取和生成圖片功能，鑑於圖片太大，存入記憶體太大，且只進行了水平梯度特徵圖對隨機數的濾波，但已足以瞭解卷積的並行操作。如需處理檔案，可將圖片資料匯入，輸出檔案已存於程式當前執行目錄的image_out.txt中。
4、總結
本文介紹了卷積，卷積與影象處理，使用Sobel運算元展示lena的水平梯度與垂直梯度的特徵圖，並通過OpenCL實現了卷積核，展示了使用Sobel水平梯度卷積核對隨機數的卷積操作結果。
如有不當，請您斧正，謝謝。
途次客
2018年7月21日