針對脊柱生物力學測量的重要意義,在LabVIEW開發環境中構建基于圖像處理技術快速獲取脊柱運動范圍(ROM)的測量系統。運用LabVIEW中的機器視覺模塊對采集到的原始脊柱圖像進行預處理后提取標記針角度信息,最終獲得脊柱的ROM。選取6具T7-T10節段人體脊柱標本,實施左/右側彎、前屈、后伸、順/逆時針扭轉6個方向上的加載,利用該系統測量在載荷為1~5 N·m時T8-T9節段的ROM。實驗結果表明,該系統能夠精確、快速地測定脊柱的ROM,為脊柱生物力學研究者提供了操作簡便可靠性強的測量工具。
引用本文： 李曉芳, 鄧林紅, 盧虎, 何斌. 基于LabVIEW圖像處理技術的脊柱運動范圍測量系統研究. 生物醫學工程學雜志, 2014, 31(4): 767-770. doi: 10.7507/1001-5515.20140143 復制
引言
對脊柱進行離體的生物力學穩定性測試可以幫助人們評價新的植入器材和手術方法[1]。脊柱的力學穩定性通常用脊柱節段的運動范圍(range of motion，ROM)來分析表示。脊柱的靈活性導致它的運動方式復雜多樣，因而怎樣精確地測量ROM仍是一個難題。目前測量ROM的方法大體分為兩種：在體和離體測量[2]。在體測量最常用的手段有X線[3]、計算機斷層成像(CT)[4-5]和磁共振成像(MRI)[6]等。離體測量脊柱標本的ROM時，最直接的方法為傳感器法，其缺點為傳感器本身會影響測量結果。間接測量大部分基于光學原理[7-8]，通過采集脊柱標本原始圖像，測量脊柱上的特定標記點在受載前后的位移來間接反映ROM值，提取和計算過程復雜。
本文基于光學原理，提出了一種簡單有效地從圖像中獲取ROM的方法。用黑色標記針插入待測脊柱節段，MTS858萬能材料試驗機結合特制夾具實施不同方向加載，COMS攝像機拍攝標本在不同載荷下的脊柱標本圖像。將該圖像輸入基于LabVIEW的圖像處理測量系統，經過ROI提取、灰度變換、平滑濾波、二值化、形態學處理等步驟獲取標記針的角度信息，最終計算出脊柱節段的ROM。
1 測試系統工作原理
1.1 提取感興趣區域
如圖 1所示，攝像頭采集的原始圖像包含夾具、脊柱節段等干擾信息，圖像背景雜亂不能直接用于后續處理。因此在采集原始圖像后，首先要提取只包含標記針信息的感興趣區域(region of interest，ROI)。LabVIEW中的IMAQ Extract.vi包含一個“可選矩形”輸入項，該函數有“X Step Size”、“Y Step Size”兩個選項，分別定義水平方向的采樣列數和垂直方向的采樣行數，通過記錄鼠標的起始位置和終止位置，即可從原始圖像中截取特定區域的矩形ROI。
圖1
提取感興趣區域
Figure1.
Acquire the region of interest (ROI)
1.2 圖像增強
ROI中的標記針下側存在陰影的干擾，將影響計算結果的準確度。因此，首先對圖像進行灰度變換。IMAQ MathLookup.vi包含對數、平方、指數、平方根等7種算子，能夠實現不同類型的灰度變換。經比較后發現，算子“Power 1/X”在本系統的處理效果最理想(取X=1.5)。
經過圖像灰度變換后，圖像中仍包含一些噪聲和冗余的細節，使用空間域平滑濾波器可有效地解決該問題，如公式(1)所示，x[m,n]代表原始圖像，而y[m,n]代表濾波后圖像。分別調用IMAQ Convolute.vi、IMAQ NthOrder.vi和IMAQ GetKernel.vi三個函數共同實現濾波功能，卷積模板H如公式(2)所示：
| $x\left[ m,n \right]=\Sigma i,j\in H\text{ }\Sigma x\left[ m-i,n-i \right]*H\left[ i,j \right]$ |
| $H=\frac{1}{9}\left| \begin{align} & 111 \\ & 111 \\ & 111 \\ \end{align} \right|$ |
1.3 標記針角度判定
將濾波后的圖像二值化，設定IMAQ Threshold.vi的輸入參數下限值為155，上限值為255。二值化后標記針與背景分離，如圖 2(上)所示。標記針中存在一些小孔或者裂痕，將會影響角度判斷的精確度，因此對二值化后圖像進行形態學處理。
| $A\cdot B=\left( A\oplus B \right)\Theta B$ |
公式(3)表示的是形態學閉運算，⊕和Θ分別代表膨脹和腐蝕,對A進行閉運算是先用結構元素B對A進行膨脹運算，再用B對膨脹的結果進行腐蝕運算[9]。閉運算能夠填充目標圖像中的細小空洞，連接斷開的鄰近目標，在不顯著改變目標面積的情況下平滑邊界。IMAQ GrayMorphology.vi包含膨脹、腐蝕、開運算、閉運算等6種運算，選定“閉運算”算子，設定結構元素大小為3×3。此外，同樣基于形態學原理的IMAQ FillHole.vi和IMAQ RejectBorder.vi也用于進一步改善圖像質量。
IMAQ Particle Analysis.vi是LabVIEW機器視覺模塊中一個功能強大的函數，能夠直接獲取二進制圖像中各孤立顆粒的形狀、面積、中心、角度等信息。運用IMAQ Particle Analysis.vi計算出各個標記針的角度，如圖 2(下)所示。最后，將上下兩個標記針的角度相減，保存得到的差值，將該差值減去中立位對應的差值，便可得到該載荷下標本的運動范圍。
圖2
二值化和形態學處理后的標記針圖像
Figure2.
The image of marking needles after binarization and morphological processing
2 精度測量
首先使用AutoCAD軟件繪制由7條直線構成且相鄰直線呈30°的半圓形量角器，并將其打印在A4紙上。將實驗用的標記針分別覆蓋于A4紙的7條直線上，并與之重合，使得相鄰標記針之間呈30°。調節攝像頭的位置，使得攝像頭與紙面保持平行，固定攝像頭。在緩慢勻速地旋轉A4紙的同時，每隔10 s拍攝一張照片。用本系統測量各個標記針的角度，結果如表 1所示，最大的絕對誤差為0.66°。
表1
系統測量結果與理論值對比
Table1.
Comparison of results between theoretical and system measurement data
3 實驗及結果
3.1 標本制備
在用X線掃描排除骨折、畸形、腫瘤等疾患后，最終選取6具成人新鮮脊柱標本用于本實驗。截取T7-T10節段，保留椎間韌帶和關節囊，剔除椎骨的椎旁肌、脂肪和結締組織。將椎體T7和T10包埋于自凝牙托粉中，通過模具凝固成圓柱形。將處理好的標本用雙層塑料袋密封，放入－20 ℃冰柜中保存，待實驗前置于常溫下自然解凍。
3.2 生物力學測試
將包埋了牙托粉的T7和T10節段分別置入圓形加載盤并用4個加壓螺釘固定，6枚標記針分別插入T8和T9椎體的前部和兩側。兩臺互成90°的攝像機在距離標本75 cm處分別拍攝標本矢狀面、冠狀面圖像。利用MTS858材料試驗機施加非破壞性載荷[10]，包括6個方向：前屈、后伸、左/右側彎和順/逆時針扭轉。在收集最終數據前，先進行兩次最大載荷為5 N·m的預加載，每次加載之間間隔60 s以消除松弛蠕變的影響[11]。正式實驗時，設置試驗機的位移速率為2 mm/min，在載荷達到1、2、3、4、5 N·m時分別拍攝標本圖像。
3.3 實驗結果
6具人體胸椎標本T7-T10功能節段分別進行左/右側彎、前屈、后伸、順/逆時針扭轉6個方向的加載，測得T8-T9節段的ROM如表 2所示。
表2
胸椎T8-T9節段的ROM(4 討論
為了驗證圖像測試系統的可行性，測量了6個人體脊柱功能單元T8-T9的ROM。由表 2可以看出，隨著載荷的增加，T8-T9節段的ROM也不斷增加，但是整體數值較小，扭轉、側彎的ROM要明顯大于前屈、后伸。從解剖學上來講，這是因為胸椎與肋骨相連，椎間盤較薄，關節突關節面呈冠狀，棘突成疊瓦狀，這些因素限制了胸椎的運動，故胸椎前后屈曲的ROM較小，而左右旋轉ROM較大[12]。2010年Deniz等[13]測量了完整狀態下的人體胸椎T8-T9節段在載荷為7.5 N·m的ROM，分別為前屈1.9°、后伸2.2°、側彎3.3°、扭轉3.9°。由此可見，Deniz等的測量結果與本實驗結果趨勢相符，只是在后伸方向略有差別，這可能是由于所使用標本的個體差異導致。
以往以圖像處理為基礎的脊柱ROM測量方法大多采用標記點法--在椎體的特定位置貼上與椎體對比度差異較大的標記物，通過計算標記點的位移間接反映脊柱節段的ROM。本文采用標記針作為標記工具、以LabVIEW為開發平臺，利用圖像處理技術獲得標記針角度信息。通過把受載狀態下的相鄰節段標記針角度差減去中立位的角度差，獲得脊柱功能節段的ROM值。實驗證明，所測人體脊柱T8-T9節段的ROM符合胸椎解剖學特點以及前人的實驗結果。因此，該系統提供了一種運算簡便、準確度高、非接觸性測量脊椎ROM的測量方法，為評價脊柱的內固定器械和手術治療方法的優劣提供了有效手段。
引言
對脊柱進行離體的生物力學穩定性測試可以幫助人們評價新的植入器材和手術方法[1]。脊柱的力學穩定性通常用脊柱節段的運動范圍(range of motion，ROM)來分析表示。脊柱的靈活性導致它的運動方式復雜多樣，因而怎樣精確地測量ROM仍是一個難題。目前測量ROM的方法大體分為兩種：在體和離體測量[2]。在體測量最常用的手段有X線[3]、計算機斷層成像(CT)[4-5]和磁共振成像(MRI)[6]等。離體測量脊柱標本的ROM時，最直接的方法為傳感器法，其缺點為傳感器本身會影響測量結果。間接測量大部分基于光學原理[7-8]，通過采集脊柱標本原始圖像，測量脊柱上的特定標記點在受載前后的位移來間接反映ROM值，提取和計算過程復雜。
本文基于光學原理，提出了一種簡單有效地從圖像中獲取ROM的方法。用黑色標記針插入待測脊柱節段，MTS858萬能材料試驗機結合特制夾具實施不同方向加載，COMS攝像機拍攝標本在不同載荷下的脊柱標本圖像。將該圖像輸入基于LabVIEW的圖像處理測量系統，經過ROI提取、灰度變換、平滑濾波、二值化、形態學處理等步驟獲取標記針的角度信息，最終計算出脊柱節段的ROM。
1 測試系統工作原理
1.1 提取感興趣區域
如圖 1所示，攝像頭采集的原始圖像包含夾具、脊柱節段等干擾信息，圖像背景雜亂不能直接用于后續處理。因此在采集原始圖像后，首先要提取只包含標記針信息的感興趣區域(region of interest，ROI)。LabVIEW中的IMAQ Extract.vi包含一個“可選矩形”輸入項，該函數有“X Step Size”、“Y Step Size”兩個選項，分別定義水平方向的采樣列數和垂直方向的采樣行數，通過記錄鼠標的起始位置和終止位置，即可從原始圖像中截取特定區域的矩形ROI。
圖1
提取感興趣區域
Figure1.
Acquire the region of interest (ROI)
1.2 圖像增強
ROI中的標記針下側存在陰影的干擾，將影響計算結果的準確度。因此，首先對圖像進行灰度變換。IMAQ MathLookup.vi包含對數、平方、指數、平方根等7種算子，能夠實現不同類型的灰度變換。經比較后發現，算子“Power 1/X”在本系統的處理效果最理想(取X=1.5)。
經過圖像灰度變換后，圖像中仍包含一些噪聲和冗余的細節，使用空間域平滑濾波器可有效地解決該問題，如公式(1)所示，x[m,n]代表原始圖像，而y[m,n]代表濾波后圖像。分別調用IMAQ Convolute.vi、IMAQ NthOrder.vi和IMAQ GetKernel.vi三個函數共同實現濾波功能，卷積模板H如公式(2)所示：
| $x\left[ m,n \right]=\Sigma i,j\in H\text{ }\Sigma x\left[ m-i,n-i \right]*H\left[ i,j \right]$ |
| $H=\frac{1}{9}\left| \begin{align} & 111 \\ & 111 \\ & 111 \\ \end{align} \right|$ |
1.3 標記針角度判定
將濾波后的圖像二值化，設定IMAQ Threshold.vi的輸入參數下限值為155，上限值為255。二值化后標記針與背景分離，如圖 2(上)所示。標記針中存在一些小孔或者裂痕，將會影響角度判斷的精確度，因此對二值化后圖像進行形態學處理。
| $A\cdot B=\left( A\oplus B \right)\Theta B$ |
公式(3)表示的是形態學閉運算，⊕和Θ分別代表膨脹和腐蝕,對A進行閉運算是先用結構元素B對A進行膨脹運算，再用B對膨脹的結果進行腐蝕運算[9]。閉運算能夠填充目標圖像中的細小空洞，連接斷開的鄰近目標，在不顯著改變目標面積的情況下平滑邊界。IMAQ GrayMorphology.vi包含膨脹、腐蝕、開運算、閉運算等6種運算，選定“閉運算”算子，設定結構元素大小為3×3。此外，同樣基于形態學原理的IMAQ FillHole.vi和IMAQ RejectBorder.vi也用于進一步改善圖像質量。
IMAQ Particle Analysis.vi是LabVIEW機器視覺模塊中一個功能強大的函數，能夠直接獲取二進制圖像中各孤立顆粒的形狀、面積、中心、角度等信息。運用IMAQ Particle Analysis.vi計算出各個標記針的角度，如圖 2(下)所示。最后，將上下兩個標記針的角度相減，保存得到的差值，將該差值減去中立位對應的差值，便可得到該載荷下標本的運動范圍。
圖2
二值化和形態學處理后的標記針圖像
Figure2.
The image of marking needles after binarization and morphological processing
2 精度測量
首先使用AutoCAD軟件繪制由7條直線構成且相鄰直線呈30°的半圓形量角器，并將其打印在A4紙上。將實驗用的標記針分別覆蓋于A4紙的7條直線上，并與之重合，使得相鄰標記針之間呈30°。調節攝像頭的位置，使得攝像頭與紙面保持平行，固定攝像頭。在緩慢勻速地旋轉A4紙的同時，每隔10 s拍攝一張照片。用本系統測量各個標記針的角度，結果如表 1所示，最大的絕對誤差為0.66°。
表1
系統測量結果與理論值對比
Table1.
Comparison of results between theoretical and system measurement data
3 實驗及結果
3.1 標本制備
在用X線掃描排除骨折、畸形、腫瘤等疾患后，最終選取6具成人新鮮脊柱標本用于本實驗。截取T7-T10節段，保留椎間韌帶和關節囊，剔除椎骨的椎旁肌、脂肪和結締組織。將椎體T7和T10包埋于自凝牙托粉中，通過模具凝固成圓柱形。將處理好的標本用雙層塑料袋密封，放入－20 ℃冰柜中保存，待實驗前置于常溫下自然解凍。
3.2 生物力學測試
將包埋了牙托粉的T7和T10節段分別置入圓形加載盤并用4個加壓螺釘固定，6枚標記針分別插入T8和T9椎體的前部和兩側。兩臺互成90°的攝像機在距離標本75 cm處分別拍攝標本矢狀面、冠狀面圖像。利用MTS858材料試驗機施加非破壞性載荷[10]，包括6個方向：前屈、后伸、左/右側彎和順/逆時針扭轉。在收集最終數據前，先進行兩次最大載荷為5 N·m的預加載，每次加載之間間隔60 s以消除松弛蠕變的影響[11]。正式實驗時，設置試驗機的位移速率為2 mm/min，在載荷達到1、2、3、4、5 N·m時分別拍攝標本圖像。
3.3 實驗結果
6具人體胸椎標本T7-T10功能節段分別進行左/右側彎、前屈、后伸、順/逆時針扭轉6個方向的加載，測得T8-T9節段的ROM如表 2所示。
表2
胸椎T8-T9節段的ROM(4 討論
為了驗證圖像測試系統的可行性，測量了6個人體脊柱功能單元T8-T9的ROM。由表 2可以看出，隨著載荷的增加，T8-T9節段的ROM也不斷增加，但是整體數值較小，扭轉、側彎的ROM要明顯大于前屈、后伸。從解剖學上來講，這是因為胸椎與肋骨相連，椎間盤較薄，關節突關節面呈冠狀，棘突成疊瓦狀，這些因素限制了胸椎的運動，故胸椎前后屈曲的ROM較小，而左右旋轉ROM較大[12]。2010年Deniz等[13]測量了完整狀態下的人體胸椎T8-T9節段在載荷為7.5 N·m的ROM，分別為前屈1.9°、后伸2.2°、側彎3.3°、扭轉3.9°。由此可見，Deniz等的測量結果與本實驗結果趨勢相符，只是在后伸方向略有差別，這可能是由于所使用標本的個體差異導致。
以往以圖像處理為基礎的脊柱ROM測量方法大多采用標記點法--在椎體的特定位置貼上與椎體對比度差異較大的標記物，通過計算標記點的位移間接反映脊柱節段的ROM。本文采用標記針作為標記工具、以LabVIEW為開發平臺，利用圖像處理技術獲得標記針角度信息。通過把受載狀態下的相鄰節段標記針角度差減去中立位的角度差，獲得脊柱功能節段的ROM值。實驗證明，所測人體脊柱T8-T9節段的ROM符合胸椎解剖學特點以及前人的實驗結果。因此，該系統提供了一種運算簡便、準確度高、非接觸性測量脊椎ROM的測量方法，為評價脊柱的內固定器械和手術治療方法的優劣提供了有效手段。
