(57)【要約】【課題】光干渉断層撮影法により取得した断層画像に含まれるノイズを除去するのに好適な光干渉断層画像処理方法を提供する。【解決手段】本発明の光干渉断層画像処理方法は、光干渉断層撮影法により取得した断層画像を処理する方法であって、(1) 断層画像の深さ方向の強度分布の特徴量を算出する特徴量算出ステップと、(2) 測定対象の表面位置および深達限界位置を算出する位置算出ステップと、(3) 前記特徴量および前記深達限界位置に基づいて、輝度低下がない場合の推定特徴量および推定深達限界位置を求める推定ステップと、(4) 前記特徴量，前記推定特徴量，前記深達限界位置，前記推定深達限界位置および前記表面位置に基づいて、深さ方向の各位置での補正量を算出する補正量算出ステップと、(5) の断層画像の深さ方向の強度分布に前記補正量を加えることで断層画像を補正する加算ステップと、を備える。【選択図】図２
2016080659.tif 000002
【特許請求の範囲】
【請求項１】
光干渉断層撮影法により取得した断層画像を処理する方法であって、
各ラインの断層画像の深さ方向の強度分布の特徴量を算出する特徴量算出ステップと、
各ラインの測定対象の表面位置および深達限界位置を算出する位置算出ステップと、
各ラインについて、前記特徴量および前記深達限界位置に基づいて、輝度低下がない場合の推定特徴量および推定深達限界位置を求める推定ステップと、
各ラインについて、前記特徴量，前記推定特徴量，前記深達限界位置，前記推定深達限界位置および前記表面位置に基づいて、深さ方向の各位置での補正量を算出する補正量算出ステップと、
各ラインの断層画像の深さ方向の強度分布に前記補正量を加えることで断層画像を補正する加算ステップと、
を備える光干渉断層画像処理方法。
【請求項２】
前記特徴量算出ステップにおいて、深さ方向の全範囲または限定された範囲での前記断層画像の強度分布の積分値を前記特徴量とする、
請求項１に記載の光干渉断層画像処理方法。
【請求項３】
前記特徴量算出ステップにおいて、深さ方向の全範囲または限定された範囲で前記断層画像の強度分布にフィッティングを適用したときのフィッティング係数を前記特徴量とする、
請求項１に記載の光干渉断層画像処理方法。
【請求項４】
処理対象の断層画像が、分光計で測定したスペクトルから構築された断層画像である、
請求項１〜３の何れか１項に記載の光干渉断層画像処理方法。
【請求項５】
処理対象の断層画像が、管状の測定対象に挿入されたプローブの回転スキャンによって得られた断層画像である、
請求項１〜４の何れか１項に記載の光干渉断層画像処理方法。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、光干渉断層撮影法により取得した断層画像を処理する方法に関するものである。
【背景技術】
【０００２】
光干渉断層撮影法（ＯＣＴ: Optical Coherence Tomography）は、測定対象の内部の構造を高分解能で観察することができる。ＯＣＴは、バイオ・医療分野で注目されており、測定対象としての生体の診断にも使われ始めている。ＯＣＴにより取得された断層画像はノイズを含む場合があるので、断層画像からノイズを除去した後に診断をすることが望まれる。
【０００３】
断層画像からノイズを除去する方法として、非特許文献１，２に記載された方法が知られている。非特許文献１に記載されたノイズ除去方法は、ガウシアンフィルタ等の平滑フィルタにより断層画像を平滑化することで、ノイズを除去するものである。非特許文献２に記載されたノイズ除去方法は、複数の断層画像の平均をとることで、スッペクルノイズを除去するものである。
【先行技術文献】
【非特許文献】
【０００４】
【非特許文献１】J. Rogowska et al., Physics inMedicine and Biology Vol.47, pp.641-655, 2002.
【非特許文献２】M. Szkulmowskiet al., Opt Exp, Vol.20, No.2, pp.1337-1359, 2012.
【発明の概要】
【発明が解決しようとする課題】
【０００５】
本発明者による知見によれば、光干渉断層撮影法により取得した断層画像に含まれる或る種のノイズは、非特許文献１，２に記載されたノイズ除去方法を用いても、除去することができない。
【０００６】
本発明は、上記問題点を解消する為になされたものであり、光干渉断層撮影法により取得した断層画像に含まれるノイズを除去するのに好適な光干渉断層画像処理方法を提供することを目的とする。
【課題を解決するための手段】
【０００７】
本発明の光干渉断層画像処理方法は、光干渉断層撮影法により取得した断層画像を処理する方法であって、(1) 各ラインの断層画像の深さ方向の強度分布の特徴量を算出する特徴量算出ステップと、(2)各ラインの測定対象の表面位置および深達限界位置を算出する位置算出ステップと、(3) 各ラインについて、前記特徴量および前記深達限界位置に基づいて、輝度低下がない場合の推定特徴量および推定深達限界位置を求める推定ステップと、(4) 各ラインについて、前記特徴量，前記推定特徴量，前記深達限界位置，前記推定深達限界位置および前記表面位置に基づいて、深さ方向の各位置での補正量を算出する補正量算出ステップと、(5) 各ラインの断層画像の深さ方向の強度分布に前記補正量を加えることで断層画像を補正する加算ステップと、を備える。
【発明の効果】
【０００８】
本発明によれば、光干渉断層撮影法により取得した断層画像に含まれるノイズを除去することができる。
【図面の簡単な説明】
【０００９】
【図１】モーションアーチファクトの影響が現れている断層画像の一例を示す図である。
【図２】本実施形態の光干渉断層画像処理方法を説明するフローチャートである。
【図３】１次元断層画像の強度分布Ｉ(i,z)の例を示す図である。
【図４】１次元断層画像の強度分布Ｉ(i,z)のうちノイズフロアより強度が大きい範囲を直線で近似して示す図である。
【図５】補正量ΔＩ(I,z)の分布を示す図である。
【図６】アーチファクト部の１次元断層画像の強度分布Ｉ(i,z)に補正量ΔＩ(I,z)を加算して得られた強度分布Ｉ'(i,z)を示す図である。
【発明を実施するための形態】
【００１０】
本発明の光干渉断層画像処理方法は、光干渉断層撮影法により取得した断層画像を処理する方法であって、(1) 各ラインの断層画像の深さ方向の強度分布の特徴量を算出する特徴量算出ステップと、(2)各ラインの測定対象の表面位置および深達限界位置を算出する位置算出ステップと、(3) 各ラインについて、前記特徴量および前記深達限界位置に基づいて、輝度低下がない場合の推定特徴量および推定深達限界位置を求める推定ステップと、(4) 各ラインについて、前記特徴量，前記推定特徴量，前記深達限界位置，前記推定深達限界位置および前記表面位置に基づいて、深さ方向の各位置での補正量を算出する補正量算出ステップと、(5) 各ラインの断層画像の深さ方向の強度分布に前記補正量を加えることで断層画像を補正する加算ステップと、を備える。
【００１１】
前記特徴量算出ステップにおいて、深さ方向の全範囲または限定された範囲での前記断層画像の強度分布の積分値を前記特徴量とするのが好適であり、或いは、深さ方向の全範囲または限定された範囲で前記断層画像の強度分布にフィッティングを適用したときのフィッティング係数を前記特徴量とするのも好適である。
【００１２】
処理対象の断層画像は、分光計で測定したスペクトルから構築された断層画像であるのが好適であり、管状の測定対象に挿入されたプローブの回転スキャンによって得られた断層画像であるのも好適である。
【００１３】
以下、添付図面を参照して、本発明を実施するための形態を詳細に説明する。なお、図面の説明において同一の要素には同一の符号を付し、重複する説明を省略する。本発明は、これらの例示に限定されるものではなく、特許請求の範囲によって示され、特許請求の範囲と均等の意味および範囲内でのすべての変更が含まれることが意図される。
【００１４】
ＯＣＴは、バイオ・医療分野で生体の診断にも使われることから、生体（患者）の負担をできるだけ軽減することが重要である。それ故、ＯＣＴは、撮影時間を短くすること（すなわち、撮影速度を速くすること）が要求されている。
【００１５】
ＯＣＴは、ＴＤ-ＯＣＴ（Time-DomainＯＣＴ）とＦＤ-ＯＣＴ（Fourier-Domain ＯＣＴ）とに大別される。ＴＤ-ＯＣＴは、撮影時に干渉計内の参照光の光路長を動的に変えることにより断層画像を取得する方式である。ＦＤ-ＯＣＴは、干渉縞のスペクトルを測定して、そのスペクトルをフーリエ変換することにより断層画像を取得する方式である。現在では、前述の撮影速度に対する要求に伴い、干渉計内にメカニカルな動的機構が必要なＴＤ-ＯＣＴより、動的機構が不要なＦＤ-ＯＣＴが主流となっている。
【００１６】
ＯＣＴは、測定対象の或る位置に測定光を照射することで、その位置における深さ方向の１次元の断層画像を取得することができる。また、ＯＣＴは、測定対象に対する測定光の照射位置を走査することで、２次元または３次元の断層画像を取得することができる。ＯＣＴにより２次元または３次元の断層画像を高速に取得する為には、測定対象に対する測定光照射の走査を速くする必要がある。
【００１７】
しかし、干渉光を受光する受光素子の露光時間に対して、走査による測定光照射位置の変位量が大きいと、当該露光中に干渉縞の位相が変化する場合がある。この場合、受光素子により得られる干渉光の強度は、露光時間に亘って累積されたものとなるので、静止時に比べて振幅が減衰され、フーリエ変換によって得られる断層画像の輝度が下がってしまう。この現象はモーションアーチファクトと呼ばれている。
【００１８】
ＦＤ-ＯＣＴの一種であるＳＤ-ＯＣＴ（Spectral-Domain ＯＣＴ）は、分光計で測定した干渉光のスペクトルから断層画像を取得する手法であり、分光された各波長の比較的弱い光を比較的長い露光時間に亘って受光する必要があることから、モーションアーチファクトの影響が顕著である。
【００１９】
ＳＤ-ＯＣＴにおいて、露光時間内で干渉縞の位相が直線的にΔφだけ変化したとすると、理論上、干渉縞の振幅の減衰率は、sin（Δφ／２）／（Δφ／２） なる式で表され、干渉縞の周波数つまり断層画像での深さに依存しない。また、実際の断層画像ではノイズフロアが存在するので、ノイズフロアよりも上に出ている範囲においてモーションアーチファクトの影響に因る減衰が生じる。
【００２０】
図１は、モーションアーチファクトの影響が現れている断層画像の一例を示す図である。測定光を９０度側方に反射させる光学部品が先端に取り付けられた光ファイバを含むＯＣＴ用カテーテルを用い、このＯＣＴ用カテーテルをＯＣＴ装置に接続した。ＯＣＴ用カテーテルをイントラリピッド（人工脂肪）溶液中に入れ、ＯＣＴ用カテーテル内部の光ファイバを高速回転させて断層画像を取得した。同図は、１回転の間に撮影した複数ラインの一次元断層画像を円形に配置してできた二次元断層画像である。同図に示される円形像において、主に５時，８時および１１時の各方向に放射状に暗く写っているのが、モーションアーチファクトに因るものである。光ファイバの高速回転で突発的に生じた振動により干渉縞の位相が変動していることが原因であると考えられる。
【００２１】
このようなモーションアーチファクトは、観察者にとって判断の妨げとなるので、断層画像から除去する必要がある。ＯＣＴにより取得した断層画像からノイズを除去する方法として、非特許文献１に記載された空間フィルタによる除去や、非特許文献２に記載されたフレーム間の平均による除去が知られている。しかし、この種のモーションアーチファクトに因るノイズは、ランダムなパターンではなく、深さ方向で略一様に減衰しているので、これらのノイズ除去方法では除去できない可能性がある。また、フレーム間の平均によるノイズ除去については、動いているものが対象である場合には、フレーム間で画像にずれが生じているので、処理によってぶれた画像になってしまう。さらに、このような処理は、断層画像を平滑化する傾向があるので、観察者に異なる印象を与える可能性もある。
【００２２】
以下に説明する本実施形態の光干渉断層画像処理方法は、ＯＣＴにより取得した断層画像からモーションアーチファクトに因るノイズを除去するのに好適である。本実施形態の光干渉断層画像処理方法は、分光計で測定したスペクトルから構築された断層画像（すなわち、ＳＤ-ＯＣＴにより取得した断層画像）を処理する場合に好適であり、また、管状の測定対象（例えば血管）に挿入されたプローブの回転スキャンによって得られた断層画像を処理する場合に好適である。
【００２３】
図２は、本実施形態の光干渉断層画像処理方法を説明するフローチャートである。本実施形態の光干渉断層画像処理方法は、ＯＣＴにより取得した断層画像を処理する方法であって、断層画像取得ステップＳ１０，特徴量算出ステップＳ１１，位置算出ステップＳ１２，推定ステップＳ１３，補正量算出ステップＳ１４および加算ステップＳ１５を備える。
【００２４】
断層画像取得ステップＳ１０では、処理対象の断層画像を取得する。この断層画像は、２次元または３次元のものであり、深さ方向の１次元の断層画像を複数ライン含む。これら複数ラインの１次元断層画像のうち第ｉラインの１次元断層画像の強度分布をＩ(i,z)と表す。ｚは深さ方向の位置を表す変数である。
【００２５】
特徴量算出ステップＳ１１では、各ラインの断層画像の深さ方向の強度分布Ｉ(i,z)の特徴量Ｓ(i)を算出する。この特徴量Ｓ(i)は、モーションアーチファクトに因る輝度低下の程度を表す。この特徴量Ｓ(i)は、深さ方向の全範囲または限定された範囲での断層画像の強度分布Ｉ(i,z)の積分値であってもよいし、また、深さ方向の全範囲または限定された範囲で断層画像の強度分布Ｉ(i,z)にフィッティングを適用したときのフィッティング係数であってもよい。
【００２６】
位置算出ステップＳ１２では、各ラインの測定対象の表面位置ｚｓ(i)および深達限界位置ｚｅ(i)を算出する。これらの位置の算出は、断層画像の強度分布Ｉ(i,z)に基づいてエッジ検出や特定関数でのフィッティング等をすることで実現可能である。
【００２７】
推定ステップＳ１３では、各ラインについて、輝度低下がない場合の推定特徴量Ｓ'(i)および推定深達限界位置ｚｅ'(i)を、その周囲の特徴量Ｓ(i)および深達限界位置ｚｅ(i)に基づいて求める。各ラインでの特徴量Ｓ(i)および深達限界ｚｅ(i)の傾向については、モーションアーチファクトが生じているラインでは共に値が局所的に変動する。そこで、ベースライン（ノイズフロア）を求めることで、モーションアーチファクトがない場合の推定特徴量Ｓ'(i)および推定深達限界位置ｚｅ'(i)を求めることができる。ベースラインは、例えばフーリエ級数でのフィッティングや前後の範囲での最大あるいは最小あるいは条件分岐を使った補間等で求めることが可能である。
【００２８】
補正量算出ステップＳ１４では、各ラインについて、特徴量Ｓ(i)，推定特徴量Ｓ'(i)，深達限界位置ｚｅ(i)，推定深達限界位置ｚｅ'(i)および表面位置ｚｓ(i)に基づいて、深さ方向の各位置での補正量ΔＩ(I,z)を算出する。
【００２９】
加算ステップＳ１５では、各ラインの断層画像の深さ方向の強度分布Ｉ(I,z)に補正量ΔＩ(I,z)を加えることで、断層画像を補正する。これにより、モーションアーチファクトに因るノイズが除去された断層画像の強度分布Ｉ'(I,z)（＝Ｉ(I,z)＋ΔＩ(I,z)）を得ることができる。
【００３０】
次に、図１に示したイントラリピッド溶液の断層画像の例を用いて、本実施形態の光干渉断層画像処理方法について更に説明する。
【００３１】
図３は、１次元断層画像の強度分布Ｉ(i,z)の例を示す図である。同図は、図１に示したイントラリピッド溶液の断層画像のうちの、モーションアーチファクトの影響があるライン（アーチファクト部）の１次元断層画像の強度分布、および、このアーチファクト部の近傍に位置しモーションアーチファクトの影響がないライン（正常部）の１次元断層画像の強度分布を示す。
【００３２】
同図に示されるように、アーチファクト部および正常部それぞれの表面位置ｚｓ(i)は同じである。アーチファクト部および正常部の何れの１次元断層画像の強度分布においても、深さ方向の或る位置より深い範囲では同程度のベースライン（ノイズフロア）となっている。ノイズフロアより強度が大きい範囲では、アーチファクト部の１次元断層画像の強度分布は、正常部の１次元断層画像の強度分布に比べて略一様に減衰している。また、ノイズフロアより強度が大きい範囲では、アーチファクト部および正常部の何れの１次元断層画像の強度分布も、深さ方向の位置変数ｚに対して直線近似し得る。
【００３３】
図４は、１次元断層画像の強度分布Ｉ(i,z)のうちノイズフロアより強度が大きい範囲を直線で近似して示す図である。同図には、ノイズフロアのレベル、表面位置ｚｓ(i)および深達限界位置ｚｅ(i)も示されている。アーチファクト部および正常部それぞれの表面位置ｚｓ(i)は同じである。しかし、アーチファクト部の深達限界位置ｚｅ(i)は正常部の深達限界位置より浅い。周囲の１次元断層画像の強度分布の積分値Ｓ(i)に基づいて求められた推定特徴量をＳ'(i)とし、アーチファクト部の１次元断層画像の強度分布の積分値をＳ(i)としたとき、両者の差ΔＳ(i)（＝Ｓ'(i)−Ｓ(i)）は、同図中で４本の実線で囲まれた範囲（破線で囲まれた範囲）の面積である。このΔＳ(i)は、補正量ΔＩ(I,z)の積分値と等しい。
【００３４】
図５は、補正量ΔＩ(I,z)の分布を示す図である。モーションアーチファクトに因る減衰量は深さ方向の位置ｚに依存しないとしてよいから、各位置での補正量ΔＩ(I,z)は、下記の式により比較的容易に求めることができる。同図に示される補正量ΔＩ(I,z)の分布は、この式に従って求めたものである。
【００３５】
【数１】
2016080659.tif 000003
【００３６】
図６は、アーチファクト部の１次元断層画像の強度分布Ｉ(i,z)に補正量ΔＩ(I,z)を加算して得られた強度分布Ｉ'(i,z)を示す図である。同図に示されるように、このようにして補正されたアーチファクト部の１次元断層画像の強度分布Ｉ'(i,z)は、正常部の１次元断層画像の強度分布と傾向が一致したものとなっている。
【００３７】
以上のように、本実施形態の光干渉断層画像処理方法は、元の断層画像に対してオフセットを与えるだけであり、空間フィルタやフレーム平均のように分布を直接変える処理ではないので、元の画質を変えることなくモーションアーチファクトに因るノイズを効果的に除去することができる。
【図２】
2016080659.tif 000004
【図３】
2016080659.tif 000005
【図４】
2016080659.tif 000006
【図５】
2016080659.tif 000007
【図６】
2016080659.tif 000008
【図１】
2016080659.tif 000009