郭雁文

(中北大學(xué) 光電儀器廠，山西 太原 030051)

0 引言

角膜地形圖參數(shù)的精確計(jì)算，對(duì)于眼科醫(yī)學(xué)應(yīng)用具有重要的意義，而角膜地形圖預(yù)處理的精度將決定參數(shù)的計(jì)算精度［1-2］.在角膜地形圖的獲取過(guò)程中，由于受光斑、睫毛等不同噪聲的特異性的影響，就會(huì)導(dǎo)致在直角坐標(biāo)系中無(wú)法獲得較高精度的預(yù)處理結(jié)果.對(duì)臨床診斷來(lái)說(shuō)，角膜地形圖的處理有三個(gè)重要的因素［3］:①角膜地形圖重建的精度;②角膜覆蓋度(表示為最大角膜直徑)理想狀況應(yīng)該是覆蓋整個(gè)角膜(即眼球中黑色虹膜的部分);③對(duì)于眼角膜實(shí)際拍攝捕捉過(guò)程中存在的各種問(wèn)題的解決能力，例如被檢測(cè)者頭部的晃動(dòng)，眼球的移動(dòng)，淚膜的質(zhì)量，粘液的存在，特別是睫毛的遮擋等造成的影響［4-5］.而要有效地解決這些問(wèn)題，首先應(yīng)先對(duì)CCD 相機(jī)捕捉到的角膜圓環(huán)圖像進(jìn)行預(yù)處理，從而找到圖像重建方案［6］.

現(xiàn)有的基于直角坐標(biāo)系的預(yù)處理算法除計(jì)算量較大之外，還會(huì)導(dǎo)致兩個(gè)極端:①導(dǎo)致有用的數(shù)據(jù)缺失嚴(yán)重;②人工擬合的數(shù)據(jù)過(guò)多，在進(jìn)行參數(shù)計(jì)算前無(wú)法保證數(shù)據(jù)的精度［7-8］.因此，根據(jù)角膜地形圖原始圖像的特性，本文提出了一種基于極坐標(biāo)的原始圖像的處理算法，在圖像坐標(biāo)轉(zhuǎn)換完成后，重點(diǎn)研究在極坐標(biāo)系下角膜地形圖的自適應(yīng)邊緣檢測(cè)及邊緣生長(zhǎng)，有效地克服了原始圖像的數(shù)據(jù)丟失和獲取的數(shù)據(jù)不完全而導(dǎo)致的誤差問(wèn)題.

1 傳統(tǒng)圖像預(yù)處理算法存在的問(wèn)題

圖1 為實(shí)驗(yàn)中所獲取的基于Placido 環(huán)的角膜地形原始圖像.由于患者眼瞼很難張開到所要求的寬度，淚膜不穩(wěn)定，角膜表面的可涵蓋粘液等，都導(dǎo)致了所獲取的圖像信息丟失嚴(yán)重［9-10］.尤其是睫毛，大量遮蓋了有用的信息［11-12］.如圖1 所示，眼瞼難以張開，睫毛又濃厚，這些都會(huì)導(dǎo)致所捕捉到的角膜地形圖質(zhì)量很差.如果利用傳統(tǒng)的基于直角坐標(biāo)系的預(yù)處理算法，那么對(duì)其直接進(jìn)行預(yù)處理所得到的圓環(huán)圖就會(huì)受到睫毛等噪聲信號(hào)的嚴(yán)重影響，難以進(jìn)行下一步的參數(shù)計(jì)算，如圖2 所示.

根據(jù)直角坐標(biāo)和極坐標(biāo)系的特點(diǎn)，直角坐標(biāo)系中的一個(gè)圓，如果以圓心作為極點(diǎn)，轉(zhuǎn)換到極坐標(biāo)系后將是一條直線［13］.從圖1中可以看到，Placido 環(huán)圖像呈現(xiàn)為準(zhǔn)圓形的同心圓，因此如果能找到這些同心圓的圓心，將其作為極點(diǎn)，將Placido 環(huán)圖像轉(zhuǎn)換為極坐標(biāo)域，那么這些亮環(huán)將成為單一方向的準(zhǔn)直線型，這樣對(duì)于去除眼瞼與眼皮的工作，以及有效檢測(cè)亮環(huán)的邊緣工作，都將會(huì)被大大地簡(jiǎn)化.基于此，本文提出了一種新的基于極坐標(biāo)的角膜圖像預(yù)處理算法.

圖1 Placido 原始圖像Fig.1 Original image of Placido

圖2 傳統(tǒng)算法預(yù)處理結(jié)果Fig.2 Traditional algorithm preprocessing results

2 基于極坐標(biāo)系的角膜地形圖像自適應(yīng)預(yù)處理算法

基于極坐標(biāo)系的角膜圖像預(yù)處理，首先要完成坐標(biāo)系的轉(zhuǎn)換，同時(shí)對(duì)轉(zhuǎn)換坐標(biāo)系之后的圖像進(jìn)行重新生成，然后在新坐標(biāo)系環(huán)境下對(duì)圖像進(jìn)行進(jìn)一步處理.改進(jìn)算法的具體流程如圖3 所示.

圖3 新算法工作步驟示意圖Fig.3 Steps for implementation of the new algorithm

2.1 極坐標(biāo)圖像生成

在實(shí)際拍攝中的Placido 圖像，圖像的中心往往不是圓環(huán)的中心，因此圓環(huán)中心的定位非常重要.利用原始圖像的準(zhǔn)圓環(huán)特性，設(shè)原角膜地形圖像為f(x，y)，nx，ny分別是該圖像在X 方向與Y 方向的像素點(diǎn)數(shù).定義一個(gè)方形子圖像fES(m，n)，對(duì)子圖像fES(m，n)實(shí)現(xiàn)Hough 變換檢測(cè)圓心，利用原始圖像近似同心圓的特征，當(dāng)檢測(cè)到兩個(gè)不同半徑的圓，并且其圓心相同時(shí)，就可以確定圓心，子圖像圓心的定位結(jié)果如圖4 所示.利用該圓心作為新極坐標(biāo)系的極坐標(biāo)原點(diǎn)，將原始圖像重新調(diào)整和生成為以極坐標(biāo)原點(diǎn)為中心的方形新圖像，目的是為了在下一步進(jìn)行極坐標(biāo)轉(zhuǎn)換后，圓環(huán)能被有效地展開為直線的形式.矯正后的圖像如圖5 所示.

圖4 子圖像圓心定位結(jié)果Fig.4 The results of sub image central location

圖5 定位后重新生成的圖像轉(zhuǎn)換至極坐標(biāo)系后如圖6 所示.圖6 中下方如山丘一樣的黑色弧形為重定位時(shí)所擴(kuò)展的部分.

圖5 矯正后的圖像Fig.5 The corrected image

圖6 轉(zhuǎn)換極坐標(biāo)后的圖像Fig.6 Image of polar coordinate conversion

2.2 基于極坐標(biāo)系的自適應(yīng)水平邊緣檢測(cè)

從圖6 中可以看到，有用的圓環(huán)呈現(xiàn)近似橫線的形式，而睫毛等噪聲信號(hào)以縱向線條或不同的角度存在，因此復(fù)雜的噪聲干擾去除就變?yōu)樘崛∷椒较虻奶荻人阕樱瑥亩蟠蠼档土嗽肼暤挠绊?因此，先對(duì)極坐標(biāo)轉(zhuǎn)換后的角膜地形圖像進(jìn)行水平邊緣檢測(cè).自適應(yīng)水平邊緣檢測(cè)步驟如圖7 所示.

圖7 自適應(yīng)水平邊緣檢測(cè)步驟Fig.7 Steps of self-adaptive horizontal edge detection

首先對(duì)圖像進(jìn)行平滑處理.選擇自適應(yīng)維納濾波器，使原始圖像和其濾波后的圖像之間的均方誤差達(dá)到最小，根據(jù)局部方差來(lái)調(diào)整濾波器效果.利用Sobel 算子對(duì)圖像進(jìn)行水平方向的檢測(cè).非極大值抑制實(shí)質(zhì)上是找到邊緣檢測(cè)強(qiáng)度數(shù)據(jù)中的最高點(diǎn)，利用上述水平邊緣檢測(cè)的結(jié)果，對(duì)邊緣方向法線兩側(cè)的兩個(gè)點(diǎn)插入邊緣強(qiáng)度［14］，如果感興趣點(diǎn)的邊緣強(qiáng)度比這兩個(gè)值大，則返回該值，否則丟棄.當(dāng)非極大值抑制定義的邊緣點(diǎn)大于上限閾值時(shí)(或小于上限閾值時(shí))，標(biāo)記為邊緣點(diǎn)，然后利用8 連通域查找它的相鄰點(diǎn)以確定它們是否也大于上限閾值(或小于上限閾值)，如果形成相鄰點(diǎn)連起來(lái)的支徑，則為檢測(cè)到的邊緣圖像，處理結(jié)果如圖8 所示.

圖8 自適應(yīng)水平邊緣檢測(cè)結(jié)果Fig.8 Results on self-adaptive horizontal edge detection

2.3 形態(tài)學(xué)后續(xù)處理

如圖8 所示，在角膜地形邊緣圖像中，因?yàn)榻廾?、光斑等影響而丟失了部分有用的信息，導(dǎo)致圖像中線段有短缺區(qū)域.而在形態(tài)學(xué)算子中，利用結(jié)構(gòu)元素對(duì)圖像做膨脹運(yùn)算，可以填充目標(biāo)內(nèi)部狹窄的裂縫和長(zhǎng)細(xì)的窄溝，消去小的孔洞，而腐蝕運(yùn)算可以消除圖像中小的成分［15］.這非常適合用來(lái)對(duì)原始角膜地形圖像進(jìn)行有效的彌補(bǔ)工作與噪聲消除，最后對(duì)圓環(huán)周圍殘余的睫毛及眼瞼進(jìn)行連通區(qū)域處理，轉(zhuǎn)換回直角坐標(biāo)系中即可.

如果選擇結(jié)構(gòu)元素SE1(x，y)為線性算子，那么膨脹可填充圖像中的細(xì)長(zhǎng)窄溝，對(duì)極坐標(biāo)下的角膜地形圖中的準(zhǔn)直線具有彌補(bǔ)作用.而腐蝕可以消除圖像中小的成分，并將圖像縮小，從而使其補(bǔ)集擴(kuò)大.因此可選用菱形結(jié)構(gòu)元素SE2(x，y)對(duì)圖像進(jìn)行殘余噪聲去除.對(duì)極坐標(biāo)邊緣圖像fE作膨脹運(yùn)算，可表示為

最后標(biāo)記連通分量，計(jì)算連通面積，面積過(guò)小的像素點(diǎn)去除掉.這里，為了避免損壞有用的信息，對(duì)于連通區(qū)域的查找采用除水平方向的連通區(qū)域查找，結(jié)果如圖9 所示.由圖9 可以看到，因?yàn)榻廾绊懰斐傻臄?shù)據(jù)損傷已基本修復(fù)，大的數(shù)據(jù)短缺不應(yīng)全部補(bǔ)齊，否則會(huì)影響到數(shù)據(jù)精度，最后將圖像轉(zhuǎn)換到直角坐標(biāo)系中.

圖9 形態(tài)學(xué)后續(xù)處理結(jié)果Fig.9 Results on morphological follow-up processing

圖10 為處理結(jié)果的比較，其中圖10(a)為原始圖像矯正后的圖像，圖10(b)為圖10(a)紅線處的灰度值.可以看到，由于光斑的影響，圖10(b)中有多處圓環(huán)漏檢.圖10(c)為處理結(jié)果圖像，圖10(d)為相同行的邊緣檢測(cè)結(jié)果圖，可以看到圓環(huán)均能有效地檢測(cè)出來(lái)，這充分說(shuō)明了自適應(yīng)預(yù)處理算法使得數(shù)據(jù)精度有了大幅度提高.

圖10 最終處理結(jié)果比較Fig.10 Compare of final processing results

3 實(shí)驗(yàn)驗(yàn)證

圖11 為實(shí)驗(yàn)中獲取的三組被測(cè)者角膜地形圖像及處理結(jié)果顯示，其中原始圖像有光斑影響，中心點(diǎn)錯(cuò)位，睫毛影響等共同存在的現(xiàn)象.

經(jīng)過(guò)自適應(yīng)預(yù)處理后，從處理結(jié)果可以看到，圖(a)中被測(cè)者甲的圖像受睫毛的影響較為嚴(yán)重，導(dǎo)致直接獲取邊緣圖像時(shí)數(shù)據(jù)損失嚴(yán)重;圖(b)的處理結(jié)果能夠去除其大量?jī)?nèi)外睫毛的影響，實(shí)現(xiàn)了如圖(c)中三維顯示的較為完整的角膜建模圖像.圖(d)為被測(cè)者乙的角膜原始圖像，其中心點(diǎn)嚴(yán)重錯(cuò)位，且有光斑影響，直接進(jìn)行參數(shù)計(jì)算誤差很大;圖(e)能準(zhǔn)確矯正圓心位置，并準(zhǔn)確獲取其有效數(shù)據(jù).圖(g)中被測(cè)者丙的角膜原始圖像光斑影響嚴(yán)重，直接檢測(cè)時(shí)圓環(huán)圖像有嚴(yán)重的斷裂現(xiàn)象;圖(h)中有效去除了光斑的影響，彌補(bǔ)了由光斑導(dǎo)致的數(shù)據(jù)確失的現(xiàn)象，檢測(cè)結(jié)果數(shù)據(jù)損失小，且精度較高，能滿足后續(xù)測(cè)量參數(shù)計(jì)算的要求.

圖11 三組原始實(shí)驗(yàn)數(shù)據(jù)及其處理結(jié)果Fig.11 Three sets of original experimental data and their processing results

4 結(jié)論

本文根據(jù)Placido 角膜原始圖像的特點(diǎn)，采用極坐標(biāo)圖像的轉(zhuǎn)換方式，把復(fù)雜的噪聲干擾去除方法轉(zhuǎn)化為極坐標(biāo)系下的水平邊緣檢測(cè)和邊緣跟蹤問(wèn)題，大大簡(jiǎn)化了角膜地形圖原始圖像的預(yù)處理難度，同時(shí)增加了算法的處理精度，獲得了圓環(huán)邊緣檢測(cè)完整的結(jié)果圖像，滿足了后續(xù)參數(shù)計(jì)算的要求.大量實(shí)驗(yàn)數(shù)據(jù)證明:該算法具有數(shù)據(jù)損失小、精度高，算法易實(shí)現(xiàn)的特點(diǎn)，滿足眼科醫(yī)學(xué)的實(shí)際使用要求.

［1］劉祖國(guó).角膜地形圖學(xué)［M］.廣州:廣東科技出版社，2001.

［2］Alkhaldi W，Iskander D R，Zoubir M.Model-order selection in Zernike polynomial expansion of corneal surfaces using the efficient detection criterion［J］.IEEE Transactions on Biomedical Engineering，2010，57(10):2429-2438.

［3］Zhu Mingxia，Collins M J，Yeo A C.Stability of corneal topography and wavefront aberrations in young Singaporeans［J］.Clinical and Experimental Optometry，2013，96(5):486-493.

［4］Wizert A，Iskander D R，Cwiklik L.Organization of lipids in the tear film:a molecular-level view［J］.Journal of Biomedical Optics，2014，9(3):e92461.

［5］Atwood D K，Andersen H E，Matthiss B，et al.Impact of topographic correction on estimation of aboveground boreal biomass using multi-temporal，l-band backscatter［J］.IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing，2014，7(8):3262-3273.

［6］Alonso-Caneiro D，Iskander D R，Collins M J.Estimating corneal surface topography in videokeratoscopy in the presence of strong signal interference［J］.IEEE Transactions on Biomedical Engineering，2008，55(10):2381-2387.

［7］Polette A，Auvinet E，Mari J L，et al.Construction of a mean surface for the variability study of the cornea［J］.IEEE Conference Publications，2014:328-335.

［8］周洪亞，沈建新，高紹雷，等.基于精確中心定位的角膜復(fù)原算法［J］.生物醫(yī)學(xué)工程學(xué)雜志，2011，28(5):872-875.Zhou Hongya，Shen Jianxin，Gao Shaolei，et al.An algorithm of corneal reconstruction based on precise location of corneal center［J］.Journal of Biomedical Engineering，2011，28(5):872-875.(in Chinese)

［9］Snellenburg J J，Braaf B，Hermans E A，et al.A forward ray tracing for image projection prediction and surface reconstruction in the evaluation of corneal topography systems［J］.Optics Express，2010，18 (18):19324-19338.

［10］Florindo J B，Soaresa S H M，de Carvalho L A V，et al.Mumford-Shah algorithm applied to videokeratography image processing and consequences to refractive power values［J］.Computer Methods and Programs in Biomedicine，2007，87(1):61-67.

［11］Carvalho L A.Absolute accuracy of placido-based videokeratographs to measure the optical aberrations of the cornea［J］.Optometry and Vision Science，2004，81(8):616-628.

［12］Alkhaldi W，Iskander D R，Zoubir A M，et al.Enhancing the standard operating range of a placido disk videokeratoscope for corneal surface estimation［J］.IEEE Transactions on Biomedical Engineering，2009，56，(3):800-809.

［13］Thalji Z，Alsmadi M.Iris recognition using robust algorithm for eyelid，eyelash and shadow avoiding［J］.World Applied Sciences Journal，2013，25(6):858-865.

［14］Nixon M S.特征提取與圖像處理［M］.第2 版.北京:電子工業(yè)出版社，2010:109-115.

［15］Burger W，Burge M J.Principles of Digital Image Processing:Advanced Methods (Undergraduate Topics in Computer Science)［M］.New York:Springer，2013.