計算 entropy 的概念

一個 kernel 中算出 entropy，利用 CPU 的

是很容易的，我們只要循序的讀取 volume data 中每個 pixel 的像素值，然後讀 取到的每個像素值所對應的存取位置就加 1，便可得出 histogram。

用 CPU 實作的演算法如下:

for(int i = 0; i < to

result[temp]= result[temp]+1;

利用 CUDA 實作 histogram，由於是利 Co

p=data[i ;) ， 但 是 在 寫 入 result 的 時 候 就 會 產 生 Write collisions (result[temp]= result[temp]+1;)，例如兩個 threads

同時對 result[temp]進行加 1，result[temp]只會加 1 一次而不是兩次，

這樣與我們想要的結果不同，這就是產生了 Write collisions。

為了要解決 Write collisions 的問題，由於 CUDA 的 shared memory 並 沒有 ardware 支援的 atomic operaion，所以我們必須用程式邏輯的方式去解 決 s

h

hared memory 的 Write collisions。由於 CUDA 的 device 在執行的 時候，會以 block 為單位，把每個的 block 分配給各別的 multiprocessor 進行 目前 CUDA 的 warp 大小都是 32，也就是 32 個 thread 會被群組成一個 warp 來一起執行，並且享有共同的 shared memory range，所以我們以一個 warp 為基礎去解決 shared memory 的 write collisions，整個流程如下:

1.給每一個不同的 thread 都有一個獨立的編號，用 threadtag 這個變數中最高位的 5 個 bits 去儲存。

2.讀取 histogram 中已經累積的數目，然後存到 count 這個變數的低位的 27 個 bits。 thread 存回 histogram 覆蓋掉的情形下，則繼續迴 圈，直到是該 thread 的 count 最後存入為止。

演算法的要點如下:

dTag = threadIdx.x << 27;// 每一個 thread 有自己的編號

addD

o{

lt [temp] & 0x07FFFFFFU;

count = threadTag | (count + 1);

，result [temp]是要存入的位置，threadTag 則 把一個 thread 的 thread id 存入到高位的五個 bits，如果有兩個以上的 threads 要存

[ unsigned int threa

ata256(result, temp , threadTag);

unsigned int count;

d

count = resu

result [temp] = count;//存取 }while(result [temp] != count);

其中 temp 是讀取到的 intensity 是

入同一個記憶體位置(result [temp])，由於每一個 thread 的 thread id 不同，count 高位的 5 個 bits 就會不同，存入 result [temp]之後高位的 5 個 bits 就會不同。接 下來利用判斷式(result [data] != count)，如果是目前的 threads 存取的 result [temp]，這個 thread 就可以跳出迴圈，如果不是則在迴圈內繼續跑，這樣就可以 確保每一個讀取到同一個 intensity 的 thread 都會對 result temp]的位置加一次，

shared memory 的 write collisions 就不會影響到最後的結果，這個方法如果像素值 越平均，threads 就不會都集中在同一個迴圈內，write collisions 就比較能在平行 處理中解決，所以執行速度就越快，如果像素值集中在同一個數值，程式就會卡 在同一個迴圈，就會從 parallel 的演算法變成 linear 的情況，這樣執行時間會大 幅上升。由於 cuda 的 shared memory 的空間只有 16KB，存每一個像素值的數目 要 4 個 bytes，如果是 1D histogram 通常可以存在 shared memory 內。而 2D 的 joint

histogram 就要動用到 global memory 來實作，由於 shared memory 的存入跟讀取 的時間較快，我們還是會將 2D 的 joint histogram 像素值比較集中的部分由 shared memory 存入跟讀取，其他的部分再由 global memory 來存入。

如果要解決 global memory 的 write collisions 問題，在 cuda 中通常 8600GT 卡會支援 atomic 函式，也就是硬體支援的 atomic operation 可以確 保一

制的，但是 一個 block 裡面的 threads 是可以透過__syncthreads 函式進行同步的，用 __sy 可以避免 Write collisions，並且達到我們想要的存取結果。

在 cuda 的架構中，一個 kernel 的 grid 是沒有 syntronization 的機 同

ncthreads 可以確保在同一個 block 裡面的所有 threads 同樣執行到程式的 那ㄧ行，同一個 block 中所有的 threads 都要做完那ㄧ行之前的指令後，才會往 下執行。如果 block 中的其中任何一個 thread 沒有完成__syncthreads 之前的指 令，那麼整個 block 中的所有 threads 都不會執行__syncthreads 後面的指令，這 是 cuda 所提供的 block level syntronization，在 cuda 架構下同一個 block 裡面的 thread 有 syntronization 的機制，但是不同的 block 之間，是沒有 syntronization 的 機制的，這也影響到我們實作上的設計。

由於果蠅腦的影像，有很廣大的背景，

問題，我們特別針對 intensity 較集中的部份，我們會給每一個 thread 自己 的 shared memory 的位置，專門存入各別 thread 讀取到較集中的像素值的數目，

__syncthreads 後再合併到同一個 shared memory 的位置，由於 shared memory 有 16KB，我們可以在 shared memory 內解決這樣的問題，這樣可以平行的解決

出現次數最多的像素值，剩下的出現次數較少的像素值，就可以透過原本 atomic 與 source volume data 相同大 小的

kern

，gpu 在處理浮點數上的精確度比 cpu 要來的稍微低，所以跟 cpu 相比較會 有一點點的誤差，由於差距非常小，不會影響到影像對位的結果。