JP5885758B2 - 情報信号生成方法 - Google Patents
情報信号生成方法 Download PDFInfo
- Publication number
- JP5885758B2 JP5885758B2 JP2013551079A JP2013551079A JP5885758B2 JP 5885758 B2 JP5885758 B2 JP 5885758B2 JP 2013551079 A JP2013551079 A JP 2013551079A JP 2013551079 A JP2013551079 A JP 2013551079A JP 5885758 B2 JP5885758 B2 JP 5885758B2
- Authority
- JP
- Japan
- Prior art keywords
- signal
- information
- extended
- extrapolated
- peak
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Expired - Fee Related
Links
- 238000000034 method Methods 0.000 title claims description 65
- 230000007274 generation of a signal involved in cell-cell signaling Effects 0.000 title claims description 54
- 230000003287 optical effect Effects 0.000 claims description 114
- 238000001228 spectrum Methods 0.000 claims description 98
- 238000005259 measurement Methods 0.000 claims description 43
- 238000012545 processing Methods 0.000 claims description 40
- 230000001131 transforming effect Effects 0.000 claims description 29
- 238000003325 tomography Methods 0.000 claims description 11
- 238000001514 detection method Methods 0.000 claims description 8
- 239000000284 extract Substances 0.000 claims description 3
- 238000003384 imaging method Methods 0.000 claims description 2
- 239000013307 optical fiber Substances 0.000 description 25
- 230000003595 spectral effect Effects 0.000 description 14
- 238000012014 optical coherence tomography Methods 0.000 description 13
- 238000004364 calculation method Methods 0.000 description 11
- 238000004458 analytical method Methods 0.000 description 8
- 238000010586 diagram Methods 0.000 description 8
- 239000010410 layer Substances 0.000 description 6
- 230000035945 sensitivity Effects 0.000 description 6
- 238000007493 shaping process Methods 0.000 description 6
- 238000006243 chemical reaction Methods 0.000 description 5
- 230000008569 process Effects 0.000 description 5
- 230000008859 change Effects 0.000 description 4
- 239000002131 composite material Substances 0.000 description 4
- 238000013213 extrapolation Methods 0.000 description 4
- 230000009467 reduction Effects 0.000 description 4
- 230000000694 effects Effects 0.000 description 3
- 230000003252 repetitive effect Effects 0.000 description 3
- 230000002194 synthesizing effect Effects 0.000 description 3
- 238000005033 Fourier transform infrared spectroscopy Methods 0.000 description 2
- 230000015556 catabolic process Effects 0.000 description 2
- 238000006731 degradation reaction Methods 0.000 description 2
- 239000000835 fiber Substances 0.000 description 2
- 230000000644 propagated effect Effects 0.000 description 2
- 239000002356 single layer Substances 0.000 description 2
- 238000010521 absorption reaction Methods 0.000 description 1
- 238000013459 approach Methods 0.000 description 1
- 230000002238 attenuated effect Effects 0.000 description 1
- 230000005540 biological transmission Effects 0.000 description 1
- 230000015572 biosynthetic process Effects 0.000 description 1
- 238000007796 conventional method Methods 0.000 description 1
- 238000012937 correction Methods 0.000 description 1
- 230000006866 deterioration Effects 0.000 description 1
- 238000000605 extraction Methods 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 230000004044 response Effects 0.000 description 1
- 230000005236 sound signal Effects 0.000 description 1
- 230000001360 synchronised effect Effects 0.000 description 1
- 238000003786 synthesis reaction Methods 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/14—Fourier, Walsh or analogous domain transformations, e.g. Laplace, Hilbert, Karhunen-Loeve, transforms
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01B—MEASURING LENGTH, THICKNESS OR SIMILAR LINEAR DIMENSIONS; MEASURING ANGLES; MEASURING AREAS; MEASURING IRREGULARITIES OF SURFACES OR CONTOURS
- G01B9/00—Measuring instruments characterised by the use of optical techniques
- G01B9/02—Interferometers
- G01B9/02083—Interferometers characterised by particular signal processing and presentation
- G01B9/02084—Processing in the Fourier or frequency domain when not imaged in the frequency domain
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01B—MEASURING LENGTH, THICKNESS OR SIMILAR LINEAR DIMENSIONS; MEASURING ANGLES; MEASURING AREAS; MEASURING IRREGULARITIES OF SURFACES OR CONTOURS
- G01B9/00—Measuring instruments characterised by the use of optical techniques
- G01B9/02—Interferometers
- G01B9/0209—Low-coherence interferometers
- G01B9/02091—Tomographic interferometers, e.g. based on optical coherence
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T5/00—Image enhancement or restoration
- G06T5/10—Image enhancement or restoration using non-spatial domain filtering
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T5/00—Image enhancement or restoration
- G06T5/70—Denoising; Smoothing
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/10—Image acquisition modality
- G06T2207/10072—Tomographic images
- G06T2207/10101—Optical tomography; Optical coherence tomography [OCT]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/20—Special algorithmic details
- G06T2207/20048—Transform domain processing
- G06T2207/20056—Discrete and fast Fourier transform, [DFT, FFT]
Landscapes
- Physics & Mathematics (AREA)
- Engineering & Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Mathematical Physics (AREA)
- Signal Processing (AREA)
- Data Mining & Analysis (AREA)
- Pure & Applied Mathematics (AREA)
- Computational Mathematics (AREA)
- Mathematical Analysis (AREA)
- Mathematical Optimization (AREA)
- General Health & Medical Sciences (AREA)
- Radiology & Medical Imaging (AREA)
- Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
- Health & Medical Sciences (AREA)
- Algebra (AREA)
- Databases & Information Systems (AREA)
- Software Systems (AREA)
- General Engineering & Computer Science (AREA)
- Investigating Or Analysing Materials By Optical Means (AREA)
- Length Measuring Devices By Optical Means (AREA)
- Optical Radar Systems And Details Thereof (AREA)
Description
本発明は、光スペクトル干渉信号を取得し、前記光スペクトル干渉信号をフーリエ変換することで測定対象の断層に関する情報信号を生成する信号生成方法等、フーリエ変換を用いて情報信号を生成する信号生成方法に関する。
フーリエ変換を用いて、各種信号を生成することは、画像処理分野、計測各種分野、音声解析分野等で広く行われている。
例として、光スペクトル干渉信号をフーリエ変換することによって被測定対象の断層情報信号を得る光干渉断層撮像装置(Fourier Domain - Optical Coherence Tomography: FD-OCT)がある。
FD-OCTでは、光源出力を二つ以上に分け、一つを参照光とし他の光を測定光として被測定対象に照射する。被測定対象から散乱光あるいは反射光が戻り、前述の参照光との光スペクトル干渉信号を取得する。該光スペクトル干渉信号は波数空間軸上で観測され、参照光路と測定光路との光路長差に応じて波数空間軸上で振動した信号が得られる。そのため、取得された光スペクトル干渉信号をフーリエ変換することで、光路長差に応じてピークを示す断層情報信号が得られる。
ここで、光スペクトル干渉信号の強度は参照光と測定対象からの戻り光の強度の積に比例するため、測定対象からの戻り光が吸収や散乱、あるいは透過によって減衰しても高感度に断層情報信号を得ることができる。
そして、光スペクトル干渉信号をフーリエ変換した断層情報信号は、光路長差に応じた周波数を有する正弦波のフーリエ変換信号とスペクトル形状をフーリエ変換した形状の畳みこみ演算となる。
そのため、スペクトル帯域(例えば、光源の波長帯域に依存)が広ければ広いほど奥行き方向に高い分解能(層構造を分解表示できる能力)を有した断層情報信号が得られる。
しかし、一般にスペクトルの帯域は有限であり、またスペクトル形状はある形状を有している。そのため、断層情報信号の形状はスペクトル形状のフーリエ変換した形状を反映する。このことにより、本来一層である断層が複数の層に分かれているように画像表示される、あるいは複数層を分解できず一層に表示される、あるいは強度の大きな断層情報信号の裾に強度の小さな断層情報信号が埋もれてしまうといった断層情報信号の劣化が生じる。
そこで、前記断層情報信号の劣化を抑制するため、これまで、断層情報信号を整形するために、取得される光スペクトル干渉信号に窓関数を乗じて波形を整形する方法が取られている。この方法によって、断層情報信号は単峰でノイズが抑制される。しかし、窓関数を乗じることで光スペクトル帯域が狭くなり、断層情報信号の幅が広がり断層像の奥行き方向の分解能が劣化してしまう。
このような信号の劣化は、光断層計測に限らず、フーリエ変換を伴う信号処理・信号解析処理には共通の課題である。
前述のような課題を解決するために、これまでレーダ等の場合において、例えば、特許文献1のようなSuper spatially variant apodizationと呼ばれる、複数回の補外を行い、回折限界を超えて分解を行う方法が提案されている。この方法では、フーリエ変換された信号のピークを中心に広帯域に波形を抽出し逆フーリエ変換することで外挿信号を生成する。そして、外挿信号を元の波形に外挿することで分解能を高める。この信号処理を繰り返し行うことで、より高い分解能で信号を取得している。
また、音声分析合成の分野でも、例えば特許文献2のような周波数分析を行い、その分析結果から信号を合成する方法が提案されている。この方法では、メインの周波数成分のみを抽出し、位相補正を行って音声を合成している。
特許文献1に記載の方法では、ピーク近傍ではなく広帯域に波形を抽出している。そのため、抽出した帯域内にノイズ成分が存在する場合には、該ノイズ成分をも信号として抽出してしまうため、ノイズレベルを十分に下げることができない。
特許文献2に記載の方法は、近接したピークを消去してしまい、本来あるべき情報をも消去してしまう怖れがある。
このようにこれまでの方法では、信号の分解能を高めることと、ノイズ成分を抑制することとを同時に達成するのは、困難であるというのが実状である。
本発明により提供される情報信号生成方法は、
信号f(a)をフーリエ変換することで情報信号F(b)を生成する信号生成方法において、
前記情報信号F(b)のピーク強度の最大値を検出する工程と、
前記ピーク強度の最大値に対応する前記信号f(a)の振幅と位相と周波数とを算出する工程と、
前記ピーク強度の最大値に対応する前記信号f(a)を前記振幅及び位相及び周波数情報に基づきa軸上に延長した信号を生成する工程と、
a1 < a2として、前記延長した信号からa1よりも小さい領域の信号と、前記延長した信号のうちa2よりも大きい領域の信号とを抽出して外挿信号を生成する工程と、
前記外挿信号と前記a1からa2までの信号f(a)とを合成して合成信号を生成する工程と、
前記合成信号をさらにフーリエ変換する工程と、
を有することを特徴とする。
信号f(a)をフーリエ変換することで情報信号F(b)を生成する信号生成方法において、
前記情報信号F(b)のピーク強度の最大値を検出する工程と、
前記ピーク強度の最大値に対応する前記信号f(a)の振幅と位相と周波数とを算出する工程と、
前記ピーク強度の最大値に対応する前記信号f(a)を前記振幅及び位相及び周波数情報に基づきa軸上に延長した信号を生成する工程と、
a1 < a2として、前記延長した信号からa1よりも小さい領域の信号と、前記延長した信号のうちa2よりも大きい領域の信号とを抽出して外挿信号を生成する工程と、
前記外挿信号と前記a1からa2までの信号f(a)とを合成して合成信号を生成する工程と、
前記合成信号をさらにフーリエ変換する工程と、
を有することを特徴とする。
本発明の信号生成方法では、情報信号F(b)のピーク強度の最大値に対応する信号f(a)に基づき外挿信号を生成し、これを用いて合成信号を生成し、さらに合成信号をフーリエ変換する。
即ち、情報信号のピーク周辺の広範囲なデータを用いるのではなく、ピークが示す振幅と位相と周波数とを用いて外挿信号を生成し、合成信号を生成することで、ピークに付随したサイドノイズを抑制できる。これにより、サイドノイズに埋もれていたピーク強度値の小さな情報信号が消去されることなく残り、より詳細な情報信号を得ることが可能となる。
本発明の信号f(a)をフーリエ変換することで情報信号F(b)を生成する情報信号生成方法は、以下の工程を備える。
即ち、
前記情報信号F(b)のピーク強度の最大値を検出する工程。
前記ピーク強度の最大値に対応する前記信号f(a)の振幅と位相と周波数とを算出する工程。
前記ピーク強度の最大値に対応する前記信号f(a)を前記振幅及び位相及び周波数情報に基づきa軸上に延長した信号f’(a)を生成する工程。
a1 < a2として、前記延長した信号からa1よりも小さい領域の信号と、前記延長した信号のうちa2よりも大きい領域の信号とを抽出して外挿信号を生成する工程。
前記外挿信号と前記a1からa2までの信号f(a)とを合成して合成信号を生成する工程。
そして前記合成信号をさらにフーリエ変換する工程である。
前記情報信号F(b)のピーク強度の最大値を検出する工程。
前記ピーク強度の最大値に対応する前記信号f(a)の振幅と位相と周波数とを算出する工程。
前記ピーク強度の最大値に対応する前記信号f(a)を前記振幅及び位相及び周波数情報に基づきa軸上に延長した信号f’(a)を生成する工程。
a1 < a2として、前記延長した信号からa1よりも小さい領域の信号と、前記延長した信号のうちa2よりも大きい領域の信号とを抽出して外挿信号を生成する工程。
前記外挿信号と前記a1からa2までの信号f(a)とを合成して合成信号を生成する工程。
そして前記合成信号をさらにフーリエ変換する工程である。
さらに、前記信号生成方法を、複数回繰り返すことも本発明は包含する。
本発明におけるf(a)及びF(b)としては、フーリエ領域光干渉断層撮像方法(Fourier Domain - Optical Coherence Tomography: FD-OCT)における光スペクトル干渉信号及び断層情報信号や、フーリエ変換赤外分光法(FT-IR)における干渉信号と分光情報信号や、光強度信号とその解析信号や、音声信号と音声解析信号等が採用可能である。
本発明の実施の形態について、以下、FD-OCTを例に詳細に説明する。
FD-OCTにおいて、取得される光スペクトル干渉信号I(k)は、以下の式(1)で表される。
したがって、前記式(1)で表される光スペクトル干渉信号をフーリエ変換して取得される断層情報信号は、
となり、光路長差に応じた周波数を有する正弦波のフーリエ変換信号とスペクトル形状をフーリエ変換した形状の畳みこみ演算となる。そのため、スペクトル帯域(例えば、光源の波長帯域に依存)が広ければ広いほど奥行き方向に高い分解能(層構造を分解表示できる能力)を有した断層情報信号が得られる。
しかし、一般にスペクトルは無限に広帯域ではなく、ある形状S(k)を有しているため、断層情報信号の形状はスペクトル形状のフーリエ変換形状を反映する。
そこで、本発明では、取得された干渉信号の正弦波成分に注目し、スペクトル形状にとらわれない断層情報信号を以下の信号生成方法によって取得可能とする。
例えば、本発明では、外挿信号を広帯域に取ることにより、窓関数を広帯域に乗じることが可能となり、外挿信号を外挿されたピーク強度の断層情報信号の幅を狭め、断層に関する情報信号の分解能を高める。
一例として、ガウシアン窓を乗じてフーリエ変換する場合、断層情報信号の分解能Δzは以下の式(3)で表される。
図5を参照して本発明の効果について説明する。図5において301は、本発明の方法で得られるスペクトルを示しており、301は、光源のスペクトル帯域を示す波数ksからkeの領域以外の領域、即ち、ks未満及びkeを超える領域について外挿信号を加えて得られたものである。
より具体的には、外挿信号により光源出力のスペクトル帯域の端(ks及びke)が1/2となるようなガウシアン窓を乗じてフーリエ変換することで得られたスペクトルである。
一方、302は、本発明を適用しない従来技術に相当するスペクトルを示しており、光源出力のスペクトル帯域の端が1/e2となるようなガウシアン窓を乗じてフーリエ変換して得られるものである。
ここで、スペクトル301と302とを対比すると、本発明によるスペクトルの半値全幅Δkは、302に比べ、およそ2.4倍広くなるため、分解能はおよそ2.4倍改善することとなる。
外挿信号を広帯域に取ることにより、窓関数を広帯域に乗じることが可能となり、外挿信号を外挿されたピーク強度の断層情報信号周辺のサイドノイズを下げることが可能となる。
次に図6を参照して、本発明の方法で得られる信号のノイズ低減について説明する。
図6において、(A)は本発明を適用しない従来技術に相当する信号スペクトル(ガウシアン関数が半値で0となるような関数)を示しており、(B)はその信号スペクトルをフーリエ変換して得られる信号を示している。
図6(C)は本発明を適用した信号スペクトルを示しており、(D)はその信号スペクトルをフーリエ変換して得られる信号を示している。
図6(C)に示すように、窓関数の裾の値が十分小さくなるまで外挿信号を外挿とすることが可能であり、図6(D)に示すようにフーリエ変換した信号ノイズレベル(サイドノイズ)をピーク強度値によっては10桁以上小さくすることができる。ここで、延長する範囲を広げれば広げるほど、ノイズを低減することができる。ガウシアン窓関数の値が1/e2となるまで延長した場合には、ノイズレベルは1桁程度改善され、ガウシアン窓関数の値が0.001となるまで延長した場合には、ノイズレベルは12桁程度改善される。
一方で、延長すればするほどデータ数が増え、信号処理に時間がかかることになる。そこで、有限の延長範囲を決定する必要がある。
OCTでは、画像の信号対雑音比(SNR)が96dB以上が望まれており、SNRが96dBを超えるノイズレベルとなるように、延長する幅を決定する。
以下、図1及び図2を参照して、本発明の信号生成方法における工程を具体的に説明する。
図1(A)は、スペクトル帯域が波数ksからkeの範囲の光源を用いて得られた光干渉スペクトル信号を示している。これを具体例で説明すると、発光波数を広帯域に有する光源あるいは、発光波数が時間的に広帯域に変化する光源から出力された光を少なくとも2つに分割し、前記光の一方を測定対象に照射することで得られる測定対象からの散乱光あるいは反射光と、前記光の他方と、を干渉させることで光スペクトル干渉信号を示している。
図1(B)は、本発明において必須の工程ではなく必要に応じて行われるものであるが、図1(A)の光スペクトル干渉信号にエンベロープを整形する窓関数を乗じて得られた信号を示している。
図1(C)は、図1(B)の光スペクトル干渉信号をフーリエ変換することで得られる測定対象の断層に関する情報信号を示している。ここでは、窓関数を乗じた信号をフーリエ変換する例を示しているが、窓関数を乗じない図1(A)の干渉信号をフーリエ変換する場合を本発明は包含する。
図1(D)に示すように情報信号からピーク強度の最大値101を検出する。
ピーク強度の最大値101に対応する光スペクトル干渉信号の振幅と位相と周波数 とを算出し、図2(E)に示すようにピーク強度の最大値101に対応する前記光スペクトル干渉信号を前記振幅及び位相及び周波数情報に基づき波数空間軸上に延長した信号を生成する。次に、本発明において必須の工程ではなく必要に応じて行われるものであるが、延長した信号(図2(E))に窓関数を乗じて、図2(F)を得る。
次いで、窓関数を乗じた信号から、光源が発光する波数帯域中の波数ksと波数keをks < keとおいて、ksよりも小さい領域の信号と延長した信号のうち波数keよりも大きい領域の信号とを抽出し外挿信号を生成する(図2(G))。
光スペクトル干渉信号から波数ksから波数keまでの信号を抽出した抽出信号と前記外挿信号とを合成して合成信号を生成する(図2(H))。
得られた合成信号をさらにフーリエ変換する(図2(I))。
図2(I)に示した情報信号は、図1(C)に示した情報信号に比べてスペクトル幅が狭く、格段にノイズが抑制されたものとなっている。
これより上述した一連の工程によって、ピーク強度の最大値を検出することで、正しく断層情報信号のみを抽出することができ、サイドノイズ成分を抑制することができることが理解される。
本発明の信号生成方法を断層に関する情報信号を生成する方法に適用する形態では、本発明は、以下の工程を有する信号生成方法を包含する。
即ち、前記信号f(a)を光スペクトル干渉信号とし、前記情報信号F(b)を断層に関する情報信号として、前記断層に関する情報信号のピーク強度の最大値を検出する工程。
前記断層に関する情報信号のピーク強度の最大値に対応する前記光スペクトル干渉信号の振幅と位相と周波数とを算出する工程。
前記振幅、位相及び周波数情報に基づき、前記a軸を波数軸として該波数軸上に延長した信号を生成する工程。
波数ks < keとして、前記延長した信号からksよりも小さい領域の信号と、前記延長した信号のうちkeよりも大きい領域の信号とを抽出して外挿信号を生成する工程。
前記外挿信号と前記ksからkeまでの光スペクトル干渉信号とを合成して合成信号を生成する工程。及び前記合成信号をさらにフーリエ変換する工程である。
本発明によると抽出信号のスペクトル帯域ksからkeまでのスペクトル形状が、理想とする窓関数とは大きく異なる形状をしている場合には、スペクトル形状を整形することで、断層情報信号を単峰化し断層画像の画質を改善することができる。
次に、測定を行う検体が複数層からなる場合に有用な繰返し外挿信号を生成し合成波形を取得することで、断層情報信号全体に渡って高感度化、高分解能化する手法について、図4と図3を参照しながら説明する。この手法によると、サイドノイズ成分に埋もれていた小さな断層情報信号を見出すことができる。
図4に示すようにまず、光スペクトル干渉信号を取得する(201)。ここでは、取得された干渉信号は等波数間隔となるように記録する、あるいは信号処理により等波数間隔となるようにデータ変換することが好適である。干渉信号のエンベロープは光源出力のスペクトル強度と一致している。そこで別に取得された光源出力のスペクトル強度を用いて、前記干渉信号のエンベロープをフラットトップ化し矩形波形状に整形することも可能である。 また、前記矩形波整形された干渉信号に窓関数を乗じることも可能である。窓関数を乗ずることは必要に応じて行うことができる。
次いで、干渉信号をフーリエ変換することで前記測定対象の断層に関する情報信号を取得する(202)。
前記取得された断層に関する情報信号からピーク強度がn番目の値のピークを検出する(203)。
本願発明では、得られた情報信号(例えば、断層に関する情報信号)のピーク強度の最大値に対応して、外挿信号を生成すること、合成信号の生成と、フーリエ変換を所定の回数Nまで繰り返すことを包含する。
ここで、第n番目(ここで、nは1以上の整数)のピーク検出が所定のNではない場合には、前記検出された第n番目のピーク強度の値が対応する振幅、位相、周波数を算出する(205)。尚、第n番目のピーク強度は、参照光路と測定光路の光路長差に応じた周波数で振動する干渉信号の振動信号を表している。算出された振幅、位相、周波数情報に基づき、外挿する波数帯域まで延長(波数空間軸上に延長)された振動波形を生成する(206)。
窓関数を用いる場合には、延長された信号に前記窓関数を乗じる工程を付加しても良い。波数空間軸上に延長された信号から(または、前記窓関数を乗じられた信号から)波数ksよりも小さい領域の信号と波数がkeよりも大きい領域の信号とを抽出し外挿信号を生成する(207)。
前記光スペクトル干渉信号から波数ksから波数keまでの信号を抽出した抽出信号と前記外挿信号とを合成して合成信号を生成する(208)。
前記合成信号をフーリエ変換し断層情報信号を取得する(209)。
ここまでの処理により、検出されたピーク強度の信号は高分解能化され、該ピーク強度周辺のサイドノイズ成分は抑制される。
上記の203〜209の工程をピーク強度の値が大きい方からn=1番目、2番目、3番目・・・N番目まで繰り返す。
上記203の工程では、ノイズとは区別された信号成分のみを抽出している。そのため、工程201〜209及び、210により光源スペクトル帯域に依存して生じるサイドノイズ成分を抑制し、高分解能、高感度な断層情報信号を得ることが可能となる。
図3の(A)は、光スペクトル干渉信号を示しており、図3の(B)は、光スペクトル干渉信号をフーリエ変換して得られた情報信号を示している。図3(B)において351は、情報信号のピーク強度の最大値である(図4における203の工程)。第n番目のピーク値でいえば1番目となる。
図3(C)は、ピーク強度の最大値351に対応する光スペクトル干渉信号の振幅、位相及び周波数に基づき波数空間軸上に延長された信号(外挿信号)を用いた合成信号(図4における208の工程)を示している。
図3(D)は、合成信号をフーリエ変換して得られる情報信号であり、352は、最初にフーリエ変換して得られた情報信号の2番目のピーク強度を示している。ここで、352は、第(n+1)番目のピーク値では、2番目のピーク値にあたる。2番目のピーク強度352について同様の工程を繰り返す。図3(F)における353、図3(H)における354、図3(J)における355は、最初にフーリエ変換して得られた情報信号の3番目、4番目、5番目のピーク強度をそれぞれ示している。図3の(B)と図3(L)との対比から理解されるように、本願発明のフーリエ変換のサイクルを繰り返すことにより格段にノイズ成分が抑制された情報信号が得られる。
また、ピーク強度の大きな情報信号から順に上記信号処理を繰返し行うことで、サイドノイズに埋もれていた断層情報信号を正しく抽出して、断層に関する情報信号全体を高分解能化、高感度化することが可能である。
前記Nは前記断層に関する情報信号のピークの数と捉えることもできる。この観点で信号生成を行うと、最小の繰返し回数で断層に関する情報信号全体を高分解能化、高感度化することが可能となる。ピーク数Nが分からない場合には、前記Nは、検出されたピーク強度が大きい方から数えて所定の強度以下となるピークまでの数としても良い。このことにより、演算を無限に繰り返すことなく、ある値以上の断層に関する信号を高分解能化、高感度化することが可能となる。
また、ピーク強度の値の大きい方からN回フーリエ変換を含む工程を繰り返すのではなく、ある閾値を設けて閾値以上のピークを検出する方法も採用することができる。
この方法は、例えば以下の工程を取り得る。
(a) 取得された干渉信号は等波数間隔となるように記録する、あるいは信号処理により等波数間隔となるようにデータ変換する。
(b) 前記干渉信号のエンベロープは光源出力のスペクトル強度と一致している。そこで別に取得された光源出力のスペクトル強度を用いて、前記干渉信号のエンベロープをフラットトップ化し矩形波形状に整形する。
(c) 前記矩形波整形された干渉信号に窓関数を乗じる。
(d) 前記窓関数を乗じられた干渉信号をフーリエ変換することで前記測定対象の断層に関する情報信号を取得する。
(e) 前記取得された断層に関する情報信号から最大ピーク値を検出し、該最大ピーク値に所定の割合をm回乗じた値をm番目の閾値とし、ピーク強度が該m番目の閾値以上のピークを全て検出する。
(f) 前記m番目の閾値以上の各々のピーク強度は、参照光路と測定光路の光路長差に応じた周波数で振動する干渉信号の振動信号を表している。そこで、前記検出された各々のピーク強度の値が対応する振幅、位相、周波数を算出する。
(g) 算出された各々の振幅、位相、周波数情報に基づき、各々の外挿する波数帯域分の振動波形を生成し、足し合わせて外挿信号とする。
(h) 前記外挿信号に前記窓関数を乗じる。
(i) 前記窓関数を乗じられた外挿信号から波数ksよりも小さい領域の信号と前記外挿信号のうち波数がkeよりも大きい領域の信号とを抽出する。
(j) 前記抽出された外挿信号を、前記窓関数を乗じられた干渉信号に外挿した合成信号を生成する。
(k) 前記合成信号をフーリエ変換し断層情報信号を取得する。ここまでの処理により、検出されたピーク強度の信号は高分解能化され、該ピーク強度周辺のサイドノイズ成分は抑制される。
(l) 上記(e)〜(k)の操作を、閾値を前記所定の割合だけ下げながらm=1番目、2番目、3番目・・・M番目まで繰り返す。
(a) 取得された干渉信号は等波数間隔となるように記録する、あるいは信号処理により等波数間隔となるようにデータ変換する。
(b) 前記干渉信号のエンベロープは光源出力のスペクトル強度と一致している。そこで別に取得された光源出力のスペクトル強度を用いて、前記干渉信号のエンベロープをフラットトップ化し矩形波形状に整形する。
(c) 前記矩形波整形された干渉信号に窓関数を乗じる。
(d) 前記窓関数を乗じられた干渉信号をフーリエ変換することで前記測定対象の断層に関する情報信号を取得する。
(e) 前記取得された断層に関する情報信号から最大ピーク値を検出し、該最大ピーク値に所定の割合をm回乗じた値をm番目の閾値とし、ピーク強度が該m番目の閾値以上のピークを全て検出する。
(f) 前記m番目の閾値以上の各々のピーク強度は、参照光路と測定光路の光路長差に応じた周波数で振動する干渉信号の振動信号を表している。そこで、前記検出された各々のピーク強度の値が対応する振幅、位相、周波数を算出する。
(g) 算出された各々の振幅、位相、周波数情報に基づき、各々の外挿する波数帯域分の振動波形を生成し、足し合わせて外挿信号とする。
(h) 前記外挿信号に前記窓関数を乗じる。
(i) 前記窓関数を乗じられた外挿信号から波数ksよりも小さい領域の信号と前記外挿信号のうち波数がkeよりも大きい領域の信号とを抽出する。
(j) 前記抽出された外挿信号を、前記窓関数を乗じられた干渉信号に外挿した合成信号を生成する。
(k) 前記合成信号をフーリエ変換し断層情報信号を取得する。ここまでの処理により、検出されたピーク強度の信号は高分解能化され、該ピーク強度周辺のサイドノイズ成分は抑制される。
(l) 上記(e)〜(k)の操作を、閾値を前記所定の割合だけ下げながらm=1番目、2番目、3番目・・・M番目まで繰り返す。
閾値以上のピークを一括で検出し、検出された全てのピークについて外挿信号を生成し、元の信号に外挿することで、演算回数を減らし高速な信号生成方法が可能となる。
また、信号生成工程(l)において、Mは予め決めておく。あるいは、M番目と決めておくのではなく、検出されたピークの数が所定の数となるまで繰り返しても良い。このことにより、演算回数を限定することが可能となる。あるいは、閾値が予め決めておいた値以下になるまで繰り返しても良い。このことにより、演算を無限に繰り返すことなく、ある値以上の断層に関する信号全体を高分解能化、高感度化することが可能となる。
また、前記閾値を下げる所定の割合は、干渉信号のエンベロープ波形を取得し、該エンベロープ波形をフーリエ変換することで得られる信号の強度の最大値と2番目のピーク値との比としても良い。このことにより、前記閾値を下げる所定の割合を最大にし、効率良く閾値を下げることが可能となる。
本発明の信号生成方法は、信号f(a)に、関数が示す波形の裾部が0となる窓関数を乗じ、延長した信号を前記窓関数が0となるaまで延長されたものである形態を包含する。
また、前記信号f’(a)を延長する範囲は、前記延長した範囲の前記窓関数をフーリエ変換した信号の強度の最大のピーク値とこれに次ぐ2番目のピーク値との比が所定の値となる範囲以上である形態を包含する。
また、前記信号f(a)は、該信号f(a)を該信号f(a)のエンベロープ情報で割ることで該信号f(a)の形状をフラットトップ形状に整形し、これに窓関数を乗じたものとすることができる。
本発明の信号生成方法では、ピーク強度値の大きな情報信号から順番に繰返し外挿信号を生成し外挿することで、ピーク強度値の大きさに比例したサイドノイズ成分が順番に抑制される。ここで抑制されるサイドノイズ成分は検出されたピークに付随しているため、サイドノイズに埋もれていた情報信号をも消去することはない。したがって、ピーク強度値の大きな情報信号から順番にサイドノイズ成分を抑制することで、正しく情報信号のみを抽出し、情報信号全体に渡ってノイズ成分を抑制することが可能となる。
また、ピーク強度値の大きな情報信号から一つずつ順番に処理を行うと演算回数が多くなり、信号生成時間が長くなってしまう。そこで、閾値を設けて、閾値以上のピークを一括で検出し、検出された全てのピークについて外挿信号を生成し、元の信号に外挿することで、演算回数を減らし高速な信号生成方法が可能となる。
前記閾値はなるべく小さくすることで、正しい情報信号のピークを一度に多く検出できる。しかし、閾値を小さくし過ぎるとサイドノイズをも情報信号として検出してしまう。一方で、サイドノイズの大きさは、元の信号のエンベロープ形状をフーリエ変換した波形で決まる。そこで、該エンベロープ形状をフーリエ変換した波形のピーク強度の最大値と、2番目の大きさを有するピーク値の比を計算し、実際に得られた情報信号のピーク強度の最大値に乗じた値を閾値とすることで、サイドノイズ成分を回避して正しく情報信号に関わるピーク強度を検出することが可能となる。
本発明により、外挿信号は理論上スペクトル帯域を無限に広げることができる。しかしながら、スペクトル帯域を無限に取ると、データ数が大きくなりすぎ、フーリエ変換などの工程で時間がかかってしまう。また情報信号の形状は、元の信号のエンベロープ形状をフーリエ変換した波形形状を反映する。そのため、一般に元の信号にcos関数やハニング関数、あるいはガウシアン関数等の、対称かつ裾が0に近づく窓関数を乗じることで、情報信号の形状を整形する。そこで、外挿信号の範囲を窓関数の裾が十分0に近づく範囲までとすることで、データ数を限定でき情報信号の高分解能化、低ノイズ化を高速に行うことが可能となる。
また、元の信号のエンベロープ形状がフラットトップでない場合、元の信号に直接窓関数を乗じても、元の信号のエンベロープ形状の影響が残り、情報信号にも影響が及ぶ。そこで、エンベロープ形状情報で光スペクトル干渉信号を割り、該信号をフラットトップ化することで、元の信号のエンベロープ形状の影響を抑制し、窓関数による整形を行うことが可能となる。
本発明は、信号f(a)をフーリエ変換することで情報信号F(b)を生成する情報信号処理装置を包含する。この情報信号処理装置は、以下を有して構成される。
即ち、前記情報信号F(b)のピーク強度の最大値を検出する処理部。前記情報信号F(b)のピーク強度の最大値に対応する前記信号f(a)の振幅と位相と周波数とを算出する処理部。前記振幅及び位相及び周波数情報に基づきa軸上に延長した信号を生成する処理部。a1 < a2として、前記延長した信号からa1よりも小さい領域の信号と、前記延長した信号のうちa2よりも大きい領域の信号とを抽出して外挿信号を生成する処理部。前記外挿信号と前記a1からa2までの信号f(a)とを合成して合成信号を生成する処理部。及び前記合成信号をさらにフーリエ変換する処理部である。これらの処理部は、所定のプログラムを実行するための演算装置で構成することができ、具体例としては、コンピュータの演算部が挙げられる。
また、本発明は光干渉断層撮像装置を包含する。
本発明の光干渉断層撮像装置は、光源部と、前記光源部からの光を検体に照射し、検体からの反射光を伝達させる検体測定部と、前記光源部からの光を参照光として伝達させる参照部と、前記反射光と前記参照光とを干渉させる干渉部と、前記干渉部からの干渉光を検出する光検出部と、前記光検出部で検出された光に基づいて前記検体の断層像を得る画像処理部と、を有する。そして、画像処理部は、本発明の信号生成装置を備えることを特徴とする。
本実施例では、本発明の信号生成方法を用いた光干渉断層撮像装置の例を示す。
本実施例で示す光干渉断層撮像装置は、FD-OCTの中でも波長掃引光源型光干渉断層撮像装置(Swept Source - Optical Coherence Tomography: SS-OCT) であり、図7に模式図を示す。図7では、光ファイバを用いて装置を構成しているが、光ファイバを用いず空間で構成しても良い。
図7に示した撮像装置(SS-OCT装置)は、発光波長が10μsの間に800nmから880nmまで変化する光源501(光源部)を用いる。光源501は例えば、広帯域光源から切出す波長を時間的に変化させる光源が用いられる。
SS-OCT装置では、光源501から出力された光を光ファイバカプラ510によって二つに分割し、前記分割された光の一方を測定光とし測定対象515に照射し、測定対象からの散乱光あるいは反射光を再度光ファイバカプラ510(干渉部)に結合させ、光ファイバカプラ518へと伝搬させる。この時、測定光は直交した二枚のガルバノミラー512、513によって構成された測定光走査光学系(検体測定部)を通る。測定光は、該測定光走査光学系により、測定対象上を走査させられる。
一方、前記分割された光の他方を参照光として伝達させるディレイライン等(参照部)を伝搬して、光ファイバカプラ518へと伝搬させる。ここでディレイラインは、射出される光ファイバ端516と入射される光ファイバ端517間の距離を変化させることが可能である空間光学系である。
前記散乱光あるいは反射光と前記参照光とを干渉させ、干渉光を検出する差動光検出器519(光検出部)によって断層に関する光スペクトル干渉信号を取得する。差動光検出器519は205MHzの応答速度で光を検出する。したがって、SS-OCT装置の一波長掃引中に2050点の断層に関する光スペクトル干渉データを出力する。出力された断層に関する光スペクトル干渉データはADボードを介して断層像を得るための画像処理部を構成する計算機520に取り込まれる。
また、波長掃引光源から出力された光は光ファイバカプラ504によって、一部、マッハチェンダー干渉計に導波される。マッハチェンダー干渉計からの光スペクトル干渉信号の差動光検出データは、等波数間隔で0と交わる為、時間的な波数の変化を示す波数クロック信号として前記ADボードの別のチャンネルを介して計算機520に取り込まれる。該波数クロック信号を用いて、前記断層に関する光スペクトル干渉データを等波数間隔の干渉信号となるように計算機内部でデータ変換される。この時、データ数が2048点となるようにデータ変換を行う。そのため、波数クロック信号が0と交わる点を2048点以上取る為に、マッハチェンダー干渉計の二つのアームの光路長差を8mm以上とする。
光ファイバカプラ502で一部分けられた光は、光検出器503によって検出され、光源の強度変化をモニタしている。強度変化情報もまたADボードによって計算機520に取り込まれ、光スペクトル干渉信号のエンベロープとして利用され、フラットトップ化に用いられる。信号をフラットトップ化することで、元の信号のエンベロープ形状の影響を抑制し、窓関数による整形を行うことが可能となる。
ADボードによるデータ取り込みは、波長掃引光源から出力されたトリガ信号521を用いて同期が取られる。
前記等波数間隔にデータ変換された干渉信号を高速フーリエ変換し、測定光が照射された光軸上の断層情報信号を取得する。ここで、フーリエ変換された信号は折り返し信号を含むので、フーリエ変換された信号の半分を抽出して断層情報信号とする。
SS-OCT装置では、参照光路(すなわち、参照光路中の最初に光を分岐する光ファイバカプラ510からディレイラインを通って干渉が起こる光ファイバカプラ518までの光路)と測定光路(すなわち、測定光路中の最初に光を分岐する光ファイバカプラ510から測定対象に照射されて干渉が起こる光ファイバカプラ518までの光路)との光路長差が長くなるに従い光スペクトル干渉信号の周波数は高くなる。そして装置の構成上、前記光路長差が0mmに比べておよそ8mmで光スペクトル干渉信号の振幅は1/2の大きさになる。そこで、断層画像を取得する前に、測定対象の最も上面の断層情報信号が0付近に来るように前記ディレイライン(射出される光ファイバ端516と入射される光ファイバ端517間の距離)を調節する。
前記測定光走査光学系を構成するガルバノミラーの1枚512を1024個の方向に角度を変化させ、各方向について、光スペクトル干渉信号を取得し、等波数間隔になるようデータ変換し、高速フーリエ変換する、という一連のデータ処理を行う。前記ガルバノミラー512の各方向での断層情報信号の強度に基づいてグレースケールで濃淡をつけ、1024個の方向の断層情報を並べることで断層画像を得る。
さらに、もう1枚のガルバノミラー513の角度を1024個の方向に変化させることで、最終的に1024×1024×1024の3次元ボリューム断層データを取得することができる。
図7において505と508とは波数クロック信号発生用のマッハツェンダー干渉計内のファイバカプラであり、506と507とはコリメータレンズである。509は、差動光検出器であり、502と504はファイバカプラである。
次に、本発明による信号生成方法を適用し、断層画像の画質の向上を行う。
本発明による信号生成方法は、取得された断層情報信号の1方向ずつについて適用される。
まず、取得されたある1方向の断層情報信号について、前記断層情報信号のピーク強度の最大値を検出する。
次に、前記ピーク強度の最大値に対応する前記光スペクトル干渉信号の振幅と位相と周波数とを算出する。該振幅情報は、検出されたピーク値の大きさから算出される。該位相情報は、光スペクトル干渉信号をフーリエ変換したデータの実部と虚部から算出される。周波数は、検出されたピークの位置から算出される。該算出された振幅及び位相及び周波数情報に基づき波数空間軸上に延長した正弦波信号を生成する。
今回、光源から出力された光スペクトル形状は、波数7496982.469[m-1](波長838.1[nm]に対応)で最大値となり、波数7139983[m-1] (波長880[nm]に対応)と波数7853982[m-1] (波長800[nm]に対応)とで光スペクトル強度が半値となるガウシアン関数形状であった。
そこで、ガウシアン関数の値が最大値に対して1/e2まで落ちる波数8460395.981[m-1] (波長742.6585376[nm]に対応)から波数6533568.956[m-1] (波長961.6773542[nm]に対応)まで、前記正弦波信号を延長する。
そして、生成された前記正弦波信号に波数7496982.469[m-1] (波長838.1[nm]に対応)で最大値となり、波数7139983[m-1] (波長880[nm]に対応)と波数7853982[m-1] (波長800[nm]に対応)とで光スペクトル強度が半値となるガウシアン関数を乗じ、外挿信号を得る。
前記外挿信号から、波数7139983[m-1] (波長880[nm]に対応)よりも小さい領域の信号と前記外挿信号のうち波数7853982[m-1] (波長800[nm]に対応)よりも大きい領域の信号とを抽出する。
前記抽出した抽出信号を前記光スペクトル干渉信号に加えさらに高速フーリエ変換する。
前記の信号生成方法を前記ガルバノミラーの全ての方向について行うことで、3次元ボリューム断層データの画質の向上を行うことが可能となる。
本例のスペクトル領域光干渉断層撮像装置(Spectral Domain - Optical Coherence Tomography: SD-OCT)の模式図を図8に示す。図8では、光ファイバを用いて装置を構成しているが、光ファイバを用いず空間で構成しても良い。
このSD-OCT装置は発光波長を800[nm]から880[nm]まで有する広帯域な光源601を用いる。ここでは、光源は広帯域LD(Super Luminescent Diode: SLD)を使用したが、超広帯域連続光(Super-continuum light source: SC光)等を用いても良い。
SD-OCTでは光源から出力された光を光ファイバカプラ602によって二つに分割し、前記分割された光の一方を測定光とし測定対象607に照射し、測定対象からの散乱光あるいは反射光を再度光ファイバカプラ602に結合させる。この時、測定光は直交した二枚のガルバノミラー604, 605によって構成された測定光走査光学系を通る。測定光は、該測定光走査光学系により、測定対象上を走査させられる。
一方、前記分割された光の他方を参照光とし、光ファイバ端608からミラー609までの距離を変化させることが可能であるディレイラインを伝搬させる。前記参照光をディレイラインのミラー609に照射し、反射した参照光も再度光ファイバカプラ602に結合させる。
前記散乱光あるいは反射光と前記参照光とを干渉させ、分光器610によって光スペクトル干渉信号を取得する。分光器は2048ピクセルのラインセンサアレイで光源出力の全帯域(=800~880[nm])を検出する。検出された光スペクトル干渉信号は画像処理部を構成する計算機611に取り込まれる。
分光器610では、隣り合うピクセルが検出する光は等波数間隔ではなく等波長間隔であるので、予め用意しておいたデータを用いて、等波数間隔の干渉信号となるように計算機611内部でデータ変換される。
等波数間隔にデータ変換された干渉信号を高速フーリエ変換し、測定光が照射された光軸上の断層情報信号を取得する。ここで、フーリエ変換された信号は折り返し信号を含むので、フーリエ変換された信号の半分を抽出して断層情報信号とする。
SD-OCT装置では、参照光路(すなわち、光ファイバカプラ602の分岐点からミラー609までの光路)と測定光路(すなわち、光ファイバカプラ602の分岐点から測定対象607までの光路)との光路長差が長くなるに従い光スペクトル干渉信号の周波数は高くなる。そして装置の構成上、光スペクトル干渉信号の振幅は、前記光路長差が0[mm]に比べておよそ8[mm]で1/2の大きさになる。そこで、断層画像を取得する前に、測定対象の最も上面の断層情報信号が0付近に来るように前記ディレイライン(光ファイバ端608からミラー609までの距離)を調節する。
前記測定光走査光学系を構成するガルバノミラーの1枚604を1024個の方向に角度を変化させ、各方向について、光スペクトル干渉信号を取得し、等波数間隔になるようデータ変換し、高速フーリエ変換する、という一連のデータ処理を行う。前記ガルバノミラー604の各方向での断層情報信号の強度に基づいてグレースケールで濃淡をつけ、1024個の方向の断層情報を並べることで断層画像を得る。
さらに、もう1枚のガルバノミラー605の角度を1024個の方向に変化させることで、最終的に1024×1024×1024の3次元ボリューム断層データを取得することができる。
次に、本発明による信号生成方法を適用し、断層画像の画質の向上を行う。
本発明による信号生成方法は、取得された断層情報信号の1方向ずつについて適用される。
まず、取得されたある1方向の断層情報信号について、前記断層情報信号のピーク強度の最大値を検出する。
次に、前記ピーク強度の最大値に対応する前記光スペクトル干渉信号の振幅と位相と周波数とを算出する。該振幅情報は、検出されたピーク値の大きさから算出される。該位相情報は、光スペクトル干渉信号をフーリエ変換したデータの実部と虚部から算出される。周波数は、検出されたピークの位置から算出される。該算出された振幅及び位相及び周波数情報に基づき波数空間軸上に延長した正弦波信号を生成する。
今回、光源から出力された光スペクトル形状は、波数7496982.469[m-1] (波長838.1[nm]に対応)で最大値となり、波数7139983[m-1] (波長880[nm]に対応)と波数7853982[m-1] (波長800[nm]に対応)とで光スペクトル強度が半値となるガウシアン関数形状であった。そこで、ガウシアン関数の値が最大値に対して1/e2まで落ちる波数8460395.981[m-1] (波長742.6585376[nm]に対応)から波数6533568.956[m-1] (波長961.6773542[nm]に対応)まで、前記正弦波信号を延長する。そして、生成された前記正弦波信号に波数7496982.469[m-1] (波長838.1[nm]に対応)で最大値となり、波数7139983[m-1] (波長880[nm]に対応)と波数7853982[m-1] (波長800[nm]に対応)とで光スペクトル強度が半値となるガウシアン関数を乗じ、外挿信号を得る。
前記外挿信号から、波数7139983[m-1] (波長880[nm]に対応)よりも小さい領域の信号と前記外挿信号のうち波数7853982[m-1] (波長800[nm]に対応)よりも大きい領域の信号とを抽出する。
前記抽出した抽出信号を前記光スペクトル干渉信号に加えさらに高速フーリエ変換する。前記の信号生成方法を前記ガルバノミラー604、605の全ての方向について行うことで、3次元ボリューム断層データの画質の向上を行うことが可能となる。
まず、実施例1記載の信号処理を行い画質の向上を図る。引き続き、ピーク強度の値が大きい方から2番目のピークについて、ピークの抽出と、抽出されたピークについての振幅と位相と周波数との算出を行う。
次に、ガウシアン関数の値が最大値に対して1/e2まで落ちる波数8460395.981[m-1] (波長742.6585376[nm]に対応)から波数6533568.956[m-1] (波長961.6773542[nm]に対応)まで、前記算出された振幅と位相と周波数に基づいて正弦波信号を生成する。
該正弦波信号に波数7496982.469[m-1] (波長838。1[nm]に対応)で最大値となり、波数7139983[m-1] (波長880[nm]に対応)と波数7853982[m-1] (波長800[nm]に対応)とで光スペクトル強度が半値となるガウシアン関数を乗じ、外挿信号を得る。
該外挿信号から波数7139983[m-1] (波長880[nm]に対応)よりも小さい領域の信号と前記外挿信号のうち波数7853982[m-1] (波長800[nm]に対応)よりも大きい領域の信号とを抽出し、抽出された外挿信号と光スペクトル干渉信号を合成する。
該合成信号を高速フーリエ変換し、改善された断層情報信号を得る。断層がN層ある場合には、これをさらに3番目、4番目...N番目のピークまで繰り返す。
上記信号処理により、実施例1よりさらに断層画質を改善することができる。
断層数Nが分からない場合には、無限に演算を繰り返さないために、N=100として信号生成を行っても良い。あるいは、今回使用したADボードが12bitであるので、検出されたピーク値が10-3.6(≒2-12)となるまで繰り返す信号生成方法としても良い。このことで、ノイズレベルを10-3.6(≒2-12)まで低減することが可能となる。
実施例1と同様にして断層情報信号の取得を行う。
取得した断層情報信号において、一つずつピークを検出せず、閾値を設けて閾値以上のピークについて、それぞれのピークについての振幅と位相と周波数算出とを行う。
検出された全てのピークについて、ガウシアン関数の値が最大値に対して1/e2まで落ちる波数8460395.981[m-1] (波長742.6585376[nm]に対応)から波数6533568.956[m-1] (波長961.6773542[nm]に対応)まで、前記算出された振幅と位相と周波数に基づいて正弦波信号を生成する。
前記全ての正弦波信号に波数7496982.469[m-1] (波長838.1[nm]に対応)で最大値となり、波数7139983[m-1] (波長880nmに対応)と波数7853982[m-1] (波長800[nm]に対応)とで光スペクトル強度が半値となるガウシアン関数を乗じ、足し合わせることで、外挿信号を得る。
該外挿信号から波数7139983[m-1] (波長880[nm]に対応)よりも小さい領域の信号と前記外挿信号のうち波数7853982[m-1] (波長800[nm]に対応)よりも大きい領域の信号とを抽出し、抽出された外挿信号と光スペクトル干渉信号を合成する。該合成信号を高速フーリエ変換し、断層情報信号を得る。
このことにより、繰返し信号生成回数を減らすことができ、断層画像の向上の高速化が可能となる。
また閾値は、光源から出力される光スペクトル強度波形を取得し、該光スペクトル強度波形をフーリエ変換することで得られる信号の強度の最大値と2番目のピーク値を検出し、前記2番目のピーク値と最大のピーク値との比を算出し、光スペクトル干渉信号をフーリエ変換して得られた断層情報信号の最大のピーク値に前記比をかけることで得られる値として決定しても良い。
実施例1または実施例2と同様にして断層情報信号の取得を行う。取得した断層情報信号において、光スペクトル干渉信号をフーリエ変換して得られた断層情報信号の最大のピーク値P0に下げ幅0.2を乗じた値を閾値Tmとし、Tm=P0×0.2mとして、m=0、1、2…Mまで閾値Tm以上のピークについて、繰返し一括で外挿信号を生成、合成し、フーリエ変換を行う信号生成方法を、検出されるピーク数が断層数Nになるまで行う。
ここで、下げ幅0.2は、矩形波をフーリエ変換して得られるsinc関数の第1ピークと第2ピークとの比0.2を参考とした。
断層数Nが分からない場合には、無限に演算を繰り返さないために、N=100として信号生成を行っても良い。あるいは、繰返し計算数M=100として信号生成を行っても良い。
また、検出されたピーク値が10-3.6となるまで繰り返す信号生成方法としても良い。あるいは、Tm閾値が10-3.6となるまで繰り返す信号生成方法としても良い。このことで、ノイズレベルを10-3.6(≒2-12)まで低減することが可能となる。
下げ幅は0.2ではなく、光源から出力される光スペクトル強度波形を取得し、該光スペクトル強度波形をフーリエ変換することで得られる信号の強度の最大値と2番目のピーク値を検出し、前記2番目のピーク値と最大のピーク値との比を下げ幅としても良い。
実施例1または実施例2と同様にして断層情報信号の取得を行う。取得した前記断層情報信号において、前記光スペクトル干渉信号は中心波数7496982[m-1] (=波長838.1[nm])、波数半値全幅713998.3[m-1] (波数7853982~ 7139983[m-1] = 波長800~880[nm])の窓関数を乗じて整形する。前記窓関数はガウシアン関数でも良いし、cos関数やハニング関数でも良い。窓関数を乗じて整形することで、フーリエ変換後の断層に関する情報信号は単峰で低サイドノイズレベルの信号となる。また、cos窓やハニング窓の場合は、窓関数の裾部が0となる波数まで外挿信号を延長する。このことで、データ数を限定でき情報信号の高分解能化、低ノイズ化を高速に行うことが可能となる。該窓関数がガウシアン関数の場合は、有限の波数範囲で光スペクトル干渉信号を延長すると、該窓関数の裾部で0とならない。そこで、延長した有限の波数範囲は、該延長した有限の波数範囲内の窓関数をフーリエ変換した波形の最大値と2番目のピーク値との比が10-3.6となる範囲とする。このことで、データ数を限定でき情報信号の高分解能化、低ノイズ化を高速に行うことが可能となる。
また、該光スペクトル干渉信号は、別に取得しておいたスペクトル強度情報で割ることでフラットトップに整形された後、窓関数を乗じて整形されても良い。
前記光スペクトル干渉信号を取得する際に窓関数を乗じている場合には、同じ窓関数を外挿信号にも乗じる。
窓関数を乗じて光スペクトル干渉信号の包絡線を整形することで、フーリエ変換後の断層情報信号のピークを単峰化し、断層画質を改善することが可能となる。
本発明は上記実施の形態に制限されるものではなく、本発明の精神及び範囲から離脱することなく、様々な変更及び変形が可能である。従って、本発明の範囲を公にするために以下の請求項を添付する。
画像処理分野、計測各種分野、音声解析分野をはじめとする種々の分野において、フーリエ変換を用いて、各種信号を生成する際に利用可能である。
201 干渉信号の取得
202 干渉信号をフーリエ変換して情報信号生成
203 情報信号のn番目のピーク強度検出
205 ピーク強度の最大値に対応する干渉信号の振幅、位相及び周波数を算出
206 振幅、位相及び周波数情報に基づき波数空間軸上に延長した信号生成
207 外挿信号の生成
208 合成信号の生成
209 合成信号をフーリエ変換
202 干渉信号をフーリエ変換して情報信号生成
203 情報信号のn番目のピーク強度検出
205 ピーク強度の最大値に対応する干渉信号の振幅、位相及び周波数を算出
206 振幅、位相及び周波数情報に基づき波数空間軸上に延長した信号生成
207 外挿信号の生成
208 合成信号の生成
209 合成信号をフーリエ変換
Claims (14)
- 信号f(a)をフーリエ変換することで情報信号F(b)を生成する情報信号生成方法において、
前記情報信号F(b)のピーク強度の最大値を検出する工程と、
前記情報信号F(b)のピーク強度の最大値に対応する前記信号f(a)の振幅と位相と周波数とを算出する工程と、
前記振幅、位相及び周波数情報に基づきa軸上に延長した信号f’(a)を生成する工程と、
a1 < a2として、前記延長した信号f’(a)からa1よりも小さい領域の信号と、前記延長した信号f’(a)のうちa2よりも大きい領域の信号とを抽出して外挿信号を生成する工程と、
前記外挿信号と前記a1からa2までの信号f(a)とを合成して合成信号を生成する工程と、
前記合成信号をさらにフーリエ変換する工程と、
を有することを特徴とする信号生成方法。 - 前記情報信号生成方法により得られた情報信号F(b)について、強度の第n番目(ここで、nは1以上の整数)のピーク値に次ぐ強度の第(n+1)番目のピーク値を検出する工程と、
前記第(n+1)番目のピーク値に対応する前記信号f(a)の振幅と位相と周波数とを算出する工程と、
前記振幅、位相及び周波数情報に基づきa軸上に延長した信号を生成する工程と、
a1 < a2として、前記延長した信号からa1よりも小さい領域の信号と、前記延長した信号のうちa2よりも大きい領域の信号とを抽出して外挿信号を生成する工程と、
前記外挿信号と前記a1からa2までの信号f(a)とを合成して合成信号を生成する工程と、
前記合成信号をさらにフーリエ変換する工程と、
をn=1からNまで繰り返す工程と、
を有することを特徴とする請求項1に記載の信号生成方法。 - 前記Nを前記情報信号F(b)のピークの数とすることを特徴とする請求項2に記載の信号生成方法。
- 前記Nは、大きい強度を示すピークから数えて所定の強度以下となるピークまでの数であることを特徴とする請求項2に記載の信号生成方法。
- 前記情報信号F(b)の強度に対して閾値を決める工程と、
前記情報信号F(b)で閾値を超える複数のピーク値を検出する工程と、
前記複数のピーク値のそれぞれに対応する前記信号f(a)の振幅と位相と周波数とを算出する工程と、
前記振幅、位相及び周波数情報に基づきa軸上に延長した信号を生成する工程と、
前記延長した信号からa1よりも小さい領域の信号と、前記延長した信号のうちa2よりも大きい領域の信号とを抽出して外挿信号を生成する工程と、
前記外挿信号と前記a1からa2までの信号f(a)とを合成して合成信号を生成する工程と、
前記合成信号をさらにフーリエ変換する工程と、
を有することを特徴とする請求項1に記載の信号生成方法。 - 前記閾値を所定の割合で下げる度に、
前記情報信号F(b)で閾値を超える複数のピーク値を検出する工程と、
前記複数のピーク値のそれぞれに対応する前記信号f(a)の振幅と位相と周波数とを算出する工程と、
前記振幅、位相及び周波数情報に基づきa軸上に延長した信号を生成する工程と、
前記延長した信号からa1よりも小さい領域の信号と、前記延長した信号のうちa2よりも大きい領域の信号とを抽出して外挿信号を生成する工程と、
前記外挿信号と前記a1からa2までの信号f(a)とを合成して合成信号を生成する工程と、
前記合成信号をさらにフーリエ変換する工程と、
を行うことを特徴とする請求項5に記載の信号生成方法。 - 前記閾値を決める工程は、f(a)のエンベロープ波形を取得し、該エンベロープ波形をフーリエ変換することで得られる信号の最大のピーク値とこれに次ぐ2番目のピーク値を検出する工程と、
前記2番目のピーク値と最大のピーク値の比を算出する工程とを有し、
前記比を前記情報信号F(b)の最大のピーク値に乗じて得られる値を、前記閾値とすることを特徴とする請求項5に記載の信号生成方法。 - 前記信号f(a)に、関数が示す波形の裾部が0となる窓関数を乗じ、前記延長した信号は前記窓関数が0となるaまで延長されたものであることを特徴とする請求項1に記載の信号生成方法。
- 前記信号f’(a)を延長する範囲は、前記延長した範囲の前記窓関数をフーリエ変換した信号の強度の最大のピーク値とこれに次ぐ2番目のピーク値との比が所定の値となる範囲以上であることを特徴とする請求項8に記載の信号生成方法。
- 前記信号f(a)は、該信号f(a)を該信号f(a)のエンベロープ情報で割ることで該信号f(a)の形状をフラットトップ形状に整形し、これに窓関数を乗じたものであることを特徴とする請求項1に記載の信号生成方法。
- 前記信号f(a)を光スペクトル干渉信号とし、前記情報信号F(b)を断層に関する情報信号として、
前記断層に関する情報信号のピーク強度の最大値を検出する工程と、
前記断層に関する情報信号のピーク強度の最大値に対応する前記光スペクトル干渉信号の振幅と位相と周波数とを算出する工程と、
前記振幅、位相及び周波数情報に基づき、前記a軸を波数軸として該波数軸上に延長した信号を生成する工程と、
波数ks < keとして、前記延長した信号からksよりも小さい領域の信号と、前記延長した信号のうちkeよりも大きい領域の信号とを抽出して外挿信号を生成する工程と、
前記外挿信号と前記ksからkeまでの光スペクトル干渉信号とを合成して合成信号を生成する工程と、
前記合成信号をさらにフーリエ変換する工程と、
を有することを特徴とする請求項1乃至10に記載の信号生成方法。 - 信号f(a)をフーリエ変換することで情報信号F(b)を生成する情報信号処理装置において、
前記情報信号F(b)のピーク強度の最大値を検出する処理部と、
前記情報信号F(b)のピーク強度の最大値に対応する前記信号f(a)の振幅と位相と周波数とを算出する処理部と、
前記振幅及び位相及び周波数情報に基づきa軸上に延長した信号を生成する処理部と、
a1 < a2として、前記延長した信号からa1よりも小さい領域の信号と、前記延長した信号のうちa2よりも大きい領域の信号とを抽出して外挿信号を生成する処理部と、
前記外挿信号と前記a1からa2までの信号f(a)とを合成して合成信号を生成する処理部と、
前記合成信号をさらにフーリエ変換する処理部と、
を有することを特徴とする信号生成装置。 - 前記信号f(a)を光スペクトル干渉信号、前記情報信号F(b)を断層に関する情報信号、
前記a軸を波数軸、とそれぞれすると共に、前記a1 < a2を、波数ks < keとした
ことを特徴とする請求項12に記載の信号生成装置。 - 光源部と、前記光源部からの光を検体に照射し、検体からの反射光を伝達させる検体測定部と、前記光源部からの光を参照光として伝達させる参照部と、前記反射光と前記参照光とを干渉させる干渉部と、前記干渉部からの干渉光を検出する光検出部と、前記光検出部で検出された光に基づいて前記検体の断層像を得る画像処理部と、を有する光干渉断層撮像装置であって、前記画像処理部は、請求項13に記載の断層に関する情報信号を生成する信号生成装置を備えることを特徴とする光干渉断層撮像装置。
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
PCT/JP2011/080223 WO2013098942A1 (ja) | 2011-12-27 | 2011-12-27 | 情報信号生成方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
JPWO2013098942A1 JPWO2013098942A1 (ja) | 2015-04-30 |
JP5885758B2 true JP5885758B2 (ja) | 2016-03-15 |
Family
ID=48655395
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
JP2013551079A Expired - Fee Related JP5885758B2 (ja) | 2011-12-27 | 2011-12-27 | 情報信号生成方法 |
Country Status (5)
Country | Link |
---|---|
US (1) | US9600444B2 (ja) |
EP (1) | EP2799838B1 (ja) |
JP (1) | JP5885758B2 (ja) |
CN (1) | CN104011497B (ja) |
WO (1) | WO2013098942A1 (ja) |
Cited By (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
KR20200105931A (ko) * | 2018-01-18 | 2020-09-09 | 로베르트 보쉬 게엠베하 | 레이더 신호의 보정 방법 및 그 장치, 그리고 레이더 장치 |
Families Citing this family (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JP6437364B2 (ja) * | 2015-03-30 | 2018-12-12 | 株式会社Screenホールディングス | 信号処理方法および画像処理装置 |
EP3109646A1 (de) * | 2015-06-23 | 2016-12-28 | Siemens Aktiengesellschaft | Verfahren zur analyse eines signals sowie vorrichtung zur durchführung des verfahrens |
US11357399B2 (en) * | 2018-08-02 | 2022-06-14 | Nidek Co., Ltd. | OCT apparatus |
CN114518162B (zh) * | 2022-01-24 | 2023-08-04 | 中国人民解放军海军工程大学 | 一种光纤水听器干涉信号强度补偿方法和系统 |
Family Cites Families (16)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JPS61111427A (ja) * | 1985-10-04 | 1986-05-29 | Hitachi Ltd | フ−リエ変換スペクトロメ−タの干渉信号処理方法 |
JPS63302307A (ja) * | 1987-06-02 | 1988-12-09 | Hitachi Ltd | 光学的膜厚測定方法 |
PL320611A1 (en) * | 1994-12-09 | 1997-10-13 | Foss Electric As | Method of receiving information |
US5686922A (en) | 1995-09-29 | 1997-11-11 | Environmental Research Institute Of Michigan | Super spatially variant apodization (Super - SVA) |
US5612784A (en) * | 1996-01-04 | 1997-03-18 | Bio-Rad Laboratories, Inc. | Digital signal processing technique for a FT-IR spectrometer using multiple modulations |
JP3928106B2 (ja) * | 1999-07-02 | 2007-06-13 | 株式会社日立メディコ | 磁気共鳴撮影装置 |
US6252668B1 (en) * | 1999-11-19 | 2001-06-26 | Zygo Corporation | Systems and methods for quantifying nonlinearities in interferometry systems |
WO2007133961A2 (en) * | 2006-05-10 | 2007-11-22 | The General Hospital Corporation | Processes, arrangements and systems for providing frequency domain imaging of a sample |
EP1870029A1 (en) * | 2006-06-23 | 2007-12-26 | OPTOPOL Technology Spolka z o.o. | Apparatus and method for frequency domain optical coherence tomography |
JP4816507B2 (ja) | 2007-02-28 | 2011-11-16 | カシオ計算機株式会社 | 音声分析合成装置、及びプログラム |
JP5448353B2 (ja) * | 2007-05-02 | 2014-03-19 | キヤノン株式会社 | 光干渉断層計を用いた画像形成方法、及び光干渉断層装置 |
JP5194963B2 (ja) * | 2008-04-03 | 2013-05-08 | 株式会社ニコン | 波形解析装置、波形解析プログラム、干渉計装置、パターン投影形状測定装置、及び波形解析方法 |
JP2010223670A (ja) * | 2009-03-23 | 2010-10-07 | Toyama Univ | 光断層画像表示システム |
GB2472059B (en) * | 2009-07-23 | 2012-09-19 | Univ Loughborough | Apparatus for the absolute measurement of two dimensional optical path distributions using interferometry |
MX338568B (es) * | 2009-08-14 | 2016-04-21 | Pyrotek Inc | Sistema de calor residual. |
TWI417519B (zh) * | 2009-12-10 | 2013-12-01 | Ind Tech Res Inst | 干涉相位差量測方法及其系統 |
-
2011
- 2011-12-27 JP JP2013551079A patent/JP5885758B2/ja not_active Expired - Fee Related
- 2011-12-27 EP EP11878685.4A patent/EP2799838B1/en not_active Not-in-force
- 2011-12-27 WO PCT/JP2011/080223 patent/WO2013098942A1/ja active Application Filing
- 2011-12-27 CN CN201180076033.9A patent/CN104011497B/zh not_active Expired - Fee Related
-
2012
- 2012-12-20 US US13/721,403 patent/US9600444B2/en active Active
Cited By (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
KR20200105931A (ko) * | 2018-01-18 | 2020-09-09 | 로베르트 보쉬 게엠베하 | 레이더 신호의 보정 방법 및 그 장치, 그리고 레이더 장치 |
KR102680054B1 (ko) | 2018-01-18 | 2024-07-02 | 로베르트 보쉬 게엠베하 | 레이더 신호의 보정 방법 및 그 장치, 그리고 레이더 장치 |
Also Published As
Publication number | Publication date |
---|---|
EP2799838A1 (en) | 2014-11-05 |
EP2799838B1 (en) | 2016-07-27 |
JPWO2013098942A1 (ja) | 2015-04-30 |
US9600444B2 (en) | 2017-03-21 |
CN104011497B (zh) | 2017-10-03 |
EP2799838A4 (en) | 2015-09-02 |
WO2013098942A1 (ja) | 2013-07-04 |
US20130166239A1 (en) | 2013-06-27 |
CN104011497A (zh) | 2014-08-27 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
JP5680826B2 (ja) | 1以上のスペクトルを符号化する内視鏡技術によるデータ生成システム | |
JP5371315B2 (ja) | 光干渉断層撮像方法および光干渉断層撮像装置 | |
JP4389032B2 (ja) | 光コヒーレンストモグラフィーの画像処理装置 | |
US8457440B1 (en) | Method and system for background subtraction in medical optical coherence tomography system | |
US8504141B2 (en) | Optical tomographic image generating apparatus and optical tomographic image generating method | |
JP5885758B2 (ja) | 情報信号生成方法 | |
JP2006047264A (ja) | オプティカル・コヒーレント・トモグラフィー装置及びこれに用いる可変波長光発生装置並びに可変波長発光光源 | |
JP5557397B2 (ja) | 半透明物質の画像化の方法および装置 | |
CN105865613A (zh) | 海洋立体监测中的水下光学检测与成像传感方法及系统 | |
US9964397B2 (en) | Multiple reference OCT system | |
JP2008151697A (ja) | 光断層画像化装置 | |
US11118894B2 (en) | Systems, methods, and media for multiple reference arm spectral domain optical coherence tomography | |
CN107407601A (zh) | 用于补偿由光谱仪系统产生的干涉图的时间周期扰动的光谱仪系统和方法 | |
JP6214020B2 (ja) | 光断層イメージング法、その装置およびプログラム | |
JP2010223670A (ja) | 光断層画像表示システム | |
JP2014016318A (ja) | 光断層画像取得方法 | |
JP5664564B2 (ja) | 光断層画像取得方法 | |
Rao et al. | Signal processing of spectral domain optical coherence tomography | |
JP5746741B2 (ja) | 画像生成装置、画像生成システム及び画像生成方法 | |
JP2014092438A (ja) | 光学的測定方法 | |
JP2013057549A (ja) | 光断層画像取得方法 |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
TRDD | Decision of grant or rejection written | ||
A01 | Written decision to grant a patent or to grant a registration (utility model) |
Free format text: JAPANESE INTERMEDIATE CODE: A01 Effective date: 20160112 |
|
A61 | First payment of annual fees (during grant procedure) |
Free format text: JAPANESE INTERMEDIATE CODE: A61 Effective date: 20160209 |
|
R151 | Written notification of patent or utility model registration |
Ref document number: 5885758 Country of ref document: JP Free format text: JAPANESE INTERMEDIATE CODE: R151 |
|
LAPS | Cancellation because of no payment of annual fees |