JP3766975B1 - パラメトリックな時間引き延ばしパルス生成装置 - Google Patents

Abstract

【課題】
従来のＴＳＰを用いたインパルス応答の測定では、特に暗騒音が多く、歪み成分と本来のＴＳＰ入力に対する応答成分が時間的に接近する低周波において、十分なＳ／Ｎ比を達成することが困難だった。また測定媒体の周波数特性を逆補正する再にもＳ／Ｎ比の低下を許していた。
【解決手段】
このため掃出線のグラフを任意に設計可能なＴＳＰ波形を測定に用い、畳み込みにより逆フィルタを施すことで、測定環境及び測定媒体の特性を考慮でき、より一層のＳ／Ｎ比向上が期待できる。
本発明ではこのＴＳＰ波形の生成手段を開発した。
【選択図】図１

Description

本発明は、系の未知の動特性、例えば音楽ホールの各場所に於ける音響特性を測定するために行われるインパルス応答の測定方法に係り、特に高精度のインパルス応答を測定するために用いる測定用入力波形とその逆フィルタの生成に関する。

インパルス応答のディジタルによる測定方法として特許文献１で示された標本化周波数の２分の１（ナイキスト周波数）付近において位相特性が滑らかに実数となるように最適化された時間引き延ばしパルス（ＯＡＴＳＰ）を測定用入力信号として用い、被測定系からの出力を、該入力信号の逆フィルタにより直線畳み込み、あるいは巡回畳み込みを行うことでインパルス応答を復元する方法が提唱されている。
また、時間引き延ばしパルスの掃出線（群遅延特性）を変更する試みとして、非特許文献１などが考案されている。
また、ＯＡＴＳＰを用いた測定では、しばしば入力波形を連続で再生し、位相が完全に一致する周期で被測定系からの出力を平均化することで、暗騒音に対して差別的に入力波形に関連した被測定系からの出力を高める同期加算という手法が平行して用いられる。

また、インパルス応答は通常ラウドスピーカーとマイクを介して測定されるため、例えば音楽ホールのインパルス応答を測定するのに従来のＴＳＰによる測定方法を単独で用いる場合、係るスピーカとマイク固有の周波数特性までもが同時に畳み込まれた結果としてのインパルス応答が得られるのみであり、音楽ホールの純粋なインパルス応答は得られない。この問題を回避あるいは軽減する措置として、予め無響室などで同じスピーカとマイク自体のインパルス応答を求めておき、この逆特性を持つフィルタにて、測定後のインパルス応答を畳み込むことで、係るスピーカとマイクの周波数特性を相殺する方法が用いられている。
特許２７２５８３８ 日本音響学会講演論文集（１９９９年９月−１０月，４３３−４３４）

ＯＡＴＳＰの掃出線は周波数が時間に対して直線的に変化する正弦波のようなものなので、被測定系においてしばしば重要な意味を持ち、より高精度な測定が求められる低周波の占める時間的割合は小さい。
例えば人間が音色や音程を判断するのに使う周波数帯域は、概ね５ＫＨｚ程度迄であり、それよりも高い周波数帯域は、高音質感を得るために必要であるが、知覚される音声の中身としての意味が薄い帯域である。
しかし、従来のＯＡＴＳＰでは、例えばサンプリング周波数を４４．１ＫＨｚとした場合、有効周波数帯域の上限であるナイキスト周波数は２２．０５ＫＨｚとなるため、５ＫＨｚから下の帯域が占める時間割合はおよそ４分の１程度でしかない。
人間の感覚として重要度の低い高周波に、多くの時間（言い換えればエネルギー）を使ってしまうため、人間の感覚としてより重要と言える低周波は、相対的にノイズの影響を受け易くなってしまう。
また、自然環境では発音源に対してあらゆる物体、空気自体もが低域濾過器として働く傾向を持つため暗騒音のスペクトラムも通常は低周波が優勢である。よって時間割合が小さい低周波は、一層ノイズの影響を受け易くなる。
さらに被測定系は、通常好ましくない非線形効果（歪み）を含んでおり、該歪みによって生じるノイズ（歪み成分）のスペクトラムは、特定時刻において正弦波と見立てたＴＳＰ波形の整数倍の周波数に集中する傾向にある。該歪み成分は、ＴＳＰの入力に伴って発生するため、同期加算を用いても影響を軽減することは適わない。
このため例えば余韻時間が非常に長い音楽ホールの測定では、図２に示すように該歪み成分の余韻が、本来のインパルス応答の成分に混入してしまう。
この現象自体は不可避であり、しかも特に低周波において、本来のインパルス応答の成分と該歪み成分とが時間的に接近してしまうため、より多くの成分が混入してしまう。
さらに、測定環境に於ける暗騒音のスペクトラムは、騒音源の種類や共鳴の具合により測定環境毎に変化するのが普通である。
また、前記したＴＳＰと測定媒体の逆フィルタを組み合た測定方法では、測定媒体の、特に再生周波数域限界付近のパワーの落ち込みを補う際、相対的な増幅率が大きくなるため、係る帯域に存在する騒音も同時に増幅されてしまい、Ｓ／Ｎ比が顕著に低下してしまう。また、ＴＳＰの入力からインパルス応答を復元するための畳み込みと、測定媒体の逆フィルタによる畳み込みの２回の畳み込みが必要となるため、計算量のロスが大きくなってしまう。

線形時間不変な連続時間の被測定系に対する非インパルスの任意波形入力信号と該入力信号の逆フィルタとを用いて上記被測定系の離散的インパルス応答を測定する方法において、振幅特性、群遅延特性、位相特性を数値積分、及び数値微分により関連付けることで、包絡線の直線性を維持しつつ、任意のパワースペクトラムを有しながらもその逆フィルタにより畳み込みを実施することでインパルス応答を厳密に再現できる時間引き延ばしパルス（ＴＳＰ）とその逆フィルタの波形を生成する装置。
及び上記測定方法において、初期時刻に周波数０を通り、単調増加し、連続な任意の時刻をパラメータとした周波数の関数によって表現される時刻−周波数グラフを掃出線として持つ時間引き延ばしパルス（ＴＳＰ）を入力信号として用いながらも、同時に生成される逆フィルタの波形により畳み込みを実施することでインパルス応答を厳密に再現できる該入力信号及び、逆フィルタの波形の生成装置。
及び上記測定方法において、該測定時とは別の時点で用意されたインパルス応答Ｂに基づいて群遅延特性及び振幅特性が補正されたＴＳＰ信号を被測定系の入力信号として用いながらも、同時に生成されるフィルタ信号により被測定系からの出力信号を畳み込みすることにより、被測定系からＢの周波数特性を取り除いたインパルス応答を復元可能な、該入力信号、及び畳み込みに用いる信号の生成装置。
及び上記のＴＳＰ生成処理において、任意設計された複数の群遅延特性を振幅特性レベルの乗算により合成する装置。
及び上記ＴＳＰ生成処理において、ＴＳＰのパワースペクトラムを、該測定場所にて事前に測定した暗騒音のパワースペクトラムに基づき、同等の特性となるように変更する装置。

時間に対して周波数が直線的に変化していたＯＡＴＳＰの掃出線に比べ、掃出線の変化グラフを自由に設計できる本発明に於けるＴＳＰは、被測定系の好まざる特性を考慮した測定用入力波形とその逆フィルタの波形を生成し、これを測定に用いることにより、一層的確なインパルス応答の測定が可能となる。例えば歪みを有する被測定系において、図３に示すように、低周波帯域において、歪み成分の混入を最小限に抑えることができ、インパルス応答の厳密な復元能も維持される。
さらに測定媒体を含む被測定系のインパルス応答の復元後に、測定媒体の逆フィルタを畳み込む場合に比べ、ＴＳＰのパワースペクトラムの変更により周波数帯域毎の補正量が平均化されるため、特定の帯域に於けるＳ／Ｎ比の低下を抑えることができる。
また測定媒体の逆フィルタ特性がＴＳＰからインパルスを復元するフィルタに組み込まれるため、測定後の畳み込み回数を２回から１回に減らすことができる。
さらに任意設計されたＴＳＰの群遅延特性及び振幅特性と上記逆フィルタ特性とを合成することにより、暗騒音や歪みの影響を考慮したＴＳＰの設計手法と併用できる。
さらに事前に測定した暗騒音に基づきＴＳＰのパワースペクトラムを同等の特性となるように変更することにより、暗騒音に対して最適なＴＳＰの設計が可能となる。

従来のＴＳＰ（時間引き延ばしパルス）はインパルスを構成する各スペクトルの位相を周波数の２乗に比例して変化させることで、インパルスのエネルギーを時間軸上に分散する技術である。その逆フィルタはＴＳＰ自身を時間軸で逆転したものとなる。
位相を周波数の２乗に比例して変化させる理由は、群遅延特性が位相特性の周波数微分として定義されているためで、実際には群遅延特性を周波数に対して直線的に変化させているものとなる。
ここで包絡線の直線性を維持しつつもパワースペクトラムが平坦ではない時間引き延ばしパルスの定義を考える。
ある帯域が他の帯域に比べ大きなパワーを持つとすると、群遅延グラフを直線とするままでは係る帯域の包絡線は高くなってしまい、直線性を維持できなるなる。
このため、相対的に大きなパワーを持つ帯域では、より長い時間幅にパワーを分散するため、群遅延グラフの傾きを大きく取る（つまり単位周波数当たりの遅延量を大きくする）必要がある。
パワーは振幅の２乗で表され、また群遅延特性は直感的にパワースペクトラムの周波数積分であることが分かる。
また、群遅延特性は位相特性の周波数微分として定義されたものであるから、ここで定義しようとしている時間引き延ばしパルスの周波数領域での設計式Ｈ（）を、振幅特性を表す関数Ａ（）を起点にして、以下のような式に表すことができる。

数１

Ｆ（ｋ）＝Σ（Ａ（ｋ）＾２） （０≦ｋ≦１）
※Σ（シグマ）は積分を表す ※＾はべき乗を表す

数２

Ｓ（ｋ）＝ΣＦ（ｋ） （０≦ｋ≦１）

数３

ｈ（ｋ）＝ｅｘｐ（−Ｓ（ｋ）Ｑ） （０≦ｋ≦１）
※ｅｘｐ（ｘ）＝ｃｏｓ（ｘ）−ｊｓｉｎ（ｘ）
Ｈ（ｋ）＝ｈ（ｋ）Ａ（ｋ） （０≦ｋ≦１）
Ｈ（ｋ）＝Ｃ（ｈ（２−ｋ））Ａ（２−ｋ） （１≦ｋ＜２）
※Ｃ（）は複素共役を表す ※Ｑは正規化のための定数
Ｆ（）は周波数に対する群遅延量を表す関数、Ｓ（）は位相回転量を表す関数である。
ｋは周波数を表すパラメータで、ｋ＝０の時０Ｈｚ、ｋ＝１の時ナイキスト周波数を表す。従来のＯＡＴＳＰの正式な定義ではＴＳＰの掃出線は高い周波数から始まり低い周波数へ下降する直線であるため、これを入力波形に用いると被測定系が歪みを有する場合、図２に示す通り歪み成分がインパルス応答の主成分の側に現れてしまう。そこで、本発明の於ける入力波形用のＴＳＰは、位相の符号を反転し、ＴＳＰの掃出線が上昇カーブを描くようにする。
ｈ（）は位相回転の結果であり、これに振幅特性Ａ（）を乗じてＨ（）を求める。
Ａ（）の値に０が含まれる場合は、逆フィルタを求めることができない。このため少なくとも（０≦ｋ≦１）の範囲でＡ（ｋ）＞０とすると、Ｆ（）のグラフは単調増加し、Ｆ（１）の時に最大の群遅延量となる。

Ｈ（）は、逆フーリエ変換し実時間領域へ展開した際、全てのエネルギーを実数へと変換するため、Ｈ（１）が実数であり、Ｈ（ｋ）とＨ（２−ｋ）が複素共役となる条件を満たす必要がある。
定数Ｑを用意しＳ（１）Ｑがπの整数倍となるようにすることで、ｈ（１）を実数にでき、実数値の振幅を乗じたＨ（１）も実数となる。

数４

Ｑ＝πＮ／Ｓ（１）
※Ｎは任意の自然数

しかし、直流からナイキスト周波数にかけての群遅延量、つまりＴＳＰの全長を表すＦ（１）は、当然Ａ（）の形状に従って変化することから、Ｆ（１）についても正規化を行う必要が生じる。よってＱの定義を以下のように修正する。

数５

Ｑ＝π・Ｉｎｔ（ｎ・Ｓ（１）／Ｆ（１））／Ｓ（１）
※Ｉｎｔは小数部切捨てを表す
※ｎは任意の実数を表す
ｎ＝０の時ＴＳＰはインパルス様の波形となり、ｎ＝１の時、ＴＳＰのナイキスト周波数の成分と直流成分が時間的に回り込んで重なってしまうため、ｎは通常（０＜ｎ＜１）の範囲から決まる。
Ｈ（）を逆フーリエ変換することで目的とするＴＳＰ波形が得られる。
逆フィルタの設計式は位相の符号を反転し、振幅特性の逆数を取れば良い。

数６

ｈｉ（ｋ）＝ｅｘｐ（Ｓ（ｋ）Ｑ） （０≦ｋ≦１）
Ｈｉ（ｋ）＝ｈｉ（ｋ）／Ａ（ｋ） （０≦ｋ≦１）
Ｈｉ（ｋ）＝Ｃ（ｈｉ（２−ｋ））／Ａ（２−ｋ） （１≦ｋ＜２）

Ｈｉ（）を逆フーリエ変換することで目的とするＴＳＰの逆フィルタ波形（以下、逆ＴＳＰと呼ぶ）が得られる。
上記の設計式は、Ｆ（）を起点に書き換えることもできる。

数７

Ａ（ｋ）＝ｓｑｒｔ（Ｆ’（ｋ））
※’（アポストロフィ）は微分を表す
※ｓｑｒｔ（）は平方根を表す
以降の式はＡ（）を起点とする場合と変わらないため省略する。
Ａ（ｋ）＞０の制約によりＦ（）は狭義の単調増加（傾きが常に０を超える）関数でなくてはならない。

ＴＳＰの設計方針として、パワースペクトラムに基づく設計を行う場合はＡ（）を起点にとり、掃出線（周波数−群遅延グラフ）の形状に基づいた設計を行う場合はＦ（）を起点にとることで設計が容易となる。

ここで、掃出線（周波数−群遅延グラフ）の形状に基づいたＴＳＰ設計手段について具体的に述べる。
まず周波数に対する時刻の関数ｆ（ｔ）を定義する。ｆ（ｔ）は、ＴＳＰの掃出線の設計を容易にするため、時刻ｔをパラメータとして、周波数ｋを出力とする関数である。
ｔの変化域を０≦ｔ≦１とし、この範囲内において、０Ｈｚからナイキスト周波数まで周波数が変化するものとする。即ちｔ＝１に於けるｋの値をナイキスト周波数とする。

数８

ｋ＝ｆ（ｔ）
ＴＳＰの掃出線は周波数０Ｈｚから、ナイキスト周波数までを変化域とするため、ｆ（ｔ）は、０Ｈｚを通る必要がある。
また、該掃出線は各周波数に対する群遅延グラフの結果であるため、１周期分のＴＳＰ波形内において同一周波数が別時刻には存在し得ない。このためｆ（ｔ）は単調に増加する必要がある。
また、ＴＳＰは全ての周波数を含む必要があるため、ｆ（ｔ）は連続していなければならない。

次に関数のパラメータを周波数とするため、ｆ（ｔ）の逆関数Ｆ（ｋ）を求める。またｆ（１）の値がナイキスト周波数（ｋ＝１）と対応するように正規化する。ここでｔとｋの変化域が逆転し、０≦ｋ≦１となる。

数９

ｔ＝Ｆ（ｋ） 但し、Ｆ（ｘ）＝Ｉｎｖ（ｆ（ｘ）／ｆ（１））
※Ｉｎｖ（）は逆関数を求める操作を表す
ｆ（ｔ）は単調増加するため、その逆関数であるＦ（ｋ）も単調増加する。以降は、数式７で示した変換を行うことにより振幅特性Ａ（ｋ）を、数式２、数式３、数式５、数式６で示した変換を行うことによりＴＳＰ及び逆ＴＳＰの周波数領域での設計式を求めることができる。

しかし、ｆ（ｋ）が急騰な増分変化を持つ場合（例えばｆ（ｔ）が任意の線分を組み合わせた関数などである場合）逆関数Ｆ（ｋ）も同様な増分変化を持つため、Ｆ（ｋ）の微分関数に相当する振幅特性Ａ（ｋ）に滑らかでない点が生じ、時間領域へ展開した場合、増分が変化している周波数付近において発振様のスペクトルとして現れてしまう。（図４）被測定系が例えばコンピュータのソフトウェア内部で完結する場合は、逆フィルタによって該発振様スペクトルを相殺することが可能であるが、通常の被測定系は、例えばスピーカやマイク、その他電気的素子を介し、少なからず歪みを有するものである。しかし、該発振様スペクトルは歪みを有する測定系では混変調や定常波をもたらし脆弱である。
このため、必要に応じて低域濾過器を表す関数Ｌ（）を用意し、Ｆ（ｋ）を以下のように定義する。

数１０

ｔ＝Ｆ（ｋ） 但し、Ｆ（ｘ）＝Ｌ（Ｉｎｖ（ｆ（ｘ）／ｆ（１）））

通常のＴＳＰ設計では低周波に対してより長い測定時間を与えたいことから、ｆ（ｔ）はｔの高次式や指数関数とすることが多い。
しかし、これらの式は、０Ｈｚ付近にて傾きが０になるか、０に近い値を持つ。このためｆ（ｔ）の逆関数の微分であるＡ（ｋ）は０Ｈｚ付近にて極を持つこととなり、包絡線に乱れが生じてしまう。この問題に対する一つの解決策としては、ｆ（ｔ）に小さな係数を持つｔの１次項を設けて、傾きが０になることを防ぐ方法が考えられる。

数１１

ｆ（ｔ）＝ｆｏ（ｔ）＋αｔ
※ｆｏ（）は単調増加であるが傾き０の点を有する任意関数
※αは小さな値の定数

またｆ（ｔ）がｔの指数関数となる場合、低周波の回り込みを許さなければ十分なＴＳＰ長を得られない場合がある。このため必要に応じ任意のロールオフ関数Ｒ（ｋ）を用意し、回り込みが起こる帯域の振幅を抑制する。

数１２

Ａ（ｋ）＝ｓｑｒｔ（Ｆ’（ｋ））Ｒ（ｋ）

以上までが掃出線の形状に基づいたＴＳＰ設計手段についての説明である。
次にＴＳＰのパワースペクトラムに基づいた設計手段について述べる。スピーカとマイクを使い音楽ホールのインパルス応答を測定するといった一般的なＴＳＰの使用方法に際し、パワースペクトラムの面で得に注意しなければならないのは、スピーカやマイクの周波数特性（いわゆる音のクセ）と、暗騒音のパワースペクトラムである。まず前者について説明する。

数式１から数式７で示したＴＳＰ、及び逆ＴＳＰの周波数領域での設計式の定義に対して、ＴＳＰ信号を逆ＴＳＰで畳み込んだ時、純粋なインパルス信号ではなく、任意信号の逆フィルタとしての特性を持つインパルス応答が復元されるような変形を試みる。
なお以下より係る畳み込みによる逆フィルタ復元の対象となる任意信号を補正信号と呼ぶ。

まず、補正信号をフーリエ変換し、振幅特性と位相特性を求める。

数１３

Ａａ（ｋ）＝ｓｑｒｔ（Ｒａ（ｋ）．Ｒｅ＾２＋Ｒａ（ｋ）．Ｉｍ＾２）
（０≦ｋ≦１）

数１４

Ｐａ（ｋ）＝ａｔａｎ（Ｒａ（ｋ）．Ｉｍ／Ｒａ（ｋ）．Ｒｅ）
（０≦ｋ≦１）
※ａｔａｎ（）は逆正接を表す
Ｒａは補正信号のフーリエ変換後の複素信号であり、Ｒｅは実部、Ｉｍは虚部の実数を表す。
Ａａ（）は補正信号の振幅特性、Ｐａ（）は位相特性である。ｋは、周波数を表すパラメータであり、０の時に直流、２の時にサンプリング周波数を表す。ＴＳＰの加工に用いるｋの範囲は、直流からナイキスト周波数までとする。

ＴＳＰの群遅延特性は、Ａａ（）を元に加工される。しかし通常スピーカやマイクのインパルス応答は、振幅特性に鋭い起伏を持ち、また再生限界域付近では振幅の平均値が全体的に減衰する傾向を持つ。
ＴＳＰの振幅特性が鋭い起伏を持つ場合、発振様のスペクトルが発生し、測定時に発生する歪みに対し脆弱となる。また、再生限界域付近の減衰が著しい場合、ＴＳＰ信号の全長の内、係る帯域に多くの時間が割り当てられることから、逆に中心となる再生域でのＳ／Ｎ比の低下を招く。
この問題を避けるため、Ａａ（）の振幅の最小値を制限し、また振幅を平均化する。

数１５

Ｅａ（ｋ）＝Ａｖｇ（Ｌｍｉｎ（Ａａ（ｋ））） （０≦ｋ≦１）
※Ｌｍｉｎ（）は最小値を制限する操作を表す
※Ａｖｇ（）は移動平均値を求める操作を表す
移動平均値を求める操作は、より滑らかな曲線を得るため窓関数による重み付き移動平均を用いる。

次にＥａ（ｋ）の逆数を取り、ＴＳＰ補正用の振幅特性Ａ１（）を得る。

数１６

Ａ１（ｋ）＝１／Ｅａ（ｋ） （０≦ｋ≦１）

また、数式７において振幅特性は制限及び平均化されてしまっているため、Ａａ（）とＥａ（）との違いを逆ＴＳＰの振幅特性へと盛り込む必要がある。以下の式によりその値Ｆａ（）を求める。

数１７

Ｆａ（ｋ）＝Ａａ（ｋ）／Ｅａ（ｋ） （０≦ｋ≦１）
逆ＴＳＰにおいてもＳ／Ｎ比の周波数帯域による偏りや発振様スペクトルといった問題が発生するため、Ａａ（）と同様の操作を行い、逆ＴＳＰ補正用の振幅特性Ａ２（）を得る。

数１８

Ａ２（ｋ）＝１／Ａｖｇ２（Ｌｍｉｎ２（Ｆａ（ｋ））） （０≦ｋ≦１）
※Ｌｍｉｎ２（）は最小値を制限する操作を表す
※Ａｖｇ２（）は移動平均値を求める操作を表す

次に、ＴＳＰの振幅特性を表すＡｐ（）、及び逆ＴＳＰの振幅特性を表すＡｎ（）を以下のように求める。

数１９

Ａｐ（ｋ）＝Ａ１（ｋ）Ａ（ｋ） （０≦ｋ≦１）

数２０

Ａｎ（ｋ）＝Ａ２（ｋ）／Ａ（ｋ） （０≦ｋ≦１）

そしてＴＳＰ及び逆ＴＳＰの群遅延特性を表すＳ（）を以下のように変形し、Ｓ２（）とする。

数２１

Ｓ２（ｋ）＝Σ（Σ（Ｆ’（ｋ）Ａｄ’（ｋ））） （０≦ｋ≦１）
但し、Ａｄ（ｋ）＝Σ（Ａ１（ｋ）＾２）
数式２１は、夫々周波数に対する群遅延特性を表す関数であるＦ（）とＡｄ（）を微分した値同士を掛け合わせ、その結果を２回積分する操作を表す。
結局のところ２つのパワースペクトラム同士を乗算した結果である。

数２２

Ｓ２（ｋ）＝Σ（Σ（Ａｐ（ｋ）＾２）） （０≦ｋ≦１）

最後にこれまで求めた式を統合してここでのＴＳＰ及び逆ＴＳＰの周波数領域に於ける設計式を得る。

数２３

ｈ（ｋ）＝ｅｘｐ（−Ｓ２（ｋ）Ｑ） （０≦ｋ≦１）
Ｈ（ｋ）＝ｈ（ｋ）Ａｐ（ｋ） （０≦ｋ≦１）
Ｈ（ｋ）＝Ｃ（ｈ（２−ｋ））Ａｐ（２−ｋ） （１≦ｋ＜２）
Ｑ＝π・Ｉｎｔ（ｎ・Ｓ２（１）／Ｓ２’（１））／Ｓ２（１）

数２４

ｈｉ（ｋ）＝ｅｘｐ（Ｓ２（ｋ）Ｑ−Ｐａ（ｋ） （０≦ｋ≦１）
Ｈｉ（ｋ）＝ｈｉ（ｋ）Ａｎ（ｋ） （０≦ｋ≦１）
Ｈｉ（ｋ）＝Ｃ（ｈｉ（２−ｋ））Ａｎ（２−ｋ） （１≦ｋ＜２）
数式２３はＴＳＰ、数式２４は逆ＴＳＰの周波数領域での設計式であり、これらを逆フーリエ変換することにより夫々実時間信号へ展開される。
数式２４においては補正信号の位相特性Ｐａ（）を差し引くことで、位相を揃えるためのフィルタ特性を盛り込む。

ところで、通常スピーカやマイク等の測定媒体の周波数特性は、各々の距離や角度、設置する床の材質、温度や湿度といった様々な条件により、大きく変動するため、厳密かつ環境を選ばずに使える逆フィルタ特性は仮定できないことが多い。
そこで、逆ＴＳＰ側に補正信号の逆フィルタ特性を盛り込まず、逆フィルタはインパルス応答の測定後にケースバイケースの対応として施すこともできる。
この場合、数式２４は次式のように定義し直される。

数２５

ｈｉ（ｋ）＝ｅｘｐ（Ｓ２（ｋ）Ｑ） （０≦ｋ≦１）
Ｈｉ（ｋ）＝ｈｉ （ｋ）／Ａｐ（ｋ） （０≦ｋ≦１）
Ｈｉ（ｋ）＝Ｃ（ｈｉ（２−ｋ））／Ａｐ（２−ｋ） （１≦ｋ＜２）

次に、暗騒音のパワースペクトラムに基づくＴＳＰの設計手段について説明する。
測定環境にて比較的騒音の少ないタイミングで一定時間サンプリングした結果を切り出し、任意の時間窓を施した後、フーリエ変換すれば暗騒音の周波数特性が得られる。ここから振幅特性を求めて数式１のＡ（）とすれば暗騒音と同じパワースペクトラムを持つＴＳＰ及び逆ＴＳＰを設計できるものと思われる。
しかし、自然環境の暗騒音は低周波帯域に集中する傾向にあり、一般的な測定環境において直流付近とナイキスト周波数付近の平均的な振幅の差は１０００倍（６０ｄＢ）を超えることも珍しくはない。
仮にこれと同じパワースペクトラムを持つようなＴＳＰを生成した場合、直流付近とナイキスト周波数付近の測定時間比は１００万倍を超え、直流に近い成分が大半の測定時間を占め、より高い周波数成分の測定時間は非常に短く、インパルス様の波形となってしまう。低周波帯域でのＳ／Ｎ比向上は重要な課題ではあるが、これほど極端な違いとなってしまうと測定時の非線形性の影響を受け易く、却って測定誤差が増加してしまうであろう。
そこで、暗騒音から求めた振幅特性の最大値を制限する操作を行い、次式を得る。

数２６

Ａｍ（ｋ）＝Ｌｍａｘ（Ａｂ（ｋ）） （０≦ｋ≦１）
※Ａｂ（）は暗騒音の振幅特性を表す
※Ｌｍａｘ（）は最大値を制限する操作を表す
ここでＡｍ（）には、Ａｂ（）に急騰な変化が元々含まれていたり、最大値を制限する操作によって急騰な変化が生じることから、このままＴＳＰに変換すると、発振様のスペクトルが発生してしまう。そこで、Ａｍ（）を平均化し、Ａｃ（）とする。

数２７

Ａｃ（ｋ）＝Ａｖｇ３（Ａｍ（ｋ）） （０≦ｋ≦１）
※Ａｖｇ３（）は移動平均値を求める操作を表す

ここでＡｃ（）に数式１９で示したように他の振幅補正を乗じることもできる。

数２８

ＡＣ（ｋ）＝Ａｃ（ｋ）Ａｘ（ｋ）
※Ａｘ（）は任意の振幅特性
以降は、数式１、数式２、数式３、数式５、数式６で示した変換を行うことによりＴＳＰ及び逆ＴＳＰの周波数領域での設計式を求めることができる。

図１に本発明に於けるＴＳＰ及び逆ＴＳＰ生成フローのブロック図を示す。
なお、本発明は、コンピュータシステムにおいて実現されるため、上記の数式で示した変換手順をアルゴリズムにより再現する必要がある。
まず、数式９に示した逆関数への変換であるが、これは二分法を用いて０≦ｔ≦１においてｔ−ｆ（ｔ）が０となるｔを逐一求め、これをメモリに格納することで実現される。ｆ（ｔ）は原点０を通り、単調増加するため、二分法に於ける開始値Ａは、一つ前に求めたｔ若しくは０とし、終了値Ｂは、開始値Ａから始まり現在値に対して必ず大きくなる任意の演算を繰り返すことにより求まる。

積分や微分は任意の数値積分及び数値微分のアルゴリズムを利用できる。簡単な例として積分は各々の関数において、一つ前にまでに求めていた値の総和であり。微分は一つ前の値との差分として求まる。
またｅｘｐ（）は実数用と虚数用のメモリバッファを各々設け、組み込み関数のｓｉｎ（）とｃｏｓ（）を用い、定義に従って計算し格納する。

本発明に於けるＴＳＰ及び逆ＴＳＰの生成フローを示した図である。 ＯＡＴＳＰを用いた場合に歪み成分がインパルス応答の主成分に混入する様子を示した図である。 本発明に於けるＴＳＰを用いた場合に歪み成分がインパルス応答の主成分に混入する現象が如何に軽減されるかを示した図である。 ｆ（ｔ）の増分が急騰に変化した際に生じる、発振様スペクトラムの様子を示した図である。

Claims (5)

1. 時間引き延ばしパルス（ＴＳＰ）を用いて線形時間不変な連続時間の被測定系の離散的インパルス応答を測定する装置において、時間引き延ばしパルスの周波数領域での設計式を、振幅特性、又は群遅延特性を起点とした数値的（非解析的）算法により生成することを特徴とする、時間引き延ばしパルスとその逆フィルタの波形を生成する装置。
2. 請求項１のＴＳＰ生成処理において、初期時刻に周波数Ｏを通り、単調増加し、連続な任意の時刻をパラメータとした周波数の関数によって表現される時刻一周波数グラフを掃出線として持つ時間引き延ばしパルスを入力信号として用いることを特徴とする、時間引き延ばしパルスとその逆フィルタの波形を生成する装置。
3. 請求項１のＴＳＰ生成処理において、測定時とは別個のインパルス応答に基づいて群遅延特性及び振幅特性を補正したＴＳＰ信号を被測定系の入力信号として用いて、被測定系から該別個のインパルス応答の周波数特性を取り除いたインパルス応答を復元することを特徴とする、時間引き延ばしパルスとその逆フィルタの波形を生成する装置。
4. 請求項１のＴＳＰ生成処理において、任意設計された複数の群遅延特性を振幅特性レベルの乗算により合成する装置。
5. 請求項１のＴＳＰ生成処理において、ＴＳＰのパワースペクトラムを、該測定場所にて事前に測定した暗騒音のパワースペクトラムに基づき、同等の特性となるように変更する装置。

Images