cuda 簡介

第二章影像對位和 cuda 簡介

2.2 cuda 簡介

GPGPU（General-Purpose computation on GPU）是以 GPU（Graphic Processing Unit）來進行通用的計算，而不是只局限在把顯示卡用在於圖形的顯示上，隨著發的 CUDA(Compute Unified Device Architecture) 是一種 GPGPU 的技術；

cuda 平台是透過 C 語言的函式庫和一些 CUDA 的函式延伸出來的，因為不用使用圖形函式庫，而不會被傳統的 render pipeline 綁住，使得在程式設計上更方便。

圖 2-1 GPU 效能和 CPU 效能對比

CUDA 的 GPGPU 是採取 SIMD 架構[2]，也就是不同的資料，要用同樣的程式來做計算。而在這種情況下，就可以透過把每一組資料的計算，都用一個 thread

去做，在 device 中要執行的 thread 又可以分成好幾個 block，其中 thread block 的型式最多可以到三維。以此平行化來加速計算。其中，演算法的平行度也要高，如此才能避免平行化的負擔大於平行化的加速。

在 CUDA 的程式架構裡，程式執行的區域會分成 Host 和 Device 其中，host 指的就是 CPU，而「device」就是 GPU 了在 CUDA 的程式架構中，主程式還是由 CPU 來執行(圖 2.2 取自官方文件) 。

圖 2.2 CUDA 執行架構

CUDA的程式是遇到了資料要平行化處理的部分，就會將要在 GPU 跑的程式編譯成 device 能執行的程式，再丟給 device 執行了，而這個程式在 CUDA 裡把他叫做kernel。 CUDA 的是標準的single instruction multiple data (SIMD)架構，

這代表了會產生許多在 device 上執行的 thread，每一個 thread 都會去執行 kernel 這個程式；雖然程式都是同一份，但是會因為其 thread index 的不同，

而取得不同的資料來計算。而在同一個 block 裡的 thread，有部分的記憶體

（shared memory）是共用的，（shared memory）是cuda加快整體效能的核心；

shared memory是內建在晶片上(built on-chip) ，而且可以非常快速的被存取，可以大幅的增加記憶體的存取速度，而由於 thread block 的最大大小是有限制的，大約是512個，所以一般來說不會把所有的 thread 都塞到同一個 block 裡。所以在同一個kernel裡面，必須用同樣維度和大小的 thread block，來組成一個 grid 做批次處理。這個 grid是指block的數目，同樣的可以是一維或二維陣列的形式。所以，實際上 CUDA 在執行一個 kernel 程式的時候，必須要指定他的 grid 和 block 的相關資訊（維度和大小）；而會產生的 thread 數目，

就會由這兩者的資訊來決定，這必須是固定的，同一個kernel不能夠變換他的 block跟threads的數目。

一個CUDA實作大致上包含了以下五個步驟:

1.在GPU上宣告出我們想要的記憶體空間，並且初始化GPU上的資料 2.把資料從CPU的記憶體上傳到GPU的記憶體空間

3.決定kernel的gridDim (要有多少個blocks)和blockDim(每個block有多少個thread) 4.執行kernel結果會存在GPU的memory內

5.把我們想要的結果從GPU的記憶體送回到CPU的記憶體

在 CUDA 中，要讓 thread 可以使用的變數，都必須要先把資料傳送到 device 的記憶體裡，而 device 的記憶體，記憶體架構可分為 DRAM 和 chip 上的記憶體兩種。在使用的時候，global memory 是存在於 DRAM，可以由整個 kernel 裡面的所有 thread 讀取跟存入，但是讀取跟存入的速度較慢。local memory 如果用的記憶體少則會存取在是在 chip 上，如果用到的記憶體多則是存在於 DRAM 上。shared memory 是內建在晶片上(built on-chip) ，shared memory 的記憶體空間只有 16KB，可以非常快速的被存取，並且大幅的增加記憶體的存取速度，

texture 則是唯獨的記憶體，也是在 DRAM 上但是讀取速度較快。(圖 2.3 取自官方文件)

圖 2.3 CUDA 記憶體架構

CUDA 的軟體架構大概分為 Library、runtime、Driver 三個部分；而我們在開發程式的時候，可以透過這三個部分，來使用 GPU 的計算能力。(圖 2.4 取自官方文件)

圖 2.4 CUDA 軟體架構

2.3 影像對位簡介

影像對位(Image Registration)是一種方法用來了解兩張影像間彼此特徵的對應關係，在許多大型的影像處理系統都使用著這樣的技術。影像對位被認為是一個基本的技術用來分析醫學影像，因為當空間上的相對應關係被清楚的定義之後，我們可以整合它所獲得的資訊，以利於很多影像的研究。

在醫學影像中，由於影像資料的成像可能由多種方式產生，這樣的分類主要來自於利用不同的儀器，產生出相對應的影像，影像模型大概可以分為二個種類:

結構型(Anatomical)與功能型(Functional)。結構型影像就例如；X-ray、電腦斷層掃描(Computed tomography)、核磁共振(Magnetic Resonance Imaging)之類的影像結構；功能型影像；正電子發射電腦斷層掃描(Positron Emission Tomography) 單光子射出電腦斷層掃描(Single Photon Emission Computed Tomography)等，包含了器官或是組織功能的內部資訊。

三維(3D)影像對位的流程是將其中一組影像資料 A(稱為 Floating Image or Transformed Image)選出一組空間轉換的參數與另一組資料 B(稱為 Reference Image)找出最相似的轉換過程[3]。這個演算法將不斷的調整空間上的參數去做影像 A、B 相似度(Similarity Measure)的測量，循序的找出最大的相似度，這篇論文中我們是以 Mutual Information 來判斷相似度。

圖 2-5 egistration 流程圖

在最近十年，越來越多研究者對於資訊理論的測量感到興趣，我們可以把一張影像當做一種資訊，而影像對位(Image Registration)可以被想成求取最大的兩張影像的共有資訊，最常被使用到的資訊量測量方法為起源於通訊理論的 Shannon Entropy，我們將在第三章會詳細介紹。

第三章 Mutual Information 簡介

3.1 前言

利用最求取最大的 Mutual Information (MI)的方法來完成 Image Registration 將這章節詳述。使用 MI 很大的好處就是不需先對影像做任何處理的動作，不必 (uncertainty)，當資訊量越大時，代表我們會有更大的不確定性。在 1948 年 Shannon 定義出一種的資訊量測量的方法[5]，成為現在 entropy 的定義，如果有一個系统 S 内存在多個事件 S = {E₁,...,E_n}，每個事件的機率分布 P = {p₁, ...,

個事件發生的機率愈平均，我們就比較不容易預期哪個事件會發生，不確定性就上升，entropy 也就比較大，如果其中幾個事件發生的機率很大，那不確定就會隨之下降，entropy 也就會比較小。

在影像處理上我們把影像當作一個隨機變數，把不同的灰階值當作不同的事

3.3 共同資訊(Mutual Information)

我們在 3.2 提到了的定義，利用 entropy 我們可以進而求出 Mutual 件，利用計算影像的 histogram 來求出每個灰階值出現的機率，如此就可以利用 Shannon Entropy 的定義去求出影像的 entropy。

entropy

Information，計算 Mutual Information 主要表示兩個隨機變數的相關聯程度，

Mutual Information 愈大表示兩個隨機變數的關聯程度愈大。 MI (Mutual ibution MPD)。

計算聯合熵(Joint E

已知，則其聯合熵可利用以下公式獲得 H(A,B)=

計算 Joint Entropy 需要的是統計 joint histogram，也就是統計兩張影像中相

air)累積的次數 A 圖的像素

值是 5，而 B 圖的像素值是 (5,10)這一個 pair 存進我們的 joint histogram。

計算共同資訊 (Mutual Information)

共同資訊可用來量測一個變數包含在另一個變數的資訊量。得到A影像與B 影像之臨界熵及聯合熵後，A與B影像之共同資訊即可依下式計算求得。

對應點的灰階值對(gray level p ，比如說在同一個位置 10，我們就把

藉由取共同資訊的最大值，而找到使得影像對位的空間轉換 最理想的參數 t：

為影像 A 和 B 使用參數 t 轉換的共同資訊， 為最大值時也就表示這

是我們要轉換的空間參數的最佳值，也就是。

3.4 MI特性

MI 與 Entropy 之間的關係，簡單來說 Entropy 可以用來測量

我們現在可以瞭解

影像包含的資訊量，也就是影像的灰階值分布情形，而 MI 則是測量兩個影像間彼此包含多少資訊量。MI 會有以下特性[6]；

1.Non-negativity 2.Symmetry

3.Independence

4.Self 5.Bounded

上述的特性我們都可以了解到，共同資訊量是不會為負的，A 和 B 的 MI 跟 B 和 A 的 MI 是一樣的，若是當 MI 為 0 的時候，也就是代表彼此的隨機變數是獨立無相關聯，自己跟自己的 MI 也就是自己的資訊量。

第四章 Cuda 實作流程

4.1 實

章節是這篇論文中最重要的部分，我們利用 CUDA 平行的算出 histogram 算出 entropy 跟 joint entropy，並且從 GPU 的記憶體

在這個部分將介紹 CUDA 實作的概念，我們將先求出 histogram，並且在同 sequential algorithm 計算 histogram

taldatanumber; i++) temp=data[i];

用 GPU 平行運算的架構，首先 ncurrent Read 是 CUDA 記憶體讀取的架構，所以平行的讀取資料是可以的

(tem ]

作簡介

這個

值，進而用 CUDA 平行的

傳回 CPU 的記憶體，求出最大的 MI 值。

4.2 計算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

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

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;

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 就比較能在平行

在文檔中利用圖形處理晶片研究共同資訊量的影像對位 (頁 9-0)

第二章 影像對位和 cuda 簡介