JP7263271B2 - 演算装置およびプログラム - Google Patents

Description

本発明は演算装置およびプログラムに関する。

周波数解析を行なう技術としては一般的にＤＦＴ（Discrete Fourier Transform：離散フーリエ変換）やそれを高速に実行するＦＦＴ（Fast Fourier Transform：高速フーリエ変換）といった方法が知られている。

特開２０１７－１１９２８５号公報 特開２００５－２０７８７０号公報 特開２００３－２３２６７４号公報 特開平７－３１１２２７号公報

ＤＦＴやＦＦＴを用いた周波数解析においては、事前に入力される波形の基本周波数が分かっていなければ結果が保証されないという問題がある。また、計算コストも高いという問題がある。

このため、ある程度長い期間の情報を入力することで誤差を最小化するであるとか、経験的な値を適用するであるとか、事前に波形を確認して基本周波数と思われる周波数を推測するといった方法によって誤差を最小化する試みが行われてきた。
１つの側面では、本発明は、入力される波の基本周期を容易に把握することを目的とする。

上記目的を達成するために、開示の周期フィルタが提供される。この周期フィルタは、入力される離散時間信号から、予め割り当てられたそれぞれの時間長が異なる基本周期候補区間に含まれる連続した離散点をパラレルに出力ポートへ出力する複数の信号出力部と、各信号出力部により出力される信号の値を、出力ポート毎にそれぞれ足し合わせることで基本周期候補区間毎の累積値を算出する加算部と、基本周期候補区間毎の累積値の頂点の軌跡により生じる傾きの差分を基礎に基本周期を導出する演算を実行する演算部と、を有している。

１態様では、入力される波の基本周期を容易に把握することができる。

第１の実施の形態の集積回路を示す図である。 第１の実施の形態の周期フィルタの一例を示す図である。 加算器の処理を説明する図である。 加算器の処理を説明する図である。 第１の実施の形態の周期フィルタの処理を説明するフローチャートである。 第１の実施の形態の周期フィルタの処理を説明するフローチャートである。 基本周期を正確に識別できることを説明する図である。 基本周期を正確に識別できることを説明する図である。 基本周期を正確に識別できることを説明する図である。 基本周期を正確に識別できることを説明する図である。 基本周期を正確に識別できることを説明する図である。 傾斜したＡＲＳスペクトルを説明する図である。 第２の実施の形態の周期フィルタの一例を示す図である。 第２の実施の形態の第１演算部および第２演算部の処理の概要を説明する図である。 第２演算を説明する図である。 第２演算を説明する図である。 第２演算を説明する図である。 第２演算の応用例を説明する図である。 第２演算の応用例を説明する図である。 第２の実施の形態の周期フィルタの処理を説明するフローチャートである。 第２の実施の形態の周期フィルタの処理を説明するフローチャートである。 第２の実施の形態の周期フィルタの処理の第１の実施例を説明する図である。 第２の実施の形態の周期フィルタの処理の第１の実施例を説明する図である。 第２の実施の形態の周期フィルタの処理の第２の実施例を説明する図である。 第２の実施の形態の周期フィルタの処理の第３の実施例を説明する図である。 第３の実施の形態の周期フィルタの処理を説明するフローチャートである。 第３の実施の形態の周期フィルタの処理を説明するフローチャートである。 第３の実施の形態の周期フィルタの処理の実施例を説明する図である。 第３の実施の形態の周期フィルタの変形例の処理を説明するフローチャートである。 第３の実施の形態の周期フィルタの変形例の処理を説明するフローチャートである。

以下、実施の形態の集積回路を、図面を参照して詳細に説明する。

以下の図面等において示す各構成の位置、大きさ、形状、範囲などは、発明の理解を容易にするため、実際の位置、大きさ、形状、範囲などを表していない場合がある。このため、本発明は、必ずしも、図面等に開示された位置、大きさ、形状、範囲等に限定されない。
実施の形態において単数形で表される要素は、文面で明らかに示されている場合を除き、複数形を含むものとする。
＜第１の実施の形態＞
図１は、第１の実施の形態の集積回路を示す図である。

集積回路１は、送信部（Transmitter）２と、受信部（Receiver）３と、ＡＤＣ（Analog Digital Converter）４と、周期フィルタ（周期決定部）５とを有している。

送信部２は、周期を測定する対象（測定対象）にアンテナ２１を介して電磁波や音声、超音波などの周期信号（電磁波等）を発信する。送信部２と受信部３との組み合わせは、図示したものに限定されず、例えば、送信側がスピーカで受信側がマイクであってもよい。
測定対象としては、例えば、生体の心拍および呼吸によって生じる体表面の偏位等が挙げられる。測定対象は１つでもよいし、図１に示すように複数でもよい。
受信部３は、測定対象周波数の偏差を検出する。測定対象周波数を検出した波形には、ノイズ（η）が乗っている。
ＡＤＣ４は、受信部３が受信したアナログ値をデジタル値（離散時間信号）に変換する。
なお、集積回路１が有する機能の一部または全部をソフトウェアで実現してもよい。
機能の一部をソフトウェアで実現する場合には、コンピュータ上でそのソフトウェアを動作させてもよい。

なお、受信部３は、実数、虚数を分離する機能を備えていてもよい。その場合、ＡＤＣ４は、実数、虚数それぞれをデジタル値に変換するようになっていてもよい。

実施の形態の周期フィルタ５は、基本周期の違いによって重なった波を分離して周期スペクトルを出力する。この周期フィルタ５は、入力される波についての事前情報がなくても波を分離することができる。

なお、周期フィルタ５は、ＦＰＧＡ（Field-Programmable gate array）等で専用回路を生成してもよいし、ＤＳＰ（Digital Signal Processor）で実現することもできる。
図２は、第１の実施の形態の周期フィルタの一例を示す図である。

周期フィルタ５は、ｋ（ｋ＝１・・・Ｖ）個の標本化回路（図面においてはｖ（ｖ＝１、２、・・・Ｖ）個のＳＰＣ：Serial Parallel Converterと表記する場合もある）５ａ１、５ａ２、・・・、５ａＶと、標本化回路５ａ１、５ａ２、・・・、５ａＶに対応するｖ（ｖ＝１、２、・・・Ｖ）個の積算器（Accumulator）として機能する加算器５ｂ１、５ｂ２、・・・、５ｂＶと、最大値検出部５ｃとを有している。

標本化回路５ａ１、５ａ２、・・・、５ａＶは、それぞれ入力される離散時間信号を標本化してパラレル出力する。より具体的には、標本化回路５ａ１、５ａ２、・・・、５ａＶは入力される離散時間信号から、予め割り当てられた、それぞれの時間長が異なる区間（基本周期候補区間）に含まれる連続した離散点を抽出してパラレルに出力する。

すなわち、基本周期候補区間において抽出する離散点の数（サンプル数）は、標本化回路５ａ１、５ａ２、・・・、５ａＶそれぞれで異なる。より具体的には、これらの標本化回路５ａ１、５ａ２、・・・、５ａＶは、予め割り当てられた周期の範囲に応じて抽出された、異なるサンプル数を出力ポートから繰り返し出力する。

例えばＡＤＣ４の出力が１４４サンプルである場合であり、標本化回路５ａ１に設定された予測される周期の範囲が３サンプルである場合、連続する３サンプルを４８回繰り返し出力する。また、標本化回路５ａ２に設定された予測される周期の範囲が４サンプルである場合、連続する４サンプルを３６回繰り返し出力する。なお、繰返し出力する回数を標本化回路５ａ１、５ａ２、・・・、５ａＶで一律にしてもよい。

本実施の形態では、標本化回路５ａ１、５ａ２、・・・、５ａＶのうち、標本化回路５ａ１のサンプル数が最も少なく、標本化回路５ａ２、５ａ３と順々に１つずつサンプル数を増やしている。

加算器５ｂ１、５ｂ２、・・・、５ｂＶは、標本化回路５ａ１、５ａ２、・・・、５ａＶの出力ポートからそれぞれ入力される各サンプルの値を、同一ポート毎に足し合わせる。
加算器５ｂ１、５ｂ２、・・・、５ｂＶは、それぞれに保存している値と新たに入力される値を加算して加算結果を更新する。

最大値検出部５ｃは、加算器５ｂ１、５ｂ２、・・・、５ｂＶの演算結果の絶対値を比較し、最も大きな加算結果の絶対値を出力した加算器５ｂ１、５ｂ２、・・・、５ｂＶに値を出力した標本化回路５ａ１、５ａ２、・・・、５ａＶに設定された区間（標本化される最初の離散点と最後の離散点の間）の長さにより決定される距離を基本周期であると判断して出力する。

なお、演算結果の絶対値の比較方法については、特に限定されないが、例えば加算器５ｂ１、５ｂ２、・・・、５ｂＶの演算結果を、加算の回数で割った平均値の絶対値を比較器で比較するようにしてもよい。他の例では、例えば加算器５ｂ１、５ｂ２、・・・、５ｂＶの全ての演算結果を正規化して比較するようにしてもよい。
図３および図４は、加算器の処理を説明する図である。

図３は、予測される周期の範囲が８サンプルである場合を例示しており、標本化回路５ａＸの出力ポート数は「８」である。このとき、パラメータｇ（予測される周期を１区間の長さとして入力波形を分割した場合に、入力波形の初めから数えた該当区間の順番）＝０、１、２、３の場合の各出力ポートそれぞれの値（出力値）はほぼ等しく、同じ波形が出現している。このため同一の出力ポート毎に出力される信号の値をそれぞれ足し合わせることで、振幅の値が大きな波形となる。従って、８サンプルが１周期であると判断することができる。

図４は、予測される周期の範囲が７サンプルである場合を例示しており、標本化回路５ａＹの出力ポート数は「７」である。このとき、パラメータｇ＝０、１、２、３の場合の各出力ポートそれぞれの値（出力値）は異なりバラバラの値となる。このため同一ポート毎に値を足し合わせると互いに打ち消しあい、振幅の値が小さな波形となる。

このように、周期が一致している波の離散値を足してくと、どんどん絶対値が大きくなっていく。周期が一致していない波の離散値を足していっても互いに打ち消し合うため絶対値は大きくならない。従って、絶対値が最も大きくなる加算器に対応する標本化回路の抽出区間の長さを検出することにより、基本周期を容易に発見することができる。

また、最大値検出部５ｃは、最も大きな積算結果を１となるように全ての積算結果の絶対値を正規化し、正規化された積算結果の絶対値と、０以上１以下の閾値とを比較して、入力される離散時間信号に含まれる信号の数と、それぞれの周期と、を判定することができる。

次に、周期フィルタ５の処理を、フローチャートを用いて説明する。なお、以下に示す処理は一例であり、一部の処理を省略してもよい。一部の処理を他の処理に置き換えてもよい。処理と処理の間に他の処理を実行してもよい。一部の処理の順序を入れ替えてもよい。
図５および図６は、周期フィルタの処理を説明するフローチャートである。
［ステップＳ１］ 周期フィルタ５は、以下の演算に使用する標本化回路の個数Ｖを設定する。その後、ステップＳ２に遷移する。
［ステップＳ２］ 周期フィルタ５は、パラメータｖ＝１に設定する。その後、ステップＳ３に遷移する。
［ステップＳ３］ 周期フィルタ５は、ｖ番目の標本化回路の出力ポート数をｕ（ｖ）に決定する。その後、ステップＳ４に遷移する。
［ステップＳ４］ 周期フィルタ５は、全ての要素が０の１行ｕ（ｖ）の配列ｗ（ｖ）を作成する。その後、ステップＳ５に遷移する。
［ステップＳ５］ 周期フィルタ５は、パラメータｇ＝１、ｋ＝１を設定する。その後、ステップＳ６に遷移する。
［ステップＳ６］ 周期フィルタ５は、累積する標本化回路出力の個数Ｇを決定する。その後、ステップＳ７に遷移する。

［ステップＳ７］ 加算器５ｂ１、５ｂ２、・・・、５ｂＶは、ｖ番目の標本化回路に入力サンプルｘ［ｋ］からｘ［ｋ＋ｕ_ｖ－１］までを入力する。得られた出力を１行ｕ（ｖ）の配列ｘ_ｖ（ｇ）とする。その後、ステップＳ８に遷移する。

［ステップＳ８］ 加算器５ｂ１、５ｂ２、・・・、５ｂＶは、以下の演算を行う。Ｗ（ｖ）＝Ｗ（ｖ）＋ｘ_ｖ（ｇ）、ｋ＝ｇ・ｕ（ｖ）。その後、ステップＳ９に遷移する。
［ステップＳ９］ 周期フィルタ５は、以下の演算を行う。ｇ＝ｇ＋１。その後、ステップＳ１０に遷移する。

［ステップＳ１０］ 周期フィルタ５は、パラメータｇが個数Ｇ以下か否かを判断する。パラメータｇが個数Ｇ以下である場合（ステップＳ１０のＹｅｓ）、ステップＳ７に遷移する。パラメータｇが個数Ｇより大きい場合（ステップＳ１０のＮｏ）、ステップＳ１１に遷移する。
［ステップＳ１１］ 周期フィルタ５は、以下の演算を行う。Ｗ（ｖ）＝Ｗ（ｖ）／ｇ。その後、ステップＳ１２に遷移する。

［ステップＳ１２］ 最大値検出部５ｃは、配列Ｗ（ｖ）の要素のうち、絶対値が最大の要素を特定し、その要素の値をｐ（ｖ）とする。その後、ステップＳ１３に遷移する。
［ステップＳ１３］ 周期フィルタ５は、以下の演算を行う。ｖ＝ｖ＋１。その後、ステップＳ１４に遷移する。

［ステップＳ１４］周期フィルタ５は、パラメータｖが個数Ｖ以下か否かを判断する。パラメータｖが個数Ｖ以下である場合（ステップＳ１４のＹｅｓ）、ステップＳ３に遷移する。パラメータｖが個数Ｖより大きい場合（ステップＳ１４のＮｏ）、ステップＳ１５に遷移する。

［ステップＳ１５］ 周期フィルタ５は、ｐ（１）からｐ（Ｖ）のうち、最大となる値を特定し、その値をｐ１とする。このｐ１が１となるようにｐ（１）からｐ（Ｖ）までを正規化し、得られた値を、それぞれＰ（１）からＰ（Ｖ）として１行Ｖ列の行列Ｐ＝［Ｐ（１），Ｐ（２），・・・，Ｐ（Ｖ）］を作る。その後、ステップＳ１６に遷移する。
［ステップＳ１６］ 周期フィルタ５は、閾値Ｑ（０＜Ｑ≦１）を設定する。その後、ステップＳ１７に遷移する。
［ステップＳ１７］ 周期フィルタ５は、行列Ｐの要素のうち、閾値Ｑを超えるものを全て選択する。選択された要素数をＮとする。
［ステップＳ１８］ 周期フィルタ５は、行列Ｐの要素のうち、閾値Ｑを超える要素の番号をＶ_１１，・・・Ｖ_Ｎ１として記録する。

［ステップＳ１９］ 周期フィルタ５は、要素数Ｎを信号数として出力する。また、周期フィルタ５は、それぞれの周期としてＶ_１１，・・・Ｖ_Ｎ１を出力する。また、周期フィルタ５は、推定基本波形としてＷ（Ｖ_１１），・・・Ｗ（Ｖ_Ｎ１）とする。その後、図５および図６の処理を終了する。
図７～図１１は、基本周期を正確に識別できることを説明する図である。
図７に示すように、１周期を３１サンプルと仮定した場合と、１つの０を付加して１周期を３２サンプルに仮定した場合を例示する。
図８に示すように、１周期を３１サンプルと仮定した場合と、３２サンプルと仮定した場合では、周期が少しずつずれていく。

図９は、図８のアナログ波形に対する周期フィルタ５のステップＳ１５における１行Ｖ列の行列Ｐ＝［Ｐ１Ｐ（２）・・・Ｐ（３９）］を示している。１周期を３１サンプルと仮定した場合の絶対値が最も大きい。従って周期フィルタ５は、３１サンプルの長さにより決定される距離を、基本周期であると判断して出力する。
図１０は、受信部３が受信したアナログ波形の一例を示している。
図１０に示すアナログ波形は、図８に示すアナログ波形に比べてノイズが加わり、一見では基本周期の特定が難しい。
図１１は、アナログ波形に対する周期フィルタ５のステップＳ１５における１行Ｖ列の行列Ｐ＝［Ｐ１Ｐ（２）・・・Ｐ（３９）］を示している。

１周期を３１サンプルと仮定した場合と、３２サンプルと仮定した場合が大きく突出しているが、１周期を３２サンプルに仮定した場合の値が最も大きい。従って周期フィルタ５は、３２サンプルの長さにより決定される距離を、基本周期であると判断して出力する。このように、周期フィルタ５によれば、僅かな周期の違いであっても基本周期を容易に検出することができる。すなわち、通常、隣り合った波を区別することや、波形のよく似た高相関の波を分離することはとても難しく、同時に、検出される波と検出されない波が隣り合っている場合でもここまで明確に区別することは困難である。

これに対し、図１１に示すように、実施の形態の周期フィルタ５によれば、図１０に示す信号を入力した場合に正しく周期が検出できることを示していることと、高相関の二つの波である３１サンプルの波と３２サンプルの波について同時に入力された場合でもそれを明確に区別できること、それ以外の波については明らかに検出されていないことなどが示されている。

以上述べたように、集積回路１によれば、標本化回路５ａ１、５ａ２、・・・、５ａＶと加算器５ｂ１、５ｂ２、・・・、５ｂＶ、および最大値検出部５ｃとを有する周期フィルタ５を備え、入力として離散時間信号を用いて、入力波形を構成する一つ又は複数の波形の周期を求めることができる。

従って、波形の情報を事前に取得しなくても、単純な標本化と足し合わせの組み合わせによる方法を用いて、１つまたは複数の入力波形を解析することができる。また、高い相関を持った異なる基本周波数の波形が同時に存在する場合に、これらを区別して出力できる。

また、周期フィルタ５は、ＦＦＴやＤＦＴ等一般的な周波数解析方法と比較して低周波域の解像度が特に高くなる。このため、より高精度に基本周期を求めることができる。
また、固定小数点演算が可能なことや、加算処理を主に用いて演算を行なうことから、高速且つ高精度、並びに、低計算資源での周期解析を実現できる。
また、異なる基本周期の波形が複数合成された入力波形を解析し、各波形の基本周期並びに１／ｎ倍周期の高調波の周期を求めることができる。
また、上記異なる基本周期の波形は、相関の高い波形であっても良い。
また、上記異なる基本周期の波形は、周波数帯域が隣接する、若しくは重なった波形であっても良い。
また。周期フィルタ５により特定された基本周期を用いてＦＦＴ解析を行なうことも出来る。

なお、周期フィルタ５が行なった処理が複数の装置によって分散処理されるようにしてもよい。例えば、一つの装置が、標本化処理を行ない、他の装置が、その出力されたパラレル値を用いて基本周期を特定するようにしても良い。
＜第２の実施の形態＞
次に、第２の実施の形態の集積回路について説明する。
以下、第２の実施の形態の集積回路について、前述した第１の実施の形態との相違点を中心に説明し、同様の事項については、その説明を省略する。

第１の実施の形態の周期フィルタ５は、比較的少ないサンプル数を抽出してパラレルに出力する標本化回路（以下、少ＳＰＣとも言う）から比較的多いサンプル数を抽出してパラレルに出力する標本化回路（以下、多ＳＰＣとも言う）を多数実装することができる。少ＳＰＣの出力ポートから出力される各サンプルの値の加算結果の累積値は、多ＳＰＣの出力ポートから出力される各サンプルの値の加算結果の累積値より大きくなる。そのため、標本化回路の出力をそのまま累積すると、少ＳＰＣの累積数が多いので、傾斜したＡＲＳスペクトル（周期に対する各標本化回路の出力ポートから出力される値の累積値）が得られる。
図１２は、傾斜したＡＲＳスペクトルを説明する図である。
縦軸は各標本化回路の出力ポートから出力される値の累積値を示している。横軸は周期を示している。
出力ポートの数が少ない標本化回路ほど、出力数が多いので、左側にいくほど標本化回路の出力の累積値が多くなる。

この問題を解決するため、例えば、標本化回路の出力を累積したあと、累積数で割って平均する方法が考えられる。これにより累積数の差がなくなり、ＡＲＳスペクトルが平らになる。

しかしながら、標本化回路の数と同じ数の割り算を実行すると、演算量がその分大幅に増えることになる。また、平均をとったとしてもＳ／Ｎ比が低いとき、今度はＡＲＳスペクトルが右肩上がりとなってしまう。このようにＡＲＳスペクトルが傾斜していると、ピーク検出がうまくいかない場合がある。第２の実施の形態の周期フィルタは、この問題を解決することを目的の１つとしている。
図１３は、第２の実施の形態の周期フィルタの一例を示す図である。

第２の実施の形態の周期フィルタ５０は、第１演算部５ｄと第２演算部５ｅとを有している。第１演算部５ｄと第２演算部５ｅは、基本周期候補区間毎の累積値の頂点の軌跡により生じる傾きの差分を基礎に基本周期を導出する演算を実行する。
図１４は、第２の実施の形態の第１演算部および第２演算部の処理の概要を説明する図である。

第１演算部５ｄは、ＡＲＳスペクトルの値の右側から左側の値の差分（以下、左差分と言う）を演算する。また、第１演算部５ｄは、ＡＲＳスペクトルの値の左側から右側の値の差分（以下、右差分と言う）を演算（第１演算）する。

図１４に示す例では、ＡＲＳスペクトルの値を識別する（１）～（６）の符号を付している。そして、左差分、右差分の演算式を示している。図１４に示す例では、隣接する値同士の差分を演算しているが、１以上の間隔を開けて左差分及び右差分を演算するようにしてもよい（この例については後述する）。

そして、第２演算部５ｅは、求めた左差分及び右差分を用いて第２演算を行う。具体的には、第２演算部５ｅは、第１演算で演算した間隔だけ左差分と右差分をずらして左差分と右差分を足し合わせる演算を行う。例えば、隣接するＡＲＳスペクトルの値同士で右差分と左差分を演算した場合は、第２演算部５ｅは、１つ分、左差分と右差分をずらして左差分と右差分を足し合わせる演算を行う。これにより、ＡＲＳスペクトルの傾斜が相殺され、ピークが浮き上がってより分かり易くなる。
図１５～図１７は、第２演算を説明する図である。

図１５に示すＡＲＳスペクトルは、ピークが１つのケースである。図１５中アスタリスクで示したＡＲＳスペクトルの値がピーク値である。左差分及び右差分を用いて第２演算を行うことでピークを検出することができる。

図１６に示すＡＲＳスペクトルは、ピークが２つのケースである。図１６中アスタリスクで示したＡＲＳスペクトルの値がピーク値である。この場合も図１５に示す場合と同様に左差分及び右差分を用いて第２演算を行うことでピークを検出することができる。

図１７に示すＡＲＳスペクトルは、ピークが３つのケースである。図１７中アスタリスクで示したＡＲＳスペクトルの値がピーク値である。この場合も図１５に示す場合と同様に左差分及び右差分を用いて第２演算を行うことでピークを検出することができる。
図１８および図１９は、第２演算の応用例を説明する図である。
図１８に示すＡＲＳスペクトルは、図１７に示すＡＲＳスペクトルと同様に、ピークが３つのケースである。

前述したように、第２演算部５ｅは、１以上の間隔を開けて左差分及び右差分を演算するようにしてもよい。図１８に示す例では、第２演算部５ｅは、１つの間隔を開けて左差分及び右差分を第１演算している。このような第１演算を行った場合、第２演算部５ｅは、２つ分、左差分と右差分をずらして左差分と右差分を足し合わせる演算を行う。
また、複数のスペクトル間で差分をとることにより、隣り合った信号も検出することができる。
図２０および図２１は、第２の実施の形態の周期フィルタの処理を説明するフローチャートである。
［ステップＳ２１］ 周期フィルタ５０は、以下の演算に使用する標本化回路の個数Ｖを設定する。その後、ステップＳ２２に遷移する。
［ステップＳ２２］ 周期フィルタ５０は、パラメータｖ＝１に設定する。その後、ステップＳ２３に遷移する。
［ステップＳ２３］ 周期フィルタ５０は、ｖ番目の標本化回路の出力ポート数をｕ（ｖ）に決定する。その後、ステップＳ２４に遷移する。
［ステップＳ２４］ 周期フィルタ５０は、全ての要素が０の１行ｕ（ｖ）の配列Ｗ（ｖ）を作成する。その後、ステップＳ２５に遷移する。
［ステップＳ２５］ 周期フィルタ５０は、パラメータｇ＝１、ｋ＝１を設定する。その後、ステップＳ２６に遷移する。
［ステップＳ２６］ 周期フィルタ５０は、累積する標本化回路出力の個数Ｇを決定する。その後、ステップＳ２７に遷移する。

［ステップＳ２７］ 加算器５ｂ１、５ｂ２、・・・、５ｂＶは、ｖ番目の標本化回路に入力サンプルｘ［ｋ］からｘ［ｋ＋ｕ_ｖ－１］までを入力する。得られた出力を１行ｕ（ｖ）の配列Ｘ_ｖ（ｇ）とする。その後、ステップＳ２８に遷移する。

［ステップＳ２８］ 加算器５ｂ１、５ｂ２、・・・、５ｂＶは、以下の演算を行う。Ｗ（ｖ）＝Ｗ（ｖ）＋Ｘ_ｖ（ｇ）、ｋ＝ｇ・ｕ（ｖ）。その後、ステップＳ２９に遷移する。
［ステップＳ２９］ 周期フィルタ５０は、以下の演算を行う。ｇ＝ｇ＋１。その後、ステップＳ３０に遷移する。

［ステップＳ３０］ 周期フィルタ５０は、パラメータｇが個数Ｇ以下か否かを判断する。パラメータｇが個数Ｇ以下である場合（ステップＳ３０のＹｅｓ）、ステップＳ２７に遷移する。パラメータｇが個数Ｇより大きい場合（ステップＳ３０のＮｏ）、ステップＳ３１に遷移する。
［ステップＳ３１］ 周期フィルタ５０は、以下の演算を行う。ｖ＝ｖ＋１。その後、ステップＳ３２に遷移する。

［ステップＳ３２］周期フィルタ５０は、パラメータｖが個数Ｖ以下か否かを判断する。パラメータｖが個数Ｖ以下である場合（ステップＳ３２のＹｅｓ）、ステップＳ２３に遷移する。パラメータｖが個数Ｖより大きい場合（ステップＳ３２のＮｏ）、ステップＳ３３に遷移する。

［ステップＳ３３］ 周期フィルタ５０は、ｐ（１）からｐ（Ｖ）のうち、最大となる値を特定し、その値をｐ１とする。このｐ１が１となるようにｐ（１）からｐ（Ｖ）までを正規化し、得られた値を、それぞれＰ（１）からＰ（Ｖ）として１行Ｖ列の行列Ｐ＝［Ｐ（１），Ｐ（２），・・・，Ｐ（Ｖ）］を作成する。その後、ステップＳ３４に遷移する。

［ステップＳ３４］ 第１演算部５ｄは、行列Ｐについて第１演算を実行する。また、第２演算部５ｅは、第２演算を実行する。そして、１行（Ｖ－２）列の行列Ｐ０＝［Ｐ０（２），Ｐ０（３），・・・，Ｐ０（Ｖ－１）］を作成する。その後、ステップＳ３５に遷移する。

［ステップＳ３５］ 周期フィルタ５０は、行列Ｐ００の要素のうち、Ｓ／Ｎ比がｓ［ｄＢ］（Ｓは０を超える任意の実数）以上の要素の番号（要素番号）をＶ_１１，・・・Ｖ_Ｎ１として記録する。その後、ステップＳ３６に遷移する。

［ステップＳ３６］ 周期フィルタ５０は、要素数Ｎを信号数として出力する。また、周期フィルタ５０は、それぞれの周期としてＶ_１１，・・・Ｖ_Ｎ１を出力する。また、周期フィルタ５０は、推定基本波形としてＷ（Ｖ_１１），・・・Ｗ（Ｖ_Ｎ１）とする。その後、図２０および図２１の処理を終了する。
以下、実際の計測結果を用いて第２の実施の形態の周期フィルタ５０の処理を説明する。
図２２および図２３は、第２の実施の形態の周期フィルタの処理の第１の実施例を説明する図である。

図２２（ａ）は、第１の実施例のＡＲＳスペクトルを示している。１周期を５１サンプル、５２サンプル、１０２サンプル、１０４サンプルと仮定した場合が大きく突出している。

図２２（ｂ）は、図２２（ａ）に示すＡＲＳスペクトルに対し周期フィルタ５０が第１演算を施した演算結果を示している。傾斜が取り除かれ、大きく突出しているサンプルがより強調されている。
図２３は、図２２（ｂ）の第１演算の演算結果に対し第２演算を実行した演算結果を示している。
図２４は、第２の実施の形態の周期フィルタの処理の第２の実施例を説明する図である。

図２４（ａ）は、第２の実施例のＡＲＳスペクトルを示している。第２の実施例では、累積値が大きく突出する複数の連続したサンプル（１０１～１０５）が存在する。また、連続したサンプルから少し離れた位置に大きく突出するサンプル（１０８）が存在する。

図２４（ｂ）は、図２４（ａ）に示すＡＲＳスペクトルに対し周期フィルタ５０が第１演算を施した演算結果を示している。第２の実施例のようなＡＲＳスペクトルに対しても傾斜が取り除かれ、大きく突出しているサンプルがより強調されている。
第２の実施の形態の集積回路によれば、さらに、比較的簡易な演算によりピークを容易かつ確実に検出することができる。
＜第３の実施の形態＞
次に、第３の実施の形態の集積回路について説明する。

以下、第３の実施の形態の集積回路について、前述した第１の実施の形態および第２の実施の形態との相違点を中心に説明し、同様の事項については、その説明を省略する。
第２の実施の形態の集積回路の処理においては、ノイズによって少し孤立した箇所に存在するサンプルは、相対的に累積値が大きくなる。
図２５は、第２の実施の形態の周期フィルタの処理の第３の実施例を説明する図である。
第３の実施例は、Ｓ／Ｎ比の低い入力（受信部３が受信したアナログ値）に対して周期フィルタ５０が処理を行った場合を示している。

図２５（ａ）は、第３の実施例のＡＲＳスペクトルを示している。第３の実施例では、累積値が大きく突出する複数の連続したサンプル（１０１～１０５）が存在する。また、連続したサンプルから少し離れた位置に大きく突出するサンプル（１０８）が存在する。
図２５（ｂ）は、図２５（ａ）に示す第３の実施例のＡＲＳスペクトルに対し周期フィルタ５０が第１演算を施した演算結果を示している。

Ｓ／Ｎ比の低い入力の場合、周期フィルタ５０が第１処理を施しても、サンプル（１０１～１０５及び１０８）の他に、本来は存在しないはずのサンプル（１２７）が検出されていることがわかる。このようなＳ／Ｎ比の低い入力の場合、検出されたサンプルが、信号であるのか、ノイズであるのかの判断が難しい場合がある。

第３の実施の形態の集積回路は、以下の処理を実行することにより、ノイズの影響を低減させて、ピークとノイズを容易に見分けられるようにすることを目的の１つとしている。

図２６および図２７は、第３の実施の形態の周期フィルタの処理を説明するフローチャートである。第３の実施の形態の周期フィルタの処理では、新たにステップＳ４２、Ｓ４３、Ｓ５７、Ｓ５８の処理が追加されている。
［ステップＳ４１］ 周期フィルタ５０は、以下の演算に使用する標本化回路の個数Ｖを設定する。その後、ステップＳ４２に遷移する。
［ステップＳ４２］ 周期フィルタ５０は、パラメータｍ＝１に設定する。その後、ステップＳ４３に遷移する。
［ステップＳ４３］ 周期フィルタ５０は、繰り返し回数Ｍの入力を受け付ける。その後、ステップＳ４４に遷移する。
［ステップＳ４４］ 周期フィルタ５０は、パラメータｖ＝１に設定する。その後、ステップＳ４５に遷移する。
［ステップＳ４５］ 周期フィルタ５０は、ｖ番目の標本化回路の出力ポート数をｕ（ｖ）に決定する。その後、ステップＳ４６に遷移する。
［ステップＳ４６］ 周期フィルタ５０は、全ての要素が０の１行ｕ（ｖ）の配列Ｗ（ｖ）を作成する。その後、ステップＳ４７に遷移する。
［ステップＳ４７］ 周期フィルタ５０は、パラメータｇ＝１、ｋ＝１を設定する。その後、ステップＳ４８に遷移する。
［ステップＳ４８］ 周期フィルタ５０は、累積する標本化回路出力の個数Ｇを決定する。その後、ステップＳ４９に遷移する。

［ステップＳ４９］ 加算器５ｂ１、５ｂ２、・・・、５ｂＶは、ｖ番目の標本化回路に入力サンプルｘ［ｋ］からｘ［ｋ＋ｕ_ｖ－１］までを入力する。得られた出力を１行ｕ（ｖ）の配列Ｘ_ｖ（ｇ）とする。その後、ステップＳ５０に遷移する。

［ステップＳ５０］ 加算器５ｂ１、５ｂ２、・・・、５ｂＶは、以下の演算を行う。Ｗ（ｖ）＝Ｗ（ｖ）＋Ｘ_ｖ（ｇ）、ｋ＝ｇ・ｕ（ｖ）。その後、ステップＳ５１に遷移する。
［ステップＳ５１］ 周期フィルタ５０は、以下の演算を行う。ｇ＝ｇ＋１。その後、ステップＳ５２に遷移する。

［ステップＳ５２］ 周期フィルタ５０は、パラメータｇが個数Ｇ以下か否かを判断する。パラメータｇが個数Ｇ以下である場合（ステップＳ５２のＹｅｓ）、ステップＳ４９に遷移する。パラメータｇが個数Ｇより大きい場合（ステップＳ５２のＮｏ）、ステップＳ５３に遷移する。
［ステップＳ５３］ 周期フィルタ５０は、以下の演算を行う。ｖ＝ｖ＋１。その後、ステップＳ５４に遷移する。

［ステップＳ５４］周期フィルタ５０は、パラメータｖが個数Ｖ以下か否かを判断する。パラメータｖが個数Ｖ以下である場合（ステップＳ５４のＹｅｓ）、ステップＳ４５に遷移する。パラメータｖが個数Ｖより大きい場合（ステップＳ５４のＮｏ）、ステップＳ５５に遷移する。

［ステップＳ５５］ 周期フィルタ５０は、ｐ（１）からｐ（Ｖ）のうち、最大となる値を特定し、その値をｐ１とする。このｐ１が１となるようにｐ（１）からｐ（Ｖ）までを正規化し、得られた値を、それぞれＰ（１）からＰ（Ｖ）として１行Ｖ列の行列Ｐ＝［Ｐ（１），Ｐ（２），・・・，Ｐ（Ｖ）］を作成する。その後、ステップＳ５６に遷移する。

［ステップＳ５６］ 第１演算部５ｄは、行列Ｐについて第１演算を実行する。また、第２演算部５ｅは、第２演算を実行する。そして、１行（Ｖ－２）列の行列Ｐ０＝［Ｐ０（２），Ｐ０（３）・・・，Ｐ０（Ｖ）］を作成する。その後、ステップＳ５７に遷移する。
［ステップＳ５７］ 周期フィルタ５０は、行列Ｐ０を積算した行列Ｐ００を作成する。その後、ステップＳ５８に遷移する。

［ステップＳ５８］ 周期フィルタ５０は、パラメータｍが繰返し回数Ｍ以下か否かを判断する。パラメータｍが繰返し回数Ｍ以下である場合（ステップＳ５８のＹｅｓ）、ステップＳ４４に遷移する。パラメータｍが繰返し回数Ｍより大きい場合（ステップＳ５８のＮｏ）、ステップＳ５９に遷移する。

［ステップＳ５９］ 周期フィルタ５０は、行列Ｐ００の要素のうち、Ｓ／Ｎ比がｓ［ｄＢ］（Ｓは０を超える任意の実数）以上の要素の番号（要素番号）をＶ_１１，・・・Ｖ_Ｎ１として記録する。その後、ステップＳ６０に遷移する。

［ステップＳ６０］ 周期フィルタ５０は、要素数Ｎを信号数として出力する。また、周期フィルタ５０は、それぞれの周期としてＶ_１１，・・・Ｖ_Ｎ１を出力する。また、周期フィルタ５０は、推定基本波形としてＷ（Ｖ_１１），・・・Ｗ（Ｖ_Ｎ１）とする。その後、図２６および図２７の処理を終了する。
この処理により、ノイズの影響を抑え、Ｓ／Ｎ比の低い入力であっても、一つまたは複数の信号を検出することができる。
図２８は、第３の実施の形態の周期フィルタの処理の実施例を説明する図である。
図２８（ａ）は、繰返し回数Ｍ＝５の場合の処理結果を示している。図２８（ｂ）は、繰返し回数Ｍ＝１０の場合の処理結果を示している。

前述したように、Ｓ／Ｎ比の低い入力の場合、第１処理を施しても、図２５（ｂ）においてはサンプル（１０１～１０５及び１０８）の他に本来は存在しないはずのサンプル（１２７）の値が検出されていたところ、第３の実施の形態の周期フィルタの処理を施すことにより、信号に対するノイズの値が相対的に小さくなり、サンプル（１２７）がノイズであることが容易に判断できる。
＜変形例＞

図２９および図３０は、第３の実施の形態の周期フィルタの変形例の処理を説明するフローチャートである。前述した第３の実施の形態の周期フィルタの処理とは、ステップＳ５０ａの処理が異なっている。

［ステップＳ５０ａ］ 加算器５ｂ１、５ｂ２、・・・、５ｂＶは、以下の演算を行う。Ｗ（ｖ）＝Ｗ（ｖ）×Ｘ_ｖ（ｇ）、ｋ＝ｇ・ｕ（ｖ）。その後、ステップＳ５１に遷移する。

この処理により、計算資源の消費が少ない足し算及び引き算のみを用いて、ノイズの影響を抑えつつ、一つまたは複数の複数の信号を検出することができる。また、強度比の大きな信号を検出したい場合であっても、信号とノイズを明確に区別することができる。

以上、本発明の演算装置およびプログラムを、図示の実施の形態に基づいて説明したが、本発明はこれに限定されるものではなく、各部の構成は、同様の機能を有する任意の構成のものに置換することができる。また、本発明に、他の任意の構成物や工程が付加されていてもよい。
また、本発明は、前述した各実施の形態のうちの、任意の２以上の構成（特徴）を組み合わせたものであってもよい。

なお、上記の処理機能は、コンピュータによって実現することができる。その場合、周期フィルタ５が有する機能の処理内容を記述したプログラムが提供される。そのプログラムをコンピュータで実行することにより、上記処理機能がコンピュータ上で実現される。処理内容を記述したプログラムは、コンピュータで読み取り可能な記録媒体に記録しておくことができる。コンピュータで読み取り可能な記録媒体としては、磁気記憶装置、光ディスク、光磁気記録媒体、半導体メモリ等が挙げられる。磁気記憶装置には、ハードディスクドライブ、フレキシブルディスク（ＦＤ）、磁気テープ等が挙げられる。光ディスクには、ＤＶＤ、ＤＶＤ－ＲＡＭ、ＣＤ－ＲＯＭ／ＲＷ等が挙げられる。光磁気記録媒体には、ＭＯ（Magneto-Optical disk）等が挙げられる。

プログラムを流通させる場合には、例えば、そのプログラムが記録されたＤＶＤ、ＣＤ－ＲＯＭ等の可搬型記録媒体が販売される。また、プログラムをサーバコンピュータの記憶装置に格納しておき、ネットワークを介して、サーバコンピュータから他のコンピュータにそのプログラムを転送することもできる。

プログラムを実行するコンピュータは、例えば、可搬型記録媒体に記録されたプログラムもしくはサーバコンピュータから転送されたプログラムを、自己の記憶装置に格納する。そして、コンピュータは、自己の記憶装置からプログラムを読み取り、プログラムに従った処理を実行する。なお、コンピュータは、可搬型記録媒体から直接プログラムを読み取り、そのプログラムに従った処理を実行することもできる。また、コンピュータは、ネットワークを介して接続されたサーバコンピュータからプログラムが転送される毎に、逐次、受け取ったプログラムに従った処理を実行することもできる。

また、上記の処理機能の少なくとも一部を、ＤＳＰ（Digital Signal Processor）、ＡＳＩＣ（Application Specific Integrated Circuit）、ＰＬＤ（Programmable Logic Device）等の電子回路で実現することもできる。

１ 集積回路
２ 送信部
３ 受信部
４ ＡＤＣ
５、５０ 周期フィルタ
５ａ１、５ａ２、・・・、５ａＶ 標本化回路
５ｂ１、５ｂ２、・・・、５ｂＶ 加算器
５ｃ 最大値検出部
５ｄ 第１演算部
５ｅ 第２演算部

Claims (2)

1. 入力される離散時間信号から、予め割り当てられたそれぞれの時間長が異なる基本周期候補区間に含まれる連続した離散点をパラレルに出力ポートへ出力する複数の信号出力部であって、前記基本周期候補区間の時間長の順番に並べたときに隣接する信号出力部の前記基本周期候補区間に含まれる連続した離散点の数の差が全て等しい複数の信号出力部と、
   前記各信号出力部により出力される信号の値を前記出力ポート毎にそれぞれ足し合わせた加算値を求め、当該加算値の前記基本周期候補区間毎の最大値を前記基本周期候補区間毎の累積値として算出する加算部と、
   前記基本周期候補区間を前記時間長の順番に並べたときの、二つの基本周期候補区間の順番の差を基本周期候補区間の間隔とし、各基本周期候補区間について、当該基本周期候補区間の累積値から、当該基本周期候補区間から一定の間隔だけ前の順番の基本周期候補区間の累積値を減算して第１差分を演算し、当該基本周期候補区間の累積値から、当該基本周期候補区間から一定の間隔だけ後の順番の基本周期候補区間の累積値を減算して第２差分を演算する第１演算部と、各基本周期候補区間の第１差分及び第２差分を足し合わせる第２演算部と、を備え、前記第２演算部の演算値のうち、Ｓ／Ｎ比が一定値以上の演算値に対応する基本周期候補区間を基本周期として出力する演算部と、
   を有することを特徴とする演算装置。
2. コンピュータに、
   入力される離散時間信号から、予め割り当てられたそれぞれの時間長が異なる基本周期候補区間に含まれる連続した離散点をパラレルに出力ポートへ出力する複数の信号出力部であって、前記基本周期候補区間の時間長の順番に並べたときに隣接する信号出力部の前記基本周期候補区間に含まれる連続した離散点の数の差が全て等しい複数の信号出力部により出力される信号の値を前記出力ポート毎にそれぞれ足し合わせた加算値を求め、当該加算値の前記基本周期候補区間毎の最大値を前記基本周期候補区間毎の累積値として算出し、
   前記基本周期候補区間を前記時間長の順番に並べたときの、二つの基本周期候補区間の順番の差を基本周期候補区間の間隔とし、各基本周期候補区間について、当該基本周期候補区間の累積値から、当該基本周期候補区間から一定の間隔だけ前の順番の基本周期候補区間の累積値を減算して第１差分を演算し、
   当該基本周期候補区間の累積値から、当該基本周期候補区間から一定の間隔だけ後の順番の基本周期候補区間の累積値を減算して第２差分を演算し、
   各基本周期候補区間の第１差分及び第２差分を足し合わせた演算値のうち、Ｓ／Ｎ比が一定値以上の演算値に対応する基本周期候補区間を基本周期として出力する、
   処理を実行させることを特徴とするプログラム。

Images