▶ 学校法人文理学園の特許一覧 ▶ 株式会社ＡＫシステムの特許一覧
(19)【発行国】日本国特許庁(JP)
(12)【公報種別】特許公報(B2)
(11)【特許番号】
(24)【登録日】2024-08-01
(45)【発行日】2024-08-09
(54)【発明の名称】伝達特性測定装置および伝達特性測定プログラム
(51)【国際特許分類】
G01S 15/08 20060101AFI20240802BHJP
G01S 7/52 20060101ALI20240802BHJP
B64C 39/02 20060101ALI20240802BHJP
B64D 45/00 20060101ALI20240802BHJP
G01S 5/30 20060101ALI20240802BHJP
【ＦＩ】
G01S15/08
G01S7/52 Z
B64C39/02
B64D45/00 Z
G01S5/30
【請求項の数】
5
(21)【出願番号】P 2020071404
(22)【出願日】2020-04-13
(65)【公開番号】P2021167778
(43)【公開日】2021-10-21
【審査請求日】2023-04-13
(73)【特許権者】
【識別番号】596056519
【氏名又は名称】学校法人文理学園
(74)【代理人】
【識別番号】100197642
【弁理士】
【氏名又は名称】南瀬 透
(74)【代理人】
【識別番号】100099508
【弁理士】
【氏名又は名称】加藤 久
(74)【代理人】
【識別番号】100219483
【弁理士】
【氏名又は名称】宇野 智也
(74)【代理人】
【識別番号】100182567
【弁理士】
【氏名又は名称】遠坂 啓太
(73)【特許権者】
【識別番号】501008462
【氏名又は名称】株式会社ＡＫシステム
(74)【代理人】
【識別番号】100197642
【弁理士】
【氏名又は名称】南瀬 透
(74)【代理人】
【識別番号】100099508
【弁理士】
【氏名又は名称】加藤 久
(74)【代理人】
【識別番号】100219483
【弁理士】
【氏名又は名称】宇野 智也
(72)【発明者】
【氏名】福島 学
(72)【発明者】
【氏名】森竹 隆広
(72)【発明者】
【氏名】石井 秀樹
【審査官】▲高▼場 正光
(56)【参考文献】
【文献】米国特許出願公開第２０１５／０１６０６５８（ＵＳ，Ａ１）
【文献】山下涼介 他，“伝達特性に基づく距離と反射特性計測に関する一検討”，日本音響学会 2020 年 春季研究発表会講演論文集［CD-ROM］，2020年03月02日，Article 1-Q-27，Pages 519-520，ISSN 1880-7658
【文献】福島学 他，“観測信号のみによる距離推定における誤推定抑制のためのＤＬＲ－ＣＳ法の活用の一検討”，日本音響学会 2010 年 春季研究発表会講演論文集［CD-ROM］，2010年03月01日，Article 2-P-11，Pages 815-816，ISSN 1880-7658
【文献】近藤善隆 他，“距離スペクトルの位相変化に着目した反射物微振動検出に関する一検討”，日本音響学会 2011 年 秋季研究発表会講演論文集［CD-ROM］，2011年09月13日，Article 1-P-10，Pages 731-732
(58)【調査した分野】(Int.Cl.，ＤＢ名)
Ｇ０１Ｓ ５／１８ － Ｇ０１Ｓ ５／３０
Ｇ０１Ｓ ７／５２ － Ｇ０１Ｓ ７／６４
Ｇ０１Ｓ １５／００ － Ｇ０１Ｓ １５／９６
ＪＳＴＰｌｕｓ（ＪＤｒｅａｍＩＩＩ）
(57)【特許請求の範囲】
【請求項1】
時間ｔ０から時間ｔ１までの音源からの音源信号である入力信号と前記音源信号が反射して観測点で観測される観測信号である出力信号とを切り出す切出手段と、前記入力信号と前記出力信号とから時間ごとの伝達特性を求める伝達特性演算手段と、
時間ｔ１より以前の時間ｔ３から、時間ｔ１以後の時間ｔ２までの入力信号に０を与え、前記伝達特性から因果性が正しい因果成分と因果性が前後した非因果成分とを分離させる歪分離手段と、
前記歪分離手段により得られた前記非因果成分からノイズレベルを推定するノイズ推定手段と、
前記ノイズ推定手段により推定されたノイズレベルを前記因果成分から除去した伝達特性を求める打ち切り手段と、
前記打ち切り手段により除去された伝達特性の減衰量を推定して増幅させた伝達特性を求める減衰補正手段と
を備えた伝達特性測定装置において、
前記伝達特性をフーリエ変換して、観測点からの距離についての距離振幅スペクトルと距離位相スペクトルとを求める距離スペクトル演算手段
を備え、
前記距離振幅スペクトルにあわらわれるピークのうち、前記距離位相スペクトルが正である場合のピークの位置で固体の反射物までの距離情報を取得する一方、前記距離位相スペクトルが負である場合のピークの位置を空気密度による反射とする伝達特性測定装置。
【請求項2】
前記ノイズ推定手段は、前記非因果成分のエネルギー量によりノイズレベルを推定する請求項１記載の伝達特性測定装置。
【請求項3】
前記減衰補正手段は、前記因果成分の減衰量をヒルベルト変換により包絡線抽出を行い、包絡線の傾斜に応じて逆数を乗算することで減衰量を増幅する請求項１または２記載の伝達特性測定装置。
【請求項4】
前記入力信号と前記出力信号は、ドローン装置の回転翼を支持するアームに配置した入力用信号用のマイクと出力信号用のマイクで観測したものであって、
前記入力用信号用のマイクと前記出力信号用のマイクとを１組として前記ドローン装置に３組配置して、
前記距離スペクトル演算手段から得られる各ピークから立体の散布図を得る座標情報取得手段を備えた請求項１から３のいずれかの項に記載の伝達特性測定装置。
【請求項5】
コンピュータを、時間ｔ０から時間ｔ１までの音源からの音源信号である入力信号と前記音源信号が反射して観測点で観測される観測信号である出力信号とを切り出す切出手段、
前記入力信号と前記出力信号とから時間ごとの伝達特性を求める伝達特性演算手段、
時間ｔ１より以前の時間ｔ３から、時間ｔ１以後の時間ｔ２までの入力信号に０を与え、前記伝達特性から因果性が正しい因果成分と因果性が前後した非因果成分とを分離させる歪分離手段、
前記歪分離手段により得られた前記非因果成分からノイズレベルを推定するノイズ推定手段、
前記ノイズ推定手段により推定されたノイズレベルを前記因果成分から除去した伝達特性を求める打ち切り手段、
前記打ち切り手段により除去された伝達特性の減衰量を推定して増幅させた伝達特性を求める減衰補正手段、
前記伝達特性をフーリエ変換して、観測点からの距離についての距離振幅スペクトルと距離位相スペクトルとを求める距離スペクトル演算手段
として機能させ、
前記距離振幅スペクトルにあわらわれるピークのうち、前記距離位相スペクトルが正である場合のピークの位置で固体の反射物までの距離情報を取得する一方、前記距離位相スペクトルが負である場合のピークの位置を空気密度による反射とする
伝達特性測定プログラム。
【発明の詳細な説明】
【技術分野】
【0001】
本発明は、信号源から送信された入力信号が伝送路に入力され、伝送路から出力された出力信号を受信して、入力信号と出力信号とから伝送路の伝達特性を演算する伝達特性測定装置および伝達特性測定プログラムに関するものである。
【背景技術】
【0002】
伝達特性は、信号源からの入力信号と、伝送路を通じて入力した出力信号とに基づいて演算することで求められる。例えば、音源からの音源信号を入力信号、反射物に反射した観測信号を出力信号として伝達特性を求めることで、反射物（観測対象）までの距離や、反射物が硬いものか、または柔らかいものか、空気密度による反射などが特定できる。このように伝達特性を求めることは、ビデオによる画像解析では特定できないものまでも特定することができるため、様々な分野に応用可能であり、有用な技術である。
伝達特性を求める際に、精度を向上させる技術として特許文献１－３および非特許文献１に記載されたものが知れている。
伝達特性を求める際に、精度を向上させる技術として特許文献１－３および非特許文献１に記載されたものが知れている。
【0003】
特許文献１に記載の応答特性測定装置は、所定の時刻ｔ0から所定の第１の時間幅Ｔ1の間の入力信号を切り出す共に、第１の時間Ｔ1よりも短い所定の第２の時間幅をＴ2としたときに時刻ｔ0＋Ｔ2から時刻ｔ0＋Ｔ1までの間の入力信号を零に置換することにより置換入力信号を生成し、また上記所定の時刻ｔ0から上記第１の時間幅Ｔ1の間の出力信号を切り出し、これら置換入力信号と切り出された出力信号とに基づいて上記計測対象の伝達関数もしくはインパルス応答関数を求める、というものである。
【0004】
特許文献２に記載の応答特性測定装置は、出力側の時間窓をｎ＝Ｌ～Ｎ－１＋Ｌにずらし、ｎ＝Ｎ～Ｎ－１＋Ｎの出力信号をｎ＝０～Ｌ－１に移動して、クロススペクトル法によりインパルス応答を求める、というものである。
【0005】
特許文献３に記載の応答特性測定装置は、インパルス応答の継続時間に比べ十分に長い時間窓を用いることができない場合に、求めた伝達関数もしくはインパルス応答関数を、時間窓長が有限であること、および時間窓として用いた関数の双方の影響を取り除くように補正する、というものである。
【0006】
特許文献１－３に記載の応答特性測定装置では、計測結果の歪について記載されたものであり、ここでは、歪を防止し、ノイズレベルの推定に利用するというものである。
【0007】
また、非特許文献１には、推定インパルス応答の振幅を変形せずに、クロススペクトル法を用いて推定伝達関数のノイズ係数を決定する方法が記載されている。
【先行技術文献】
【特許文献】
【0008】
【文献】特開平７－８３９７７号公報
【文献】特開平７－１２０５１６号公報
【文献】特開平７－１２０５１７号公報
【非特許文献】
【0009】
【文献】福島 学，外４名，“Ｆｒ２，Ｅ．３”，［online］，平成１６年４月４日，ＩＣＡ２００４，［令和２年３月２日検索］，インターネット＜URL：https://www.icacommission.org/Proceedings/ICA2004Kyoto/pdf/Fr2.E.3.pdf＞
【発明の概要】
【発明が解決しようとする課題】
【0010】
特許文献１－３と非特許文献１などにより求める伝送路の伝達特性の精度を高めるものではあるが、伝達特性を求める際に歪やノイズが含まれてしまうことから、更に精度を噛める技術が求められている。
【0011】
そこで本発明は、歪やノイズを推定して除去することで、更に精度を高めた伝達特性を求めることが可能な伝達特性測定装置および伝達特性測定プログラム並びにこれらを用いたドローン装置を提供することを目的とする。
【課題を解決するための手段】
【0012】
本発明の伝達特性測定装置は、時間ｔ０から時間ｔ１までの音源からの音源信号である入力信号と前記音源信号が反射して観測点で観測される観測信号である出力信号とを切り出す切出手段と、前記入力信号と前記出力信号とから時間ごとの伝達特性を求める伝達特性演算手段と、時間ｔ１より以前の時間ｔ３から、時間ｔ１以後の時間ｔ２までの入力信号に０を与え、前記伝達特性から因果性が正しい因果成分と因果性が前後した非因果成分とを分離させる歪分離手段と、前記歪分離手段により得られた前記非因果成分からノイズレベルを推定するノイズ推定手段と、前記ノイズ推定手段により推定されたノイズレベルを前記因果成分から除去した伝達特性を求める打ち切り手段と、前記打ち切り手段により除去された伝達特性の減衰量を推定して増幅させた伝達特性を求める減衰補正手段と
を備えた伝達特性測定装置において、前記伝達特性をフーリエ変換して、観測点からの距離についての距離振幅スペクトルと距離位相スペクトルとを求める距離スペクトル演算手段を備え、前記距離振幅スペクトルにあわらわれるピークのうち、前記距離位相スペクトルが正である場合のピークの位置で固体の反射物までの距離情報を取得する一方、前記距離位相スペクトルが負である場合のピークの位置を空気密度による反射とすることを特徴としたものである。
を備えた伝達特性測定装置において、前記伝達特性をフーリエ変換して、観測点からの距離についての距離振幅スペクトルと距離位相スペクトルとを求める距離スペクトル演算手段を備え、前記距離振幅スペクトルにあわらわれるピークのうち、前記距離位相スペクトルが正である場合のピークの位置で固体の反射物までの距離情報を取得する一方、前記距離位相スペクトルが負である場合のピークの位置を空気密度による反射とすることを特徴としたものである。
【0013】
また、本発明の伝達特性測定プログラムは、コンピュータを、時間ｔ０から時間ｔ１までの音源からの音源信号である入力信号と前記音源信号が反射して観測点で観測される観測信号である出力信号とを切り出す切出手段、前記入力信号と前記出力信号とから時間ごとの伝達特性を求める伝達特性演算手段、時間ｔ１より以前の時間ｔ３から、時間ｔ１以後の時間ｔ２までの入力信号に０を与え、前記伝達特性から因果性が正しい因果成分と因果性が前後した非因果成分とを分離させる歪分離手段、前記歪分離手段により得られた前記非因果成分からノイズレベルを推定するノイズ推定手段、前記ノイズ推定手段により推定されたノイズレベルを前記因果成分から除去した伝達特性を求める打ち切り手段、前記打ち切り手段により除去された伝達特性の減衰量を推定して増幅させた伝達特性を求める減衰補正手段、前記伝達特性をフーリエ変換して、観測点からの距離についての距離振幅スペクトルと距離位相スペクトルとを求める距離スペクトル演算手段として機能させ、前記距離振幅スペクトルにあわらわれるピークのうち、前記距離位相スペクトルが正である場合のピークの位置で固体の反射物までの距離情報を取得する一方、前記距離位相スペクトルが負である場合のピークの位置を空気密度による反射とすることを特徴とする。
【0014】
本発明によれば、歪分離手段が、時間ｔ３より大きい時間からｔ２までには０を与え、ｔ０より時間ｔ３までの非因果成分を、時間ｔ２より大きい時間に分離させることができる。減衰補正手段が、伝達特性の減衰量を推定して増幅させた伝達特性を求めるときに、ノイズ推定手段が非因果成分からノイズレベルを推定し、打ち切り手段が、ノイズレベルを因果成分から除去した伝達特性を求めるので、補正してもノイズが拡大せず、精度のよい伝達特性を得ることができる。
【0015】
前記ノイズ推定手段は、前記非因果成分のエネルギー量によりノイズレベルを推定することができる。そうすることで、非因果成分によるノイズレベルを推定することができる。
【0016】
前記減衰補正手段は、前記因果成分の減衰量をヒルベルト変換により包絡線抽出を行い、包絡線の傾斜に応じて逆数を乗算することで減衰量を増幅するものとすることができる。そうすることで、減衰するノイズを抑えた伝達特性を示す信号を増幅することができる。
【0018】
前記入力信号と前記出力信号は、ドローン装置の回転翼を支持するアームに配置した入力用信号用のマイクと出力信号用のマイクで観測したものであって、前記入力用信号用のマイクと前記出力信号用のマイクとを１組として前記ドローン装置に３組配置して、前記距離スペクトル演算手段から得られる各ピークから立体の散布図を得る座標情報取得手段を備えたものとすることができる。そうすることで、三次元空間の反射物の位置や材質、状況などを観測することができるので、光学的に計測不可能な空気密度による反射も観測することができる。
【発明の効果】
【0021】
本発明は、非因果成分からノイズレベルを推定し、ノイズレベルを因果成分から除去した伝達特性を求めるので、補正してもノイズが拡大せず、精度のよい伝達特性を得ることができる。よって、本発明は、歪やノイズを推定して除去することで、更に精度を高めた伝達特性を求めることが可能である。
【図面の簡単な説明】
【0022】
【図1】本発明の実施の形態に係る伝達特性測定装置の構成を示す図である。
【図2】音源からの音源信号ｘ（ｔ）が反射物に反射して観測信号ｙ（ｔ）を観測点で観測することを説明するための図である。
【図3】（Ａ）は入力信号と出力信号と伝達特性を示すモデル図、（Ｂ）は入力信号ｘ（ｎ）と、ｘ（ｎ）が伝達特性ｈ（ｎ）により生じる出力信号ｙ（ｎ）の関係を示す図である。
【図4】（Ａ）はｘ（ｎ）とｙ（ｎ）とからγxy（τ）を求めることを説明するための図、（Ｂ）は観測用時間窓として時間ｔ１×２－１までを観測対象としたことを説明するための図、（Ｃ）は非因果成分を除去した状態を説明するための図である。
【図5】因果律を満たさない歪が分離された状態を説明するための波形図である。
【図6】音源近傍信号と反射物側信号とが計測できるよう２個のマイクを設置した波形であり、（Ａ）は振幅をリニア表示した波形図、（Ｂ）は振幅を対数表示（ｄＢ）した波形図、（Ｃ）は打ち切り処理した状態を説明するための波形図である。
【図7】（Ａ）は打ち切り処理前の歪なし区間における包絡線を示す波形図、（Ｂ）は打ち切り処理後の歪なし区間に逆数補正を行った状態の波形図である。
【図8】（Ａ）はパワースペクトルを示す波形図、（Ｂ）は距離振幅スペクトルを示す波形図、（Ｃ）は距離位相スペクトルを示す波形図、（Ｄ）はアンラップ処理後の距離位相スペクトルを示す波形図である。
【図9】距離振幅スペクトルにより、反射物の位置で観測点からの距離を半径として、同心円（同心球）を描くことを説明するための図である。
【図10】本発明の実施の形態２に係るドローン装置を示す図である。
【発明を実施するための形態】
【0023】
（実施の形態１）
本発明の実施の形態１に係る伝達特性測定装置を図面に基づいて説明する。なお、本明細書における伝達特性とは、インパルス応答（時間）、伝達関数（周波数）を含む概念である。
本発明の実施の形態１に係る伝達特性測定装置を図面に基づいて説明する。なお、本明細書における伝達特性とは、インパルス応答（時間）、伝達関数（周波数）を含む概念である。
【0024】
図１に示す本実施の形態１に係る伝達特性測定装置１０は、音響測距法により軸上の距離を求め、それを空間情報に展開して、３次元空間の反射物の位置を特定するものである。
伝達特性測定装置１０は、切出手段１０１と、伝達特性演算手段１０２と、歪分離手段１０３と、ノイズ推定手段１０４と、打ち切り手段１０５と、減衰補正手段１０６と、パワースペクトル演算手段１０７と、距離スペクトル演算手段１０８と、座標情報取得手段１０９とを備えている。
伝達特性測定装置１０は、切出手段１０１と、伝達特性演算手段１０２と、歪分離手段１０３と、ノイズ推定手段１０４と、打ち切り手段１０５と、減衰補正手段１０６と、パワースペクトル演算手段１０７と、距離スペクトル演算手段１０８と、座標情報取得手段１０９とを備えている。
【0025】
切出手段１０１は、所定の時間ｔ０から時間ｔ１までの入力信号と出力信号とを切り出す。
伝達特性演算手段１０２は、切り出された入力信号と出力信号とから時間ごとの伝達特性を求める。
歪分離手段１０３は、伝達特性演算手段１０２から得られた伝達特性であり、時間ｔ０から時間ｔ１までをフーリエ変換により伝達特性を求めるときに、時間ｔ１以後の時間ｔ２であり、時間ｔ１より以前の時間ｔ３から時間ｔ２までの入力信号に０を与え、因果性が正しい因果成分と因果性が前後した非因果成分とを分離させる伝達特性を求める。本実施の形態では、時間ｔ２を時間ｔ１×２－１となる時間としている。
伝達特性演算手段１０２は、切り出された入力信号と出力信号とから時間ごとの伝達特性を求める。
歪分離手段１０３は、伝達特性演算手段１０２から得られた伝達特性であり、時間ｔ０から時間ｔ１までをフーリエ変換により伝達特性を求めるときに、時間ｔ１以後の時間ｔ２であり、時間ｔ１より以前の時間ｔ３から時間ｔ２までの入力信号に０を与え、因果性が正しい因果成分と因果性が前後した非因果成分とを分離させる伝達特性を求める。本実施の形態では、時間ｔ２を時間ｔ１×２－１となる時間としている。
【0026】
ノイズ推定手段１０４は、歪分離手段１０３により得られた非因果成分からノイズレベルを推定する。
打ち切り手段１０５は、ノイズ推定手段１０４により推定されたノイズレベルを因果成分から除去した伝達特性を求める。
減衰補正手段１０６は、打ち切り手段１０５により除去された伝達特性の減衰量を推定して増幅させた伝達特性を求める。
打ち切り手段１０５は、ノイズ推定手段１０４により推定されたノイズレベルを因果成分から除去した伝達特性を求める。
減衰補正手段１０６は、打ち切り手段１０５により除去された伝達特性の減衰量を推定して増幅させた伝達特性を求める。
【0027】
パワースペクトル演算手段１０７は、減衰補正手段１０６により補正された伝達特性に基づいてパワースペクトルを求める。
距離スペクトル演算手段１０８は、パワースペクトルをフーリエ変換して、振幅スペクトルと位相スペクトルとを求める。
座標情報取得手段１０９が、振幅スペクトルの各ピークの位置で半径を取得して同心球を描き、振幅スペクトルの大きさで対象物の材質情報を取得し、位相スペクトルの符号で対象物の状況情報を取得することを三方向で行い、立体の散布図を得る。
距離スペクトル演算手段１０８は、パワースペクトルをフーリエ変換して、振幅スペクトルと位相スペクトルとを求める。
座標情報取得手段１０９が、振幅スペクトルの各ピークの位置で半径を取得して同心球を描き、振幅スペクトルの大きさで対象物の材質情報を取得し、位相スペクトルの符号で対象物の状況情報を取得することを三方向で行い、立体の散布図を得る。
【0028】
本実施の形態１に係る伝達特性測定装置１０では、例えば、図２に示すように、音源からの音源信号ｘ（ｔ）が反射物に反射して観測信号ｙ（ｔ）を観測点で２個のマイクにて観測することができる。但し、ｃは光速である。このときのパワースペクトルは、以下の式（１）により表すことができる。但し、ＦＴはフーリエ変換である。
【0029】
【0030】
本実施の形態１に係る伝達特性測定装置１０では、音源信号ｘ（ｔ）と観測信号ｙ（ｔ）とから伝達特性ｈ（ｔ）／Ｈ（ｆ）を導出し、そのパワースペクトルに生じる変調周期からｔ１－ｔ０を計測する。
このため，「音源信号ｘ（ｔ）」がパルス信号の必要が無く，例えばドローン飛翔音のような信号を音源として使用可能である。信号源は、空気中を伝搬する音や、物体を伝搬する振動、液体を伝搬する流体、気体の流れなど、任意のものでよく、信号を媒介する伝送路も任意のものでよいので、光学的に写らない対象を可視化することができる。
このため，「音源信号ｘ（ｔ）」がパルス信号の必要が無く，例えばドローン飛翔音のような信号を音源として使用可能である。信号源は、空気中を伝搬する音や、物体を伝搬する振動、液体を伝搬する流体、気体の流れなど、任意のものでよく、信号を媒介する伝送路も任意のものでよいので、光学的に写らない対象を可視化することができる。
【0031】
次に、伝達特性測定装置１０の動作を説明するが、まず、モデルを簡単にするために図３（Ａ）に示すように、入力信号をｘ（ｎ）、出力信号をｙ（ｎ）、インパルス応答（伝達特性）をｈ（ｎ）とする。なお、ｎは時間ｔと同じ離散値である。
ここで、入力信号ｘ（ｎ）と、ｘ（ｎ）が伝達特性ｈ（ｎ）により生じる出力信号ｙ（ｎ）の関係を、図３（Ｂ）に示す。図３（Ｂ）では、縦軸を入力信号、横軸を出力信号としている。出力信号は、同離散時間のレスポンスの総和（縦の和）である。
ここで、入力信号ｘ（ｎ）と、ｘ（ｎ）が伝達特性ｈ（ｎ）により生じる出力信号ｙ（ｎ）の関係を、図３（Ｂ）に示す。図３（Ｂ）では、縦軸を入力信号、横軸を出力信号としている。出力信号は、同離散時間のレスポンスの総和（縦の和）である。
【0032】
切出手段１０１（図１参照）は、これからｘ（ｎ）における時間ｔとしてｔ＝０からｔ＝７まで（入力用時間窓）（ｙ（ｔ）の時間ｔをｔ＝０からｔ＝７まで（観測用時間窓））を切り出す。そうすることにより図３（Ｂ）に示す表では、網掛け部分が除去される。
従って、ｔ＝０からｔ＝７までとした入力用時間窓以外の時間ｔを含む信号を歪として除去することができる。
従って、ｔ＝０からｔ＝７までとした入力用時間窓以外の時間ｔを含む信号を歪として除去することができる。
【0033】
次に、図３（Ｂ）による入力信号ｘ（ｎ）とｙ（ｎ）とから、伝達特性演算手段１０２（図１参照）は、ｈ（ｎ）を求める。
図３（Ｂ）に示す例からｙ（ｎ）を求める式は、ｎ＝７までとすると、次の式（２）による連立方程式で表すことができる。
図３（Ｂ）に示す例からｙ（ｎ）を求める式は、ｎ＝７までとすると、次の式（２）による連立方程式で表すことができる。
【0034】
ここで、ｙ（０）について、演算処理する際には、x(-1)h(1)+x(-2)h(2)+x(-3)h(3)+x(-4)h(4)+x(-5)h(5)+x(-6)h(6)+x(-7)h(7)が補われる。また、ｙ（１）についても、+x(-1)h(2)+x(-2)h(3)+x(-3)h(4)+x(-4)h(5)+x(-5)h(6)+x(-6)h(7)が補われる。以下、ｙ（５）、ｙ（６）、ｙ（７）も同様である。そうすると、式（１）は、以下の式（３）に示すようになる。
【0035】
従って、式（３）をベクトルで表すと以下の式（４）のように書き直すことができる。
ここで、ｘは入力信号ベクトル、ｙは出力信号ベクトルである。
【0036】
式（４）による伝達特性ｈ（ｔ）をベクトル化すると式（５）に示す式とすることができる。
【0037】
式（５）を数列化すると式（６）に示す行列式となる。
但し、ｙは出力信号ベクトルであり、ｘは入力信号行列であり、ｈはインパルス応答ベクトルである。
【0038】
式（６）からインパルス応答ベクトル（伝達特性）ｈ（ベクトル）を求めると式（７）のようになる。
【0039】
このＸTｙをγxy（τ）として、ｘ（ｎ）とｙ（ｎ）とからγxy（τ）を求めた表にしたものを図４（Ａ）に示す。γは、ハットｈである。
この図４（Ａ）に示す表から、まず、因果成分と非因果成分とを分離することで歪みの無い期間を抽出する。因果成分とは、入力信号と出力信号とで事象の時間的順序が正常な状態の要素（因果要素）であり、非因果的成分とは、事象の時間的順序が逆となった状態の要素（非因果要素）である。図４（Ａ）では、白抜き部分が因果成分であり、網掛け部分が非因果成分である。非因果成分はフーリエ変換した際に巡回性の問題から発生する。この非因果成分が歪となる。
この図４（Ａ）に示す表から、まず、因果成分と非因果成分とを分離することで歪みの無い期間を抽出する。因果成分とは、入力信号と出力信号とで事象の時間的順序が正常な状態の要素（因果要素）であり、非因果的成分とは、事象の時間的順序が逆となった状態の要素（非因果要素）である。図４（Ａ）では、白抜き部分が因果成分であり、網掛け部分が非因果成分である。非因果成分はフーリエ変換した際に巡回性の問題から発生する。この非因果成分が歪となる。
【0040】
そこで、歪分離手段１０３（図１参照）では、時間ｔ０から時間ｔ１までの信号をフーリエ変換により伝達特性を求めるときに、観測用時間窓として時間ｔ１×２－１となる時間ｔ２までの伝達特性を観測対象とする。
図４（Ａ）に示す例では、時間ｔ０＝０とする時間ｔから、時間ｔ１＝７とする時間ｔまでを観測対象としている。
図４（Ａ）に示す例では、時間ｔ０＝０とする時間ｔから、時間ｔ１＝７とする時間ｔまでを観測対象としている。
【0041】
図４（Ｂ）に示すように、観測用時間窓として時間ｔ１×２－１＝１５となる時間ｔ２までを観測対象とすると、γ（０）からγ（７）までの伝達特性が、γ（８）からγ（１５）までに繰り返しあらわれる。
そのため、歪分離手段１０３は、入力用時間窓を、時間ｔ０から、時間ｔ１より以前の時間ｔ３までとすることで、時間ｔ３から時間ｔ２までの間の入力信号ｘ（ｎ）に０を与える。
そのため、歪分離手段１０３は、入力用時間窓を、時間ｔ０から、時間ｔ１より以前の時間ｔ３までとすることで、時間ｔ３から時間ｔ２までの間の入力信号ｘ（ｎ）に０を与える。
【0042】
次に、時間ｔ０から時間ｔ１までの非因果成分に０を与えると共に、時間ｔ１の次から時間ｔ２までの因果成分に０を与える。
そうすることで、図４（Ｃ）に示すように、観測対象とした時間ｔ０（ｔ＝０）から時間ｔ１（ｔ＝７）までに含まれた非因果成分が無くなり、γ（１３）からγ（１５）に移動したように見える。
そうすることで、図４（Ｃ）に示すように、観測対象とした時間ｔ０（ｔ＝０）から時間ｔ１（ｔ＝７）までに含まれた非因果成分が無くなり、γ（１３）からγ（１５）に移動したように見える。
【0043】
このように因果律を満たさない歪が分離できれば、時間ｔ０からの時間ｔ３までの伝達特性であるγ（０）からγ（４）を歪なし期間とすることできる。
例えば、図５に示す波形からも、因果律を満たさない歪が分離できていることがわかる（矩形枠部分にて示す）。
例えば、図５に示す波形からも、因果律を満たさない歪が分離できていることがわかる（矩形枠部分にて示す）。
【0044】
例えば、音源近傍信号（入力信号）と反射物側信号（出力信号）とが計測できるよう２個のマイクを設置した波形を図６に示す。図６（Ａ）は、振幅をリニア表示したものであり、同図（Ｂ）は振幅を対数表示（ｄＢ）したものである。
【0045】
ここで、音源近傍信号をｘ、反射物側信号をｙとして、ｘとｙから伝達特性ｈを推定する（図４（Ｃ）に示す処理）。
そうすると、図６（Ｂ）に示すように、区間Ｓ１から、因果律を満たさない歪を区間Ｓ２に分離することで、歪なし区間とすることができる。
そうすると、図６（Ｂ）に示すように、区間Ｓ１から、因果律を満たさない歪を区間Ｓ２に分離することで、歪なし区間とすることができる。
【0046】
次に、区間Ｓ２が特定できると、ノイズ推定手段１０４が、切り離した非因果成分からノイズレベルを推定する。ノイズレベルの推定は、非因果成分のエネルギー量の平均により求めることができる。
図６（Ｂ）に示す例では、区間Ｓ２レベルがほぼ平均的なレベルであることから、この区間Ｓ２のノイズレベルは、区間Ｓ２の波形の波高値付近の値となり、レベルを示す水平な直線Ｌ１にて示すことができる。
図６（Ｂ）に示す例では、区間Ｓ２レベルがほぼ平均的なレベルであることから、この区間Ｓ２のノイズレベルは、区間Ｓ２の波形の波高値付近の値となり、レベルを示す水平な直線Ｌ１にて示すことができる。
【0047】
ノイズレベルが推定できると、次に、打ち切り手段１０５（図１参照）により打ち切り処理を行う。打ち切り処理とは、歪なし区間の伝達特性ｈに対して、ノイズ推定手段１０４（図１参照）が推定したノイズレベルによりノイズを除去する処理を指す。
図６（Ｃ）に示す例では、区間Ｓ２の任意の区間を切り出して、ＳＮ比が大きい部分（直線Ｌ１より下の部分）を除去している。
図６（Ｃ）に示す例では、区間Ｓ２の任意の区間を切り出して、ＳＮ比が大きい部分（直線Ｌ１より下の部分）を除去している。
【0048】
打ち切り処理後には、減衰補正手段１０６が、伝達特性（インパルス応答（ｈ））の減衰量を推定して、伝達特性ｈに対する補正を行う。
補正は、まず、打ち切り処理前の歪なし区間Ｓ２の伝達特性ｈからヒルベルト変換により包絡線抽出を行い、得られた包絡線である近似直線の傾きから減衰定数を求める。次に、減衰定数に基づいて打ち切り処理後のｈに対して補正を行う。
補正は、まず、打ち切り処理前の歪なし区間Ｓ２の伝達特性ｈからヒルベルト変換により包絡線抽出を行い、得られた包絡線である近似直線の傾きから減衰定数を求める。次に、減衰定数に基づいて打ち切り処理後のｈに対して補正を行う。
【0049】
図７（Ａ）にて、打ち切り処理前の歪なし区間Ｓ２（図６（Ｂ）参照）における包絡線Ｌ２を示す。この包絡線Ｌ２の傾きを求めることで、歪なし区間Ｓ２における減衰定数を推定することができる。
そして、この減衰定数が推定できると、打ち切り処理後のｈに逆数補正をする。逆数補正は、負の傾きに対して符号を反転させた正の傾きの直線に合わせることで波高値を平坦化する（図７（Ｂ）参照）。このように、予め打ち切り処理を行っているので、逆数補正を行ってもノイズが拡大してしまうことを抑止することができる。
そして、この減衰定数が推定できると、打ち切り処理後のｈに逆数補正をする。逆数補正は、負の傾きに対して符号を反転させた正の傾きの直線に合わせることで波高値を平坦化する（図７（Ｂ）参照）。このように、予め打ち切り処理を行っているので、逆数補正を行ってもノイズが拡大してしまうことを抑止することができる。
【0050】
このようにして、補正された伝達特性ｈが得られることにより、精度がよい伝達特性ｈを得ることができる。
【0051】
次に、パワースペクトル演算手段１０７が、パワースペクトルを求める。
パワースペクトルは、上記式（１）の以下の式（８）に基づいてパワースペクトルを求める。
Ｈ*（ｆ）Ｈ（ｆ）により図８（Ａ）に示すパワースペクトルが得られ、これをフーリエ変換すると最終式のｃｏｓ成分が得られる。これを距離スペクトル演算手段１０８（図１参照）がフーリエ変換する。この結果を「距離スペクトル」と称することにする。
パワースペクトルは、上記式（１）の以下の式（８）に基づいてパワースペクトルを求める。
【0052】
距離スペクトルでは、図８（Ｂ）に示す距離に応じた振幅を示す距離振幅スペクトル（振幅スペクトル）と、図８（Ｃ）に示す距離に応じた位相を示す距離位相スペクトル（位相スペクトル）とを求めることができる。また、図８（Ｄ）にアンラップ処理後の距離位相スペクトルを示す。
図８（Ｂ）に示す距離振幅スペクトルは、縦軸が振幅であり、横軸が時間である。
また、図８（Ｃ）および同図（Ｄ）に示す距離位相スペクトルは、縦軸が角度（－π～＋π）であり、横軸が時間である。
図８（Ｂ）に示す距離振幅スペクトルは、縦軸が振幅であり、横軸が時間である。
また、図８（Ｃ）および同図（Ｄ）に示す距離位相スペクトルは、縦軸が角度（－π～＋π）であり、横軸が時間である。
【0053】
距離振幅スペクトルは、反射物がピーク位置にあらわれる。
距離位相スペクトルは、位相が正となる値では、固体の反射物の反射率で材質を予測することができる。また、位相が負となる値では、空気密度による反射（例えば、開口端反射など。）であることが予測できる。
従って、距離スペクトル演算手段１０８（図１参照）は、距離振幅スペクトルに基づいて、反射物の位置により観測点からの距離を取得し、距離振幅スペクトルの大きさで反射物の材質情報を取得し、更に、距離位相スペクトルの符号で反射物の状況情報を取得する。
距離位相スペクトルは、位相が正となる値では、固体の反射物の反射率で材質を予測することができる。また、位相が負となる値では、空気密度による反射（例えば、開口端反射など。）であることが予測できる。
従って、距離スペクトル演算手段１０８（図１参照）は、距離振幅スペクトルに基づいて、反射物の位置により観測点からの距離を取得し、距離振幅スペクトルの大きさで反射物の材質情報を取得し、更に、距離位相スペクトルの符号で反射物の状況情報を取得する。
【0054】
このように固体での反射だけでなく、光学的に計測不可能な空気密度による反射も観測することができる。また、打ち切り処理を行ってから逆数補正を行っているため、ノイズが抑止された伝達特性であるからこそ、正確に反射物の位置や材質、状況など正確に特定することができる。
【0055】
後述する３組の入力用信号用のマイクと出力信号用のマイクを配置したドローン装置では、距離スペクトルは、いずれか一つの組でのある特定の方向に向かう直線上（一方向目）の情報である。そこで、座標情報取得手段１０９は、距離振幅スペクトルに基づいて、反射物の位置により観測点からの距離を半径として、図９に示すように同心円（同心球）を描く。
次に、他の組での他の方向（二方向目）の距離振幅スペクトルから描く同心円（同心球）により交点が求められれば平面（２次元）の散布図が得られる。更に、他の組でのもう一方向（三方向目）の距離振幅スペクトルから描く円により交点が求められれば立体（３次元）の散布図が得られる。
このようにして、反射物までの軸上の距離を求め、それを空間情報に展開することで、３次元空間の反射物の位置や材質、状況などを特定することができる。
次に、他の組での他の方向（二方向目）の距離振幅スペクトルから描く同心円（同心球）により交点が求められれば平面（２次元）の散布図が得られる。更に、他の組でのもう一方向（三方向目）の距離振幅スペクトルから描く円により交点が求められれば立体（３次元）の散布図が得られる。
このようにして、反射物までの軸上の距離を求め、それを空間情報に展開することで、３次元空間の反射物の位置や材質、状況などを特定することができる。
【0056】
（実施の形態２）
本発明の実施の形態２に係る伝達特性測定装置を搭載してドローン装置を図面に基づいて説明する。
図１０に示すドローン装置２０は、障害物を検知して回避運動を行う自律型のドローン装置である。
ドローン装置２０は、８本の回転翼２１を支持するアーム２２を備えている。アーム２２の回転翼２１側には、入力用信号用のマイク２３ａが配置されており、マイク２３ａが配置されたアーム２２の裏側には、出力信号用のマイク２３ｂが配置されている。
ドローン装置２０には、入力用信号用のマイク２３ａと出力信号用のマイク２３ｂとを一組として３組配置されている。
本発明の実施の形態２に係る伝達特性測定装置を搭載してドローン装置を図面に基づいて説明する。
図１０に示すドローン装置２０は、障害物を検知して回避運動を行う自律型のドローン装置である。
ドローン装置２０は、８本の回転翼２１を支持するアーム２２を備えている。アーム２２の回転翼２１側には、入力用信号用のマイク２３ａが配置されており、マイク２３ａが配置されたアーム２２の裏側には、出力信号用のマイク２３ｂが配置されている。
ドローン装置２０には、入力用信号用のマイク２３ａと出力信号用のマイク２３ｂとを一組として３組配置されている。
【0057】
また、ドローン装置２０には、制御ボックス３０が配置されている。制御ボックス３０には、本実施の形態２に係る伝達特性測定装置が内蔵されている。
ここで、本実施の形態２に係る伝達特性測定装置を図９に基づいて説明する。なお、図１１においては、図１と同じ構成のものは同符号を付して説明を省略する。また、図１１においては、回転翼２１（図１０参照）を駆動するモーターなどは図示していない。
ここで、本実施の形態２に係る伝達特性測定装置を図９に基づいて説明する。なお、図１１においては、図１と同じ構成のものは同符号を付して説明を省略する。また、図１１においては、回転翼２１（図１０参照）を駆動するモーターなどは図示していない。
【0058】
図１１に示す伝達特性測定装置１０ｘでは、座標情報取得手段１０９により位置や材質、状況などが特定された３次元空間の反射物を回避するための操舵手段１１０を備えている。
【0059】
座標情報取得手段１０９により得た立体の散布図から、衝突すると飛行に支障がある反射物が進行方向に検出されると、操舵手段１１０が、回避運動をとるように飛行手段であるモーターに指示する。
例えば、反射物が固体であり、反射物までの距離が、飛行速度に基づいて回避可能な距離の限界値から決定された所定閾値未満である反射物が特定された場合に、左方への移動や右方への移動、上昇、下降するなどの回避運動をモーターに指示する。
そうすることで、図１０に示すドローン装置２０は、正確に反射物（障害物）を避けながら、飛行することができる。
例えば、反射物が固体であり、反射物までの距離が、飛行速度に基づいて回避可能な距離の限界値から決定された所定閾値未満である反射物が特定された場合に、左方への移動や右方への移動、上昇、下降するなどの回避運動をモーターに指示する。
そうすることで、図１０に示すドローン装置２０は、正確に反射物（障害物）を避けながら、飛行することができる。
【産業上の利用可能性】
【0060】
本発明は、信号源として、空気中を伝搬する音や、物体を伝搬する振動、液体を伝搬する流体、気体の流れなど、任意のものでよく、光学的に写らない対象を可視化することができるので、様々な分野の測定に応用することができ、特に、自律飛行するドローン装置に最適である。
【符号の説明】
【0061】
１０，１０ｘ 伝達特性測定装置
１０１ 切出手段
１０２ 伝達特性演算手段
１０３ 歪分離手段
１０４ ノイズ推定手段
１０５ 打ち切り手段
１０６ 減衰補正手段
１０７ パワースペクトル演算手段
１０８ 距離スペクトル演算手段
１０９ 座標情報取得手段
１１０ 操舵手段
２０ ドローン装置
２１ 回転翼
２２ アーム
２３ａ，２３ｂ マイク
３０ 制御ボックス
Ｓ１，Ｓ２ 区間
Ｌ１ 直線
Ｌ２ 包絡線
１０１ 切出手段
１０２ 伝達特性演算手段
１０３ 歪分離手段
１０４ ノイズ推定手段
１０５ 打ち切り手段
１０６ 減衰補正手段
１０７ パワースペクトル演算手段
１０８ 距離スペクトル演算手段
１０９ 座標情報取得手段
１１０ 操舵手段
２０ ドローン装置
２１ 回転翼
２２ アーム
２３ａ，２３ｂ マイク
３０ 制御ボックス
Ｓ１，Ｓ２ 区間
Ｌ１ 直線
Ｌ２ 包絡線
【図1】
【図2】
【図3】
【図4】
【図5】
【図6】
【図7】
【図8】
【図9】
【図10】
【図11】