"エンジン音における、所定の周波数の周波数信号を所定の時間ごとに算出する周波数分析手段と、時間の経過に伴って前記周波数信号の位相が加速度的に増加するか又は加速度的に減少するかを判定することにより、エンジン回転数の増加又は減少を判定する回転数判定手段とを備える回転数増減判定装置。前記回転数判定手段は、時間の経過に伴って前記位相が加速度的に増加する場合にはエンジン回転数が増加していると判定し、時間の経過に伴って前記位相が加速度的に減少する場合にはエンジン回転数が減少していると判定する請求項1記載の回転数増減判定装置。さらに、前記周波数信号の位相の時間変化を近似する位相曲線を算出する位相曲線算出部を備え、前記回転数判定手段は、前記位相曲線の形状に基づいて前記周波数信号の位相が加速度的に増加するか又は加速度的に減少するかを判定することにより、エンジン回転数の増加又は減少を判定する請求項1記載の回転数増減判定装置。前記回転数判定手段は、前記位相曲線が下に凸の場合に前記周波数信号の位相が加速度的に増加していると判定することにより、エンジン回転数が増加していると判定する請求項3記載の回転数増減判定装置。前記回転数判定手段は、前記位相曲線が上に凸の場合に前記周波数信号の位相が加速度的に減少していると判定することにより、エンジン回転数が減少していると判定する請求項3記載の回転数増減判定装置。前記回転数判定手段は、時間の経過に伴う位相の変化値が所定の閾値以下の場合にのみ前記エンジン回転数の増加又は減少を判定する請求項3記載の回転数増減判定装置。前記位相曲線は、2次多項式で表される曲線である請求項3記載の回転数増減判定装置。さらに、所定数の前記位相との差が小さくなるように、前記所定数の前記位相とは異なる他の前記位相に±2π×mラジアンmは自然数を加算することにより、他の前記位相を補正する位相補正部を備える請求項3記載の回転数増減判定装置。さらに、前記位相曲線と前記周波数信号の位相との誤差を算出する誤差算出手段と、互いに異なる角度範囲ごとに、当該角度範囲内に収まるように前記位相に±2π×mラジアンmは自然数を加算することにより、前記位相を補正する位相補正部とを備え、前記位相曲線算出部は、前記角度範囲ごとに、前記位相曲線を算出し、前記誤差算出手段は、前記角度範囲ごとに、前記誤差を算出し、前記位相補正部は、さらに、前記位相曲線と、前記周波数信号の位相との誤差が最小となるときの角度範囲を選択し、前記回転数判定手段は、選択された前記角度範囲における前記位相曲線の形状に基づいて前記周波数信号の位相が加速度的に増加するか又は加速度的に減少するかを判定することにより、エンジン回転数の増加又は減少を判定する請求項3記載の回転数増減判定装置。前記周波数分析手段は、雑音とエンジン音とを含む混合音における、前記所定の周波数の周波数信号を前記所定の時間ごとに算出し、前記位相曲線算出部は、前記混合音の周波数信号の位相の時間変化を近似する位相曲線を算出し、前記回転数増減判定装置は、さらに、前記位相曲線と、前記混合音の周波数信号の位相との誤差を算出する誤差算出手段と、前記誤差に基づいて、前記混合音がエンジン音か否かを識別する音響信号識別手段とを備え、前記回転数判定手段は、前記音響信号識別手段でエンジン音であると識別された前記混合音の位相について、エンジン回転数の増加又は減少を判定する請求項3記載の回転数増減判定装置。前記周波数分析手段は、各々がエンジン音の入力を受け付ける互いに離間して配置された複数のマイクロホンで受け付けられた複数のエンジン音のそれぞれについて周波数信号を算出し、前記回転数増減判定装置は、さらに、前記複数のマイクロホンで受け付けられた複数の前記エンジン音の到達時間差に基づいて、前記エンジン音の音源方向を検知し、前記回転数判定手段によりエンジン回転数が増加していると判定された場合にのみ、前記音源方向の検知結果を出力する方向検知部を備える請求項1記載の回転数増減判定装置。前記回転数判定手段は、さらに、エンジン回転数が増加している場合に前記エンジン音を発する車両が加速していると判定し、エンジン回転数が減少している場合に前記エンジン音を発する車両が減速していると判定する請求項1記載の回転数増減判定装置。エンジン音における、所定の周波数の周波数信号を所定の時間ごとに算出する周波数分析ステップと、時間の経過に伴って前記周波数信号の位相が加速度的に増加するか又は加速度的に減少するかを判定することにより、エンジン回転数の増加又は減少を判定する回転数判定ステップとを含む回転数増減判定方法。エンジン音における、所定の周波数の周波数信号を所定の時間ごとに算出する周波数分析ステップと、時間の経過に伴って前記周波数信号の位相が加速度的に増加するか又は加速度的に減少するかを判定することにより、エンジン回転数の増加又は減少を判定する回転数判定ステップとをコンピュータに実行させるためのプログラム。第1の従来技術では、周囲音を音圧レベル信号に変換し、音圧レベル信号の特定周波数帯域における絶対量と判断レベルとを比較することで、周囲車両の有無を判断するとともに、音圧レベル信号の時間的変化から周囲車両が接近するか否かを判断している。したがって、第1の従来技術によっては、さらに詳細な接近の状況である周囲車両のエンジン回転数の増減の状況又は周囲車両の加減速の状況を判定することができないという問題がある。また、周囲車両のエンジン回転数の増減又は周囲車両の接近や加速などの判定には、一般的に、エンジン音の周波数の変化や音圧の変化が観測できる十分に長い時間数秒の音声信号を要する。このため、周囲車両のエンジン回転数の増減の状況又は周囲車両の加減速の状況を運転者に短時間に知らせる必要がある安全運転支援などのアプリケーションに、従来の技術を利用することが困難である。本発明は、上述の課題を解決するためになされたものであり、自車の周囲に存在する周囲車両のエンジン回転数の増減をリアルタイムで判断することができる回転数増減判定装置等を提供することを目的とする。上記目的を達成するために、本発明のある局面に係る回転数増減判定装置は、エンジン音における、所定の周波数の周波数信号を所定の時間ごとに算出する周波数分析手段と、時間の経過に伴って前記周波数信号の位相が加速度的に増加するか又は加速度的に減少するかを判定することにより、エンジン回転数の増加又は減少を判定する回転数判定手段とを備える。具体的には、前記回転数判定手段は、時間の経過に伴って前記位相が加速度的に増加する場合にはエンジン回転数が増加していると判定し、時間の経過に伴って前記位相が加速度的に減少する場合にはエンジン回転数が減少していると判定する。エンジン回転数が増加する場合には、エンジン音の周波数が時間の経過とともに増加し、エンジン音の周波数信号の位相は加速度的に増加する。一方、エンジン回転数が減少する場合には、エンジン音の周波数が時間の経過とともに減少し、エンジン音の周波数信号の位相は加速度的に減少する。位相が加速度的に増加しているか又は加速度的に減少しているかは、短い時間範囲に含まれる位相から判断することが可能である。このため、この構成によると、自車の周囲に存在する周囲車両のエンジン回転数の増減をリアルタイムで判断することができる。好ましくは、上述の回転数増減判定装置は、さらに、前記周波数信号の位相の時間変化を近似する位相曲線を算出する位相曲線算出部を備え、前記回転数判定手段は、前記位相曲線の形状に基づいて前記周波数信号の位相が加速度的に増加するか又は加速度的に減少するかを判定することにより、エンジン回転数の増加又は減少を判定する。具体的には、前記回転数判定手段は、前記位相曲線が下に凸の場合に前記周波数信号の位相が加速度的に増加していると判定することにより、エンジン回転数が増加していると判定する。また、前記回転数判定手段は、前記位相曲線が上に凸の場合に前記周波数信号の位相が加速度的に減少していると判定することにより、エンジン回転数が減少していると判定する。位相が加速度的に増加する場合には、位相曲線は下に凸の形状を有し、位相が加速度的に減少する場合には、位相曲線は上に凸の形状を有するという性質がある。この性質を利用することにより、精度良く、位相が加速度的に増加しているか又は加速度的に減少しているかを判定することができ、これにより、エンジン回転数が増加しているか減少しているかを判定することができる。好ましくは、前記回転数判定手段は、時間の経過に伴う位相の変化値が所定の閾値以下の場合にのみ前記エンジン回転数の増加又は減少を判定する。周囲車両がギアチェンジしたような場合には、位相が急激に変化する。このため、そのような場合を除外して、上記判定を行うことができるようになる。好ましくは、上記回転数増減判定装置は、さらに、所定数の前記位相との差が小さくなるように、前記所定数の前記位相とは異なる他の前記位相に±2π×mラジアンmは自然数を加算することにより、他の前記位相を補正する位相補正部を備える。これにより、他の時刻の位相と大きく外れた位相を補正することができ、エンジン回転数の増減を精度良く判定することが可能となる。また、上記回転数増減判定装置は、さらに、前記位相曲線と前記周波数信号の位相との誤差を算出する誤差算出手段と、互いに異なる角度範囲ごとに、当該角度範囲内に収まるように前記位相に±2π×mラジアンmは自然数を加算することにより、前記位相を補正する位相補正部とを備え、前記位相曲線算出部は、前記角度範囲ごとに、前記位相曲線を算出し、前記誤差算出手段は、前記角度範囲ごとに、前記誤差を算出し、前記位相補正部は、さらに、前記位相曲線と、前記周波数信号の位相との誤差が最小となるときの角度範囲を選択し、前記回転数判定手段は、選択された前記角度範囲における前記位相曲線の形状に基づいて前記周波数信号の位相が加速度的に増加するか又は加速度的に減少するかを判定することにより、エンジン回転数の増加又は減少を判定してもよい。これにより、他の時刻の位相と大きく外れた位相を補正することができ、エンジン回転数の増減を精度良く判定することが可能となる。好ましくは、前記周波数分析手段は、雑音とエンジン音とを含む混合音における、前記所定の周波数の周波数信号を前記所定の時間ごとに算出し、前記位相曲線算出部は、前記混合音の周波数信号の位相の時間変化を近似する位相曲線を算出し、前記回転数増減判定装置は、さらに、前記位相曲線と、前記混合音の周波数信号の位相との誤差を算出する誤差算出手段と、前記誤差に基づいて、前記混合音がエンジン音か否かを識別する音響信号識別手段とを備え、前記回転数判定手段は、前記音響信号識別手段でエンジン音であると識別された前記混合音の位相について、エンジン回転数の増加又は減少を判定する。この構成によると、雑音の影響を除去して、エンジン音についてのみエンジン回転数の増加又は減少を判定することができる。このため、判定の精度を向上させることができる。さらに好ましくは、前記周波数分析手段は、各々がエンジン音の入力を受け付ける互いに離間して配置された複数のマイクロホンで受け付けられた複数のエンジン音のそれぞれについて周波数信号を算出し、前記回転数増減判定装置は、さらに、前記複数のマイクロホンで受け付けられた複数の前記エンジン音の到達時間差に基づいて、前記エンジン音の音源方向を検知し、前記回転数判定手段によりエンジン回転数が増加していると判定された場合にのみ、前記音源方向の検知結果を出力する方向検知部を備える。エンジン回転数が増加していると判定された場合にのみ、音源方向の検知結果を出力することができる。このため、周囲車両が加速しながら接近してくるなどのように、特に危険な場合にのみ、周囲車両が接近してくる方向を運転者に提示することが可能となる。なお、本発明は、このような特徴的な手段を備える回転数増減判定装置として実現することができるだけでなく、回転数増減判定装置に含まれる特徴的な手段をステップとする回転数増減判定方法として実現したり、回転数増減判定方法に含まれる特徴的なステップをコンピュータに実行させるプログラムとして実現したりすることもできる。そして、そのようなプログラムは、CD-ROMCompactDiscReadOnlyMemory等の不揮発性の記録媒体やインターネット等の通信ネットワークを介して流通させることができるのは言うまでもない。本発明によると、自車の周囲に存在する周囲車両のエンジン回転数の増減をリアルタイムで判断することができる。本発明の特徴は、エンジン音などの周期音であって、周波数が時間にともなって変化する音の位相の時間変化に着目して、車両の加減速を判断することである。なお、本願発明における周期音とは、位相が一定の音又は位相の変化が連続的な音を示す。ここで、図1を用いて本発明で用いる位相の定義を行う。図1aには、入力したエンジン音の例が模式的に示されている。横軸は時間を表しており、縦軸は振幅を表している。ここではエンジンの回転数が時刻に対して一定であり、エンジン音の周波数が変化しない場合の例が示されている。また、図1bには、フーリエ変換を用いて周波数分析を行う場合の基底波形である周波数fの正弦波ここではエンジン音の周波数と同じ値を所定の周波数fとしているが示されている。横軸と縦軸は図1aと同じである。この基底波形と入力した混合音との畳み込み処理を行うことで周波数信号位相を求める。この例では、基底波形を時間軸方向に移動させずに固定し、入力したエンジン音と畳み込み処理を行うことで、時刻ごとの周波数信号位相を求めている。この処理で求めた結果を図1cに示す。横軸は時間を表しており、縦軸は位相を表している。この例では、エンジンの回転数が時刻に対して一定であり、入力したエンジン音の周波数が時刻に対して一定である。このため、所定の周波数fでの位相は加速度的に増加又は加速度的に減少はしていない。この例では回転数が一定であるエンジン音の周波数と同じ値を所定の周波数fとしたが、エンジン音の周波数よりも小さい値を所定の周波数fとした場合には位相は一次関数的に増加する。また、エンジン音の周波数よりも大きい値を所定の周波数fとした場合には、位相は一次関数的に減少する。いずれの場合も、所定の周波数fでの位相は、加速度的に増加又は加速度的に減少はしていない。なお、音声信号分野や高速フーリエ変換FFTなどでは基底波形を時間軸方向にずらしながら畳み込みを行うのが一般的である。この基底波形を時間軸方向にずらしながら畳み込みを行う場合は後に位相を補正することで本発明で定義する位相へと変換することが可能である。以下、図を用いて説明する。図2は位相を説明する図である。図2aには、入力したエンジン音の例が模式的に示されている。横軸は時間を表しており、縦軸は振幅を表している。また、図2bには、フーリエ変換を用いて周波数分析を行う場合の基底波形である周波数fの正弦波ここではエンジン音の周波数と同じ値を所定の周波数fとしているが示されている。横軸と縦軸は図2aと同じである。この基底波形と入力した混合音との畳み込み処理を行うことで周波数信号位相を求める。この例では、基底波形を時間軸方向に移動させながら、入力したエンジン音と畳み込み処理を行うことで、時刻ごとの周波数信号位相を求めている。この処理で求めた結果を図2cに示す。横軸は時間を表しており、縦軸は位相を表している。入力したエンジン音は周波数fであるため、周波数fでの位相のパターンは、1/fの時刻の周期で規則的に繰り返されることとなる。そこで算出された位相ψtから規則的に繰り返される位相を補正ψ ́t=mod2πψt-2πftfは分析周波数することで図2dに示すような位相が得られる。つまり、位相補正を行うことにより、図1cに示す、本発明で定義される位相へと変換することが可能となる。次に、エンジン回転数に伴うエンジン音の時刻に対する周波数の変化について説明する。図3は自動車のエンジン音を後述するDFT分析部において分析したスペクトログラムである。縦軸は周波数を、横軸は時間をそれぞれ示しており、色の濃度は周波数信号のパワーの大きさを示している。濃い色黒い色はパワーが大きいことを示している。図3は、風などの雑音をなるべく除去したデータであり、色の濃い部分黒っぽい部分がおおむねエンジン音を示している。一般的にこのようにエンジン音は、時間とともに回転数が変化したデータであり、スペクトログラムから時間の経過と共に周波数が変化していることが分かる。エンジンは、所定数のシリンダーがピストン運動を行うことで駆動系を回転させている。そして、車両から発せられるエンジン音は、このエンジンの回転に依存した音と、エンジンの回転には依存しない固定振動音や非周期音とからなる。特に車両の外部から検知できる主な音は、エンジンの回転に依存した周期音であり、本実施の形態では、このエンジンの回転に依存する周期音に着目し、加減速の判定を行う。図3の点線の円501、502及び503に示すように、エンジン音は、回転数が変化することで周波数が部分的に時刻に応じて変化していることが分かる。ここで、周波数の変化に着目すると、周波数がランダムに変化したり、離散的に飛んだりすることはほとんどなく、所定の時間間隔でみると、所定の増減を示していることが分かる。例えば、区間Aでは右肩下がりに周波数が減少していることが分かる。この区間ではエンジン回転数は減少しており車両は減速している。区間Bでは右肩上がりに周波数が増加していることが分かる。この区間ではエンジン回転数は増加しており車両は加速している。また、区間Cではほぼ一定の周波数で推移していることが分かる。この区間ではエンジン回転数は一定であり車両は定常走行している。ここで、エンジン回転数の増減とエンジン音の位相との関係について分析する。図4aは、区間Cにおける、エンジン回転数が一定のときのエンジン音を模式的に示した図である。ここではエンジン音の周波数をfとする。図4bは基底波形を示す図である。ここでは基底波形の周波数をエンジン音の周波数fと同じ値にしている。図4cは、基底波形に対する位相を示す図である。図4cに示すように、エンジン回転数が一定であるエンジン音は図1に示す正弦波のように一定の周期を有する。このため、所定の周波数fでの位相は、時間変化に対して加速度的に増加又は加速度的に減少はしない。なお、対象とする音が一定の周波数であり、基底波形の周波数が低い場合、位相は徐々に遅れることとなる。しかし、減少量は一定となるため、位相の形状は線形的に減少することとなる。一方、対象とする音が一定の周波数であり、基底波形の周波数が高い場合、位相は徐々に早くなる。しかし、その増加量は一定となるため、位相の形状は線形的に増加することとなる。図5aは、区間Bにおける、エンジン回転数が増加して車両が加速するときのエンジン音を模式的に示した図である。このときエンジン音の周波数は時間とともに増加する。図5bは基底波形を示す図である。例えば、基底波形の周波数はfとする。図5cは、基底波形に対する位相を示す図である。エンジン音は正弦波のように周期性を有しつつ、徐々に周期が高くなる波形を有することから、図5cに示すように、基底波形に対する位相は、時間変化に対して加速度的に増加する。図6aは、区間Aにおける、エンジン回転数が減少して車両が減速するときのエンジン音を模式的に示した図である。このときエンジン音の周波数は時間とともに減少する。図6bは基底波形を示す図である。例えば、基底波形の周波数はfとする。図6cは、基底波形に対する位相を示す図である。エンジン音は正弦波のように周期性を有しつつ、徐々に周期が低くなる波形を有することから、図6cに示すように、基底波形に対する位相は、時間変化に対して加速度的に減少する。したがって、図5c又は図6cに示すように、基底波形に対する位相を用いて、位相の時間変化に対する加速度的な増減を求めることで、エンジン回転数の増減、すなわち車両の加減速を判定することができる。また、本実施の形態では短時間で大きく変化する位相の性質を利用することで、スペクトルのパワーの変化によって加減速を求める従来技術と比較して、短時間のデータで瞬時に加減速を判定することが可能となる。よって、周囲車両の加減速の状況を短時間で運転者に知らせることができる。例えば、こちらの車両が走行する道路が優先道路で、相手の車両が走行する道路に一旦停止線がある死角交差点において、相手の車両が加速又は定常走行で交差点を通過しようとしているのか、一旦停止で止まろうとしているのかを運転者に知らせることができる。以下、本発明の実施の形態について、図面を参照しながら説明する。実施の形態1実施の形態1に係る加減速判定装置について説明する。この加減速判定装置は、請求の範囲の回転数増減判定装置に対応する。図7は、本発明の実施の形態1における雑音除去装置の構成を示すブロック図である。図7において、加減速判定装置3000は、DFT分析部3002と、位相補正部3003jj=1~Mと、周波数信号選択部3004jj=1~Mと、位相曲線算出部3005jj=1~Mと、加減速判定部3006jj=1~Mとを含む。位相補正部3003jj=1~Mは、M個の位相補正部から構成され、j番目の位相補正部3003jは、後述する周波数帯域jについての処理を実行する。本明細書中において、同様の参照符号の記載を行う処理部についても、同様である。DFT分析部3002は、請求の範囲の周波数分析手段に対応する。加減速判定部3006jは、請求の範囲の回転数判定手段に対応する。DFT分析部3002は、入力されたエンジン音3001に対してフーリエ変換処理を施し、エンジン音3001の位相情報を含む周波数信号を複数の周波数帯域のそれぞれについて求める。なお、DFT分析部3002は、高速フーリエ変換や離散コサイン変換やウェーブレット変換などの別の周波数変換方法で周波数変換を行ってもよい。以下では、DFT分析部3002から求められた周波数帯域の個数をMとして、それらの周波数帯域を指定する番号を記号jj=1~Mで表すこととする。位相補正部3003jj=1~Mは、DFT分析部3002が求めた周波数帯域jの周波数信号に対して、時刻tの周波数信号の位相をψtラジアンとするときに、ψ ́t=mod2πψt-2πftfは分析周波数に位相を補正する。周波数信号選択部3004jj=1~Mは、所定の時間幅において、位相補正部3003jj=1~Mが位相補正した周波数信号の中から、位相曲線の算出に用いる周波数信号を選択する。位相曲線算出部3005jj=1~Mは、周波数信号選択部3004jj=1~Mが選択した周波数信号の補正された位相ψ ́tを用いて、時間経過に伴い位相が変化する位相形状を二次曲線として計算する。加減速判定部3006jj=1~Mは、位相曲線算出部3005jj=1~Mが算出した位相曲線から、位相の増加量をもとに、エンジン回転数の増減、つまり車両の加減速を判定する。時間の経過に伴い、エンジン回転数が増加しているときは、車両が加速しているときであり、エンジン回転数が減少しているときは、車両が減速しているときである。これらの処理を、所定の時間幅を時間方向に移動させながら行うこととなる。なお、本発明の必須の構成要件は、図7に示したDFT分析部3002と、加減速判定部3006jとである。DFT分析部3002が、図1cに示した本発明で定義する位相を直接導き出すことができるのであれば、位相補正部3003jは不要である。次に、以上のように構成された加減速判定装置3000の動作について説明する。以下では、j番目の周波数帯域について説明を行う。ここでは、周波数帯域の中心周波数と基底波形の周波数とが一致する場合を例にして説明を行う。つまり、分析周波数fに対して位相ψ ́t=mod2πψt-2πftにおける周波数fが増加するか否かを判定することになる。なお、本実施の形態において、DFT分析部3002は、いわゆる基底波形を時間軸方向にずらしながら行う一般的な周波数分析とし、得られる位相はψtとなる。そこで前述で定義した位相ψ ́へと補正する処理ψ ́t=mod2πψt-2πftを行うこととする。図8は、加減速判定装置3000の動作手順を示すフローチャートである。初めに、DFT分析部3002は、エンジン音3001を受付けて、エンジン音3001に対してフーリエ変換処理を施し、周波数信号を周波数帯域jごとに求めるステップS101。次に、位相補正部3003jは、DFT分析部3002が求めた周波数帯域jの周波数信号に対して、時刻tの周波数信号の位相をψtラジアンとするときに、ψtをψ ́t=mod2πψt-2πftfは分析周波数に変換することで位相補正を行うステップS102j。ここで、本発明において位相を用いる理由及び位相補正を行う方法例について図を用いて説明する。図3は自動車のエンジン音をDFT分析部3002において分析したスペクトログラムである。縦軸は周波数、横軸は時間をそれぞれ示しており、色の濃度は周波数信号のパワーの大きさを示している。濃い色はパワーが大きいことを示している。図3は、風などの雑音をなるべく除去したデータであり、色の濃い部分がおおむねエンジン音を示している。一般的にこのようにエンジン音は、時間とともに回転数が変化したデータであり、スペクトログラムから時間の経過と共に周波数が変化していることが分かる。エンジンは、所定数のシリンダーがピストン運動を行うことで駆動系を回転させている。そして、車両から発せられるエンジン音は、このエンジンの回転に依存した音と、エンジンの回転には依存しない固定振動音又は非周期音とからなる。特に、車両の外部から検知できる主な音は、エンジンの回転に依存した周期音である。本実施の形態では、周期音がこのエンジンの回転に依存する周期音である点に着目し、位相の時間変化をもとに加減速の判定を行う。図3の点線の円501、502、503に示すように、エンジン音は、回転数が変化することで周波数が時刻に応じて変化していることが分かる。ここで周波数の変化に着目すると、周波数がランダムに変化したり、離散的に飛んだりすることはほとんどなく、所定の時間間隔でみると、所定の増減を示していることが分かる。例えば区間Aでは右肩下がりに周波数が減少していることが分かる。この区間ではエンジン回転数は減少しており車両は減速している。区間Bでは右肩上がりに周波数が増加していることが分かる。この区間ではエンジン回転数は増加しており車両は加速している。また、区間Cではほぼ一定の周波数で推移していることが分かる。この区間ではエンジン回転数は一定であり車両は定常走行している。図9はDFT分析におけるパワーと位相について説明する図である。図9aは図3と同様に、自動車のエンジン音をDFT分析したスペクトログラムである。図9bはDFT分析の概念を示す図である。例えばエンジン回転数が増加して加速している区間である時刻t1から所定の時間窓幅の所定の窓関数ハニング窓を用いて複素空間上に周波数信号601を表したものである。周波数f1、f2、f3等、各周波数の振幅と位相が算出される。周波数信号601の長さが振幅の大きさパワーを示し、周波数信号601と実軸とのなす角が位相を示している。そして、時間シフトを行いながら各時刻における周波数信号を求めることとなる。ここで、一般的にスペクトログラムは各時刻における各周波数のパワーを示すのみであり、位相については省略されている。図3や図9aに示すスペクトログラムも同様に、DFT分析したパワーの大きさのみを表示したものである。周波数信号の位相ψt及び大きさパワーPtは、周波数信号の実部をxtと表すこととして、周波数信号の虚部をytと表すこととすると、及びである。ここでの記号tは周波数信号の時刻を表している。図9cには、図9aにおいて、エンジン回転数が増加して加速している区間である周波数例えば周波数f4のパワーの時間変動が示されている。横軸は時間軸である。縦軸は周波数信号の大きさパワーを表している。図9cより、パワーの変動はランダムであり、増加や減少を観測することはできない。図9cに示すように、一般的にスペクトログラムは位相情報を省略し、パワーのみで信号の変化を表す。このため、エンジン音の音圧の変化を観測するためには、十分に長い時間数秒の音声信号を必要とする。さらに、風などの雑音を含む場合、音圧の変化はノイズに埋もれてしまうため、観測が困難となる。このため、周囲車両の加減速の状況を運転者に短時間に知らせる必要がある安全運転支援などのアプリケーションには利用することが、従来困難であった。図9dには、図9aにおいて、エンジン回転数が増加して加速している区間の所定の周波数間f4からf5へと回転数が増加しているとするでの時間変動が示されている。横軸は時間軸である。縦軸は周波であり、斜線で塗りつぶした部分902を一定のパワーを有する区間として表している。図9dより、周波数の変動はランダムであり、エンジン回転数の増加や減少を観測することはできないことが分かる。図9cに示すように、一般的にスペクトログラムでは位相情報を省略し、パワーのみで信号の変化を表すため、エンジン音の周波数の変化を観測するためには、十分に長い時間数秒の音声信号を必要とする。さらに、風などの雑音を含む場合、周波数変化はさらにノイズに埋もれてしまうため、観測が困難となる。例えばエンジン音が周波数f4から周波数f5へと変化していても、その間、雑音があれば変化を周波数情報から観測することはできない。このため、周囲車両の加減速の状況を運転者に短時間に知らせる必要がある安全運転支援などのアプリケーションには利用することが困難であった。そこで本実施の形態では位相に着目し、位相の時間変化をもとに加減速を判定することとする。上記エンジン音の回転数の増減と、位相の時間変化との関係を数式で表すと以下の関係式で表すことができる。図3等に示すようにエンジン音の周波数の変化は、周波数がランダムに変化したり、離散的に飛んだりすることはほとんどなく、所定の時間間隔でみると、所定の増減を示していることが分かる。したがって、この増減を例えば下記の式4で示すような、一次の区分線形で近似することとする。具体的には所定の時間区間で見た場合、時刻tにおける周波数fは、初期値f0から時刻tに比例比例係数Aして増減する線分で線形近似できると考えられる。そして、周波数fを上記式4で表した場合、時刻tにおける位相ψは、とあらわせる。ここで右辺の第3項のψ0は初期位相であり、第2項2πf0tは時刻tに比例して角周波数2πf0tだけ位相が進むことを示している。そして第1項πAt2から、位相は二次曲線で近似できることを示している。次に位相の時間変化の近似処理を容易にするための位相補正処理について説明を行う。一般的にFFTやDFTで得られる位相は、基底波形を時間軸にずらしつつ算出しているため、図2cと図2dに示すように、位相ψtを位相ψ ́t=mod2πψt-2πftfは分析周波数に変換することで位相補正を行う必要がある。以下、詳細を説明する。初めに、位相補正部3003jは、基準の時刻を決定する。図10aは、図9aにおける時刻t1からの所定時間区間における位相を示す図であって、図10aの黒丸印の時刻t0を基準の時刻に決定している。次に、位相補正部3003jは、位相を補正する周波数信号の複数の時刻を決定する。この例では、図10aの5個の白丸印の時刻t1、t2、t3、t4、t5を、位相を補正する周波数信号の時刻に決定している。ここで、基準の時刻t0における周波数信号の位相をと表すこととして、位相を補正する5個の時刻における周波数信号の位相をと表すことにする。これらの補正する前の位相を図10aにおいて×印で示してある。また、対応する時刻の周波数信号の大きさはで表すことができる。次に、図11に、時刻t2における周波数信号の位相を補正する方法を示す。図11aと図10aとは同じ内容のものである。また、図11bは、1/ffは分析周波数の時間間隔で、等角速度で0~2πラジアンまで規則的に変化する位相を実線で表している。ここで、補正したあとの位相をと表すことにする。図11bにおいて、基準の時刻t0と時刻t2との位相を比較すると、時刻t2の位相は時刻t0の位相よりだけ大きい。そこで、図11aにおいて、基準の時刻t0の位相ψt0との時間差に起因する位相のずれを補正するために、時刻t2の位相ψt2からΔψを差し引いてψ ́t2を求める。これが位相補正後の時刻t2の位相である。このとき、時刻t0の位相は基準の時刻における位相であるので位相補正後も同じ値となる。具体的には、位相補正後の位相をにより求める。位相補正したあとの周波数信号の位相を図10bに×印で示す。図10bの表示の方法は図10aと同様であるため説明を省略する。次に、位相曲線算出部3005jは、位相補正部3003jが求めた補正後の位相情報を用いて、位相の時間変化を曲線として算出する。再度、図8を参照して、周波数信号選択部3004jは、位相補正部3003jが求めた所定の時間幅における位相補正された周波数信号から、位相曲線算出部3005jが位相の形状を計算する際に用いる周波数信号を選択するステップS103j。ここでは、分析の対象とする時刻をt0として、時刻t0と時刻t1、t2、t3、t4、t5とにおける周波数信号の位相から位相の形状を算出する。このとき、位相曲線を求めるときに用いた周波数信号時刻t0~t5の6個の周波数信号は所定の値以上の数から構成されている。このことは、位相距離を求めるために選択された周波数信号の数が少ない場合に、位相の時間変化の規則性を判定することが困難になるからである。ここでの所定の時間幅の時間長は、抽出音の位相の時間変化の性質に基づいて決定することとしてもよい。次に、位相曲線算出部3005jは、位相曲線を算出するステップS104j。位相曲線は例えば以下の二次多項式式12で近似して算出することとする。図12は位相曲線の算出処理を説明する図である。図12に示すように、所定の数の点から二次曲線を算出することができる。本実施の形態では、二次曲線を重回帰曲線として算出する。具体的には、各時刻tii=0,1,2,3,4,5おける補正後の位相をψ ́tiとした場合、二次曲線Ψtの各係数A2、A1、A0はそれぞれ、と表せる。なお、各係数はである。再度、図8を参照して、加減速判定部3006jj=1~Mは、位相曲線算出部3005jj=1~Mが算出した位相曲線から、位相の増加量をもとに、エンジン回転数の増減、つまり車両の加減速を判定するステップS105j。つまり、加減速判定部3006jは、位相曲線算出部3005jが算出した曲線からの加減速を判定する。具体的には、位相曲線算出部3005jが算出した二次曲線の凸とつの向きで、加減速を判定することとなる。式12で得られた係数A2が正、すなわち下に凸とつの場合は、エンジン回転数が増加している、つまり車両が加速していると判定する。一方、係数A2が負、すなわち上に凸とつの場合は、エンジン回転数が減少している、つまり車両が減速していると判定する。なお、本実施の形態では、分析対象である時刻t0に対し、時刻t1、t2、t3、t4、t5の位相から位相の形状を算出している。例えば時刻t2を分析対象とする場合つまり時刻t2を時刻t0 ́とした場合、あらたに時刻t1 ́、t2 ́、t3 ́、t4 ́、t5 ́の位相から位相曲線を算出して加減速を判定してもよいし、既に算出されたt0、t1、t2、t3、t4、t5の位相から算出された位相曲線から加減速を判定することとしてもよい。後者の判定方法を行うことにより、計算量の削減の効果を奏する。さらに、時刻ごとに加減速を判定するのではなく、分析対象を所定の区間として、所定の区間ごとに加減速を判定することとしてもよい。なお、位相補正部3003jは、位相補正において以下に説明する位相補正の処理を、さらに行うこととしてもよい。以下に説明する位相補正を行う場合には、位相曲線の算出や位相曲線との誤差算出等の処理が付随して行われる。このため、位相補正部3003jは、位相曲線算出部3005jでの計算結果を随時参照しながら処理を行うものとする。図13は、さらに実施される位相補正を説明する図である。図13のグラフはいずれもエンジン音の一部を周波数分析したグラフであり、横軸は時間を、縦軸は位相をそれぞれ示す。各白い丸印は位相補正部3003iで位相補正された周波数信号である。図13aにおいて、白い丸印で示される周波数信号の位相を用いて位相曲線を算出すると太い破線で示す曲線が算出される。細い破線は誤差閾値である。細い破線は、エンジン音と雑音との疆界を示す線であり、位相が2つの細い破線の内側にあれば、エンジン音の位相を示し、外側にあれば、雑音の位相を示す。算出された位相曲線との誤差を算出すると、各周波数信号と曲線との誤差は大きく、閾値から大きく外れた点が多いことが分かる。ここで時刻t6、t7、t8、t9の周波数信号の位相に着目すると、他の時刻の位相と大きく外れていることが分かる。これは、位相が0~2πの周期でトーラス状になっていることに起因する。そこで、トーラス状に起因する現象を考慮して位相曲線を算出することとしてもよい。これにより、他の時刻の位相と大きく外れた位相を補正することができ、位相の時間変化を精度良く曲線近似することが可能となる。例えば、前、後、又は前後N個の位相を用いて位相を補正することとしてもよい。図13bにおいて、例えば、時刻t1からt5例えば、N=5の位相の平均を算出する。平均位相がψ=2π×10/360であったとする。次に時刻t6の位相がψ6=2π×170/360であったとする。しかし位相はトーラスになっているため、ψ6=2π×170/360±2πの可能性がある。なお、実際は±2π×mmは自然数の可能性があるが、ここではm=1の場合のみを考慮することとする。なお、周波数が大きく変化する場合、位相変化も大きくなるため、分析する音に応じて考慮するmを可変としてもよい。また、平均を算出するための位相の選択時刻は、時刻t1~t5に限定されるものではなく、任意の時刻を用いることが可能である。次に、時刻t6の位相ψ6を平均位相ψとの誤差が小さい値に補正する。図13bの場合、ψ6=2π×170/360-2πとなる。同様に時刻t7の位相を時刻t2~t5の位相と、補正後の時刻t6の位相とを用いて補正する。本例の場合ψ7=ψ7-2πと補正される。同様の処理を時刻t8、t9と行っていく。図13cに補正後の位相を示す。時刻t6、t7、t8、t9の位相が補正されていることが分かる。この補正後の位相情報を用いて位相曲線を算出した場合、太い破線で示される曲線が算出される。図13cの場合、曲線とその閾値内に各周波数信号が含まれるため、エンジン音として適切に抽出されることとなる。なお、位相の補正手法はこれに限ったものではない。例えば、まず位相曲線を算出し、算出された形状との誤差が大きい各点に±2πの位相補正を行うこととしてもよい。また、位相がとり得る角度範囲を補正することとしてもよい。以下、図を用いて説明する。図14は位相の補正処理を説明する図である。図14のグラフはいずれも縦軸が位相を、横軸が時間を示している。白い丸印は各時刻における周波数信号の位相を示す。図14aは0から2πを角度範囲とした場合の周波数信号の位相を示す。各位相を元に位相曲線が算出され黒い曲線で示している。図14cは曲線との誤差を元に位相を補正している。具体的には時刻t1の位相に2πを加算する補正がされている。また時刻t8の位相に-2πを加算する補正がされている。一方、図14bは-πからπを角度範囲とした場合の周波数信号の位相を示す。図14aと同様に、各位相を元に位相曲線が算出され黒い曲線で示している。図14dは曲線との誤差を元に位相を補正している。具体的には時刻t10の位相に-2πを加算する補正がされている。図14cの角度範囲の場合の曲線との誤差と、図14dの角度範囲の場合の曲線との誤差を比較した場合、図14cの角度範囲の場合の曲線との誤差が小さくなる。そこで図14cの角度範囲を用いた位相曲線を用いる。このように、角度範囲の制御を行い位相曲線を算出することとしてもよい。これにより、他の時刻の位相と大きく外れた位相を補正することができ、より加減速の判定を精度良く行うことが可能となる。以上説明したように、エンジン回転数が増加する場合には、エンジン音の周波数が時間の経過とともに増加し、エンジン音の周波数信号の位相は加速度的に増加する。一方、エンジン回転数が減少する場合には、エンジン音の周波数が時間の経過とともに減少し、エンジン音の周波数信号の位相は加速度的に減少する。位相が加速度的に増加しているか加速度的に減少しているかは、短い時間範囲に含まれる位相から判断することが可能である。このため、実施の形態1によると、自車の周囲に存在する周囲車両のエンジン回転数の増減をリアルタイムで判断することができる。これにより、周囲車両が加速しているのか減速しているのかをリアルタイムで判断することができる。実施の形態2次に、実施の形態2に係る雑音除去装置について説明する。この雑音除去装置は、請求の範囲の回転数増減判定装置に対応する。実施の形態1では、エンジン音を受け付け、位相の時間変化をもとに加減速の判定を行う手法について説明を行った。本実施の形態では、エンジン音と風などの雑音との混合音を受け付け、混合音からエンジン音を抽出し、位相の時間変化をもとに加減速の判定を行う手法について説明を行う。図15及び図16は、本発明の実施の形態2における雑音除去装置の構成を示すブロック図である。図15において、雑音除去装置1500は、マイクロホン2400と、DFT分析部2402と、雑音除去処理部1504と、加減速判定部3006jとを含む。DFT分析部2402は、図7に示したDFT分析部3002と同様の処理を行う。このため、その詳細な説明はここでは繰り返さない。以下では、DFT分析部2402によって求められた周波数帯域の個数をMとして、それらの周波数帯域を指定する番号を記号jj=1~Mで表すこととする。雑音除去処理部1504は、位相補正部1501jj=1~Mと、抽出音判定部1502jj=1~Mと、音抽出部1503jj=1~Mとを含む。音抽出部1503jは、請求の範囲の音響信号識別手段に対応する。位相補正部1501jj=1~Mは、DFT分析部2402が求めた周波数帯域jの周波数信号に対して、時刻tの周波数信号の位相をψtラジアンとするときに、ψ ́t=mod2πψt-2πftfは分析周波数に位相を補正する。抽出音判定部1502jj=1~Mは、所定の時間幅において、分析の対象とする時刻の位相補正された周波数信号から位相の時間変化を近似する位相曲線近似曲線を算出し、算出された位相曲線と分析対象とする時刻の位相との誤差を算出する。このとき、位相距離位相曲線と分析対象とする時刻の位相との誤差を求めるときに用いた周波数信号の数は第1のしきい値以上の数から構成されている。このとき位相距離はψ ́tを用いて計算する。音抽出部1503jj=1~Mは、抽出音判定部1502jj=1~Mが計算した誤差位相距離をもとに、誤差が第2の閾値以下になる周波数信号を抽出音として抽出する。加減速判定部3006jj=1~Mは、音抽出部1503jj=1~Mにより抽出されたエンジン音についてのみ、位相曲線算出部3005jj=1~Mが算出した位相曲線から、位相の増加量をもとに、エンジン回転数の増減、つまり車両の加減速を判定する。これらの処理を、所定の時間幅を時間方向に移動させながら行うことにより、時間‐周波数領域ごとに抽出音の周波数信号2408を取り出すことができる。そして、抽出されたエンジン音における位相曲線の形状具体的には凸の向きをもとに加減速判定部3006jにおいて加減速を判定する。つまり、加減速判定部3006jj=1~Mは、音抽出部1503jj=1~Mにより抽出されたエンジン音についてのみ、位相曲線算出部3005jj=1~Mが算出した位相曲線から、位相の増加量をもとに加減速を判定する。図16に、抽出音判定部1502jj=1~Mの構成を示すブロック図を示す。抽出音判定部1502jj=1~Mは、周波数信号選択部1600jj=1~Mと、位相距離判定部1601jj=1~Mと、位相曲線算出部1602jj=1~Mとから構成される。位相距離判定部1601jは、請求の範囲の誤差算出手段に対応する。周波数信号選択部1600jj=1~Mは、所定の時間幅において、位相補正部1501jj=1~Mが位相補正した周波数信号の中から位相曲線の算出及び位相距離を計算するのに用いる周波数信号を選択する。位相曲線算出部1602jj=1~Mは、周波数信号選択部1600jj=1~Mが選択した周波数信号の補正された位相ψ ́tを用いて時間経過に伴い位相が変化する位相形状を二次曲線として計算する。そして位相距離判定部1601jj=1~Mは、位相曲線算出部1602jj=1~Mが算出した位相曲線と分析対象とする時刻の補正後の位相ψ ́tとの位相距離を判定する。次に、以上のように構成された雑音除去装置1500の動作について説明する。以下では、j番目の周波数帯域について説明を行うが、他の周波数帯域についても同様の処理が行われる。ここでは、周波数帯域の中心周波数と分析周波数位相距離を求めるψ ́t=mod2πψt-2πftにおける周波数fであって、周波数fに抽出音が存在するか否かを判定することになるとが一致する場合を例にして説明を行う。他の方法として、周波数帯域を含む複数の周波数を分析周波数として抽出音の判定を行ってもよい。この場合は、中心周波数の周辺の周波数に抽出音が存在するか否かを判定することができる。図17及び図18は、雑音除去装置1500の動作手順を示すフローチャートである。初めに、マイクロホン2400は、外部からの混合音2401を集音し、集音した混合音を、DFT分析部2402に出力するS200。DFT分析部2402は、混合音2401を受け付けて、混合音2401に対してフーリエ変換処理を施し、混合音2401の周波数信号を周波数帯域jごとに求めるステップS300。次に、位相補正部1501jは、DFT分析部2402が求めた周波数帯域jの周波数信号に対して、時刻tの周波数信号の位相をψtラジアンとするときに、位相ψtを、位相ψ ́t=mod2πψt-2πftfは分析周波数に変換することで位相補正を行うステップS1700j。ここで、本発明において位相を用いる理由について図を用いて説明する。図19はDFT分析におけるパワーと位相について説明する図である。図19aは図3と同様に、自動車のエンジン音をDFT分析したスペクトログラムである。図19bは時刻t1から所定の時間窓幅のハニング窓を用いて複素空間上に周波数信号601を表したものである。周波数f1、f2、f3等、各周波数のパワーと位相が算出される。周波数信号601の長さがパワーを示し、周波数信号601と実軸とのなす角が位相を示している。そして、図19aにt1、t2、t3と示すように、時間シフトを行いながら各時刻における周波数信号を求めることとなる。なお、一般的にスペクトログラムは各時刻における各周波数のパワーを示すのみであり、位相については省略されている。図3や図19aに示すスペクトログラムも同様に、DFT分析したパワーの大きさのみを表示したものである。図19cには、図19aにおける所定の周波数例えば周波数f4での時間方向の位相の変動を示す。横軸は時間を表している。縦軸は周波数信号の位相を表しており、0~2πラジアンの間の値で示される。図19dには、図19aにおける、所定の周波数例えば周波数f4のパワーの時間変動を示す。横軸は時間軸である。縦軸は周波数信号の大きさパワーを表している。図20は、風などの雑音がある場合の自動車のエンジン音を説明する図である。図20aは図3と同様、自動車のエンジン音をDFT分析部したスペクトログラムである。縦軸は周波数を、横軸は時間をそれぞれ示しており、色の濃度は周波数信号のパワーの大きさを示している。しかし、図3と異なり、風などの雑音が含まれているため、エンジン音以外の周波数にも色の濃い部分が存在し、エンジン音なのか風雑音なのかパワーのみではまったく分からない状態となっている。図20bは、時刻t2におけるエンジン音部分のある周波数f4の所定時間のパワーの推移を示すグラフである。風雑音の影響でパワーが乱れていることが分かる。図20cは、時刻t3におけるエンジン音がない部分である周波数f4の所定時間のパワーの推移を示すグラフである。非定常のパワーが存在することが分かる。また、図20b及び図20cを比較しても、パワーのみでは風雑音なのか、エンジン音が存在するのかまったく区別がつかないことが分かる。そこで、本発明では、エンジン音を抽出すべく、位相の時間変化を用いる。まず、エンジン音の位相特性について説明する。エンジンは、所定数のシリンダーがピストン運動を行うことで駆動系を回転させている。そして、車両から発せられるエンジン音は、このエンジンの回転に依存した音と、エンジンの回転には依存しない固定振動音又は非周期音とからなる。特に、車両の外部から検知できる主な音は、エンジンの回転に依存した周期音であり、本発明ではこのエンジンの回転に依存する周期音をエンジン音として抽出する。図3に示すようにエンジン音は、回転数が変化することで周波数が変化していることが分かる。ここで周波数の変化に着目すると、周波数がランダムに変化したり、離散的に飛ぶことはほとんどなく、所定の時間間隔でみるとほぼ時刻に応じて周波数が変化していることが分かる。したがって、エンジン音は上記式4で示すような区分線形で近似できる。具体的には所定の時間区間で見た場合、時刻tにおける周波数fは、初期値f0から時刻tに比例比例係数Aして増減する線分で線形近似できると考えられる。そして周波数fを上記式4で表した場合、時刻tにおける位相ψは、上記式5で表すことができる。位相補正部1501jは、位相の時間変化の近似処理を容易にするための位相補正処理を行う。つまり、位相補正部1501jは、図19cに示されている周波数信号の位相ψtを位相ψ ́t=mod2πψt-2πftfは分析周波数に変換することで位相補正を行う。この位相補正処理の詳細は、図10及び図11を参照して説明した実施の形態1に係る位相補正部3003jが実行する位相補正処理と同様である。このため、その詳細な説明はここでは繰り返さない。再度、図17を参照して、抽出音判定部1502jは、位相補正部1501jが求めた補正後の位相情報を用いて、位相の形状を算出する。そして分析の対象とする時刻の周波数信号と、分析の対象とする時刻とは異なる複数の時刻における周波数信号との位相距離誤差を求めるステップS1701j。図18は、抽出音の周波数信号を判定する処理ステップS1701jの動作手順を示すフローチャートである。周波数信号選択処理S1800j及び位相曲線算出処理S1801jは、実施の形態1で説明した周波数信号選択処理図8のS103j及び位相曲線算出処理S104jとそれぞれ同様である。このため、その詳細な説明は、ここでは繰り返さない。図18を参照して、位相距離判定部1601jは、位相曲線算出部1602jが算出した形状からの位相距離を計算するステップS1802j。この例では、位相距離誤差E0は位相の差分誤差であり、で求める。なお、分析対象とする点を除いて形状を算出し、算出した形状と分析対象とする点との位相の差を計算することとしてもよい。この計算方法によると、分析対象とする点に算出した形状から著しく外れるノイズが含まれる場合、より精度良く形状を近似することが可能となる。なお、本例では分析対象である時刻t0に対し、時刻t1、t2、t3、t4、t5の位相から位相の形状を算出している。例えば時刻t2を分析対象とする場合つまり時刻t0 ́とした場合、あらたに時刻t1 ́、t2 ́、t3 ́、t4 ́、t5 ́の位相から位相曲線を算出して誤差を算出してもよいし、既に算出されたt0、t1、t2、t3、t4、t5による位相曲線から誤差を算出することとしてもよい。すなわち、すでに算出された位相曲線を用いた誤差は、となる。この方法によると、位相曲線の算出回数が減るため、計算量の削減の効果を奏する。さらに、分析対象を所定の区間とし、誤差の平均によって、分析対象区間中の全ての周波数信号が誤差か否かを区別することとしてもよい。例えば誤差の平均は以下の式23であらわすことができる。再度、図17を参照して、音抽出部1503jは、位相距離誤差が閾値以下である分析の対象とする周波数信号の各々を抽出音として抽出するステップS1702j。そして、加減速判定部3006jは、抽出されたエンジン音部分の、位相曲線の形状凸の向きをもとに加減速を判定するステップS105j。図21は、位相距離を求める所定の時間幅96msにおける、混合音の周波数信号の位相補正された位相ψ ́tを模式的に示した図である。横軸は時間tを示し、縦軸は位相補正された位相ψ ́tを示す。黒丸印は分析の対象とする周波数信号の位相を示し、白丸印は位相曲線を求めるのに用いられた周波数信号の位相を示す。太い破線1101は算出された位相曲線である。位相補正された各点をもとに位相曲線として二次曲線が算出されていることが分かる。細い破線1102は誤差の閾値例えば20度とするを示す。つまり、上側の破線1102は、破線1101を上に閾値分だけシフトさせたものであり、下側の破線1102は、破線1101を下に閾値分だけシフトさせたものである。分析の対象とする周波数信号の位相が2つの破線1102内に収まっていれば、その周波数信号は抽出音周期音の周波数信号であると判定され、2つの破線1102内に収まっていなければ、その周波数信号は雑音の周波数信号であると判定される。図21aにおいて、黒丸印で示される分析対象とされる周波数信号の位相は、位相の二次曲線との誤差が閾値未満である。このため、音抽出部1503jは、その周波数信号を抽出音の周波数信号として抽出する。図21bにおいて、黒丸印で示される分析対象とされる周波数信号の位相のそれぞれは、位相の二次曲線との誤差が閾値以上である。このため、音抽出部1503jは、これらの周波数信号を抽出音の周波数信号として抽出せずに、雑音として除去する。図22は、本実施の形態に示す手法によるエンジン音抽出処理を説明する図である。式3に示すようにエンジン音を区分線形で近似した場合、位相は式12に示すように二次曲線で近似できる。図22aは、図19aに示したのと同じスペクトログラムである。図22b~図22eは、図22aに四角印で示す4つの領域における周波数信号を示すグラフである。4つの領域のそれぞれは、1つの周波数帯域を有する領域である。図22b~図22eに示すグラフにおいて、横軸は時間を、縦軸は位相をそれぞれ示す。白い丸印は実際の分析された周波数信号を示し、太い破線は算出された近似曲線を示す。また、細い破線は、抽出音と雑音との閾値を示す。図22bは、エンジンの回転数が下がっている、すなわち時間周波数空間で周波数の時間変化が負の傾きの一次式で近似できるエンジン音部分の補正後の位相を示したグラフである。位相曲線は上に凸の形状をしていることが分かる。そして分析された各周波数信号はほぼ閾値以内に収まっていることが分かる。図22cは、エンジンの回転数が上がっている、すなわち時間周波数空間で周波数の時間変化が正の傾きの一次式で近似できるエンジン音部分の補正後の位相を示したグラフである。位相曲線は下に凸の形状をしていることが分かる。そして分析された各周波数信号はほぼ閾値以内に収まっていることが分かる。図22dは、エンジンの回転数が一定、すなわち時間周波数空間で周波数が変化しない二次の係数が0で近似できるエンジン音部分の補正後の位相を示したグラフである。位相曲線は二次の項が0となり、一次直線の形状をしていることが分かる。そして分析された各周波数信号はほぼ閾値以内に収まっていることが分かる。このことから、二次曲線による表現は周波数が変化しないエンジン音も含めて識別できることがわかる。図22eは、風雑音部分の補正後の位相を示したグラフである。風雑音の周波数信号の位相は、ばらばらであるため、二次の近似曲線を算出したとしても、その曲線からの誤差が大きく、閾値以内である信号部分はほとんどないことが分かる。このように、算出された曲線と曲線からの誤差をもとに風雑音とエンジン音を区別することができる。図23は、位相曲線からの誤差を説明する図である。横軸はエンジン音、雨音、風雑音の各音響信号を示している。縦軸は、本手法による位相曲線からの誤差の平均及び分散を示している。つまり、縦軸における線分の幅がとり得る誤差の範囲を示し、ひし形が平均値を示している。例えば、エンジン音の場合、誤差の範囲は、1度から18度の間であり、誤差の平均は10度である。分析条件は以下となる。8kHzでサンプリングされた各音声に対し、256点32msで周波数分析を行い、768点96msを区間として位相曲線の算出を行った。そして位相曲線からの誤差の平均と分散を算出した。図23より、エンジン音は平均値が10度と位相曲線からの誤差が小さいのに対し、雨音は68度、風雑音は48度と位相の位相曲線との誤差が大きいことわかる。このようにエンジン音のような周期音と、風雑音のような非周期音とは位相曲線からの誤差に大きな違いがあることが分かる。本例の場合、例えば閾値を20度等に設定し、閾値以下をエンジン音として適切に抽出できることとなる。図24は音識別を説明する図である。各グラフの横軸は時間を示しており、縦軸は周波数を示している。図24aは風雑音とエンジン音が混合した音を周波数分析したスペクトログラムである。色の濃さはパワーの大きさを表しており、色が濃いほどパワーが大きいことを示している。分析条件は以下となる。8kHzでサンプリングされた音声に対し、512点で周波数分析を行い、1536点を区間として位相曲線の算出を行った。そして位相曲線からの誤差の閾値を20度としてエンジン音の抽出を行った。図24bは本実施の形態における手法で風雑音とエンジン音を識別したグラフである。黒い部分がエンジン音として抽出された部分である。図24aでは風などの影響により雑音が混合されているため、いったいどの部分がエンジン音であるか抽出することは困難である。しかしながら、本実施の形態における方法でエンジン音を抽出した場合、適切にエンジン音が抽出できていることを示している。特にエンジンの回転数が急激に増加する部分や、減少する部分、定常音ともに抽出できていることが分かる。以上説明したように、本実施の形態によると、時間-周波数領域ごとにエンジン音と、風雑音、雨音、暗騒音などとを区別することができる。このため、雑音の影響を除去して、エンジン音についてのみエンジン回転数の増加又は減少周囲車両の加速度の増加又は減少を判定することができる。このため、判定の精度を向上させることができる。実施の形態3次に、実施の形態3に係る車両検知装置について説明する。この車両検知装置は、請求の範囲の回転数増減判定装置に対応する。実施の形態3に係る車両検知装置は、複数のマイクロホンから入力される各々の混合音から、エンジン音抽出音の周波数信号を判定し、音の到達時間差より車両の到達方向を算出し、運転者に接近車両の方向及び存在を知らせるものである。その際、加速している接近車両のみについて方向及び存在を知らせ、減速又は等速走行している接近車両については方向及び存在を知らせない。図25及び図26は、本発明の実施の形態3における車両検知装置の構成を示すブロック図である。図25において、車両検知装置4100は、マイクロホン41071と、マイクロホン41072と、DFT分析部1100と、車両検知処理部4101と、加減速判定部3006jj=1~Mと、方向検知部4108とを備える。車両検知処理部4101は、位相補正部4102jj=1~Mと、抽出音判定部4103jj=1~Mと、音抽出部4104jj=1~Mと、方向検知部4108と、提示部4106とを含む。また、図26において、抽出音判定部4103jj=1~Mは、位相距離判定部4200jj=1~Mと、位相曲線算出部4201jj=1~Mと、周波数信号選択部4202jj=1~Mとを含む。位相距離判定部4200jは、請求の範囲の誤差算出手段に対応する。図25において、マイクロホン41071は、外部からの混合音24011を集音する。マイクロホン41072は、外部からの混合音24012を集音する。この例では、マイクロホン41071とマイクロホン41072はそれぞれ自車両の左前と右前のバンパーに設置されている。これらの混合音の各々は、例えば8kHzでサンプリングされた車両のエンジン音と風雑音とから構成されている。なお、サンプリング周波数は8kHzに限定されるものではない。DFT分析部1100は、入力された混合音24011と混合音24012の各々に対して離散フーリエ変換処理を施し、混合音24011と混合音2の周波数信号を求める。ここでのDFTの時間窓幅は256点38msである。以下では、DFT分析部1100から求められた周波数帯域の個数をMとして、それらの周波数帯域を指定する番号を記号jj=1~Mで表すこととする。この例では、車両のエンジン音が存在する10Hz~500Hzの周波数帯域を10Hz間隔ごとに分割してM=50、周波数信号を求める。位相補正部4102jj=1~Mは、DFT分析部1100が求めた周波数帯域jj=1~Mの周波数信号に対して、時刻tの周波数信号の位相をψtラジアンとするときに、ψ ́ ́t=mod2πψt-2πf ́tf ́は周波数帯域の周波数に位相を補正する。なお、この例は、ψtを分析周波数で補正するのではなく、周波数信号を求めた周波数帯域の周波数f ́で補正を行う。抽出音判定部4103jj=1~Mは、所定の時間幅において、分析の対象とする時刻の位相補正された周波数信号から位相曲線を算出し、算出された位相曲線をもとに抽出音の判定を行う。このとき、位相距離を求めるときに用いた周波数信号の数は第1のしきい値以上の数から構成されている。ここでは所定の時間幅を96msとする。なお、このとき位相距離は補正後の位相ψ ́ ́tを用いて計算する。抽出音判定部4103jj=1~Mが実施する処理は、実施の形態2に示した抽出音判定部1502jj=1~Mが実施する処理と同様である。このため、その詳細な説明は繰り返さない。図26に、抽出音判定部4103jj=1~Mの構成を示すブロック図を示す。抽出音判定部4103jj=1~Mは、位相距離判定部4200jj=1~Mと、位相曲線算出部4201jj=1~Mと、周波数信号選択部1600jj=1~Mとから構成される。周波数信号選択部4202jj=1~Mは、所定の時間幅において、位相補正部4102jj=1~Mが位相補正した周波数信号の中から位相曲線の算出及び位相距離を計算するのに用いる周波数信号を選択する。周波数信号選択部4202jj=1~Mが実施する処理は、実施の形態2に示した周波数信号選択部1600jj=1~Mが実施する処理と同様である。このため、その詳細な説明は繰り返さない。位相曲線算出部4201jj=1~Mは、周波数信号の補正された位相ψ ́ ́tを用いて時間経過に伴い位相が変化する位相形状を曲線として計算する。位相曲線算出部4201jj=1~Mが実施する処理は、実施の形態2に示した位相曲線算出部1602jj=1~Mが実施する処理と同様である。このため、その詳細な説明は繰り返さない。そして位相距離判定部4200jj=1~Mは、位相曲線算出部4201jj=1~Mが算出した位相曲線との位相距離が第2のしきい値以下か否かを判定する。具体的には位相曲線を算出する区間を768点96msとして位相曲線を算出し、位相距離を算出する。位相距離判定部4200jj=1~Mによる位相曲線の算出方法と位相距離誤差の算出方法とは、実施の形態2に示した位相距離判定部1601jj=1~Mのそれらと同様である。このため、その詳細な説明は繰り返さない。次に、音抽出部4104jj=1~Mは、抽出音判定部4103jj=1~Mが判定した位相距離をもとにエンジン音を抽出する。具体的には、誤差の閾値を20度とし、閾値以下をエンジン音として抽出する。音抽出部4104jj=1~Mが実施する処理は、実施の形態2に示した音抽出部1503jj=1~Mが実施する処理と同様である。このため、その詳細な説明は繰り返さない。なお、音抽出部4104jj=1~Mは、さらに、エンジン音が抽出されたときに、抽出音検知フラグ4105を出力する。再度、図25を参照して、加減速判定部3006jは、抽出音検知フラグ4105の有無に基づいて、音抽出部4104jにより抽出されたエンジン音についてのみ、位相曲線算出部4201jが算出した位相曲線から、位相の増加量をもとに、エンジン回転数の増減、つまり車両の加減速を判定する。方向検知部4108は、抽出されたエンジン音の時間-周波数領域に対し、車両の存在する方向を特定する。例えば到達時間差をもとに車両の方向を検知する。例えばいずれか1つのマイクロホンにおいてエンジン音が抽出された場合、両マイクロホンを用いて車両の存在する方向を特定する。風雑音は両マイクロホンに対して均一ではなく、一方のマイクロホンにのみ風雑音が存在し、他方には存在しない場合もあるからである。なお、両マイクロホンにおいてエンジン音が抽出された場合に方向を特定することとしてもよい。また、方向検知部4108は、加減速判定部3006jによりエンジン回転数が増加していると判定された場合車両が加速していると判定された場合にのみ、車両の方向の検知結果を出力する。マイクロホン41071及びマイクロホン41072の間隔をdmとする。エンジン音が自車両に対して方角θラジアンから検出されるとする。マイク間での到達時間差をΔtsとし、音速をcm/sとすると、方向θラジアンは以下の式24で表すことができる。最後に、車両検知装置4100に接続された提示部4106は、方向検知部4108で検知された車両の方向を運転者に知らせる。例えば、提示部4106は、どちらの方向から車両が来ているのかをディスプレイに表示するようにしても良い。なお、方向検知部4108からは、エンジン回転数が増加していると判定された車両の方向のみが出力されるため、提示部4106は、加速している車両の方向のみを運転者に知らせることができる。車両検知装置4100及び提示部4106は、これらの処理を、所定の時間幅を時間方向に移動させながら行う。次に、以上のように構成された車両検知装置4100の動作について説明する。以下では、j番目の周波数帯域周波数帯域の周波数はf ́について説明を行う。図27、図28は、車両検知装置4100の動作手順を示すフローチャートである。初めに、マイクロホン41071及び41072は、外部からの混合音2401をそれぞれ集音し、集音した混合音を、DFT分析部2402に出力するステップS201。DFT分析部1100は、混合音24011と混合音24012を受け付けて、混合音24011と混合音24012のそれぞれに対して離散フーリエ変換処理を施し、混合音24011と混合音24012の周波数信号を求めるステップS300。次に、位相補正部4102jは、DFT分析部1100が求めた周波数帯域j周波数f ́の周波数信号に対して、時刻tの周波数信号の位相をψtラジアンとするときに、位相ψtを位相ψ ́ ́t=mod2πψt-2πf ́tf ́は周波数帯域の周波数に変換することで位相補正を行うステップS4300j。次に、抽出音判定部4103j位相距離判定部4200jは、混合音混合音24011、混合音24012ごとに、所定の時間幅における第1のしきい値以上の数から構成されている位相補正された周波数信号第1のしきい値は、所定の時間幅における時刻の周波数信号の80の数であるの位相ψ ́ ́tを用いて、分析周波数fを設定して、設定された分析周波数fを用いて位相距離を求めるステップS4301j。図28を用いて、ステップ4301jの処理について詳細に説明する。初めに周波数信号選択部4202jは、位相補正部4102jが求めた所定の時間幅における位相補正された周波数信号から、位相曲線算出部4201jが位相の形状の計算に用いる周波数信号を選択するステップS1800j。そして、位相曲線算出部4201jが、位相曲線を算出するステップS1801j。次に、位相距離判定部4200jは、位相曲線算出部4201jが算出した形状と分析対象とする時刻の補正後の位相との位相距離を計算するステップS1802j。再度、図27を参照して、音抽出部4104jは、位相距離が第2のしきい値以下になる所定の時間幅における周波数信号をエンジン音の周波数信号に判定するステップS4302j。なお、音抽出部4104jj=1~Mは、さらに、エンジン音が抽出されたときに、抽出音検知フラグ4105を出力する。加減速判定部3006jは、抽出音検知フラグ4105の有無に基づいて、音抽出部4104jにより抽出されたエンジン音についてのみ、位相曲線算出部4201jが算出した位相曲線から、位相の増加量をもとに加減速を判定するS4303j。方向検知部4108は、音抽出部4104jで抽出されたエンジン音の時間-周波数領域に対し、車両の存在する方向を特定し、車両のエンジン回転数が増加していると判定された場合車両が加速していると判定された場合にのみ、車両の方向の検知結果を提示部4106に出力する。提示部4106は、方向検知部4108で検知された車両の方向を運転者に知らせるステップS4304。以上説明したように、実施の形態3に係る車両検知装置によると、エンジン回転数が増加していると判定された場合にのみ、音源方向の検知結果を出力することができる。このため、周囲車両が加速しながら接近してくるなどのように、特に危険な場合にのみ、周囲車両が接近してくる方向を運転者に提示することが可能となる。以上、本発明の実施の形態に係る加減速判定装置、雑音除去装置及び車両検知装置について説明したが、本発明は、これらの実施の形態に限定されるものではない。例えば、上述の実施の形態ではエンジン音の抽出を例に挙げて説明を行ったが、本発明が抽出対象とする音はエンジン音に限定されるものではなく、例えば、人間もしくは動物の音声又はモーター音などのように周期音であれば本発明を適用可能である。また、音抽出部は、周波数信号ごとに周期音か雑音かを判断したが、所定の時間幅ごとにその時間幅に含まれる周波数信号が周期音か雑音かを判断するようにしても良い。例えば、図21を参照して、音抽出部は、所定の時間幅ごとに、その時間幅に含まれる周波数信号の位相に対する、位相曲線算出部が求めた二次曲線との誤差が閾値未満となる位相の割合が所定の割合以上となる場合に、その時間幅に含まれる周波数信号の全てが周期音と判断し、上記割合が所定の割合未満となる場合に、その時間幅に含まれる周波数信号の全てが雑音と判断しても良い。また、加速度判定部は、時間の経過に伴う位相の変化値が所定の閾値以下の場合にのみ、エンジン回転数の増加又は減少周囲車両の加減速を判定するようにしても良い。例えば、連続する時刻間での位相の差の絶対値が、所定の閾値以下の場合にのみ、上記判定をするようにしても良い。周囲車両がギアチェンジしたような場合には、位相が急激に変化する。このため、そのような場合を除外して、上記判定を行うことができるようになる。また、実施の形態3において、加速している接近車両のみについてその方向を提示することとしたが、加速している接近車両と等速走行している接近車両について方向を提示し、減速している接近車両については方向を提示しないようにしてもよい。また、上記の各装置は、具体的には、マイクロプロセッサ、ROM、RAM、ハードディスクドライブ、ディスプレイユニット、キーボード、マウスなどから構成されるコンピュータシステムとして構成されても良い。RAM又はハードディスクドライブには、コンピュータプログラムが記憶されている。マイクロプロセッサが、コンピュータプログラムに従って動作することにより、各装置は、その機能を達成する。ここでコンピュータプログラムは、所定の機能を達成するために、コンピュータに対する指令を示す命令コードが複数個組み合わされて構成されたものである。さらに、上記の各装置を構成する構成要素の一部又は全部は、1個のシステムLSILargeScaleIntegration:大規模集積回路から構成されているとしても良い。システムLSIは、複数の構成部を1個のチップ上に集積して製造された超多機能LSIであり、具体的には、マイクロプロセッサ、ROM、RAMなどを含んで構成されるコンピュータシステムである。RAMには、コンピュータプログラムが記憶されている。マイクロプロセッサが、コンピュータプログラムに従って動作することにより、システムLSIは、その機能を達成する。さらにまた、上記の各装置を構成する構成要素の一部又は全部は、各装置に脱着可能なICカード又は単体のモジュールから構成されているとしても良い。ICカード又はモジュールは、マイクロプロセッサ、ROM、RAMなどから構成されるコンピュータシステムである。ICカード又はモジュールは、上記の超多機能LSIを含むとしても良い。マイクロプロセッサが、コンピュータプログラムに従って動作することにより、ICカード又はモジュールは、その機能を達成する。このICカード又はこのモジュールは、耐タンパ性を有するとしても良い。また、本発明は、上記に示す方法であるとしても良い。また、これらの方法をコンピュータにより実現するコンピュータプログラムであるとしても良いし、前記コンピュータプログラムからなるデジタル信号であるとしても良い。さらに、本発明は、上記コンピュータプログラム又は上記デジタル信号をコンピュータ読み取り可能な不揮発性の記録媒体、例えば、フレキシブルディスク、ハードディスク、CD-ROM、MO、DVD、DVD-ROM、DVD-RAM、BDBlurayDisc登録商標、半導体メモリなどに記録したものとしても良い。また、これらの不揮発性の記録媒体に記録されている上記デジタル信号であるとしても良い。また、本発明は、上記コンピュータプログラム又は上記デジタル信号を、電気通信回線、無線又は有線通信回線、インターネットを代表とするネットワーク、データ放送等を経由して伝送するものとしても良い。また、本発明は、マイクロプロセッサとメモリを備えたコンピュータシステムであって、上記メモリは、上記コンピュータプログラムを記憶しており、上記マイクロプロセッサは、上記コンピュータプログラムに従って動作するとしても良い。また、上記プログラム又は上記デジタル信号を上記不揮発性の記録媒体に記録して移送することにより、又は上記プログラム又は上記デジタル信号を上記ネットワーク等を経由して移送することにより、独立した他のコンピュータシステムにより実施するとしても良い。さらに、上記実施の形態及び上記変形例をそれぞれ組み合わせるとしても良い。今回開示された実施の形態はすべての点で例示であって制限的なものではないと考えられるべきである。本発明の範囲は上記した説明ではなくて請求の範囲によって示され、請求の範囲と均等の意味及び範囲内でのすべての変更が含まれることが意図される。"