JP7298706B2 - 光パルス試験方法及び光パルス試験装置 - Google Patents

Description

本開示は、光ファイバ振動センシングにおける周波数多重によるサンプリングレート向上とフェーディング対策が可能な光パルス試験方法及び光パルス試験装置に関する。

光ファイバに加わった物理的な振動を、光ファイバ長手方向に分布的に計測する手段として、被測定光ファイバにパルス試験光を入射し、レイリー散乱による後方散乱光を検出するＤＡＳ（Ｄｉｓｔｒｉｂｕｔｅｄ Ａｃｏｕｓｔｉｃ Ｓｅｎｓｉｎｇ）と呼ばれる手法が知られている（例えば、非特許文献１を参照。）。

ＤＡＳでは、光ファイバに加わった物理的な振動による光ファイバの光路長変化を捉え、振動のセンシングを行う。振動を検出することで、被測定光ファイバ周辺での、物体の動き等を検出することが可能である。

ＤＡＳにおける後方散乱光の検出方法として、被測定光ファイバの各地点からの散乱光強度を測定し、散乱光強度の時間変化を観測する手法があり、ＤＡＳ－Ｉ（ＤＡＳ－Iｎｔｅｎｓｉｔｙ）と呼ばれている。ＤＡＳ－Ｉは装置構成が簡便にできる特徴があるが、散乱光強度から振動によるファイバの光路長変化を定量的に計算することができないため、定性的な測定手法である（例えば、非特許文献２を参照。）。

一方で、被測定光ファイバの各地点からの散乱光の位相を測定し、位相の時間変化を観測する手法であるＤＡＳ－Ｐ（ＤＡＳ－Pｈａｓｅ）も研究開発されている。ＤＡＳ－Ｐは、装置構成や信号処理がＤＡＳ－Ｉより複雑となるが、振動によるファイバの光路長変化に対して位相が線形に変化し、その変化率も光ファイバ長手方向で同一となるため、振動の定量的な測定が可能となり、被測定光ファイバに加わった振動を忠実に再現することができる（例えば、非特許文献２を参照。）。

ＤＡＳ－Ｐによる測定では、パルス光を被測定光ファイバに入射し、パルス光を入射した時刻ｔでの、散乱された光の位相を、光ファイバの長手方向に分布的に計測する。つまり、光ファイバの入射端からの距離ｌとして、散乱光の位相θ（ｌ、ｔ）を測定する。パルス光を、時間間隔Ｔで、繰り返し被測定光ファイバに入射することで、整数ｎとして時刻ｔ＝ｎＴにおける散乱された光の位相の時間変化θ（ｌ、ｔ）を、被測定光ファイバの長手方向の各点について測定する。ただし実際は、入射端から距離ｔまでパルス光が伝搬する時間だけ、距離ｌの地点を測定する時刻は、パルスを入射した時刻より遅れる。さらに、散乱光が入射端まで戻ってくるのに要する時間だけ、測定器で測定する時刻は遅れることに注意する。距離ｌから距離ｌ＋δｌまでの区間に加わった物理的な振動の各時刻ｎＴでの大きさは、距離ｌ＋δｌでの位相θ（ｌ＋δｌ，ｎＴ）と、距離ｌでの位相θ（ｌ，ｎＴ）との差分δθ（ｌ，ｎＴ）に比例することが知られている。つまり、時刻ゼロを基準とすれば、下式を満たす。

散乱光の位相を検出するための装置構成としては、被測定光ファイバからの後方散乱光を直接フォトダイオードなどで検波する直接検波の構成や、別途用意した参照光と合波させて検出するコヒーレント検波を使用した構成がある（例えば、非特許文献１を参照。）。

コヒーレント検波を行い、位相を計算する機構では、ヒルベルト変換を用いてソフトウェアベースで処理する機構と、９０度光ハイブリッドを用いてハードウェアベースで処理する機構の二つに細分されるが、どちらの手法においても、散乱光の同相成分Ｉ（ｌ，ｎＴ）と直交成分Ｑ（ｌ，ｎＴ）を取得し、下式により位相を計算する。

ただし、４象限逆正接演算子Ａｒｃｔａｎによる出力値はラジアン単位で（－π，π］の範囲にあり、ｍを任意の整数として、２ｍπ＋θ（ｌ，ｎＴ）はｘｙ平面上で全て同じベクトル方向となるため、２ｍπだけの不確定性が上記で計算したθ_ｃａｌ（ｌ，ｎＴ）には存在する。

したがって、θ（ｌ，ｎＴ）のより正確な評価方法として、位相アンラップ等の信号処理がさらに行われる。一般的な位相アンラップでは、アンラップ後の位相を

とすれば、任意の整数をｐとして、

がπラジアンより大きくなる場合に、

がπラジアン以下になるような適切な整数ｑを選択して、アンラップ後の位相を

として、次式を計算する。

と計算する。なお、上添え字ｕｎｗｒａｐはアンラップ後の位相であることを表す。

ＤＡＳによる測定においては、光を検出するためのＰＤの熱雑音や、その後の電気段での雑音、光によるショット雑音などの、測定器の雑音が存在する。したがって、測定する散乱光の強度や位相にも、測定器の雑音による影響が現れる。

特に、散乱光の位相を測定する場合、測定器の雑音の影響が大きくなってしまうと、単に位相の不確かさが増加するだけでなく、雑音がない場合の理想的な位相値と比較して、大きく異なる測定値をとる確率が大きくなる。

例えば、コヒーレント検波の場合に、同相成分を横軸に、直交成分を縦軸にした時の、測定された散乱光のベクトルについて、雑音がない時のベクトルの向きが測定したい位相に対応するが、雑音の影響が大きいと、ベクトルの向きが反対の方向を向き、雑音がない場合の理想的な位相値と比較して、実際に測定される位相値がπラジアン程度異なる値をとる確率が大きくなる。このような点においては、式（１）から振動の大きさを計算する際に、大きな物理的な力が光ファイバに加わったとする誤認識につながる。また、雑音の影響が大きくなると、式（３）で示したアンラップ処理において、整数ｑの選択を誤る点が増加し、選択を誤った点の前後で２π以上の実際には存在しない位相値の違いが生じてしまう。このような位相値の違いも、式（１）から振動の大きさを計算する際に、大きな物理的な力が光ファイバに加わったとする誤認識につながる。

正確に位相を測定するためには、測定器の雑音の影響を低減する必要がある。測定器の雑音の影響が大きくなるのは、測定器の雑音が各地点・各時刻について同程度あるとみなせる際には、散乱光の強度そのものが小さくなる場合であるから、散乱光の強度を各地点及び各時刻で大きくすることができれば、測定器の雑音の影響を低減することが可能となる。

散乱光の強度そのものが小さくなる原因となっているのは、プローブとなるパルス光が被測定光ファイバを伝搬するのに従って発生する吸収や散乱による損失だけではない。有限な時間幅を持ったパルス光を被測定光ファイバに入射して、パルス光の散乱を検出しているため、被測定光ファイバ上の非常に細かく分布している多数の散乱体からの散乱光の干渉が起きる。干渉の結果として、各時刻における散乱体の被測定光ファイバの長手方向での分布に応じて、散乱光の強度が小さくなる地点が発生する。この現象はフェーディングと呼ばれる（例えば、非特許文献３を参照。）。

したがって、ＤＡＳ－Ｐにおける散乱光の位相を測定する場合、測定器の雑音の影響を低減するために、フェーディングによって、各時刻で散乱光の強度が小さくなる地点が発生することの防止が必要という課題がある。

当該課題を解決する手段として、単純に入射する光パルスのピーク強度を大きくする方法がある。しかし、ピーク強度を大きくすると、非線形効果が発生し、パルス光の特性が被測定光ファイバの伝搬に伴い変化する。このため、入射可能な光パルスのピーク強度は制限され、上記課題を十分に解決できない場合がある。

上記課題を解決するために、ＤＡＳ－Ｐにおける散乱光の位相を測定するときに、入射する光パルスのピーク強度を大きくせずに測定器の雑音の影響を低減できる位相測定方法及び信号処理装置が提案されている（例えば、非特許文献４を参照。）。

非特許文献４では、上記課題を解決するために、振動によるファイバ状態の変化が無視できる時間間隔で、異なる光周波数成分のパルスを並べて波長多重したパルス光を被測定光ファイバに入射し、被測定光ファイバからの各波長における散乱光を、同相成分を横軸に直交成分を縦軸にした２次元平面上にプロットして得られる散乱光ベクトルを作成し、作成した散乱光ベクトルを被測定光ファイバ上の各地点で波長ごとに回転させることで向きを一致させ、向きを一致させたベクトル同士を加算平均することで新たなベクトルを生成し、生成した新たなベクトルの同相成分と直交成分の値を用いて位相を計算している。

ＤＡＳ－Ｐにおける測定では、測定距離と測定可能な振動周波数の上限の間にトレードオフが生じる課題も存在する。単一周波数の光パルスを用いる場合、測定距離が長くなると、遠端からの散乱光が戻ってくる時刻が、パルス入射時刻に対して遅れる。したがって、遠端からの散乱光と、次の光パルスを入射した際の入射端付近からの散乱光が合波しないために、光パルスを入射する繰り返し周波数には上限が生じる。したがって、サンプリング定理から、繰り返し周波数の１／２倍のナイキスト周波数より大きい振動周波数の振動については、エイリアシングのため、正しく測定することができないという課題がある。

上記課題の解決方法として非特許文献５が提案されている。非特許文献５では、上記課題を解決するために、異なる光周波数成分のパルスを時間的に等間隔で並べて波長多重したパルス光を被測定光ファイバに入射し、被測定光ファイバからの各波長における散乱光を、同相成分を横軸に直交成分を縦軸にした２次元平面上にプロットして得られる散乱光ベクトルを作成する。異なる光周波数に対応する散乱光ベクトルの角度を、異なる時刻における光位相として連結させることにより、単一光周波数の場合に、測定距離から決まる測定可能な振動周波数ｆｖに対して、Ｎ波多重により、測定可能な振動周波数をＮ×ｆｖとすることができる。

しかし、測定距離と測定可能な振動周波数の上限とのトレードオフには、位相アンラップが正しく行われる必要があることから、さらに厳しい条件が加わる。隣り合う光パルスでサンプリングした際の位相変化の大きさの絶対値がπより大きく変化する場合には、位相アンラップを一意に行うことができなくなるため、位相アンラップの失敗につながってしまう（例えば、非特許文献６を参照。）。

したがって、隣り合うサンプリング点での位相変化の大きさの絶対値の上限はπという制約が生じる。つまり、ナイキスト周波数以下の範囲内であっても、振動周波数が高くなるほど、振動振幅が同じでも、隣り合うサンプリング点での位相変化量は大きくなるため、測定可能な振動周波数の上限にさらなる条件が生じる。

尚、位相変化の大きさと、振動によってファイバに加わった歪量との関係は、例えば非特許文献７で説明されている。

非特許文献７によれば、全長ｌのファイバが歪量εによってΔｌだけ伸びた時、Δｌだけ伸びた分による光が通過する際の位相変化の増加量Δφは下式となる。

ここで、ｋ＝２πｎ／λは伝搬定数、ｎはファイバの実効屈折率、μ_ｐはポアソン比、ｐ_１１とｐ_１２はストレイン－オプティックテンソル成分である。歪量εは、ファイバの元の長さに対する変化量の割合であるΔｌ／ｌとして定義されることに注意する。例えば、通常の通信波長帯付近のλ＝１５５５ｎｍの場合を考えると、ｎ＝１．４７、μ_ｐ＝０．１７、ｐ_１１＝０．１２１、ｐ_１２＝０．２７１の値となるため、

となることが知られている（例えば、特許文献８を参照。）。ただし、Ｋ＝４．６×１０^６ｍ^－１である。この関係式を使用すれば、位相変化の大きさの条件を歪量の条件に置き替えることが可能である。

Ａ． Ｍａｓｏｕｄｉ， Ｔ． Ｐ． Ｎｅｗｓｏｎ， "Ｃｏｎｔｒｉｂｕｔｅｄ Ｒｖｉｅｗ： Ｄｉｓｔｒｉｂｕｔｅｄ ｏｐｔｉｃａｌ ｆｉｂｒｅ ｄｙｎａｍｉｃ ｓｔｒａｉｎ ｓｅｎｓｉｎｇ，" Ｒｅｖｉｅｗ ｏｆ Ｓｃｉｅｎｔｉｆｉｃ Ｉｎｓｔｒｕｍｅｎｔｓ， ｖｏｌ． ８７， ｐ． ０１１５０１ （２０１６） 西口憲一、李哲賢、グジクアーター、横山光徳、増田欣増、「光ファイバによる分布型音波センサの試作とその信号処理」、信学技報、１１５ （２０２）， ｐｐ． ２９－３４ （２０１５） Ｇ． Ｙａｎｇ ｅｔ ａｌ．， "Ｌｏｎｇ－Ｒａｎｇｅ Ｄｉｓｔｒｉｂｕｔｅｄ Ｖｉｂｒａｔｉｏｎ Ｓｅｎｓｉｎｇ Ｂａｓｅｄ ｏｎ Ｐｈａｓｅ Ｅｘｔｒａｃｔｉｏｎ ｆｒｏｍ Ｐｈａｓｅ－Ｓｅｎｓｉｔｉｖｅ ＯＴＤＲ，" ＩＥＥＥ Ｐｈｏｔｏｎｉｃｓ Ｊｏｕｒｎａｌ， ｖｏｌ． ８， ｎｏ． ３， ２０１６． 脇坂佳史、飯田大輔、岡本圭司、押田博之、「周波数多重位相ＯＴＤＲを用いた分布振動計測方法」、電子情報通信学会２０１９年ソサイエティ大会、２０１９年８月２７日 Ｄ． Ｉｉｄａ， Ｋ． Ｔｏｇｅ， Ｔ． Ｍａｎａｂｅ， ‘Ｄｉｓｔｒｉｂｕｔｅｄ ｍｅａｓｕｒｅｍｅｎｔ ｏｆ ａｃｏｕｓｔｉｃ ｖｉｂｒａｔｉｏｎ ｌｏｃａｔｉｏｎ ｗｉｔｈ ｆｒｅｑｕｅｎｃｙ ｍｕｌｔｉｐｌｅｘｅｄ ｐｈａｓｅ－ＯＴＤＲ’， Ｏｐｔ． Ｆｉｂｅｒ Ｔｅｃｈｎｏｌ．， ３６ （２０１７） ｐｐ． １９－２５， ＤＯＩ： １０．１０１６／ｊ．ｙｏｆｔｅ．２０１７．０２．００５ Ｍａｒｉａ Ｒｏｓａｒｉｏ Ｆｅｒｎａｎｄｅｚ－Ｒｕｉｚ， Ｈｕｇｏ Ｆ． Ｍａｒｔｉｎｓ， "Ｓｔｅａｄｙ－Ｓｅｎｓｉｔｉｖｉｔｙ Ｄｉｓｔｒｉｂｕｔｅｄ Ａｃｏｕｓｔｉｃ Ｓｅｎｓｏｒｓ，" Ｊ． Ｌｉｇｈｔｗａｖｅ Ｔｅｃｈｎｏｌ． ３６， ５６９０－５６９６ （２０１８） Ｃ． Ｄ． Ｂｕｔｔｅｒ ａｎｄ Ｇ． Ｂ． Ｈｏｃｋｅｒ， " Ｆｉｂｅｒ ｏｐｔｉｃｓ ｓｔｒａｉｎ ｇａｕｇｅ，" Ａｐｐｌ． Ｏｐｔ． １７， ２８６７－２８６９ （１９７８） Ａ． Ｅ． Ａｌｅｋｓｅｅｖ ｅｔ ａｌ．， Ｌａｓｅｒ Ｐｈｙｓ．， ２９ （２０１９） ０５５１０６

（課題１）上記の、測定距離と測定可能な振動周波数の上限とのトレードオフを考慮した上で、目的とする振動周波数を測定するために、非特許文献５に記載のような周波数多重の方法を実施した際に必要となる、周波数多重数についての設計方法は明らかになっていない。光周波数の多重数を増加させるには、送信系と受信系の帯域を増加させる必要があり、コストを要するため、必要最小限の周波数多重数を設計する方法が必要となる。

（課題２）非特許文献５に記載のような周波数多重の方法を行う際に、各光周波数間の角度差を補正しないと、測定される位相変化が、実際の位相変化に対して、歪んでしまうという課題が生じる。この課題に対する解決方法は今までに提案されていない。

（課題３）測定距離と測定可能な振動周波数の上限とのトレードオフを解決するための非特許文献５に記載のような周波数多重の方法と、フェーディング対策のための非特許文献４に記載した周波数多重の方法を、同時に実施する際の、光周波数パルスの構成方法と受信信号処理方法については今までに提案されていない。

そこで、本発明は、上記課題を解決するために、次の３つを目的とする。
（目的１）上記課題１を解決すること。すなわち、上記の、測定距離と測定可能な振動周波数の上限とのトレードオフを考慮した上で、目的とする振動周波数を測定するために、非特許文献５に記載のような周波数多重の方法を実施した際に必要となる、周波数多重数の最小値の設計方法を提案すること。

（目的２）上記課題２を解決すること。すなわち、非特許文献５に記載のような周波数多重の方法を行う際に、各光周波数間の角度差を補正するための、光周波数パルスの構成方法と受信信号処理方法について提案すること。

（目的３）上記課題３を解決すること。すなわち、測定距離と測定可能な振動周波数の上限とのトレードオフを解決するための非特許文献５に記載のような周波数多重の方法と、フェーディング対策のための非特許文献４に記載した周波数多重の方法を、同時に実施する際の、光周波数パルスの構成方法と受信信号処理方法について提案すること。

上記目的１を達成するために、本発明に係る光パルス試験方法は、
異なる光周波数の光パルスを時間的に等間隔で並べた波長多重数Ｎの光パルス列を被測定光ファイバの一端に入射すること、
前記被測定光ファイバの前記一端に戻ってきた各波長の散乱光を受光すること、
前記被測定光ファイバの振動を前記散乱光の位相成分の時間変化として観測すること、及び
前記波長多重数Ｎを数Ｃ１を満たす最小値に決定すること、
を特徴とする。

なお、Ｆ_ν（ｔ）が単一の振動周波数Ａ×ｓｉｎ（２πｆ_ｖｔ）である場合、波長多重数Ｎを数Ｃ１ａを満たす最小値に決定する。

ただし、ｆ_ｖは振動周波数である。

また、本発明に係る光パルス試験装置は、
異なる光周波数の光パルスを時間的に等間隔で並べた波長多重数Ｎの光パルス列を被測定光ファイバの一端に入射する光源と、
前記被測定光ファイバの前記一端に戻ってきた各波長の散乱光を受光する受光器と、
前記被測定光ファイバの振動を前記散乱光の位相成分の時間変化として観測する信号処理部と、
前記波長多重数Ｎを数Ｃ１を満たす最小値に決定する計算部と、
を備える。

上記目的２を達成するために、本発明に係る光パルス試験方法は、
前記光パルス列を前記被測定光ファイバの一端に入射する前に、前記光パルス列の光パルス間隔より短い間隔で異なる光周波数の前記光パルスを並べた光パルス対を並べた補正用光パルス列を前記被測定光ファイバの一端に入射すること、及び
前記散乱光の位相成分の補正値を数Ｃ２で計算すること
を特徴とする。

ただし、φ（ｉ（ｍ’），ｉ（ｍ））は光周波数が異なる２つの前記光パルス間の角度差、
Ｐは前記補正用光パルス列に含まれる前記光パルス対の数、
ｐは前記光パルス対の番号（１からＰまでの整数）、
ｔ_ｐはｐ番目の前記光パルス対が前記被測定光ファイバの一端に入射される時間、
ｚは前記被測定光ファイバの前記一端からの距離、
ｒ_ｉ（ｍ）（ｚ，ｔ_ｐ＋２ｚ／ν）及びｒ_{ｉ（ｍ’）}（ｚ，ｔ_ｐ＋２ｚ／ν）は距離zからの散乱光の複素ベクトル、
ａｒｇ関数は引数の複素ベクトルの偏角を－πからπの範囲で計算し、実数を出力する関数、
Ｒ（＊）は複素平面上で角度＊だけ時計まわりに複素ベクトルｒを回転させる演算子である。

また、本発明に係る光パルス試験装置の前記光源は、前記光パルス列を前記被測定光ファイバの一端に入射する前に、前記光パルス列の光パルス間隔より短い間隔で異なる光周波数の前記光パルスを並べた光パルス対を並べた補正用光パルス列を前記被測定光ファイバの一端に入射し、前記信号処理部は、前記散乱光の位相成分の補正値を数Ｃ２で計算することを特徴とする。

上記目的３を達成するために、本発明に係る光パルス試験方法は、光パルス列を形成する光パルスを、前記光パルス列の光パルス間隔より短い間隔で異なる光周波数の微小光パルスをＭ個並べて形成すること、及び雑音レベルの１／√Ｍ倍が所定値以下となるように前記微小光パルスの数Ｍを設定することを特徴とする。

上記目的３を達成するために、次のような方法もある。
（方法Ａ）
本発明に係る光パルス試験方法は、
Ｎ×Ｍ＋１個（ＮとＭは自然数）の異なる光周波数の微小光パルスを並べた集団をＭ個形成すること、
Ｍ個の前記集団を並べて、先頭の前記集団からＭ個づつ前記微小光パルスを切り出して光パルス対をＮ×Ｍ＋１個形成すること、
Ｎ×Ｍ＋１個の前記光パルス対を時間的に等間隔で並べた光パルス列を被測定光ファイバの一端に入射すること、
前記被測定光ファイバの前記一端に戻ってきた各波長の散乱光を受光すること、
前記被測定光ファイバの振動を前記散乱光の位相成分の時間変化として観測すること、
波長多重数Ｎを数Ｃ１を満たす最小値に決定すること、及び
雑音レベルの１／√Ｍ倍が所定値以下となるように波長多重数Ｍを設定すること
を特徴とする。

（方法Ｂ）
本発明に係る光パルス試験方法は、
Ｎ×Ｍ個（ＮとＭは自然数）の異なる光周波数の微小光パルスを並べた集団をＮ＋１個形成すること、
Ｎ＋１個の前記集団を並べて、先頭の前記集団からＭ個づつ前記微小光パルスを切り出して光パルス対をＮ（Ｎ＋１）個形成すること、
Ｎ（Ｎ＋１）個の前記光パルス対のうち、Ｎ＋１毎の前記光パルス対にいずれの前記微小光パルスの光周波数とも異なる光周波数の追加微小光パルスを追加すること、
前記追加微小光パルスが追加された前記光パルス対を含むＮ（Ｎ＋１）個の前記光パルス対を時間的に等間隔で並べた光パルス列を被測定光ファイバの一端に入射すること、
前記被測定光ファイバの前記一端に戻ってきた各波長の散乱光を受光すること、
前記被測定光ファイバの振動を前記散乱光の位相成分の時間変化として観測すること、
波長多重数Ｎを数Ｃ１を満たす最小値に決定すること、及び
雑音レベルの１／√Ｍ倍が所定値以下となるように波長多重数Ｍを設定すること
を特徴とする。

（方法Ｃ）
本発明に係る光パルス試験方法は、
Ｎ×Ｍ＋１個（ＮとＭは自然数）の異なる光周波数の微小光パルスを並べた集団を（Ｎ＋１）Ｍ＋１個形成すること、
（Ｎ＋１）Ｍ＋１個の前記集団を並べ、１＋ｋ（Ｎ＋１）番目のパルス対にＭ＋１個の前記微小光パルスが含まれ、他の前記パルス対にＭ個の前記微小光パルスが含まれるように、先頭の前記集団から前記微小光パルスを切り出して前記光パルス対を（Ｎ×Ｍ＋１）（Ｎ＋１）個形成すること、
（Ｎ×Ｍ＋１）（Ｎ＋１）個の前記光パルス対を時間的に等間隔で並べた光パルス列を被測定光ファイバの一端に入射すること、
前記被測定光ファイバの前記一端に戻ってきた各波長の散乱光を受光すること、
前記被測定光ファイバの振動を前記散乱光の位相成分の時間変化として観測すること、
波長多重数Ｎを数Ｃ１を満たす最小値に決定すること、及び
雑音レベルの１／√Ｍ倍が所定値以下となるように波長多重数Ｍを設定すること
を特徴とする。

なお、上記各発明は、可能な限り組み合わせることができる。

本発明は、次の光パルス試験方法及び光パルス試験装置を提供することができる。
（１）測定距離と測定可能な振動周波数の上限とのトレードオフを考慮した上で、目的とする振動周波数を測定するために、非特許文献５に記載のような周波数多重の方法を実施した際に必要となる、周波数多重数の最小値を設計できる。
（２）非特許文献５に記載のような周波数多重の方法を行う際に、各光周波数間の角度差を補正することができる。
（３）測定距離と測定可能な振動周波数の上限とのトレードオフを解決するための非特許文献５に記載のような周波数多重の方法と、フェーディング対策のための非特許文献４に記載した周波数多重の方法を、同時に実施することができる。

本発明に係る光パルス試験装置を説明する図である。 本発明に係る光パルス試験装置から被測定光ファイバに入射する光パルス列を説明する図である。 本発明に係る光パルス試験装置で測定した被測定光ファイバの位置ｚにおける位相変化を説明する図である。（ａ）は実測値である。（ｂ）は異なる光周波数の角度差を補正した結果である。プロットの記号（○、×、△）は、光パルスの光周波数に対応している。つまり、同一の記号のプロットは光周波数が同じ光パルスで測定された位相である。 本発明に係る光パルス試験装置から被測定光ファイバに入射する光パルス列を説明する図である。 異なる光周波数の角度差を説明する図である。 本発明に係る光パルス試験装置から被測定光ファイバに入射する光パルス列を説明する図である。 本発明に係る光パルス試験装置から被測定光ファイバに入射する光パルス列を説明する図である。 本発明に係る光パルス試験装置から被測定光ファイバに入射する光パルス列を説明する図である。 本発明に係る光パルス試験装置から被測定光ファイバに入射する光パルス列を説明する図である。 本発明に係る光パルス試験装置から被測定光ファイバに入射する光パルス列を説明する図である。 本発明に関連する信号処理装置を備える振動検出装置を説明する図である。 本発明に関連する信号処理装置の構造を説明する図である。 本発明に関連する位相測定方法を説明する図である。 雑音によるベクトルの指す位置の不確かさを説明する図である。 雑音の有無によるベクトルの違いを説明する図である。 本発明に関連する位相測定方法の原理を説明する図である。 本発明に関連する位相測定方法の原理を説明する図である。 本発明に関連する位相測定方法の効果を説明する図である。 本発明に関連する信号処理装置の構造を説明する図である。 本発明に関連する位相測定方法を説明する図である。 本発明に係る光パルス試験方法を説明するフローチャートである。 本発明に係る光パルス試験方法を説明するフローチャートである。

添付の図面を参照して本発明の実施形態を説明する。以下に説明する実施形態は本発明の実施形態であり、本発明は、以下の実施形態に制限されるものではない。なお、本明細書及び図面において符号が同じ構成要素は、相互に同一のものを示すものとする。

（実施形態１）
図１は、本実施形態のＤＡＳ－Ｐで振動検出を行う振動検出装置を説明する図である。
本振動検出装置は、光パルス試験装置であって、
異なる光周波数の光パルスを時間的に等間隔で並べた波長多重数Ｎの光パルス列を被測定光ファイバの一端に入射する光源と、
前記被測定光ファイバの前記一端に戻ってきた各波長の散乱光を受光する受光器と、
前記被測定光ファイバの振動を前記散乱光の位相成分の時間変化として観測する信号処理部と、
前記波長多重数Ｎを数Ｃ１を満たす最小値に決定する計算部１８と、
を備える。

ＣＷ光源１、カプラ２、及び光変調器３が前記光源に相当する。９０度光ハイブリッド７及びバランス検出器（１３、１４）が前記受光器に相当する。前記受光器は、９０度光ハイブリッド７を用いてコヒーレント検波を行う。信号処理装置１７が前記信号処理部に相当する。

図２１は、本振動検出装置が行う光パルス試験方法を説明するフローチャートである。本光パルス試験方法は、
波長多重数Ｎを数Ｃ１を満たす最小値に決定すること（ステップＳ２１）、
異なる光周波数の光パルスを時間的に等間隔で並べた波長多重数Ｎの光パルス列を被測定光ファイバ６の一端に入射すること（ステップＳ２２）、
被測定光ファイバ６の前記一端に戻ってきた各波長の散乱光を受光すること（ステップＳ２３）、及び
被測定光ファイバ６の振動を前記散乱光の位相成分の時間変化として観測すること（ステップＳ２４）
を行う。

測定器３１は、次のように被測定光ファイバ６からの散乱光を測定する。ＣＷ光源１から周波数がｆ_０の単一波長の連続光が射出され、カプラ２により参照光とプローブ光に分岐される。プローブ光は、光変調器３によって、光パルス４のように波長多重の光パルスに整形される。光パルス４は、周波数がｆ_０＋ｆ_ｉ（ｉは整数）かつパルス幅が光ファイバ長手方向での測定の空間分解能に対応する値Ｗに設定された微小パルス４ａが、ｉ＝１、２、・・・、Ｎ（Ｎは整数）だけ並んだ構成である。隣り合う光周波数のパルスは等間隔Ｔ_Ｎで並べられ、ｉ＝Ｎのパルスの次には、ｉ＝１のパルスを同じ間隔Ｔ_Ｎで並べる。光周波数ｆ_ｉの周波数設定は、各光周波数の散乱光の帯域が重ならないように設計する。目安として、どの光周波数の差も２／Ｗ以上となるようにする。

光変調器３の種類は光パルス４を生成できるならば具体的な指定はなく、数が複数の場合もある。例えば、ＳＳＢ変調器や周波数可変なＡＯ変調器などを用いても良いし、パルス化における消光比を大きくするためにさらにＳＯＡなどによる強度変調を行っても良い。

光パルス４は、サーキュレータ５を介して、被測定光ファイバ６に入射される。光ファイバ６の長手方向の各点で散乱された光が、後方散乱光としてサーキュレータ５に戻り、９０度光ハイブリッド７の一方の入力部に入射される。カプラ２により分岐された参照光は、９０度光ハイブリッド７のもう一方の入力部に入射される。

９０度光ハイブリッド７の内部構成は、９０度光ハイブリッドの機能さえ備えていれば、なんでもよい。構成例を図１に示す。後方散乱光は、５０：５０の分岐比のカプラ８に入射され、２分岐された散乱光が、５０：５０の分岐比のカプラ１２と、５０：５０のカプラ１１の入力部に入射される。参照光は、５０：５０の分岐比のカプラ９に入射され、２分岐された参照光の一方が、カプラ１１の入力部に入射され、他方が、位相シフタ１０で位相をπ／２だけシフトされてカプラ１２の入力部に入射される。

カプラ１１の２つの出力がバランス検出器１３によって検出され、アナログの同相成分Ｉ^{ａｎａｌｏｇ}である電気信号１５が出力される。カプラ１２の２つの出力がバランス検出器１４によって検出され、アナログの直交成分Ｑ^{ａｎａｌｏｇ}である電気信号１６が出力される。

電気信号１５と電気信号１６は、信号の周波数帯域をエイリアシングなくサンプリングが可能なＡＤ変換機能素子１７ａとＡＤ変換機能素子１７ｂを備えた信号処理装置１７に送られる。信号処置装置１７では、ＡＤ変換機能素子１７ａとＡＤ変換機能素子１７ｂから出力されたデジタル化された同相成分Ｉ^{ｄｉｇｉｔａｌ}と直交成分Ｑ^{ｄｉｇｉｔａｌ}の信号に対して、信号処理部１７ｃによって光パルス４を構成する各周波数ｆ_０＋ｆ_ｉ（ｉ＝１、２、・・・、Ｎ）のパルスによる散乱光による信号を分離する。

つまり、信号処理部１７ｃは、各周波数ｆ_０＋ｆ_ｉ成分のパルスを単独で入射した場合に得られる同相成分Ｉ_ｉ ^{ｍｅａｓｕｒｅ}と直交成分Ｑ_ｉ ^{ｍｅａｓｕｒｅ}を、全てのｉに関する同相成分の重ね合わせとなっているＩ^{ｄｉｇｉｔａｌ}と、全てのｉに関する直交成分の重ね合わせとなっているＱ^{ｄｉｇｉｔａｌ}に対して信号処理を行うことで分離する。具体的な信号処理の方法は、Ｉ^{ｄｉｇｉｔａｌ}とＱ^{ｄｉｇｉｔａｌ}から、Ｉ_ｉ ^{ｍｅａｓｕｒｅ}とＱ_ｉ ^{ｍｅａｓｕｒｅ}を正確に分離できるならどんな手法を用いても良い。例えば、Ｉ^{ｄｉｇｉｔａｌ}とＱ^{ｄｉｇｉｔａｌ}を、中心周波数がｆ_０＋ｆ_ｉであり通過帯域が２／Ｗであるデジタルバンドパスフィルタにそれぞれ通した上で位相遅延を保証することで、Ｉ_ｉ ^{ｍｅａｓｕｒｅ}とＱ_ｉ ^{ｍｅａｓｕｒｅ}を計算する方法などが考えられる。

また、前記方法では、アナログの電気信号の状態にある同相成分と直交成分を、ＡＤ変換してデジタル化した後に、各周波数成分への分離を行っているが、アナログの電気信号の状態にある同相成分と直交成分をアナログ電気フィルタによって各周波数成分へ分離した後にＡＤ変換するなどしても良い。

信号処理部１７ｃによって取得されたＩ_ｉ ^{ｍｅａｓｕｒｅ}とＱ_ｉ ^{ｍｅａｓｕｒｅ}を元に、信号処理部１７ｄで位相の計算を行う。まず、同相成分をｘ軸（実数軸）、直交成分をｙ軸（虚数軸）としたｘｙ平面上における複素ベクトルｒ_ｉを作成する。

また複素ベクトルの角度θ_ｉを以下のように計算する。

ここで、各光周波数ｆｉのパルスは時刻ｉ×Ｔ_Ｎ＋ｎ×Ｎ×Ｔ_Ｎ（ｎは任意の整数）に入射していると考えると、光ファイバ上の長手方向の入射端から距離ｚの位置での光ファイバの状態は、光パルスの伝搬時間を考慮して時刻ｉ×Ｔ_Ｎ＋ｎ×Ｎ×Ｔ_Ｎ＋z／ν（ｎは任意の整数）で測定している。ここで、νは光ファイバ中での光速である。さらに、散乱された散乱光が伝搬して入射端まで戻る時間を考慮すると、測定器での測定時刻は、ｉ×Ｔ_Ｎ＋ｎ×Ｎ×Ｔ_Ｎ＋２z／ν（ｎは任意の整数）となる。そこで、距離ｚの地点での散乱光の複素ベクトルｒ_ｉを、測定器の測定時刻を陽に表して、
ｒ_i （z，ｉＴ_Ｎ＋ｎＮＴ_Ｎ＋２z／ν）
と記述する。

本実施形態では、測定時刻ｍＴ_Ｎ＋２z／ν（ｍは整数）における位相θ（ｚ，ｍＴ_Ｎ＋２z／ν）を、ｍＴ_Ｎ＋２z／ν＝ｉＴ_Ｎ＋ｎＮＴ_Ｎ＋２z／νを満たすｉとｎを用いて、以下のように計算する。

そして、光ファイバ上での距離ｚ_１から距離ｚ_２の区間に加わった振動による位相変化を、数１－３ａと数１－３ｂとの差分を数１－３ｃとして計算する。

尚、光ファイバの状態を測定した瞬間の時刻は、上述のように散乱光が入射端に戻るのに要する時間は含めないので、距離ｚ_１の地点では時刻ｍＴ_Ｎ＋z_１／ν、距離ｚ_２の地点では時刻ｍＴ_Ｎ＋ｚ_２／ν、となり、時間差(z_１－z_２)／νだけ違いがある。しかし、ｚ_１とｚ_２との距離の差は空間分解能と同等程度で、通常は数ｍから数十ｍ程度に設定するため、時間差(z_１－z_２)／νは数十から数百ｎｓとなり、測定対象となる通常の振動の時間変化のスケールに対して非常に短いため、光ファイバの状態を測定した時刻の差は無視できる。そのため、該当区間に加わった振動を正しく測定可能である。

ここで、各周波数ｆ_ｉ単体のパルス列に注目すると、ＮＴ_Ｎの間隔だけ離れている。測定距離をＺとすれば、光パルスの往復時間は２Ｚ／νとなるため、ＮＴ_Ｎ≧２Ｚ／νでなければならない。つまり、Ｔ_Ｎ≧２Ｚ／（νＮ）である。測定距離Ｚが決まっても、多重数Ｎを増やせば、サンプリングの時間間隔Ｔ_Ｎを細かくすることができる。

本実施形態では、測定距離と測定可能な振動周波数の上限とのトレードオフを考慮した上で、必要となる周波数多重数Ｎについての設計方法を示す。下記する設計方法に従い必要な多重数Ｎを計算部１８が計算し、その情報に基づき光変調器３を動作させることで、必要な多重数Ｎを備えた光パルス列４を作成する。

まず、単一の振動周波数の振動を測定する場合を考える。振動周波数ｆ_νがナイキスト周波数以下になるための条件は、

である。また、ある時刻ｔから時刻ｔ＋Ｔ_Ｎの間の位相変化の絶対値Ｃ（ｔ）は、光ファイバ上で目的とする区間での位相変化を、時間に依存しない係数Ａを用いてＡ×ｓｉｎ（２πｆ_νｔ）とすれば、

である。ただし、位相の時間変化をＡ×ｓｉｎ（２πｆ_νｔ）として、初期位相はゼロと仮定しているが、この仮定をしても、以降の議論の一般性は失われない。時刻ｔに依存して、Ｃ（ｔ）の大きさも変化するが、最も大きい値は、

である。位相アンラップが失敗しない条件は、

である。これを変形すると、

となる。π／２Ａが１より大きい場合、つまりＡがπ／２より小さい場合には、常に上式は成立するので、振動周波数ｆ_νがナイキスト周波数以下になるための条件のみを考えればよい。Ａがπ／２以上の場合は、πｆ_ν Ｔ_Ｎ＜π／２が成立するとすれば、以下のように変形可能である。

ところで、πｆ_ν Ｔ_Ｎ＜π／２が成立するという仮定は、振動周波数ｆ_νがナイキスト周波数以下になるための式（１－４）の条件と同一であり、上で得られた式（１－９）の条件はＴ_Ｎについて、より厳しい制約である。したがって、以上の議論をまとめると、Ａがπ／２より小さい場合には不等式（１－４）が成立すればよく、Ａがπ／２以上の場合には不等式（１－９）が成立すればよい。ところで、多重数Ｎについては、測定距離Ｚに対して、ＮＴ_Ｎ≧２Z／νという条件があった。

つまり、Ａがπ／２より小さい場合には、

を満たすような最小のＮを、使用する多重数に設定する。設定したＮを使用して、Ｔ_ＮはＮＴ_Ｎ≧２Z／νと不等式（１－４）を満たす範囲で設定することが可能である。測定距離を、ＮＴ_Ｎ≧２Z／νを満たす範囲で長くしたい場合には、Ｔ_Ｎはなるべく大きい値に設定し、サンプリングを細かく設定したい場合には、Ｔ_Ｎはなるべく小さい値に設定する。

また、Ａがπ／２以上の場合には、

を満たすような最小のＮを、使用する多重数に設定する。設定したＮを使用して、Ｔ_ＮはＮＴ_Ｎ≧２Z／νと不等式（１－９）を満たす範囲で設定することが可能である。測定距離を、ＮＴ_Ｎ≧２Z／νを満たす範囲で長くしたい場合には、Ｔ_Ｎはなるべく大きい値に設定し、サンプリングを細かく設定したい場合には、Ｔ_Ｎはなるべく小さい値に設定する。

次に、任意の振動波形について考える。目的とする光ファイバ上での区間での位相変化の時間波形をＦ_ν（ｔ）とする。そのフーリエ変換を

とする。数１－１１ａが占める周波数帯域の最大値を

とする。数１－１１ｂの振動周波数がナイキスト周波数以下になるための条件は、

である。また、ある時刻ｔから時刻ｔ＋Ｔ_Ｎの間の位相変化の絶対値Ｃ（ｔ）は、

となる。時刻ｔを任意に変化させた時の、Ｃ（ｔ）の最大値をＣ_ｍａｘ（Ｔ_Ｎ）とすると、位相アンラップに失敗しないための条件は、

となる。任意のＴ_Ｎに対して数１－１４が成立する場合は、つまりＦ_ν（ｔ）の最大値と最小値の差の大きさがπより小さい場合は、式（１－１２）の条件のみを考えれば良い。それ以外の場合には、式（１－１２）が成立する条件では、Ｃ_ｍａｘはＴ_Ｎに対して増加関数である。そこで、Ｃ_ｍａｘの逆関数Ｃ_ｍａｘ ^－１を式（１－１４）に演算して、

の条件が得られる。式（１－１５）の条件は式（１－１２）の条件よりもＴ_Ｎについて、より厳しい制約となる。そこで、Ｆ_ν（ｔ）の最大値と最小値の差の大きさがπより小さい場合は式（１－１２）の条件を考え、π以上の場合は式（１－１５）の条件を考えればよい。

尚、逆関数の計算の方法は、解析的に計算が可能な場合は数式を使用して導く。また、解析的に計算が不可能な場合には、Ｔ_Ｎとｔを十分細かく変化させながら、Ｃ（ｔ）を数値計算する。数値計算の方法としては、様々な公知の手法を適用可能である。例えば、Ｔ_Ｎを最初は粗く変化させて、式（１－１４）を満たさない最小のＴ_Ｎと、式（１－１４）を満たす最大のＴ_Ｎを求める。その範囲内で、Ｔ_Ｎをより細かく変化させて、式（１－１４）を満たさない最小のＴ_Ｎと、式（１－１４）を満たす最大のＴ_Ｎをアップデートする。これを繰り返して、十分な精度で式（１－１５）を満たす最大のＴ_Ｎを計算することが可能である。

多重数Ｎについては、単一の振動周波数の振動を考えた時と同様にして、Ｆ_ν（ｔ）の最大値と最小値の差の大きさがπより小さい場合は、式（１－１０）においてＦ_νを数１－１１ｂに置き換えた条件を満たす最小のＮを、周波数多重数として設定する。Ｆ_ν（ｔ）の最大値と最小値の差の大きさがπ以上の場合は、

を満たす最小のＮを、周波数多重数として設定する。

Ｔ_Ｎの設定も単一の振動周波数の振動を考えた時と同様にすればよい。

計算部１８は、上述のように設定したＮに基づきパルス列信号ＰＳを作成し、光変調器３に与える。

（実施形態２）
実施形態１の方法では、入射端から距離ｚの地点から散乱された散乱光ベクトルの、測定時刻ｍＴ_Ｎ＋２z／νでの位相値θ（ｚ，ｍＴ_Ｎ＋２z／ν）を、ｍＴ_Ｎ＋２z／ν＝ｉＴ_Ｎ＋ｎＮＴ_Ｎ＋２z／νを満たすｉとｎを用いて式（１－３）のように計算した。式（１－３）の右辺において、ｍ＝ｉ＋ｎＮであるため、ｉとｎを、ｍに依存することが分かるように、ｉ（ｍ）とｎ（ｍ）と以降では記述する。時刻ｍＴ_Ｎ＋２z／νから時刻ｍ’Ｔ_Ｎ＋２z／νでの位相変化量は、式（１－３）を用いれば、

となる。ｉ（ｍ）とｉ（ｍ’）が表す光周波数が異なる場合、式（２－１）は、実際の位相変化量に、異なる光周波数間の角度差が加算されたものとなっている。つまり、実際の位相変化量をΔθとして、異なる光周波数間の角度差をφ（ｉ（ｍ’），ｉ（ｍ））とすれば、

となる。

式（２－２）から分かるように、測定される位相変化の波形は、異なる周波数間の角度差φ（ｉ（ｍ’），ｉ（ｍ））があるため、実際の位相変化を表す波形とは、異なってしまう。この様子を図２と図３に示す。図２と図３では多重数３で正弦波の振動を測定した場合の例である。図２は、実施形態１に記載の方法により周波数多重した入射光パルス列を表している。図３（ａ）は、入射端から距離ｚの位置で測定される位相である。

本実施形態では、異なる光周波数の角度差φ（ｉ（ｍ’），ｉ（ｍ））を補正するための、周波数多重パルスの並べ方と信号処理方法について記述する。異なる光周波数の角度差を補正した場合の位相は、図３（ｂ）となり、実際の位相変化を正しく測定することが可能となる。

異なる光周波数の角度差を補正するために、前記光源は、前記光パルス列を前記被測定光ファイバの一端に入射する前に、前記光パルス列の光パルス間隔より短い間隔で異なる光周波数の前記光パルスを並べた光パルス対を並べた補正用光パルス列を前記被測定光ファイバの一端に入射し、前記信号処理部は、前記散乱光の位相成分の補正値を数Ｃ２で計算する。

図２２は、本振動検出装置が行う光パルス試験方法を説明する図である。本光パルス試験方法は、図１の光パルス試験方法に、前記光パルス列を被測定光ファイバ６の一端に入射する前に、前記光パルス列の光パルス間隔より短い間隔で異なる光周波数の前記光パルスを並べた光パルス対を並べた補正用光パルス列を被測定光ファイバ６の一端に入射すること（ステップＳ２２ａ）と、前記散乱光の位相成分の補正値を数Ｃ２で計算すること（ステップＳ２５）が付加される。

以下に詳細を説明する。
異なる光周波数の角度差は、異なる光周波数を振動による光ファイバの状態変化が無視できる時間間隔で並べて被測定光ファイバに入射し、その信号を処理することにより計算することが可能である。光ファイバの状態変化が無視できる微小パルス間の時間間隔は、振動周波数ｆνと振動の大きさに依存するが、通常は数ｎｓ程度としておけば十分である。そのため、図４に示すように、実施形態１で記載した測定用のメインパルス列（２０など）とは異なる、使用する全ての光周波数のパルスを振動によるファイバの状態変化が無視できる時間間隔で並べた光周波数多重パルスを、ＮＴ_Ｎの繰り返し間隔で並べた異なる光周波数の角度差補正用のパルス列２１を入射する。パルス列２１で繰り返し間隔をＮＴ_Ｎとしているのは、同じ光周波数のパルスの時間間隔は２z／ν以上でなければならないためである。

パルス列２１を構成する周波数多重パルスの数がＰ個であるとする。第ｐ番目（ｐは１からＰの整数）の周波数多重パルスの先頭を入射した時刻をｔ_ｐとすれば、各光周波数成分ｆ_ｉのパルスを入射した時刻はｔ_ｐ＋δｔ_ｉとなる。δｔ_ｉは周波数多重パルスの先頭と各光周波数成分のパルスの先頭との時間差であり設計時に既知である。入射端から距離ｚからの散乱光の複素ベクトルは、各光周波数成分の光パルスの往復時間を考えると測定時刻ｔ_ｐ＋２ｚ／ν＋δｔ_ｉで観測される。各光周波数パルスの測定時刻ｔ_ｐ＋２ｚ／ν＋δｔ_ｉでの散乱光ベクトルを単にｒ_ｉ（ｚ，ｔ_ｐ＋２ｚ／ν）と記述する。異なる光周波数の角度差φ（ｉ（ｍ’），ｉ（ｍ））を以下のように計算する。

ここで、ａｒｇ関数は引数の複素ベクトルの偏角を－πからπの範囲で計算する関数であり、出力は実数である。また、Ｒ（＊）は、複素ベクトルｒに対して作用する演算子であり、複素平面上で角度＊だけ時計まわりに複素ベクトルを回転させるものである。

式（２－３）は、「付録」に記載の「第１ベクトル回転手順」と同様の原理を用いており、「第１ベクトル回転手順」における「ｉ（ｍ）を基準波長とみなしたｉ（ｍ’）の波長の回転角度」として、異なる光周波数の角度差を計算する（ただし符号は正負逆である）。式（２－３）は、「付録」の式（１１）～（１３）に対応している。

角度差φ（ｉ（ｍ’），ｉ（ｍ））は、ｚには依存するが、レーザの周波数ドリフト等の影響を除けば、測定時刻には依存しないため、ｉ（ｍ’）とｉ（ｍ）に対応する周波数が同じであれば、ｍ’やｍごとに個別に角度差φ（ｉ（ｍ’），ｉ（ｍ））を計算する必要はない。例えば、ｉ（２）とｉ（１）に対応する周波数がｆ_１とｆ_２の場合に、式（２－３）で角度差φ(ｉ（２），ｉ（１）)を計算すれば、ｉ（５）とｉ（４）に対応する周波数がｆ_１とｆ_２の場合に、角度差φ(ｉ（５），ｉ（４）)としてφ(ｉ（２），ｉ（１）)の値を使用することが可能である。つまり、周波数のみが重要であるため、φ（ｉ（ｍ’），ｉ（ｍ））＝φ（ｆ_{ｉ（ｍ’）}，ｆ_ｉ（ｍ））とできる。

式（２－３）では、各光周波数の複素ベクトルの角度の差自体は、雑音を除いて一定であること利用し、複数時刻の複素ベクトルの平均をとることで、雑音の影響を低減し、φ（ｉ（ｍ’），ｉ（ｍ））を計算している。特に、複素ベクトルの回転演算子Ｒ（＊）を用いて複素ベクトルの回転を行うことで、効率よく雑音の低減が可能である。

式（２－３）で計算した異なる光周波数の角度差φ（ｉ（ｍ’），ｉ（ｍ））を用いて、式（２－２）より、実際の位相変化Δθは、

で計算できる。つまり、信号処理の手順としては、測定用メインパルス列（２０など）で取得した信号に対して実施形態１に記載の方法によりθ（ｚ，ｍＴ_Ｎ＋２z／ν）を計算し、異なる光周波数の角度差補正用のパルス列で取得した信号から式（２－３）に基づいて異なる光周波数の角度差φ（ｉ（ｍ’），ｉ（ｍ））を計算し、それらの値を用いて式（２－４）から実際の位相変化量Δθを計算する。

尚、パルス列２１を構成する周波数多重パルスの数Ｐは、式（２－３）で計算する角度差φ（ｉ（ｍ’），ｉ（ｍ））に必要な精度から設計する。数Ｐが増加するほど、測定対象となる振動をモニタするために必要なサンプリング間隔Ｔ_Ｎで測定できない時間帯が増加するため、数Ｐは必要最小限の値に設定する。

事前に、測定対象である振動を測定するために必要な角度差φ（ｉ（ｍ’），ｉ（ｍ））の精度±ξを決定する。例えば、角度差φ（ｉ（ｍ’），ｉ（ｍ））が真の値から＋ξだけ異なる値で計算された場合、ｉ（ｍ’）とｉ（ｍ）の光周波数と同じ組み合わせの光周波数に対応する、他の時刻間の位相の時間変化も全て、＋ξだけ実際の位相変化に対して歪んで、測定される。そのため、測定対象である振動の大きさや波形などの事前情報に基づき、許容できる計算の精度±ξを見積もる。

一方で、周波数多重パルスの数Ｐと式（２－３）で計算する角度差φ（ｉ（ｍ’），ｉ（ｍ））の精度±δφの関係式は以下のように導出できる。まず、フェーディングのため、入射端からの距離ｚの具体的な値に依存して、散乱光の振幅の大きさに対応する各光周波数の散乱光ベクトルの大きさである｜ｒ_{ｉ（ｍ’）}｜と｜ｒ_ｉ（ｍ）｜にバラツキが生じるため、精度±δφも距離ｚによって変化してしまう。そのため、予め数Ｐを設定する目的のために、｜ｒ_{ｉ（ｍ’）}｜と｜ｒ_ｉ（ｍ）｜のそれぞれの統計的な平均値を用いて、数Ｐと精度±δφの統計的な平均値との関係を導出する。また、式（２－３）中では、異なる時刻ｔ_ｐ＋２z／ν（ｐ＝１，・・・，Ｐ）の平均を演算しているため、異なる時刻の｜ｒ_{ｉ（ｍ’）}｜と｜ｒ_ｉ（ｍ）｜が計算結果に影響を与えるが、それらの値はパルス列２１を入射している間に光ファイバに振動が生じることで時間的に変化してしまうため、計算結果としての角度差φ（ｉ（ｍ’），ｉ（ｍ））の精度±δφも変化する。そのため、角度差φ（ｉ（ｍ’），ｉ（ｍ））の精度±δφは、振動による｜ｒ_{ｉ（ｍ’）}｜と｜ｒ_ｉ（ｍ）｜の時間変化に依存するが、予め数Ｐを設定する際に振動による変化を考慮するのは困難なため、ここでは予め数Ｐを設定する目的のため、｜ｒ_{ｉ（ｍ’）}｜と｜ｒ_ｉ（ｍ）｜の時間変化はないものとする。｜ｒ_{ｉ（ｍ’）}｜と｜ｒ_ｉ（ｍ）｜のそれぞれの統計的な平均値をＩ_{ｉ（ｍ’）}（ｚ）とＩ_ｉ（ｍ）（ｚ）とする。これらの値は、例えば、入射光強度と光ファイバの伝搬損失値から計算して評価する、あるいは、光強度損失を測定可能なハンディタイプのコヒーレンスの良くないレーザ光源を用いたＯＴＤＲ装置を用いて被測定光ファイバの光強度損失分布を予め測定する、などの方法で取得することが可能である。後者の方法では、敢えてコヒーレンスの良くないレーザ光源を用いることで、光強度損失の統計的な平均値を測定することが可能となる。得られたＩ_{ｉ（ｍ’）}（ｚ）とＩ_ｉ（ｍ）（ｚ）を用いて、図５に示すＩＱ平面上でのベクトル２３とベクトル２５を作成する。ベクトル２３は、Ｉ軸に平行なベクトルとし、その長さはＩ_ｉ（ｍ）（ｚ）とする。ベクトル２４は、式（２－３）で計算した角度差φ（ｉ（ｍ’），ｉ（ｍ））だけＩ軸に対して半時計回りに角度を有しており、長さをＩ_{ｉ（ｍ’）}（ｚ）とする。

ベクトル２３とベクトル２４には、測定におけるショット雑音やＰＤの熱雑音、その後の電気回路等での雑音などが合わさった雑音に起因する、不確かさ２５と不確かさ２６があると考える。シングルショット計測時での不確かさの大きさを、不確かさの円の半径として、“Ｎｏｉｓｅ”と表記する。Ｎｏｉｓｅ値は、測定器の性能のみから決まるため、事前に測定可能である。周波数多重パルスの時間間隔はＮＴ_Ｎであり、ＰＤのインパルス応答等よりも十分に長いため、異なる時刻ｔ_ｐ＋２z／ν同士の雑音には相関がなくランダムであるとみなせるため、不確かさ２５と不確かさ２６の大きさは、不確かさの円の半径として、周波数多重パルスの数Ｐを用いて、Ｎｏｉｓｅ／√Ｐと表せる。したがって、±δφの統計平均値は次のように近似できる。

式（２－５）はＰについて減少関数であるから、δφ≦ξを満たす最小のＰを設定すればよい。

具体的手順をまとめると、測定対象である振動を測定するために必要な角度差φ（ｉ（ｍ’），ｉ（ｍ））の精度±ξを測定対象である振動の大きさや波形などの事前情報に基づき決定し、散乱光の振幅の大きさの平均であるＩ_{ｉ（ｍ’）}（ｚ）とＩ_ｉ（ｍ）（ｚ）、および雑音の大きさＮｏｉｓｅを評価／測定し、式（２－５）に代入してδφとＰとの関係式を求め、最終的に、δφ≦ξを満たす最小のＰを計算し、図４のパルス列２１の周波数多重パルスの数として設定する。

尚、上記計算では散乱光の振幅の大きさの平均に比例するＩ_{ｉ（ｍ’）}（ｚ）とＩ_ｉ（ｍ）（ｚ）を用いて計算を行っているため、δφ≦ξを満たすのは統計的に全地点のうち５０％の地点となる。この条件をより厳しくし、より多くの地点でδφ≦ξを満たす必要がある場合には、散乱光の振幅の大きさの平均よりも小さい数値に基づいて計算を実施することも可能である。また、入射端からの距離ｚについて、特に計測をしたい箇所が特定可能な場合は、事前にその箇所における｜Ｉ_{ｉ（ｍ’）}（ｚ）｜と｜Ｉ_ｉ（ｍ）（ｚ）｜を測定し、それらの値をＩ_{ｉ（ｍ’）}（ｚ）とＩ_ｉ（ｍ）（ｚ）に替えて式（２－５）と同様の式を導出し、数Ｐを決定することも可能である。

また、パルス列２１を入射している間に光ファイバに振動が生じることで｜Ｉ_{ｉ（ｍ’）}｜と｜Ｉ_ｉ（ｍ）｜の時間変化はないものとして式（２－５）を導出したが、パルス列２１を入射している間に光ファイバに振動が生じることによる｜Ｉ_{ｉ（ｍ’）}｜と｜Ｉ_ｉ（ｍ）｜の時間変化に関する情報が事前に得られる場合は、その情報を取り入れることも可能である。

尚、レーザの周波数ドリフトにより、各光周波数の絶対値が変化すると、角度差φ（ｉ（ｍ’），ｉ（ｍ））にも変化が生じる。この対策として、図４のパルス列２１や２２で示しているように、異なる光周波数の角度差補正用のパルス列を定期的に挿入し、角度差φ（ｉ（ｍ’），ｉ（ｍ））の値をアップデートする。図４の例の場合、パルス列２１で異なる光周波数の角度差を取得し、その値を用いてパルス列２０で測定した位相を補正し、パルス列２２で改めて異なる光周波数間の角度差を取得し、その値を用いてパルス列２７で測定した位相を補正する。補正用のパルス列を挿入する時間間隔（例えば、図４の例では、パルス列２１とパルス列２２の間隔）は、レーザの周波数ドリフトの特性から決定される。定性的には、周波数ドリフトの程度が大きい場合には、異なる光周波数の角度差補正用のパルス列を挿入する時間間隔を短く設定する。

本実施形態の方法を使用すれば、実施形態１で発生する異なる光周波数間の角度差による波形歪みを補正することが可能である。一方で、本実施形態の方法では、時間間隔がＮＴ_Ｎで測定しなければならないパルス列２１や２２の時間帯が発生し、常に時間間隔Ｔ_Ｎで測定可能な実施形態１と比べて、欠点となる。そのため、常に時間間隔Ｔ_Ｎで測定することを重視する場合は実施形態１を使用し、波形歪みがなるべくないように振動計測を行うことを重視する場合は本実施形態を使用する，という使い分けを行う。

（実施形態３）
実施形態１と実施形態２は、周波数多重によりサンプリング点数を向上させることで、測定距離と測定可能な振動周波数の上限とのトレードオフを解決する方法に関する。本実施形態は、測定距離と測定可能な振動周波数の上限とのトレードオフの解決と、フェーディングによる位相検出ができない地点の発生の解決を同時に行うための、周波数多重の方法について記載する。具体的には、周波数多重のパルス列の構成方法と、信号処理の方法に関する。

測定対象となる振動の特性等から、測定距離と測定可能な振動周波数の上限とのトレードオフを解決するために必要な多重数Ｎを決定する。多重数Ｎの決定の方法は実施形態１と同様である。また、フェーディング対策に必要な多重数Ｍを決定する。多重数Ｍの決定の方法は「付録」の記載の通りである。つまり、雑音レベルの１／√Ｍ倍が所定値（被測定光ファイバの種類、測定対象の振動形態、あるいは光パルス試験装置の性能等に応じて設定される値）以下となるように前記微小光パルスの数Ｍを設定する。

決定した多重数Ｎと多重数Ｍに基づき、周波数多重のパルス列として図６のパルス列３１を被測定光ファイバに入射する。入射の方法などは実施形態１に記載した図１と同様で、変調器３を図６のパルスを生成するように動作させる。

パルス列３１は、Ｎ個の異なるパルス対からなるパターンを測定時間分だけ繰り返したものである。隣り合うパルス対の間隔はＴ_Ｎに設定し、それぞれのパルス対を入射して得られる信号から、異なる時刻のファイバの状態を測定する。各パルス対はＭ個の異なる光周波数パルスから構成する。Ｍ個の異なる光周波数パルスの時間間隔は、振動による光ファイバの状態変化が無視できるほど小さく設定する。したがって、Ｍ個の異なる光周波数パルスによる信号で、同じ時刻の光ファイバの状態を測定しているとみなせるため、それらＭ個の信号を用いてフェーディング対策をする。Ｎが２で、Ｍが２の例を、パルス列３２に示した。

具体的な信号処理の手順を記す。パルス対ｉの先頭を入射した時刻をｉ×Ｔ_Ｎ＋ｎＮＴ_Ｎ（ｎは任意の整数）とする。それぞれのパルス対の先頭の光周波数を基準波長にとり、「付録」の記載の方法に従い、パルス対を構成するＭ個の異なる光周波数の信号を平均することで、入射端から距離ｚの位置で散乱された散乱光ベクトルの位相θ_ｉ（ｉＴ_Ｎ＋ｎＴ_Ｎ＋２z／ν）を計算する。そして、実施形態１の式（１－３）のｉを周波数番号ではなくパルス対番号にかえることで、実施形態１と同様にして、θ（ｚ，ｍＴ_Ｎ＋２z／ν）を計算する。

（実施形態４）
実施形態３は、実施形態１でフェーディング対策を実施する場合への一般化となっている。実施形態３でＭ＝１とすると実施形態１となる。そのため、実施形態１で異なる光周波数の角度差φ（ｉ（ｍ’），ｉ（ｍ））が補正されておらず波形が歪んでしまう問題は、実施形態３でも同様に生じる。本実施形態では、実施形態２に記載した異なる光周波数の角度差補正の方法を一般化することで、実施形態３で生じる波形歪みを補正する。具体的には、周波数多重のパルス列の構成方法と、信号処理の方法に関する。

周波数多重のパルス列として図７（ａ）のパルス列４１を被測定光ファイバに入射する。入射の方法などは実施形態１に記載した図１と同様で、変調器３を図７（ａ）のパルスを生成するように動作させる。

パルス列４１は実施形態３で述べた測定用メインパルス列４３や４５の間に定期的に、異なる光周波数の角度差補正用のパルス列４４や４６を挿入して作成する。異なる光周波数の角度差補正用のパルス列４４や４６は、補正用パルス対４２を繰り返し間隔ＮＴ_Ｎで数Ｐだけ並べた物とする。補正用パルス対４２は、測定用メインパルス列４３（図６のパルス列３１と同じ）の各パルス対の先頭の光周波数ｆ_１，ｆ_{（Ｍ＋１）}，・・・，ｆ_{（（Ｎ－１）Ｍ＋１）}のパルスを、振動によるファイバの状態変化が無視可能な時間間隔で並べた構成とする。

測定用メインパルス列４３を入射して、実施形態３に記載の方法で、位相θ（ｚ，ｍＴ_Ｎ＋２z／ν）を得る。しかし、実施形態２の式（２－２）と同様に、測定される位相変化の波形は、異なる光周波数間の角度差φ（ｉ（ｍ’），ｉ（ｍ））があるため、実際の位相変化を表す波形とは、異なってしまう。ただし、角度差φ（ｉ（ｍ’），ｉ（ｍ））は本実施形態の場合には、ｉ（ｍ’）に対応するパルス対の先頭の周波数と、ｉ（ｍ）に対応するパルス対の先頭の周波数との角度差に対応する。この角度差φ（ｉ（ｍ’），ｉ（ｍ））は、パルス列４４の測定データを用いて式（２－３）を用いて計算することが可能である。計算した位相θ（ｚ，ｍＴ_Ｎ＋２z／ν）と角度差φ（ｉ（ｍ’），ｉ（ｍ））とを式（２－４）に代入して、異なる光周波数間の角度差を補正した位相を計算する。

尚、補正用パルス対４２に替えて、図７（ｂ）の補正用パルス対４７を使用してもよい。補正用パルス対４７は、全ての光周波数のパルスを、振動による被測定光ファイバの状態変化が無視できる時間間隔で並べて構成される。このパルス対４７を用いても上記記載の方法と同様にして角度差φ（ｉ（ｍ’），ｉ（ｍ））が計算可能である。補正用パルス対４２の使用との違いは、補正用パルス対４７を構成する周波数の数が、補正用パルス対４２と比べてＭ倍となるため、パルス列４４で測定した位相の時間変化については、計算時間は増加するが、フェーディング対策に割く光周波数の数が増えているため、感度がより良い測定が可能な点である。補正用パルス対４２を使うか、補正用パルス対４７を使うかは、パルス列４４で測定する位相の時間変化について計算時間を要しても感度よく計算したい場合は補正用パルス対４７を選択し、計算時間を少しでも減らしたい場合やパルス列４４を入射している時間帯では位相変化の測定は必要ない場合等は補正用パルス対４２を使用する。

尚、異なる光周波数の角度差補正用のパルス列４４を構成するパルス対４２（あるいはパルス対４７）の数Ｐについては、実施形態２と同様にして設定する。

また、レーザの周波数ドリフトの対策に関しても実施形態２と同様にして、異なる光周波数の角度差補正用のパルス列を定期的に挿入して、角度差φ（ｉ（ｍ’），ｉ（ｍ））を定期的に更新して使用して対応する。図７（ａ）のパルス列４１の例の場合、パルス列４３で計算した位相の補正にはパルス列４４の計算結果を用い、パルス列４５で計算した位相の補正にはパルス列４６の計算結果を用いる。異なる光周波数の角度差補正用のパルス列を挿入する時間間隔（例えば４４と４６の間の時間）は、実施形態２と同様にして決定する。

（実施形態５）
実施形態４では、異なる光周波数の角度差補正用のパルス列４４やパルス列４７を入射している時間帯は、位相の時間変化を測定する時間間隔がＮＴ_Ｎとなってしまう。本実施形態では、使用する光周波数を実施形態４の場合より、多重数ＮやＭの値にかかわらず、１つだけ増やすことにより、全ての時間帯で位相の時間変化を時間間隔Ｔ_Ｎで測定する方法を述べる。具体的には、周波数多重のパルス列の構成方法と、信号処理の方法に関する。

使用するパルス列の具体的な構成方法にはいくつもの種類が考えられるが、本実施形態では大きく３通りの構成方法を示す。３通りの構成方法を構成方法Ａ、構成方法Ｂ、そして構成方法Ｃとして区別する。

［構成方法Ａ］
まず、最も簡単な構成方法Ａについて、信号処理の方法と合わせて述べる。

周波数番号の１からＮＭ＋１をこの順に左から並べる。並べ終わったら、その右側に、再度、周波数番号の１からＮＭ＋１をこの順に左から並べる。この動作をＭ回繰り返して、図８に示すような周波数番号列５１を作成する。周波数番号列５１の最も左からＭ個をパルス対１で使用する周波数の組み合わせに選ぶ。その右からＭ個をパルス対２で使用する周波数の組み合わせに選ぶ。これを繰り返して、パルス対組み合わせ５２にあるような、ＮＭ＋１個のパルス対で使用する周波数の組み合わせを決定する。これらＮＭ＋１個のパルス対を用いて、パルス列５３のような入射パルスを作成する。パルス列５３は、ＮＭ＋１個のパルス対を間隔Ｔ_Ｎで並べた一つのパターンを測定時間分だけ繰り返して構成する。

入射パルス列５３で得られた信号を次のように処理することで、異なる光周波数の角度差を補正済みの、全ての時間帯において時間間隔Ｔ_Ｎで測定した位相を計算する。

パルス対ｉの先頭を入射した時刻をｉ×Ｔ_Ｎ＋ｎ×Ｎ×Ｔ_Ｎ （ｎは任意の整数）とする。実施形態３と同様に、それぞれのパルス対の先頭の周波数を基準波長にとり、「付録」に記載の方法に従い、パルス対を構成するＭ個の異なる光周波数の信号を平均することで、入射端から距離ｚの位置での位相θ_ｉ（ｉＴ_Ｎ＋ｎＮＴ_Ｎ＋２z／ν）を計算する。そして、実施形態１の式（１－３）のｉを周波数番号ではなくパルス対番号にかえることで、実施形態１と同様にして、位相θ（ｚ，ｍＴ_Ｎ＋２z／ν）を計算する。

位相θ（ｚ，ｍＴ_Ｎ＋２z／ν）には、実施形態４の冒頭で述べたのと同様で、異なる光周波数間の角度差が含まれる。異なる光周波数間の角度差の補正を漏れなく行うためには、任意の二つのパルス対の先頭の光周波数の角度差補正を行う必要があるが、これには隣り合うパルス対の先頭の光周波数の角度差補正ができれば十分である。この理由は、任意の光周波数ｆ_ｉとｆ_ｊの角度差φ（ｆ_ｉ，ｆ_ｊ）は、別の任意の光周波数ｆ_ｋを利用して、以下のように記述できることによる。

ただし、ｉ，ｊ及びｋは任意の正の整数である。

この式を用いると、ｉ＜ｊを満たす正の整数ｉとｊを任意に選んだ時に、パルス対ｊの先頭の光周波数をｆ_ｊ ^ｐｆとし、パルス対ｉの先頭の光周波数をｆ_ｉ ^ｐｆとすれば、角度差φ（ｆ_ｊ ^ｐｆ，ｆ_ｉ ^ｐｆ）は以下のように隣り合うパルス対の先頭の光周波数の角度差の和として記述可能である。

ただし、ｉとｊは任意の正の整数であり、ｉ＜ｊである。

そして、任意の隣り合うパルス対の先頭の光周波数の角度差φ（ｆ_ｉ＋１ ^ｐｆ，ｆ_ｉ ^ｐｆ）（ｉは任意の正の整数）は、パルス対ｉの最後尾の光周波数ｆ_ｉ ^ｐｒを使用して以下のように記述できる。

式（５－３）中の角度差φ（ｆ_ｉ＋１ ^ｐｆ，ｆ_ｉ ^ｐｒ）と角度差φ（ｆ_ｉ ^ｐｒ，ｆ_ｉ ^ｐｆ）は次のようにして求める。まず、角度差φ（ｆ_ｉ ^ｐｒ，ｆ_ｉ ^ｐｆ）についてであるが、パルス対ｉは繰り返し入射しているため、パルス対ｉに含まれるｆ_ｉ ^ｐｒとｆ_ｉ ^ｐｆの角度差φ（ｆ_ｉ ^ｐｒ，ｆ_ｉ ^ｐｆ）は式（２－３）を用いて計算できる。式（２－３）中のｐは、測定時間内でｐ回目に登場したパルス対ｉに対応させる。例えば、ｐ＝１は測定時間内の最初のパルス対ｉで、ｐ＝Ｐは測定時間内の最後のパルス対ｉに対応させる。式（２－３）式中のｔ_ｐは、パルス対ｉが測定時間内でｐ回目に入射された時の時刻とする。また、角度差φ（ｆ_ｉ＋１ ^ｐｆ，ｆ_ｉ ^ｐｒ）についてであるが、ｆ_ｉ＋１ ^ｐｆとｆ_ｉ ^ｐｒは周波数番号が１つだけ異なるため、ｆ_ｉ＋１ ^ｐｆとｆ_ｉ ^ｐｒを含むパルス対が必ず一つ以上存在する。そのパルス対について、式（２－３）を同様に用いれば角度差φ（ｆ_ｉ＋１ ^ｐｆ，ｆ_ｉ ^ｐｒ）を計算できる。このようにして異なる光周波数間の角度差の補正を漏れなく行うことが可能なことが示された。

まとめると、式（５－２）と式（５－３）を用いれば、次の式が得られる。

具体的計算手順としては、まず式（５－４）の各項を、適切なパルス対を用いて式（２－３）式を使用して計算し、次に式（５－４）を用いて角度差φ（ｆ_ｊ ^ｐｆ，ｆ_ｉ ^ｐｆ）を計算し、位相θ（ｚ，ｍＴ_Ｎ＋２z／ν）と角度差φ（ｆ_ｊ ^ｐｆ，ｆ_ｉ ^ｐｆ）を用いて式（２－４）と同様の要領で、最終的な位相変化を計算する。

尚、上記の説明では、それぞれのパルス対の先頭の光周波数を基準として位相θ（ｚ，ｍＴ_Ｎ＋２z／ν）を計算し、それぞれのパルス対の先頭と最後尾の光周波数を利用して式（５－４）のように角度差を計算した。しかし、パルス対の先頭や最後尾の周波数を必ず使用する必要はなく、それぞれのパルス対の先頭からｘ番目の光周波数を基準として位相θ（ｚ，ｍＴ_Ｎ＋２z／ν）を計算可能であるし、それぞれのパルス対の先頭からｘ番目と最後尾からｙ番目の光周波数を利用して式（５－４）と同様な式を作成した上で、右辺の各項の角度差を計算し、角度差補正を漏れなく行うことができる（ただし、同一パルス対において先頭からｘ番目と最後尾からｙ番目は別の光周波数に対応するようにｘとｙは選ぶ必要はある）。このため、装置の不具合等により、特定の光周波数の測定が困難となった場合でも、柔軟な対応が可能となる。

［構成方法Ｂ］
次に構成方法Ｂについて、信号処理方法と合わせて記述する。

まず、パルス列の構成方法について記述する。周波数番号の１からＮＭをこの順に左から並べる。並べ終わったら、その右側に、再度、周波数番号の１からＮＭをこの順に左から並べる。この動作をＮ＋１回繰り返して、図９に示すような周波数番号列５０１を作成する。周波数番号列５０１の最も左からＭ個をパルス対１で使用する周波数の組み合わせに選ぶ。その右からＭ個をパルス対２で使用する周波数の組み合わせに選ぶ。これを繰り返して、パルス対組み合わせ５０２にあるような、Ｎ（Ｎ＋１）個のパルス対で使用する周波数の組み合わせを決定する。次に、パルス対の周波数組み合わせ５０３に示すように、周波数ｆ_ＮＭ＋１をパルス対番号が１＋ｋ（Ｎ＋１）のパルス対に追加する（ｋ＝０，１，・・・，（Ｎ－１））。これにより、パルス対１＋ｋ（Ｎ＋１）を構成する周波数の数がＭ＋１個となる。これらＮ（Ｎ＋１）個のパルス対を用いて、パルス列５０４のような入射パルスを作成する。パルス列５０４は、Ｎ（Ｎ＋１）個のパルス対を先頭の時間間隔がＴ_Ｎとなるように並べた一つのパターンを測定時間分だけ繰り返して構成する。

入射パルス列５０４で得られた信号を次のように処理することで、異なる光周波数の角度差を補正済みの、全ての時間帯において時間間隔Ｔ_Ｎで測定した位相を計算する。

パルス対ｉの先頭を入射した時刻をｉ×Ｔ_Ｎ＋ｎ×Ｎ×Ｔ_Ｎ （ｎは任意の整数）とする。実施形態３と同様に、それぞれのパルス対の先頭の周波数を基準波長にとり、「付録」に記載の方法に従い、パルス対を構成するＭ個の異なる光周波数の信号を平均することで、入射端から距離ｚの位置での位相θ_ｉ（ｉＴ_Ｎ＋ｎＮＴ_Ｎ＋２z／ν）を計算する。そして、実施形態１の式（１－３）のｉを周波数番号ではなくパルス対番号にかえることで、実施形態１と同様にして、位相θ（ｚ，ｍＴ_Ｎ＋２z／ν）を計算する。

位相θ（ｚ，ｍＴ_Ｎ＋２z／ν）には、実施形態４の冒頭で述べたのと同様で、異なる光周波数間の角度差が含まれる。異なる光周波数間の角度差の補正を漏れなく行うためには、任意の二つのパルス対の先頭の光周波数の角度差補正を行う必要がある。ｉ＜ｊを満たす正の整数ｉとｊを任意に選んだ時に、パルス対ｊの先頭の光周波数をｆ_ｊ ^ｐｆとし、パルス対ｉの光周波数をｆ_ｉ ^ｐｆとすれば、角度差φ（ｆ_ｊ ^ｐｆ，ｆ_ｉ ^ｐｆ）は以下のようにｆ_ＮＭ＋１を用いて展開できる。

ただし、ｉとｊは任意の正の整数であり、ｉ＜ｊである。

本実施形態で使用するパルス対の光周波数の組み合わせ５０３では、周波数ｆ_ＮＭ＋１をパルス対番号が１＋ｋ（Ｎ＋１）のパルス対（ｋ＝０，１，・・・，（Ｎ－１））に追加しているため、周波数ｆ_ＮＭ＋１と他の周波数とは、必ず１回、同一のパルス内に存在している。そのため、式（５－５）の右辺の各項を、対応するパルス対について式（２－３）を用いることで計算可能である。得られたφ（ｆ_ｊ ^ｐｆ，ｆ_ｉ ^ｐｆ）の値を用いて、式（２－４）と同様の要領で、位相θ（ｚ，ｍＴ_Ｎ＋２z／ν）から最終的な位相を計算する。

尚、上記手順では、それぞれのパルス対の先頭の周波数を基準として位相θ（ｚ，ｍＴ_Ｎ＋２z／ν）を計算し、各パルス対の先頭の周波数とｆ_ＮＭ＋１を用いてφ（ｆ_ｊ ^ｐｆ，ｆ_ｉ ^ｐｆ）を計算し、最終的な位相を取得した。しかし、各パルス対の先頭の周波数を必ずしも使用する必要はなく、先頭からｘ番目（ｘは１≦ｘ≦Ｍを満たす整数）を任意に選んで使用可能である。ただし、パルス列の構成の仕様からｆ_ＮＭ＋１は必ず使用する必要がある。

構成方法Ｂと構成方法Ａとの違いは、構成方法Ｂでは（Ｎ＋１）個のパルス対に１つのパルス対は、含まれる多重数の数がＭ＋１となるため、そのパルス対で測定した時刻の位相については、フェーディング対策に使用する多重数が１つ増えており、より感度良く計算が可能となる利点がある。一方で、構成方法Ｂは、ＮＭ＋１番目の光周波数ｆ_ＮＭ＋１を必ず用いて、異なる光周波数間の角度差補正を行う必要があるため、装置の不具合等によりＮＭ＋１番目の光周波数の測定が困難となった場合には、異なる光周波数の角度差補正ができなくなる欠点がある。また、多重数Ｎと多重数Ｍに応じて、構成方法Ａと構成方法Ｂのパターン長にも違いが生じる。

［構成方法Ｃ］
最後に構成法Ｃについて、信号処理方法と合わせて記述する。

まず、パルス列の構成方法について記述する。周波数番号の１からＮＭ＋１をこの順に左から並べる。並べ終わったら、その右側に、再度、周波数番号の１からＮＭ＋１をこの順に左から並べる。この動作を（Ｎ＋１）Ｍ＋１回繰り返して、図１０に示すような周波数番号列５００１を作成する。次に、パルス対１＋ｋ（Ｎ＋１）に含まれる光周波数の数がＭ＋１個、それ以外のパルス対に含まれる光周波数の数がＭ個、となるように、左から順番にパルス対を作成していく（ｋ＝０，１，・・・）。結果として合計（ＮＭ＋１）（Ｎ＋１）個のパルス対ができる。これら（ＮＭ＋１）（Ｎ＋１）個のパルス対を用いて、パルス列５００３のような入射パルスを作成する。パルス列５００３は、（ＮＭ＋１）（Ｎ＋１）個のパルス対を先頭の時間間隔がＴ_Ｎとなるように並べた一つのパターンを測定時間分だけ繰り返して構成する。

パルス対ｉの先頭を入射した時刻をｉ×Ｔ_Ｎ＋ｎ×Ｎ×Ｔ_Ｎ （ｎは任意の整数）とする。実施形態３と同様に、それぞれのパルス対の先頭の周波数を基準波長にとり、「付録」に記載の方法に従い、パルス対を構成するＭ個の異なる光周波数の信号を平均することで、入射端から距離ｚの位置での位相θ_ｉ（ｉＴ_Ｎ＋ｎＮＴ_Ｎ＋２z／ν）を計算する。そして、実施形態１の式（１－３）のｉを周波数番号ではなくパルス対番号にかえることで、実施形態１と同様にして、位相θ（ｚ，ｍＴ_Ｎ＋２z／ν）を計算する。

位相θ（ｚ，ｍＴ_Ｎ＋２z／ν）には、実施形態４の冒頭で述べたのと同様で、異なる光周波数間の角度差が含まれる。異なる光周波数間の角度差の補正を漏れなく行うためには、任意の二つのパルス対の先頭の光周波数の角度差補正を行う必要がある。構成方法Ｃの並べ方では、一つのパターンの中で、全ての光周波数がＮ＋１個ずつパルス対の先頭となる。一方で、隣り合う周波数番号の光周波数が含まれるパルス対の数も、全ての周波数番号についてＮ＋１個である。したがって、特定の周波数番号に偏ることなく、異なる光周波数間の角度差の補正を漏れなく行うことが可能である。

特に、特定の周波数番号の光周波数の測定が困難となった場合でも、その光周波数を先頭とするパルス対については別の光周波数を基準にして位相θ_ｉ（ｉＴ_Ｎ＋ｎＮＴ_Ｎ＋２z／ν）が計算でき、その周波数を使用せずに異なる光周波数間の角度差の補正を依然として行うことが可能である。

まとめると、構成方法Ｃは、（Ｎ＋１）個のパルス対に１つのパルス対に含まれる多重数の数がＭ＋１となるため、そのパルス対で測定した時刻の位相については、より感度良く計算が可能となるという構成方法Ｂの利点を残しつつ、ＮＭ＋１番目の光周波数を必ず用いて、異なる光周波数間の角度差補正を行う必要があるという条件を排除した方法となる。ただし、パターン長がＡやＢに比べて長くなるという欠点がある。

実際の測定では、測定対象を考慮し、最も適切なパルス構成方法を選択し、そのパルス列に対応する信号ＰＳを図１の計算部１８から変調器３に送信し動作させる。信号処理方法は上記記載の通りである。

尚、本発明は、上記実施形態例そのままに限定されるものではなく、実施段階ではその要旨を逸脱しない範囲で構成要素を変形して具体化可能である。

［付録］
（例１）
図１１は、本実施形態のＤＡＳ－Ｐで振動検出を行う振動検出装置を説明する図である。本振動検出装置は、受信系に９０度光ハイブリッドを用いてコヒーレント検波を行う測定器３１と、信号処理装置１７とを備える。

測定器３１は、次のように被測定光ファイバ６からの散乱光を測定する。ＣＷ光源１から周波数がｆ_０の単一波長の連続光が射出され、カプラ２により参照光とプローブ光に分岐される。プローブ光は、光変調器３によって、光パルス４のように波長多重の光パルスに整形される。光パルス４は、周波数がｆ_０＋ｆ_ｉ（ｉは整数）かつパルス幅が光ファイバ長手方向での測定の空間分解能に対応する値Ｗに設定された微小パルス４ａが、ｉ＝１，２，・・・，Ｎ（Ｎは整数）だけ並んだ構成である。ｆ_ｉは、各時刻及び各地点における散乱光の強度が、異なるｉ同士で無相関とみなせる程度まで十分に離れているように選択をする。

光変調器３の種類は光パルス４を生成できるならば具体的な指定はなく、数が複数の場合もある。例えば、ＳＳＢ変調器や周波数可変なＡＯ変調器などを用いても良いし、パルス化における消光比を大きくするためにさらにＳＯＡなどによる強度変調を行っても良い。

光パルス４は、サーキュレータ５を介して、被測定光ファイバ６に入射される。光ファイバ６の長手方向の各点で散乱された光が、後方散乱光としてサーキュレータ５に戻り、９０度光ハイブリッド７の一方の入力部に入射される。カプラ２により分岐された参照光は、９０度光ハイブリッド７のもう一方の入力部に入射される。

９０度光ハイブリッド７の内部構成は、９０度光ハイブリッドの機能さえ備えていれば、なんでもよい。構成例を図１１に示す。後方散乱光は、５０：５０の分岐比のカプラ８に入射され、２分岐された散乱光が、５０：５０の分岐比のカプラ１２と、５０：５０のカプラ１１の入力部に入射される。参照光は、５０：５０の分岐比のカプラ９に入射され、２分岐された参照光の一方が、カプラ１１の入力部に入射され、他方が、位相シフタ１０で位相をπ／２だけシフトされてカプラ１２の入力部に入射される。

カプラ１１の２つの出力がバランス検出器１３によって検出され、アナログの同相成分I^{ａｎａｌｏｇ}である電気信号１５が出力される。カプラ１２の２つの出力がバランス検出器１４によって検出され、アナログの直交成分Q^{ａｎａｌｏｇ}である電気信号１６が出力される。

電気信号１５と電気信号１６は、信号の周波数帯域をエイリアシングなくサンプリングが可能なＡＤ変換機能素子１７ａとＡＤ変換機能素子１７ｂを備えた信号処理装置１７に送られる。信号処置装置１７では、ＡＤ変換機能素子１７ａとＡＤ変換機能素子１７ｂから出力されたデジタル化された同相成分Ｉ^{ｄｉｇｉｔａｌ}と直交成分Q^{ｄｉｇｉｔａｌ}の信号に対して、信号処理部１７ｃによって光パルス４を構成する各周波数ｆ_０＋ｆ_ｉ（ｉ＝１，２，・・・，Ｎ）のパルスによる散乱光による信号を分離する。

つまり、信号処理部１７ｃは、各周波数ｆ_０＋ｆ_ｉ成分のパルスを単独で入射した場合に得られる同相成分Ｉ_ｉ ^{ｍｅａｓｕｒｅ}と直交成分Ｑ_ｉ ^{ｍｅａｓｕｒｅ}を、全てのｉに関する同相成分の重ね合わせとなっているＩ^{ｄｉｇｉｔａｌ}と、全てのｉに関する直交成分の重ね合わせとなっているQ^{ｄｉｇｉｔａｌ}に対して信号処理を行うことで分離する。具体的な信号処理の方法は、Ｉ^{ｄｉｇｉｔａｌ}とQ^{ｄｉｇｉｔａｌ}から、Ｉ_ｉ ^{ｍｅａｓｕｒｅ}とＱ_ｉ ^{ｍｅａｓｕｒｅ}を正確に分離できるならどんな手法を用いても良い。例えば、Ｉ^{ｄｉｇｉｔａｌ}とQ^{ｄｉｇｉｔａｌ}を、中心周波数がｆ_０＋ｆ_ｉであり通過帯域が２／Ｗであるデジタルバンドパスフィルタにそれぞれ通した上で位相遅延を保証することで、Ｉ_ｉ ^{ｍｅａｓｕｒｅ}とＱ_ｉ ^{ｍｅａｓｕｒｅ}を計算する方法などが考えられる。

また、前記方法では、アナログの電気信号の状態にある同相成分と直交成分を、ＡＤ変換してデジタル化した後に、各周波数成分への分離を行っているが、アナログの電気信号の状態にある同相成分と直交成分をアナログ電気フィルタによって各周波数成分へ分離した後にＡＤ変換するなどしても良い。

信号処理部１７ｃによって取得されたＩ_ｉ ^{ｍｅａｓｕｒｅ}とＱ_ｉ ^{ｍｅａｓｕｒｅ}を元に、信号処理部１７ｄで位相の計算を行う。図１２は、信号処理部１７ｄの構造を説明する図である。信号処理部１７ｄは、
測定器３１で測定された、被測定光ファイバ６に入射した波長多重の光パルス４で発生した散乱光の同相成分と直交成分が入力される入力部２１と、
入力部２１に入力された前記散乱光の同相成分と直交成分のうち、任意の時刻且つ被測定光ファイバ６の任意位置における同相成分と直交成分により構成される２次元ベクトルを、光パルス４に多重される波長毎に取得するベクトル取得回路２２、
ベクトル取得回路２２が取得した前記２次元ベクトルのうち、基準時刻の各波長の前記２次元ベクトルを、該２次元ベクトルそれぞれが基準方向を向くように波長毎の基準回転量だけ回転し、ベクトル取得回路２２が取得した前記２次元ベクトルのうち、前記基準時刻と異なる他時刻の各波長の前記２次元ベクトルを、前記波長毎の基準回転量だけそれぞれ回転するベクトル回転回路２３、
ベクトル回転回路２３が回転した前記基準時刻の各波長の前記２次元ベクトルを加算平均して合成基準ベクトルを計算し、ベクトル回転回路２３が回転した前記他時刻の各波長の前記２次元ベクトルを加算平均して合成ベクトルを計算し、前記合成基準ベクトルと前記合成ベクトルとが成す角度から前記散乱光の位相変化量を計算する演算回路２４と、
を備える。

図１３は、本振動検出装置が行う位相測定方法を説明する図である。当該位相測定方法は、
被測定光ファイバ６に入射した波長多重の光パルス４で発生した散乱光の同相成分と直交成分を測定する測定手順Ｓ０１と、
測定手順Ｓ０１で測定した前記散乱光の同相成分と直交成分のうち、任意の時刻且つ被測定光ファイバ６の任意位置における同相成分と直交成分により構成される２次元ベクトルを、光パルス４に多重される波長毎に取得するベクトル取得手順Ｓ０２と、
ベクトル取得手順Ｓ０２で取得した前記２次元ベクトルのうち、基準時刻の各波長の前記２次元ベクトルを、該２次元ベクトルそれぞれが基準方向を向くように波長毎の基準回転量だけ回転し、ベクトル取得手順Ｓ０２で取得した前記２次元ベクトルのうち、前記基準時刻と異なる他時刻の各波長の前記２次元ベクトルを、前記波長毎の基準回転量だけそれぞれ回転するベクトル回転手順Ｓ０３と、
ベクトル回転手順Ｓ０３で回転した前記基準時刻の各波長の前記２次元ベクトルを加算平均して合成基準ベクトルを計算し、ベクトル回転手順Ｓ０３で回転した前記他時刻の各波長の前記２次元ベクトルを加算平均して合成ベクトルを計算し、前記合成基準ベクトルと前記合成ベクトルとが成す角度から前記散乱光の位相変化量を計算する演算手順Ｓ０４と、
を行う。

まず、単一周波数のみの同相成分Ｉ_ｉ ^{ｍｅａｓｕｒｅ}と直交成分Ｑ_ｉ ^{ｍｅａｓｕｒｅ}を用いて位相を計算する方法について説明し、位相の計算時のフェーディングの影響について説明する。
雑音がない時の同相成分Ｉ_ｉ（ｌ，ｎＴ）に雑音が加わった測定値がＩ_ｉ ^{ｍｅａｓｕｒｅ}（ｌ，ｎＴ）であり、雑音がない時の直交成分Ｑ_ｉ（ｌ，ｎＴ）に雑音が加わった測定値がＱ_ｉ ^{ｍｅａｓｕｒｅ}（ｌ，ｎＴ）である。つまり、同相成分と直交成分のそれぞれに重ね合わさる雑音を、Ｎ_ＩとＮ_Ｑとすれば、それらは次式で表せる。

位相は、測定した同相成分と直交成分から、信号処理装置１７により、

と計算される。雑音Ｎ_ＩとＮ_Ｑが存在するため、同相成分をｘ軸、直交成分をｙ軸としたｘｙ平面上における、ベクトル
（ｘ、ｙ）＝（Ｉ_ｉ ^{ｍｅａｓｕｒｅ}（ｌ，ｎＴ），Ｑ_ｉ ^{ｍｅａｓｕｒｅ}（ｌ，ｎＴ））
の指す位置に不確かさが生じ、ベクトルの指す向きである位相にも不確かさが生じる。例えば、雑音がない場合には、光ファイバに振動が加わっていない状態で、位相の計算値θ_ｉ ^ｃａｌ（ｌ，ｎＴ）は、各ｌについて時間変化せず一定値をとる。ところが、雑音があると、光ファイバに振動が加わっていない状態でも、位相の計算値θ_ｉ ^ｃａｌ（ｌ，ｎＴ）は、各ｌについて時間変化する。

図１４は、この様子を示した図である。ある位置ｌ_０でのベクトル
（Ｉ_ｉ ^{ｍｅａｓｕｒｅ}（ｌ_０，ｎＴ），Ｑ_ｉ ^{ｍｅａｓｕｒｅ}（ｌ_０，ｎＴ））
をｘｙ平面にプロットすると、雑音がない場合には、ベクトル２０１で示した
（ｘ，ｙ）＝（Ｉ_ｉ（ｌ，ｎＴ），Ｑ_ｉ（ｌ，ｎＴ））
で常に一定であり、その角度２０２も時間変化しない。しかし、実際には雑音が存在するため、各時刻での同相成分と直交成分から構成されるベクトルは、ベクトル２０３のように、ベクトル２０１とは異なり、その角度２０４も角度２０２とは異なる。そのため、実際の測定値から構成されるベクトルは、ベクトル２０１を中心としてばらつく。ばらつきの程度は、各軸方向での測定値の標準偏差を用いて評価できる。例えば、ｘ軸方向については、測定値のｘ成分のばらつきから、２０５で表した雑音Ｎ_Ｉの標準偏差σ（Ｎ_Ｉ）だけ不確かさがある。

コヒーレント検波の場合、参照光の強度を十分に大きくするため、ショット雑音が支配的となり、雑音の分布も正規分布で近似できる。また、図１１における二つのバランス検出器１３と１４に入射する光強度は同程度とみなせるため、雑音Ｎ_ＩとＮ_Ｑの標準偏差も同じ大きさとみなすことができ、不確かさはベクトル２０１を中心とした円形となる。

参照光のショット雑音以外の、ＰＤの熱雑音等の雑音が無視できない、例えばコヒーレント検波ではなくて直接検波を実施するような装置構成の場合にも、二つのバランス検出器１３と１４の雑音特性が同じであるとみなせるため、不確かさはベクトル２０１を中心とした円形と考えてよい。

ただし、直接検波を実施するような装置構成において、散乱光のショット雑音を考慮する必要がある場合などは、不確かさの程度は散乱光強度に依存して各地点で異なってくるが、散乱光強度が小さい点での不確かさは、ＰＤの熱雑音等の電気段以降での測定器の雑音が支配的となるため、以下で説明するフェーディングの現象による影響については、直接検波を実施するような装置構成においても成立する。

フェーディングの現象により、散乱光強度が小さくなる地点が発生する。そのような地点については、位相計算時の不確かさが増加してしまうため、小さな振動の検出が困難となる。特に、ＳＮ比が１未満となるような、雑音がない時の散乱光の振幅が図１５のベクトル２０６で示したように小さい場合には、測定されるベクトルは、ベクトル２０７で示したように、雑音がない時のベクトル２０６とは、大きく異なる値をとる確率が増加してしまい、振動の誤検知に繋がる。また、そのような点では、引き続く式（３）で示したアンラップ処理を行う場合にも、整数ｑの選択を誤る確率が増加するため、特に大きな振動が加わったとする誤検知に繋がってしまう。

単一波長で実験を行った際のフェーディングによる散乱光強度Ｐのばらつきの分布Ｄ（Ｐ）は、散乱光の強度の平均値を＜Ｐ＞として、以下の式を満たすことが知られている。

この式のように散乱光強度Ｐのばらつきの分布Ｄ（Ｐ）は、散乱光強度Ｐが小さい程大きくなる。したがって、単一波長のパルス強度を大きくすることで、散乱光強度の小さい点を削減しようとする場合、非常に大きなピーク強度が必要となるため、非線形効果などのパルス歪を考えると、限界がある。

そこで、信号処理部１７ｄは、以下に説明するように、ｉ＝１，２，・・・，Ｎの異なるＮ個の周波数における同相成分Ｉ_ｉ ^{ｍｅａｓｕｒｅ}と直交成分Ｑ_ｉ ^{ｍｅａｓｕｒｅ}を用いて位相計算を行い、フェーディングによる散乱光強度が小さい点での位相の不確かさの増加を防止する。

［測定手順Ｓ０１］
図１１で説明した測定系を利用して被測定光ファイバ６に入射した波長多重の光パルス４で発生した散乱光の同相成分と直交成分を測定する。

［ベクトル取得手順Ｓ０２］
測定手順Ｓ０１で測定した前記散乱光の同相成分と直交成分のうち、任意の時刻且つ被測定光ファイバ６の任意位置における同相成分と直交成分により構成される２次元ベクトルを、光パルス４に多重される波長毎に取得する。

［ベクトル回転手順Ｓ０３］
まず、時刻ゼロの時の測定値のベクトル
（Ｉ_ｉ ^{ｍｅａｓｕｒｅ}（ｌ，０），Ｑ_ｉ ^{ｍｅａｓｕｒｅ}（ｌ，０））
から位相θ_ｉ ^ｃａｌ（ｌ，０）を計算する。続いて、計算した位相値θ_ｉ ^ｃａｌ（ｌ，０）とは逆向きの回転量で、各時刻のベクトル
（Ｉ_ｉ ^{ｍｅａｓｕｒｅ}（ｌ，ｎＴ），Ｑ_ｉ ^{ｍｅａｓｕｒｅ}（ｌ，ｎＴ））
を回転させることで、各時刻及び各地点における新しいベクトルを式（Ｆ８）のように計算する。

［演算手順Ｓ０４］
そして、各波長に関する新しく計算したベクトルを式（Ｆ９）のように加算平均して、位相計算に直接用いるベクトルを計算する。

最後に、ベクトル（Ｉ^ｎｅｗ（ｌ，ｎＴ），Ｑ^ｎｅｗ（ｌ，ｎＴ））から式（Ｆ１０）のように位相θ^ｃａｌ（ｌ，ｎＴ）を計算する。

ベクトル（Ｉ^ｎｅｗ（ｌ，ｎＴ），Ｑ^ｎｅｗ（ｌ，ｎＴ））を使用してθ^ｃａｌ（ｌ，ｎＴ）を計算することで、フェーディングにより散乱光強度が減少する地点を減らすことが可能になる。以下に原理を述べる。

ｉ＝１，２，・・・，Ｎの異なるＮ個の周波数におけるθ_ｉ ^ｃａｌ（ｌ，ｎＴ）は互いに異なる値をとっている。例えば、Ｎ＝２を例にとれば、ｉ＝１とｉ＝２の時刻ゼロでの雑音がないときのベクトル（Ｉ_ｉ（ｌ，０），Ｑ_ｉ（ｌ，０））は、図１６のベクトル３０１とベクトル３０２のように、向きも大きさも異なる。振動によって、時刻ｎＴにおける入射端から距離ｌの地点よりも手前におけるファイバの正味の伸縮量が、時刻ゼロと比較して変化しているならば、ｉ＝１とｉ＝２の時刻ｎＴでのベクトル（Ｉ_ｉ（ｌ，ｎＴ），Ｑ_ｉ（ｌ，ｎＴ））は、図１６のベクトル３０３とベクトル３０４のようにそれぞれ変化している。ベクトル３０３とベクトル３０４の長さは、ベクトル３０１とべクトル３０２に対してそれぞれ変化し、その変化量はｉ＝１とｉ＝２で異なるが、ベクトル３０３とベクトル３０４の向きは、ベクトル３０１とベクトル３０２に対して、それぞれ同じ量だけ変化する。

つまり、角度３０５と角度３０６は同一であり、この量が式（１）のθ（ｌ，ｎＴ）に対応する。雑音がない場合には、式（Ｆ８）により、ベクトル３０１、ベクトル３０２、ベクトル３０３及びベクトル３０４は、図１７のようにそれぞれベクトル３０７、ベクトル３０８、ベクトル３０９及びベクトル３１０に移される。つまり、時刻ごとに全波長のベクトルの方向を揃える。

ベクトル３０７とベクトル３０８の平均として（Ｉ^ｎｅｗ（ｌ，０），Ｑ^ｎｅｗ（ｌ，０））がベクトル３１１として得られ、ベクトル３０９とベクトル３１０の平均として（Ｉ^ｎｅｗ（ｌ，ｎＴ），Ｑ^ｎｅｗ（ｌ，ｎＴ））がベクトル３１２として得られる。時刻０からｎＴの間の位相の変化量は角度３１３であるが、これは角度３０５や角度３０６と同一である。

実際の測定では、雑音によりベクトル３０１、ベクトル３０２、ベクトル３０３、ベクトル３０４、ベクトル３１１、及びベクトル３１２の向きに不確かさが存在し、結果として、角度３０５や角度３０６にも不確かさが伴う。しかし、加算平均したベクトル３１１とベクトル３１２を用いて角度３１３を計算することで、不確かさを減少させることができる。その理由は２つある。

１つの理由は、ベクトル３０７とベクトル３０８の振幅の２乗に対応する強度の確率分布が、独立に式（Ｆ７）に従うからである。このため、ベクトル３１１の振幅の２乗に対応する強度の確率分布は、理論的には図１８の分布４０２のようになり、平均値自体は同じでも、単一波長の場合にフェーディングにより散乱光強度が著しく小さくなる点を取り除くことが可能となる。比較のために、図１８には、一つの波長のみの場合の確率分布を分布４０１として示している。

他の理由は、ベクトル３１１では、ベクトル３０７とベクトル３０８を平均化することにより、雑音のレベルはベクトル３０１やベクトル３０２に対して１／√２となるからである。このため、ベクトル３１１の長さの平均値自体は、ベクトル３０１やベクトル３０２と変わらずとも、雑音レベルが小さくなることで、位相計算時の不確かさを減少させることが可能となる。ベクトル３１２についても同様である。

ここではＮが２つの周波数の場合について、本提案の効果を具体的に記述したが、一般化することは可能である。まず、多重数Ｎが大きくなる程、散乱光強度がゼロに近くなる地点の数は減少する。この様子を、Ｎ＝５の場合を分布４０３に、Ｎ＝１０の場合を分布４０４に示した。また、雑音レベルの大きさも１／√Ｎ倍となるため、同じ平均強度でも、Ｎが大きくなる程、位相計算時の不確かさは減少する。

なお、本実施形態で説明した信号処理方法は、単純にθ_ｉ ^ｃａｌ（ｌ，ｎＴ）を異なるｉについて平均化する方法、例えばＮ＝２の場合で位相３０５と位相３０６の平均を計算する方法、とは異なる。単純にθ_ｉ ^ｃａｌ（ｌ，ｎＴ）を異なるｉについて平均化する方法では、θ_ｉ ^ｃａｌ（ｌ，ｎＴ）自体は、単一波長で計算されるため、フェーディングによる散乱光強度の小さい所では、雑音がない場合と比較すると大きく測定値が異なる。このため振動の誤検知の発生を減少させることはできない。位相の平均化で雑音がない場合の理想的な位相値と測定値との相違を減少させることはできるが、当該相違の発生頻度自体は、散乱光強度自体が小さくなる地点が波長ごとに違うため増加することになる。つまり、単純にθ_ｉ ^ｃａｌ（ｌ，ｎＴ）を異なるｉについて平均化してもフェーディングによる散乱光強度の小さい点の除去にはならない。

尚、説明のため式（Ｆ９）では右辺に１／Ｎを乗じているが、これを乗じなくても式（Ｆ１０）で計算される位相値は変化しないため、実際の計算では１／Ｎを乗じなくても良い。

（例２）
例１では、式（Ｆ８）のベクトル回転の回転角度を時刻ゼロでのθ_ｉ ^ｃａｌ（ｌ，０）とした場合の例を説明した。θ_ｉ ^ｃａｌ（ｌ，０）を算出する際に使用する同相成分と直交成分は、（Ｉ_ｉ ^{ｍｅａｓｕｒｅ}（ｌ，０），Ｑ_ｉ ^{ｍｅａｓｕｒｅ}（ｌ，０））であり雑音の影響を含んだベクトルとなっている。このため、θ_ｉ ^ｃａｌ（ｌ，０）も雑音の影響を受けている。もし、時刻ゼロにおいて、θ_ｉ ^ｃａｌ（ｌ，０）の値が雑音がない場合のθ_ｉ（ｌ，０）の値とが大きく異なっている場合、式（Ｆ８）でのベクトル回転の効果が得られなくなることがある。

図１５で説明する。例えば、ファイバ長手方向のある地点、且つある周波数における、時刻ゼロでの散乱光のベクトルが、雑音がない場合はベクトル２０６となるものが、雑音によりベクトル２０７となっていたとする。式（Ｆ８）でベクトルを回転させた上で式（Ｆ９）でベクトルを加算平均する演算は、ベクトル長さがゼロに近い点（フェーディングの影響を受けている点）を減らす効果がある。しかし、この例の場合、ベクトル２０７の角度を基準としてその周波数の回転角度を決定するため、その周波数については当該効果を得られなくなる。

また、回転角度の誤りは、測定する振動の大きさを正しく評価できない場合も発生する。当該場合を図１６で説明する。例えば、ベクトル３０１と３０３（波長１）及びベクトル３０２と３０４（波長２）が雑音の影響がないベクトルであったとすると、波長の異なるベクトルを一致させる回転角度はゼロとはならない。しかし、基準とした時刻ゼロにおいて、雑音の影響により波長１と波長２のベクトルの向きがたまたま同一であった場合、実施形態１の手法では波長２のベクトルを回転角度ゼロ、つまり無回転のままで、波長１のベクトルと加算平均することになる。つまり、ベクトル３０１とベクトル３０２（無回転）を加算平均したベクトルを時刻ゼロにおける平均ベクトルとし、ベクトル３０３とベクトル３０４（無回転）を加算平均したベクトルを時刻ｎＴにおける平均ベクトルとする。このため、時刻ゼロから時刻ｎＴへの平均ベクトルの角度変化は、もはやベクトル３０１からベクトル３０３への角度変化や、ベクトル３０２からベクトル３０４への角度変化とは一致せず、正しく位相変化を捉えることができない。

本実施形態では、上記のような不具合の発生頻度を低減する手法を説明する。
図１９は、本実施形態の振動検出装置の信号処理部１７ｄの構造を説明する図である。本実施形態の信号処理部１７ｄは、
測定器で測定された、被測定光ファイバに入射した波長多重の光パルスで発生した散乱光の同相成分と直交成分が入力される入力部２１と、
入力部２１に入力された前記散乱光の同相成分と直交成分のうち、任意の時刻且つ被測定光ファイバ６の任意位置における同相成分と直交成分により構成される２次元ベクトルを、前記光パルスに多重される波長毎に取得するベクトル取得回路２２と、
ベクトル取得回路２２が取得した前記２次元ベクトルのうち、基準波長の各時刻の前記２次元ベクトルを、該２次元ベクトルそれぞれが基準方向を向くように時刻毎の基準回転量だけ回転し、ベクトル取得回路２２が取得した前記２次元ベクトルのうち、前記基準波長と異なる他波長の各時刻の前記２次元ベクトルを、前記時刻毎の基準回転量だけそれぞれ回転する第１ベクトル回転回路２３－１と、
第１ベクトル回転回路２３－１が回転した前記基準波長の各時刻の前記２次元ベクトルを加算平均して第１合成基準ベクトルを計算し、第１ベクトル回転回路２３－１が回転した前記他波長の各時刻の前記２次元ベクトルを加算平均して波長毎の第１合成ベクトルを計算し、前記第１合成基準ベクトルと前記第１合成ベクトルとが成す角度から波長毎の基準回転量を計算する第１演算回路２４－１と、
ベクトル取得回路２２が取得した前記２次元ベクトルのうち、基準時刻の各波長の前記２次元ベクトルを、第１演算回路２４－１が計算した前記波長毎の基準回転量だけ回転し、ベクトル取得回路２２が取得した前記２次元ベクトルのうち、前記基準時刻と異なる他時刻の各波長の前記２次元ベクトルを、第１演算回路２４－１が計算した前記波長毎の基準回転量だけそれぞれ回転する第２ベクトル回転回路２３－２と、
第２ベクトル回転回路２３－２が回転した前記基準時刻の各波長の前記２次元ベクトルを加算平均して第２合成基準ベクトルを計算し、第２ベクトル回転回路２３－２が回転した前記他時刻の各波長の前記２次元ベクトルを加算平均して第２合成ベクトルを計算し、前記第２合成基準ベクトルと前記第２合成ベクトルとが成す角度から前記散乱光の位相変化量を計算する第２演算回路２４－２と、
を備える。

図２０は、本振動検出装置が行う位相測定方法を説明する図である。当該位相測定方法は、
被測定光ファイバに入射した波長多重の光パルスで発生した散乱光の同相成分と直交成分を測定する測定手順Ｓ０１と、
測定手順Ｓ０１で測定した前記散乱光の同相成分と直交成分のうち、任意の時刻且つ前記被測定光ファイバの任意位置における同相成分と直交成分により構成される２次元ベクトルを、前記光パルスに多重される波長毎に取得するベクトル取得手順Ｓ０２と、
ベクトル取得手順Ｓ０２で取得した前記２次元ベクトルのうち、基準波長の各時刻の前記２次元ベクトルを、該２次元ベクトルそれぞれが基準方向を向くように時刻毎の基準回転量だけ回転し、ベクトル取得手順Ｓ０２で取得した前記２次元ベクトルのうち、前記基準波長と異なる他波長の各時刻の前記２次元ベクトルを、前記時刻毎の基準回転量だけそれぞれ回転する第１ベクトル回転手順Ｓ１３と、
第１ベクトル回転手順Ｓ１３で回転した前記基準波長の各時刻の前記２次元ベクトルを加算平均して第１合成基準ベクトルを計算し、第１ベクトル回転手順Ｓ１３で回転した前記他波長の各時刻の前記２次元ベクトルを加算平均して波長毎の第１合成ベクトルを計算し、前記第１合成基準ベクトルと前記第１合成ベクトルとが成す角度から波長毎の基準回転量を計算する第１演算手順Ｓ１４と、
ベクトル取得手順Ｓ０２で取得した前記２次元ベクトルのうち、基準時刻の各波長の前記２次元ベクトルを、第１演算手順Ｓ１４で計算した前記波長毎の基準回転量だけ回転し、ベクトル取得手順Ｓ０２で取得した前記２次元ベクトルのうち、前記基準時刻と異なる他時刻の各波長の前記２次元ベクトルを、第１演算手順Ｓ１４で計算した前記波長毎の基準回転量だけそれぞれ回転する第２ベクトル回転手順Ｓ１５と、
第２ベクトル回転手順Ｓ１５で回転した前記基準時刻の各波長の前記２次元ベクトルを加算平均して第２合成基準ベクトルを計算し、第２ベクトル回転手順Ｓ１５で回転した前記他時刻の各波長の前記２次元ベクトルを加算平均して第２合成ベクトルを計算し、前記第２合成基準ベクトルと前記第２合成ベクトルとが成す角度から前記散乱光の位相変化量を計算する第２演算手順Ｓ１６と、
を行う。

入力部２１、ベクトル取得回路２２、測定手順Ｓ０１及びベクトル取得手順Ｓ０２は、実施形態１の説明と同じである。
まず、第１ベクトル回転回路２３－１は次の第１ベクトル回転手順Ｓ１３を行う。
ｉ＝１を基準（基準波長）にして、全てのｉについて下式を計算する。

続いて、第１演算回路２３－１は次の第１演算手順Ｓ１４を行う。
式（Ｆ１１）で得られたｒ’（ｌ，ｎＴ）の全ての時刻について加算平均を下式で計算する。

Ｍは時刻方向でのサンプル数である。式（Ｆ１２）で得られたｒ_ｉ’’（ｌ）を用いて、下式を計算する。

簡潔に説明すると次のようになる。
まず、基準波長ｉ＝１について、各時刻でのベクトルがＩ軸を向くように回転させる。この回転角をθ_ｉ＝１（ｔ）とする。θ_ｉ＝１（ｔ）は時刻で変化する。また、回転後のベクトルをＡｔとする。
次に、他の波長について、各時刻でのベクトルをθ_ｉ＝１（ｔ）で回転する。回転後の波長毎のベクトルをＢｔ_ｉとする。
続いて、波長毎に回転後のベクトルを時間平均する。基準波長についての時間平均後のベクトルをΣＡｔとする。他の波長についての時間平均後のベクトルをΣＢｔ_ｉとする。
そして、ΣＡｔとΣＢｔ_ｉとの成す角度を求める。この角度が式（Ｆ１３）式のθ_ｉ’（ｌ）である。θ_ｉ’（ｌ）は他の波長毎に存在する。

第２ベクトル回転回路２３－２が行う第２ベクトル回転手順Ｓ１５、及び第２演算回路２４－２が行う第２演算手順Ｓ１６は、次の点を除いて実施形態１で説明したベクトル回転手順Ｓ０３及び演算手順Ｓ０４とそれぞれ同じである。式（Ｆ８）の計算時に、θ_ｉ ^ｃａｌ（ｌ，０）に替えて、式（Ｆ１３）式のθ_ｉ’（ｌ）を用いてベクトル（Ｉ^ｎｅｗ（ｌ，ｎＴ），Ｑ^ｎｅｗ（ｌ，ｎＴ））を計算し、式（Ｆ９）で位相の計算を行う。

本手法の意義を説明する。振動が起きても異なる周波数同士のベクトルのなす角度は、雑音がない場合には変化しない。例えば、２つの周波数の場合を例にとれば、基準時刻のベクトル３０１とベクトル３０２のなす角度と、時間ｎＴ後のベクトル３０３とベクトル３０４のなす角度は、雑音の影響がないならば同一となる。そこで、式（Ｆ１１）によるベクトルの回転を行った上で、式（Ｆ１２）でベクトルの加算平均を行い、基準時刻の加算平均ベクトルと時間ｎＴ後の加算平均ベクトルとが成す角度を回転角度θ_ｉ’（ｌ）とする。これにより、各ベクトルに付随する雑音の大きさを１／√Ｍにすることができる。したがって、Ｍを十分に大きくとることで雑音の影響を減らすことができる。

なお、本実施形態の手法は、異なる周波数同士のベクトルの位相差を各時刻で求めた上で、全ての時刻についての位相差を平均することで回転角度を求める手法（以下、比較手法と記載）と異なる。本実施形態の手法は、予め全ての波長のベクトルを波長１の回転角度で回転させている点で異なる。例えば、周波数が２つの場合に、比較手法は、単純にベクトル３０２とベクトル３０１の角度差、及びベクトル３０３とベクトル３０２の角度差を計算し、それらの平均値として回転角度を求める。比較手法は、各時刻での位相差の計算値が雑音がない場合と大きく異なる点（被測定ファイバの計測位置）の発生確率は変化しないので、式（Ｆ８）による回転の効果が十分でない。

尚、説明のため式（Ｆ１２）では１／Ｍを右辺に乗じているが、これを乗じなくても式（Ｆ１３）で計算される回転角度は変わらないので、実際の計算においては１／Ｍを乗じなくても良い。

例２の手法は、例１の手法と比べて、計算時間は増加するが、最終的な位相計算における不確かさを、例１より減少させることができる。

１：光源
２：カプラ
３：光変調器
４：光パルス
４ａ：微小パルス
５：サーキュレータ
６：被測定光ファイバ
７：９０度光ハイブリッド
８、９：カプラ
１０：位相シフタ
１１、１２：カプラ
１３、１４：バランス検出器
１５：アナログの同相成分の電気信号
１６：アナログの直交成分の電気信号
１７：信号処理装置
１７ａ、１７ｂ：ＡＤ変換素子
１７ｃ、１７ｄ：信号処理部
１８：計算部
２１：入力部
２２：ベクトル取得回路
２３：ベクトル回転回路
２４：演算回路
２３－１：第１ベクトル回転回路
２３－２：第２ベクトル回転回路
２４－１：第１演算回路
２４－２：第２演算回路
３１：測定器

Claims (8)

1. 光パルス試験方法であって、
   異なる光周波数の光パルスを時間的に等間隔で並べた波長多重数Ｎの光パルス列を被測定光ファイバの一端に入射すること、
   前記被測定光ファイバの前記一端に戻ってきた各波長の散乱光を受光すること、
   前記被測定光ファイバの振動を前記散乱光の位相成分の時間変化として観測すること、及び
   前記波長多重数Ｎを数Ｃ１を満たす最小値に決定すること、
   を特徴とする光パルス試験方法。
   

   
2. 前記光パルスを、前記光パルス列の光パルス間隔より短い間隔で異なる光周波数の微小光パルスをＭ個並べて形成すること、及び
   雑音レベルの１／√Ｍ倍が所定値以下となるように前記微小光パルスの数Ｍを設定すること
   を特徴とする請求項１に記載の光パルス試験方法。
3. 光パルス試験方法であって、
   Ｎ×Ｍ＋１個（ＮとＭは自然数）の異なる光周波数の微小光パルスを並べた集団をＭ個形成すること、
   Ｍ個の前記集団を並べて、先頭の前記集団からＭ個づつ前記微小光パルスを切り出して光パルス対をＮ×Ｍ＋１個形成すること、
   Ｎ×Ｍ＋１個の前記光パルス対を時間的に等間隔で並べた光パルス列を被測定光ファイバの一端に入射すること、
   前記被測定光ファイバの前記一端に戻ってきた各波長の散乱光を受光すること、
   前記被測定光ファイバの振動を前記散乱光の位相成分の時間変化として観測すること、
   波長多重数Ｎを数Ｃ１を満たす最小値に決定すること、及び
   雑音レベルの１／√Ｍ倍が所定値以下となるように波長多重数Ｍを設定すること
   を特徴とする光パルス試験方法。
   

   
4. 光パルス試験方法であって、
   Ｎ×Ｍ個（ＮとＭは自然数）の異なる光周波数の微小光パルスを並べた集団をＮ＋１個形成すること、
   Ｎ＋１個の前記集団を並べて、先頭の前記集団からＭ個づつ前記微小光パルスを切り出して光パルス対をＮ（Ｎ＋１）個形成すること、
   Ｎ（Ｎ＋１）個の前記光パルス対のうち、Ｎ＋１毎の前記光パルス対にいずれの前記微小光パルスの光周波数とも異なる光周波数の追加微小光パルスを追加すること、
   前記追加微小光パルスが追加された前記光パルス対を含むＮ（Ｎ＋１）個の前記光パルス対を時間的に等間隔で並べた光パルス列を被測定光ファイバの一端に入射すること、
   前記被測定光ファイバの前記一端に戻ってきた各波長の散乱光を受光すること、
   前記被測定光ファイバの振動を前記散乱光の位相成分の時間変化として観測すること、
   波長多重数Ｎを数Ｃ１を満たす最小値に決定すること、及び
   雑音レベルの１／√Ｍ倍が所定値以下となるように波長多重数Ｍを設定すること
   を特徴とする光パルス試験方法。
   

   
5. 光パルス試験方法であって、
   Ｎ×Ｍ＋１個（ＮとＭは自然数）の異なる光周波数の微小光パルスを並べた集団を（Ｎ＋１）Ｍ＋１個形成すること、
   （Ｎ＋１）Ｍ＋１個の前記集団を並べ、１＋ｋ（Ｎ＋１）番目のパルス対にＭ＋１個の前記微小光パルスが含まれ、他の前記パルス対にＭ個の前記微小光パルスが含まれるように、先頭の前記集団から前記微小光パルスを切り出して前記光パルス対を（Ｎ×Ｍ＋１）（Ｎ＋１）個形成すること、
   （Ｎ×Ｍ＋１）（Ｎ＋１）個の前記光パルス対を時間的に等間隔で並べた光パルス列を被測定光ファイバの一端に入射すること、
   前記被測定光ファイバの前記一端に戻ってきた各波長の散乱光を受光すること、
   前記被測定光ファイバの振動を前記散乱光の位相成分の時間変化として観測すること、
   波長多重数Ｎを数Ｃ１を満たす最小値に決定すること、及び
   雑音レベルの１／√Ｍ倍が所定値以下となるように波長多重数Ｍを設定すること
   を特徴とする光パルス試験方法。
   

   
6. 前記光パルス列を前記被測定光ファイバの一端に入射する前に、前記光パルス列の光パルス間隔より短い間隔で異なる光周波数の前記光パルスを並べた光パルス対を並べた補正用光パルス列を前記被測定光ファイバの一端に入射すること、及び
   前記散乱光の位相成分の補正値を数Ｃ２で計算すること
   を特徴とする請求項１から５のいずれかに記載する光パルス試験方法。
   

   

   ただし、φ（ｉ（ｍ’），ｉ（ｍ））は光周波数が異なる２つの前記光パルス間の角度差、
   Ｐは前記補正用光パルス列に含まれる前記光パルス対の数、
   ｐは前記光パルス対の番号（１からＰまでの整数）、
   ｔ_ｐはｐ番目の前記光パルス対が前記被測定光ファイバの一端に入射される時間、
   ｚは前記被測定光ファイバの前記一端からの距離、
   ｒ_ｉ（ｍ）（ｚ，ｔ_ｐ＋２ｚ／ν）及びｒ_{ｉ（ｍ’）}（ｚ，ｔ_ｐ＋２ｚ／ν）は距離zからの散乱光の複素ベクトル、
   ａｒｇ関数は引数の複素ベクトルの偏角を－πからπの範囲で計算し、実数を出力する関数、
   Ｒ（＊）は複素平面上で角度＊だけ時計まわりに複素ベクトルｒを回転させる演算子である。
7. 光パルス試験装置であって、
   異なる光周波数の光パルスを時間的に等間隔で並べた波長多重数Ｎの光パルス列を被測定光ファイバの一端に入射する光源と、
   前記被測定光ファイバの前記一端に戻ってきた各波長の散乱光を受光する受光器と、
   前記被測定光ファイバの振動を前記散乱光の位相成分の時間変化として観測する信号処理部と、
   前記波長多重数Ｎを数Ｃ１を満たす最小値に決定する計算部と、
   を備えることを特徴とする光パルス試験装置。
   

   
8. 前記光源は、前記光パルス列を前記被測定光ファイバの一端に入射する前に、前記光パルス列の光パルス間隔より短い間隔で異なる光周波数の前記光パルスを並べた光パルス対を並べた補正用光パルス列を前記被測定光ファイバの一端に入射し、
   前記信号処理部は、前記散乱光の位相成分の補正値を数Ｃ２で計算すること
   を特徴とする請求項７に記載する光パルス試験装置。
   

   

   ただし、φ（ｉ（ｍ’），ｉ（ｍ））は光周波数が異なる２つの前記光パルス間の角度差、
   Ｐは前記補正用光パルス列に含まれる前記光パルス対の数、
   ｐは前記光パルス対の番号（１からＰまでの整数）、
   ｔ_ｐはｐ番目の前記光パルス対が前記被測定光ファイバの一端に入射される時間、
   ｚは前記被測定光ファイバの前記一端からの距離、
   ｒ_ｉ（ｍ）（ｚ，ｔ_ｐ＋２ｚ／ν）及びｒ_{ｉ（ｍ’）}（ｚ，ｔ_ｐ＋２ｚ／ν）は距離zからの散乱光の複素ベクトル、
   ａｒｇ関数は引数の複素ベクトルの偏角を－πからπの範囲で計算し、実数を出力する関数、
   Ｒ（＊）は複素平面上で角度＊だけ時計まわりに複素ベクトルｒを回転させる演算子である。

Images