以FPGA實現腹腔鏡手術訓練系統之3D器械座標偵測

以 FPGA 實現腹腔鏡手術訓練系統之 3D 器械座標偵測

FPGA-based Instrument 3D Position Tracking for Laparoscopic Surgery Training System

林淵翔 許銘峰 許博凱

Lin, Yuan-Hsiang Shiu, Ming-Feng Syu, Bo-Kai

國立台灣科技大學 電子工程系

Department of Electronic and Computer Engineering National Taiwan University of Science and Technology

摘要

具有虛擬實境功能的腹腔鏡手術訓練 系統是可以幫助初學醫生練習和評估基礎 手術技能的有效工具，系統為了要準確且 即時的控制虛擬物件與回饋，需要準確且 即時的追蹤實際操作中器械的 3D 座標位 置，本專題以 FPGA 實現即時的 3D 器械 座標偵測的演算法，並加以驗證準確度，

可以有效降低電腦端的運算量。

本專題以立體視覺和影像處理方法來 偵測實際器械的 3D 座標，在專題製作的 過程中，使用 Matlab 開發相關演算法，再 加上 Matlab 團隊在 2012 年推出了 HDL coder 和 FPGA in the loop 技術，可以大幅 減少演算法的開發，因此本專題藉由引進 這個新技術到本專題中，一方面學習新的 技巧，同時也加速 3D 定位演算法的開發。

圖 1為系統架構圖，包括一手術訓練 箱、FPGA平台(DE2-115)和PC端三個部分，

手術訓練箱包括兩組器械和兩組攝影機 (D5M)。FPGA平台上開發了影像前置處理 單元和空間座標轉換單元，用來偵測箱內 手術器械頂端的3D座標，再將座標資料透 過RS-232傳送至電腦，供應用軟體整合使 用。PC端包含虛擬實境設計，將三維座標 應用於不同的訓練場景。

經過實驗驗證，本系統3D座標偵測的 平均誤差在x, y, z軸分別為0.028 ± 0.084 cm，0.125 ± 0.082 cm和-0.119 ±0.335 cm。

在FPGA上的計算速度可以達到 48 Frame Per Second (FPS)，較在Matlab環境的 0.1 FPS 和 FPGA in the loop 環境的5 FPS 高，可以達到準確且即時的運算效果。

圖 1. 系統架構圖

一、簡介 1.1 微創手術

微創手術(minimally invasive surgery) 是全球外科手術發展的主流，已經廣泛使 用在多種手術治療上。其中，腹腔鏡手術 是最早發展的一種微創手術。對病人有許 多好處，包括傷口小、可以降低疼痛、減 少感染和縮短復原時間等。但和傳統手術 相比，執行微創手術的醫生需要更複雜的 手術技巧，也比較困難。因為醫生必須透 過電視螢幕觀看手術位置並操作長長的手 術器械，會有手眼不協調的問題。因此，

需要一個訓練工具以降低醫生的學習曲線，

增加病患的安全。

訓練的方法包括實體的訓練箱方法 [1]-[3]，和最新的虛擬實境方法[4]-[8]。在 本系統中，我們希望在虛擬實境的方法中 使用真正的器械去操控電腦上的虛擬物件，

因此需要偵測訓練箱內實際的器械位置，

以回饋到電腦軟體中去操控虛擬物件。

1.2 FPGA in the Loop技術簡介

FPGA in the Loop是Matlab [9]-[11]針 對FPGA設計的一種功能，設計者只需將 系統所需的各個方塊完成接好，Matlab就 能 將 設 計 的 方 塊 轉 成 HDL code ， 配 合 Simulink就可在內部進行驗證，而FPGA in the Loop提供了實驗板回授驗證的功能，

Matlab會將資料透過網路線送入FPGA的 實驗板，經過計算後，再將輸出回傳到 Matlab，以供驗證，還可以依照需求自行 生成testbench碼，據公司估計可節省約 48%的時間。

圖 2. FPGA in the Loop [12]

二、研究方法 2.1 2D座標轉3D座標

本系統使用立體視覺[13]和影像處理 的方法，在器械上安裝兩個Marker，採用 兩個攝影機(640480 pixels)去偵測這兩個 Marker的3D座標，進而計算出器械頂端的 3D座標。首先說明空間中的點如何投影到 影像平面上。圖3，假設空間中任一點P，

其相對於CCD中心的座標值為(X, Y, Z)，經 過投射後成像於影像平面中的P’點，其相 對於影像中心點的座標值為(x, y)，而影像 中心點相對於CCD的座標值為(0, 0, f )，其 中f為CCD中心點到影像面的距離。三維空 間中的P點成像於二維影像平面的P’點，

其關係式為(1)：

(1) 其中(x, y)的單位是像素(pixel)。

圖 3. 空間中的點成像於 CCD 的影像面

圖 3為單顆CCD攝影機擷取影像的 架構，並無法計算P點的深度。所以必須 透過兩台CCD攝影機來抓取深度資訊。以 兩部CCD攝影機同時拍攝空間中任一目 標點為例，圖 4以簡單的視覺立體幾何說 明目標點與立體影像之關係。

(x,y,z)

(X

,Y

) (X

,Y

)

Object

BaseLine b Origin

Image Plane 1 Image Plane 2

f f

圖 4. 目標物成像於兩部 CCD 的影像面

關係式如下：

這邊的 和 為Camera CMOS感光元件 單位Pixel上X方向與Y方向的長度。

_{( )}

^{( )}

( )

( )

不過要計算Z座標，必須要先確定鏡頭焦 距f和CMOS感光元件單位Pixel方向上的 長 度 、 。 我 們 可 以 透 過 Matlab Calibration Toolbox [14]的棋盤格校正來 修正因鏡頭導致影像扭曲的現象，並得到 焦距f與 和 的關係比率，因此公式可簡 化如下:

_{( )}

₍ ^{( )} ₎

_{( )}

(5)

( )為物件之3D座標，( )為 左鏡頭計算出的2D座標，( )為右鏡頭 計算出的2D座標，b為兩鏡頭之距離， 為 軸之單位像距， 為 軸之單位像距。

2.2 Matlab開發環境

圖 5為Matlab的開發環境，可以不用 直接從HDL code著手，直接從影像處理演 算法著手，處理方法包含色彩空間轉換、

Maker顏色偵測、膨脹與侵蝕、2D座標計 算、3D座標計算及校正等，當演算法開發 完成後，可用Matlab HDL Coder將已開發 完的演算法轉成Verilog HDL，更進一步使 用FPGA in the Loop功能，將影像輸入到 DE2-115，再從DE2-115輸出到Matlab完成 驗證，透過Matlab的協助，可有效縮短開 發時程。

當演算法開發完成後，可視為演算法 元件(Intellectual Property)，如圖 6所示，

直接將兩組攝影機的影像輸入至FPGA中，

得到Tip的座標，再由RS232將三維座標傳 至PC端，以增加處理的速度。

圖 5. Matlab 開發環境

圖 6. 實際硬體方塊圖

2.3 色彩空間轉換(RGB轉HSV)

由於RGB轉HSV是浮點數運算，並不 適合FPGA計算，因此我們特別找了定點 數的轉換公式，使用Simulink的block拉成。

轉換公式如下：

{

{

(6)

r, g, b代表紅色、綠色、藍色，範圍從 0到255，max = max(r, g, b)，min = min(r, g, b)，h為色相( ~ )，s為飽和度(0~100)，

v為亮度(0~100)。

2.4 Marker顏色偵測

如圖 7 所示，我們在器械上貼上兩個 顏色的 Marker，在利用影像處理方法找出 每個 Marker 的中心座標，最後可以用 (7)~(9)計算出器械頂端(Tip)的 3D 座標。

圖 7. 器械 Marker

𝑇 𝑝

_𝐿 ^𝐿

𝑇 𝑝

_𝐿 ^𝐿

_{𝑇 𝑝}

_𝐿 ^𝐿

如表 1，設定各顏色的 HSV 閥值，就 能將 Marker 顏色挑出來。並設定輸出為布 林的資料型態。

表 1. 顏色偵測閥值

Red Yellow

H 345~359 and 0~15 45~75

S >30 >30

V >30 >30

如圖 8所示，左圖為攝影機偵測黃色；

右圖為攝影機偵測紅色，把顏色偵測到的 部分設為1；未偵測到的部分設為0，最後 轉成布林輸出，以縮減資料傳輸量。

圖 8. Maker 顏色偵測

2.5 膨脹與侵蝕

膨脹侵蝕為影像處理中常用的功能，

所以就不再贅述，在我們的影像中，經常 出現雜點，造成計算重心時的誤判，為了 解決這個問題，就得需要用到膨脹侵蝕，

我們先將影像侵蝕，消除雜訊，再把所需 的標記膨脹放大，已確保資料正確。

圖 9. 膨脹與侵蝕的 block

如圖 10所示，左下圖為原圖；左上圖 為二值化之後的圖；右上圖為經過一次膨 脹的圖；右下圖為經過一次侵蝕的圖，若 未經過處理，圖片上的雜點就會影響到計 算重心時的誤判，所以我們將原圖侵蝕後 再膨脹，可以消除雜點。

2.6 Marker標籤化

在本次專題中，有兩組器械，一組器 械上有兩個標示貼紙，所以一張圖上會有 四個標示物，為了在計算座標時的方便，

需要把每個單點都區分出來，Marker 標籤 化就是在完成這件事，把所有的單點都標 上 1 到 4 的編號，以便之後的重心運算。

但最後我們使用一個簡單且直觀的方 法，將四個 Marker 用四個顏色去代表，拿 掉標籤化，省去了標籤化的時間。

圖 10.膨脹與侵蝕

2.7 重心計算

重心的計算主要是將膨脹侵蝕後的標 記轉換為座標軸上的點，再利用 2D 轉 3D 的公式，算出實際上的 3D 座標。

在這次的專題中，由於執行速度的需要，

我們並沒有使用最準確的算法，而解決方 法會在下述討論中說明，雖然不是用標準 算法，但與精確值比較，結果十分相近。

為了加速我們的系統，我們只取 X, Y 的極 大、極小值做平均，嚴格來說，計算出來 的值是中心而不是重心，但是因為我們的 標記面積不大，且形狀非常單純，所以重 心及中心的值是十分相近的，最後測出來 的結果，與 Simulink 內建的 block 算出的 重心座標相比，吻合度也非常之高，所以 我們能確定，若圖形形狀較單純時，如下 圖 11，利用這種算法並不會影響到我們最 後的成效。

圖 11. 重心與中心座標比較

一般的重心算法：

( )

^{∑ ( )}

本次使用中心算法：

( )

( ) ( )

(8)

2.8 準確度驗證方法

如圖 12所示，為了要確認器械上標籤 的3D座標偵測準確度，我們在訓練箱工作 區域內使用兩種Marker的高度做測試，分 別是y = 1公分及y = 3公分，x為每2.5公分 做一測試點，z為每1公分做一測試點，x, z 平面共有80個測試點，總共測試160個點。

Stereo Camera 6cm

z = 24.1cm

Origin

z

x -x

· · ·

· ·· · · ·

· · · · ·

· · · · ·

· · · · ·

· · · · ·

· · · · ·

· · · · ·

· · · · ·

· · · · ·

· · · · ·

· · · · ·

· · · · ·

· · · · ·

· · · · ·

· · · · ·

·

z = 39.1cm

2.5cm 5cm -5cm -2.5cm

圖 12. 測試物體放置圖

三、 結果

如表2所示，Marker 3D座標的量測平 均誤差(與實際值比較)，總共有160個測試 點，在座標校準前的平均誤差x：0.011 ± 0.124 cm，y：-0.337 ± 0.415 cm，z：-0.078

±0.426 cm，我們使用線性回歸方法去校準 測量資料，校準後的平均誤差為x：0.028 ± 0.084 cm，y：0.125 ± 0.082 cm，z：-0.119

±0.335 cm。

表 2. 量測結果

x y z

Before calibration

0.011 ± 0.124 cm

-0.337 ± 0.415 cm

-0.078

±0.426 cm After

calibration

0.028 ± 0.084 cm

0.125 ± 0.082 cm

-0.119 ± 0.335 cm

另外，在執行速度上，如下表 3 所示，

在電腦端使用 Matlab 執行時，計算速度只 有 0.1 FPS，而在 FPGA in the loop 執行時，

已提高 50 倍到達 5 FPS，最後在 FPGA DE2-115 平台執行時已可以達到 48 FPS。

表 3. 不同平台的執行速度

Matlab FPGA in the loop FPGA Speed

(FPS)

0.1 5 48

四、討論

4.1 Simulink在本專題中的優與劣

Simulink確實內建非常多便利的block，

在我們模擬演算法時幫助很大，但是大多 數並不支援轉HDL code，所以幾乎所有的 block在實際轉code時，都需要重新建製出 相同功能的block，但與直接打code比起來，

所需時程還是相對短的。

4.2 顏色偵測的比較(RGB與HSV)

因 為 我 們 發 現 就 同 一 種 顏 色 來 說 ， RGB的變動範圍非常的不穩定，即使在同 一個光源下，器械的遠近變化也會造成誤 判，若把偵測範圍調寬又會把雜物也算進 去，而HSV就無此類的困擾，以H(色調) 來判斷顏色，S(飽和度)及V(明度)代表明 暗強弱的變化，便能較容易達到顏色偵測 的目的。

4.3 FPGA與Simulink中的膨脹與侵蝕 – 串列處理 v.s. 陣列處理

在Simulink裡軟體是採陣列處理，而 FPGA是採串列處理，單就3*3的遮罩來說，

軟體因為是矩陣處理，所以執行起來非常 容易，但就串列處理的FPGA來說，光是 要把一張圖給儲存下來就得消耗不少的記 憶體，為此，我們找到了一個解決方法。

如圖 13所示，假設原圖為320*240的 圖片，掃描方向由左至右、由上至下，輸 入依序為1、2、3……319、320、321、

322……639、640……，其中line buffer的 長 度 取 決 於 圖 片 row 的 長 度 ， 也 就 是 320pixels，而Matrix3x1合併成為3x1的矩 陣，再buffer兩次，再併起來成為9x1的矩 陣。而最後的(OR)閘和及(AND)閘代表膨 脹的計算和侵蝕的計算。

五、 結論

使用FPGA in the Loop開發環境能有 效地減少開發時程，但需對硬體處理的資 料流有概念，不能把它當成軟體，而Matlab 軟體的操作及Simulink block熟悉，都是相 當重要。目前本系統已經可以即時計算出 器械的3D座標，且經過驗證，量測的平均

誤差x、y、及z的誤差都小於0.2公分。並 且可以使用RS232把座標值傳到電腦，未 來結合嵌入式系統(Linux,Android…)，即 可增加可攜性、方便性、以及降低成本。

圖 13. 膨脹與侵蝕的 block

六、參考文獻

[1] M. Alfa-Wali and A. Antoniou,

“Eco-Friendly Laparoscopic Home Trainer,” Journal of the Society for Simulation in Healthcare, vol. 6, no. 3, pp. 176-179, Jun., 2011.

[2] A. M. Martinez and D. L. Espinoza,

“Novel Laparoscopic Home Trainer,”

Surgical Laparoscopy Endoscopy &

Percutaneous Techniques, vol. 17, no. 4, pp. 300-302, Aug., 2007.

[3] J. Hwang et al., “A novel laparoscopic ventral herniorrhaphy training system,”

Surgical Laparoscopy Endoscopy &

Percutaneous Techniques, vol. 20, no. 1, pp. e16-e18, Feb., 2010.

[4] R. Aggarwal, P. Crochet, A. Dias, A.

Misra, P. Ziprin and A. Darzi,

“Development of a virtual reality

training curriculum for laparoscopic cholecystectomy,” British Journal of Surgery, vol. 96, no. 9, pp. 1086-1093, Sep., 2009.

[5] C. R. Larsen et al., “Effect of virtual reality training on laparoscopic surgery:

randomized control trial,” BMJ, vol.

338:b1802, pp. 1-6, May, 2009.

[6] P. Kanumuri et al., “Virtual reality and computer-enhanced training devices equally improve laparoscopic surgical skills in novices,” JSLS, vol. 12, no. 3, pp. 219-226, Jul-Sep, 2008.

[7] K. S. Gurusamy, R. Aggarwal, L.

Palanivelu, and B. R. Davidson, “Virtual reality training for surgical trainees in laparoscopic surgery,” Cochrane Database Systemic Review, vol. 1, CD006575, Jan. 21, 2009.

[8] E. M. McDougall et al., “Preliminary study of virtual reality and model simulation for learning laparoscopic suturing skills,” The Journal of Urology, vol. 182, no. 3, pp. 1018-1025, Jul., 2009.

[9] MathWorks, Simulink - Simulation and Model-Based Design [Online]. Available:

http://www.mathworks.com/products/si mulink/

[10] MathWorks, VHDL code - HDL Coder - MATLAB & Simulink [Online].

Available:

http://www.mathworks.com/products/hdl -coder/index.html

[11] MathWorks, VHDL Test Bench - HDL Verifier - MATLAB & Simulink [Online].

Available:

http://www.mathworks.com/products/hdl -verifier/index.html

[12] MathWorks, Introduction to FPGA Design Using MATLAB and Simulink - MathWorks Webinar [Online]. Available:

http://www.mathworks.com/company/ev ents/webinars/wbnr60358.html?id=6035 8&p1=942727305&p2=942727310 [13] A. D. Marshall, Vision Systems

[Online]. Available:

http://www.cs.cf.ac.uk/Dave/Vision_lect ure/

[14] J. Y. Bouguet (2010, July, 9), Camera Calibration Toolbox for Matlab [Online].

Available:

http://www.vision.caltech.edu/bouguetj/

calib_doc/