JP2018023618A - 振動周波数計測装置 - Google Patents

振動周波数計測装置 Download PDF

Info

Publication number
JP2018023618A
JP2018023618A JP2016157671A JP2016157671A JP2018023618A JP 2018023618 A JP2018023618 A JP 2018023618A JP 2016157671 A JP2016157671 A JP 2016157671A JP 2016157671 A JP2016157671 A JP 2016157671A JP 2018023618 A JP2018023618 A JP 2018023618A
Authority
JP
Japan
Prior art keywords
frequency
unit
value
data
related value
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
JP2016157671A
Other languages
English (en)
Other versions
JP6783091B2 (ja
Inventor
裕和 山本
Hirokazu Yamamoto
裕和 山本
哲好 柴田
Tetsuyoshi Shibata
哲好 柴田
勝 村山
Masaru Murayama
勝 村山
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Sumitomo Riko Co Ltd
Original Assignee
Sumitomo Riko Co Ltd
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Sumitomo Riko Co Ltd filed Critical Sumitomo Riko Co Ltd
Priority to JP2016157671A priority Critical patent/JP6783091B2/ja
Publication of JP2018023618A publication Critical patent/JP2018023618A/ja
Application granted granted Critical
Publication of JP6783091B2 publication Critical patent/JP6783091B2/ja
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Abstract

【課題】簡易な方法により高精度に振動周波数を計測することができる振動周波数計測装置を提供する。【解決手段】振動周波数計測装置60は、振動センサ10aの検出信号から所定周波数範囲を抽出し、フィルタ処理された処理信号に対して、複数の対象周波数fについて、基準時刻T0nの第一データP(T0n)及び基準時刻T0nから半周期前の第二データP(T0n−0.5Tf)に基づいて単位関連値Q1(n,Tf),Q2(n,Tf)を演算し、複数の対象周波数fについて、基準時刻T0nを所定時間分とする複数の単位関連値Q1(n,Tf),Q2(n,Tf)に基づいて、総合関連値R1(Tf),R2(Tf)を演算し、複数の対象周波数fの総合関連値R1(Tf),R2(Tf)に基づいて、検出信号の振動周波数を特定する。【選択図】図9

Description

本発明は、振動周波数計測装置に関するものである。
振動センサが検出した信号から振動周波数を計測するためには、フーリエ変換処理を行ったり自己相関関数の演算処理を行ったりすることが一般的であった。例えば、特許文献1,2には、生体から得られた生体信号の自己相関関数を求め、自己相関関数からピーク(極大値)を検出し、心拍や呼吸などの生体情報の周波数を計測することが記載されている。
特公昭62−4971号公報 特開2015−188603号公報
しかし、フーリエ変換処理では、バタフライ演算を行うなど、多大な演算量を必要とする。また、自己相関関数の演算処理も同様である。そのため、フーリエ変換処理や自己相関関数の演算処理を行うためには、高速な演算処理装置が必要であった。
本発明は、簡易な方法により高精度に振動周波数を計測することができる振動周波数計測装置を提供することを目的とする。
本発明に係る振動周波数計測装置は、振動センサの検出信号に含まれる振動の周波数を計測する振動周波数計測装置であって、前記検出信号から所定周波数範囲を抽出するフィルタ処理部と、前記フィルタ処理部で処理された処理信号に対して、複数の対象周波数について、基準時刻の第一データ及び前記基準時刻から半周期前の第二データに基づいて単位関連値を演算する単位関連値演算部と、前記複数の対象周波数について、前記基準時刻を所定時間分とする複数の前記単位関連値に基づいて、総合関連値を演算する総合関連値演算部と、前記複数の対象周波数の前記総合関連値に基づいて、前記検出信号の振動周波数を特定する周波数特定部とを備える。
上記振動周波数計測装置によれば、単位関連値は、基準時刻の第一データと基準時刻から半周期前の第二データとに基づいて演算される。ここで、第一データと第二データとの関係について説明する。特定周波数の正弦波において、位相(π/4)における振幅が最大値を示し、位相(3π/4)における振幅が最小値を示す。従って、位相(π/4)における振幅と位相(3π/4)における振幅との差が最大となる。また、位相(π/4)における振幅の絶対値及び位相(3π/4)における振幅の絶対値が、最大値を示す。位相(π/4)と位相(3π/4)とは、半周期(π/2)ずれている。
従って、正弦波の特定周波数に一致する周波数f1においては、半周期(π/2)ずれた第一データ及び第二データを用いることにより、単位関連値の一つは、大きな値となる。一方、正弦波の特定周波数に一致しない周波数f2においては、半周期(π/2)ずれた第一データ及び第二データを用いた場合に、単位関連値の一つは、上記の場合に比べて小さな値となる。
ただし、半周期(π/2)ずれた位相であっても、位相0°、π/2における振幅はゼロとなるため、両者の差及び両者の絶対値は、最小となる。そこで、基準時刻を所定時間分とする複数の単位関連値に基づいた総合関連値が、演算されている。総合関連値は、単位関連値と同様に、周波数f1,f2に対して同様の関係を有する。つまり、正弦波の特定周波数に一致する周波数f1においては、総合関連値は、大きな値となる。一方、正弦波の特定周波数に一致しない周波数f2においては、総合関連値は、小さな値となる。従って、周波数特定部は、総合関連値に基づいて振動周波数を特定することができる。
つまり、相互に半周期ずれた第一データと第二データとの関係を利用することにより、振動周波数を特定している。従って、単位関連値の演算及び総合関連値の演算は、フーリエ変換処理や自己相関関数の演算処理に比べると、非常に簡易な処理となる。そのため、高速な演算処理装置を用いずに、振動周波数を計測することができる。
さらに、振動周波数計測装置は、上述したように、相互に半周期ずれた第一データと第二データとの関係を利用することにより、高精度に振動周波数を特定することができる。ただし、特定周波数に一致しない周波数f2における第一データと第二データとが上記関係を有するためには、ある程度ノイズが除去された信号を用いる必要がある。そこで、単位関連値演算部は、振動センサの検出信号そのものを用いるのではなく、フィルタ処理部で処理された処理信号を用いることとしている。このように、振動周波数計測装置は、検出信号をフィルタ処理することにより、確実に振動周波数を特定することができる。
第一実施形態の生体情報計測システムの全体構成図である。 図1のセンサユニットの分解斜視図である。 図1のセンサユニットの取付位置の説明図である。 振動周波数計測装置の構成を示す図である。 圧力センサユニットの計測順序及び計測周期を示す図である。 37サイクル加算部により演算される二次信号の演算周期を示す図である。 1サイクル加算部により演算される一次信号Dのグラフである。 37サイクル加算部により演算される二次信号Sum(n)のグラフである。 第一差分演算部により演算される三次信号Sig3_1のグラフである。 第二差分演算部により演算される三次信号Sig3_2のグラフである。 単位関連値Q1(a)及び総合関連値R1(a)の基本説明、及び、これらと振動周波数との関係を示す図である。 対象周波数fの周期Tが正弦波の周期Tf0に一致する場合に、第一データから第六データの位置を示す図である。 対象周波数fの周期Tが正弦波の周期Tf0に一致しない場合に、第一データから第六データの位置を示す図である。 呼吸数計測部及び心拍数計測部による第一周波数計測処理を示すフローチャートである。 第一周波数計測処理が実行された場合の各回の演算対象及び演算結果を示す。 呼吸数計測部及び心拍数計測部による第二周波数計測処理を示すフローチャートである。 第二周波数計測処理が実行された場合の各回の演算対象及び演算結果を示す。 呼吸数計測部及び心拍数計測部による第三周波数計測処理を示すフローチャートである。 第三周波数計測処理が実行された場合の各回の演算対象及び演算結果を示す。 呼吸数計測部及び心拍数計測部による第四周波数計測処理を示すフローチャートである。 第四周波数計測処理が実行された場合の各回の演算対象及び演算結果を示す。
<1.第一実施形態>
(1−1.生体情報計測システム1の構成)
生体情報計測システム1(以下、計測システムと称する)の構成について、図1−図3を参照して説明する。計測システム1は、面状に形成されたセンサユニット10に付与された身体の生体情報を計測する。第一実施形態における計測システム1は、生体情報として、呼吸数及び心拍数を計測する。計測システム1は、センサユニット10、電源装置20、スイッチ回路40、切替制御装置50、及び、振動周波数計測装置60を備える。
センサユニット10は、可撓性を有し、面状に形成される。センサユニット10は、面法線方向に、圧縮変形可能である。例えば、センサユニット10は、図2に示すように、1列の第一電極11と、4列の第二電極12と、誘電層13とを備える。なお、第一電極11と第二電極12の列数は、適宜変更可能である。
第二電極12は、センサユニット10の面法線方向に、第一電極11に対して距離を隔てて配置される。第二電極12は、帯状に形成され、相互に平行に配置される。第二電極12の延在方向は、第一電極11の延在方向(図1,2の左右方向)に対して直交する方向(図1の上下方向)である。誘電層13は、弾性変形可能な面状に形成され、第一電極11と第二電極12との間に配置される。
第一電極11及び第二電極12は、エラストマー中に導電性フィラーを配合させることにより成形される。第一電極11及び第二電極12は、可撓性を有し、伸縮自在な性質を有する。第一電極11及び第二電極12を構成するエラストマーには、例えば、シリコーンゴム、エチレン−プロピレン共重合ゴム、天然ゴム、スチレン−ブタジエン共重合ゴム、アクリロニトリル−ブタジエン共重合ゴム、アクリルゴム、エピクロロヒドリンゴム、クロロスルホン化ポリエチレン、塩素化ポリエチレン、ウレタンゴムなどが適用される。また、第一電極11及び第二電極12に配合される導電性フィラーには、導電性を有する粒子であればよく、例えば、炭素材料や金属等の微粒子が適用される。
誘電層13は、エラストマーにより成形され、可撓性を有し且つ伸縮自在な性質を有する。誘電層13を構成するエラストマーには、例えば、シリコーンゴム、アクリロニトリル−ブタジエン共重合ゴム、アクリルゴム、エピクロロヒドリンゴム、クロロスルホン化ポリエチレン、塩素化ポリエチレン、ウレタンゴムなどが適用される。
従って、第一電極11とそれぞれの第二電極12との対向位置A,B,C,Dが列状に位置する。センサユニット10は、配列された複数の対向位置A,B,C,Dにおいて、静電容量型センサとして機能する圧力センサセル10a(振動センサ)を備える。このように、第一実施形態においては、センサユニット10は、横4列に配列された4個の圧力センサセル10aを備える。そして、4個の圧力センサセル10aが、面状に配列されている。なお、第一電極11を複数列にする場合には、第一電極11と第二電極12との対向位置がマトリックス状となる。
そして、センサユニット10が面法線方向に圧縮する力を受けた場合には、誘電層13が圧縮変形することにより、第一電極11と第二電極12の離間距離が短くなる。つまり、第一電極11と第二電極12との間の静電容量が大きくなる。
ここで、センサユニット10は、図3に示すように、例えば、シート70の座面71前方の内部に配置される。詳細には、センサユニット10は、座面71前方の表皮の裏面側に配置される。第一実施形態においては、センサユニット10は、第二電極12の延在方向がシート70の前後方向に一致するように、座面71に配置される。また、センサユニット10は、左右の大腿部に相当する範囲に配置される。つまり、センサユニット10は、着座する人の左右の大腿部によって、体圧を受ける。そして、第二電極12の延在方向は、大腿部の延在方向、さらには、大腿動脈の延在方向に一致する。つまり、センサユニット10は、大腿動脈の脈波や呼吸成分の影響を受ける。
なお、センサユニット10は、シート70の座面71前方の他に、座面71後方、背面72やヘッドレスト73に配置してもよい。センサユニット10が座面後方に配置される場合は、センサユニット10は、人の臀部により体圧を受け、人の臀部における動脈の脈波や呼吸成分の影響を受ける。また、背面72に配置される場合には、センサユニット10は、人の背部により体圧を受け、人の背部における動脈の脈波や呼吸成分の影響を受ける。また、センサユニット10がヘッドレスト73に配置される場合には、センサユニット10は、人の頭部により体圧を受け、例えば首部における動脈の脈波や呼吸成分の影響を受ける。
電源装置20は、所定の電圧を発生し、センサユニット10の第一電極11に対して所定電圧を印加する。スイッチ回路40は、複数のスイッチにより構成される。各スイッチの一端は、対応する第二電極12に接続され、各スイッチの他端は、後述する振動周波数計測装置60に接続される。図1においては、左側から1列目の第二電極12に対応するスイッチがONされ、他はOFFされている。切替制御装置50は、スイッチ回路40の各スイッチのON/OFFの切替を実行する。そして、切替制御装置50は、計測対象とする圧力センサセル10aを、電源装置20及び振動周波数計測装置60に接続させる。
振動周波数計測装置60は、計測対象の圧力センサセル10aによる検出値を取得して、当該検出値に基づいて圧力センサセル10aに面する身体の脈波及び呼吸成分を抽出して、心拍数及び呼吸数を演算する。具体的には、振動周波数計測装置60は、脈波と呼吸成分の計測対象の圧力センサセル10aの静電容量の変化を取得して、脈波の影響により受ける力の変化を抽出して心拍数を演算する。さらに、振動周波数計測装置60は、当該圧力センサセル10aの静電容量の変化を取得して、呼吸の影響により受ける力の変化を抽出して呼吸数を演算する。
(1−2.振動周波数計測装置60の構成)
(1−2−1.振動周波数計測装置60の基本構成)
次に、振動周波数計測装置60(以下、計測装置と称する)の基本構成について、図4を参照して説明する。図4に示すように、計測装置60は、圧力センサセル10aによる検出信号に対してフィルタ処理を行うフィルタ処理部61と、処理信号に基づいて呼吸数を計測する呼吸数計測部62と、処理信号に基づいて心拍数を計測する心拍数計測部63とを備える。
(1−2−2.フィルタ処理部61の構成)
フィルタ処理部61について、図4−図6、図7A及び図7Bを参照して説明する。フィルタ処理部61は、圧力センサセル10aによる検出信号から所定周波数範囲を抽出する。フィルタ処理部61は、ローパスフィルタ、バンドパスフィルタ、及び、これらのフィルタ機能に相当する演算処理を含む。ここで、計測装置60が計測可能な1分間当たりの心拍数を48〜180とし、計測装置60が計測可能な1分間当たりの呼吸数を6〜27とする。この場合、心拍数48〜180は、0.8〜3.0Hzとなる。また、呼吸数6〜27は、0.10〜0.45Hzとなる。そこで、フィルタ処理部61は、少なくとも0.10〜3.00Hzを含む周波数範囲を抽出する。
まず、1サイクルの周期Ta及び計測周期Tbについて図5を参照して説明する。センサユニット10は、4個の圧力センサセル10aのそれぞれによる計測を行う。具体的には、スイッチ回路40を切り替えて、計測対象の圧力センサセル10aをA→B→C→Dの順に切り替える。従って、同一の圧力センサセル10aによる計測周期Tbは、0.25msとなる。すなわち、同一の圧力センサセル10aによる計測周波数は、4kHzとなる。
また、1サイクルの周期Taは、4msとする。つまり、1サイクルの周波数は、250Hzとなる。従って、1サイクルにおいて、同一の圧力センサセル10aによる計測は、16回ずつ行われる。なお、1サイクルの周期、同一の圧力センサセル10aによる計測周期は、適宜変更することができる。
そして、フィルタ処理部61は、1サイクル加算部611、37サイクル加算部612、第一差分演算部613、及び、第二差分演算部614を備える。1サイクル加算部611は、計測周期Tb(0.25ms)で、各圧力センサセル10aによる検出信号を一次信号Dとして取得して、記憶する。一次信号Dの挙動は、図7Aに示す。つまり、一次信号Dは、非常に多くのノイズを含んでいる。なお、図7Aの横軸は、一次信号Dの取得回数、すなわち一次信号Dのサンプリング回数である。
さらに、1サイクル加算部611は、「A」〜「D」のそれぞれの一次信号Dを1サイクル分加算して、記憶する。上述したように、各圧力センサセル10aによる計測は、1サイクルにおいて、16回実施される。従って、1サイクル加算部611は、16回分の「A」の一次信号Dを加算すると共に、「B」〜「D」についても同様に処理する。
37サイクル加算部612は、1サイクル加算部611により得られた1サイクル分のデータを、連続する37サイクル分加算する。37サイクル分加算されたデータは、二次信号Sum(n)として記憶される。つまり、37サイクル加算部612は、148msの期間(第一時間長N)に含まれる「A」の圧力センサセル10aによる592回分の検出値(一次信号D)を加算して、二次信号Sum(n)を生成する。
さらに、37サイクル加算部612は、第二時間間隔a(a≦N)ずつずらしながら、二次信号Sum(n)の生成を順次実行する。第一実施形態においては、第二時間間隔aは、1サイクル、すなわち4msとするが、1サイクルとは異なる時間としてもよい。従って、図6に示すように、1サイクル分ずつずらしながら、二次信号Sum(1)、Sum(2)、Sum(3)、・・・、Sum(n)が生成される。そして、37サイクル加算部612は、「B」〜「D」の圧力センサセル10aについても同様に、それぞれ二次信号Sum(n)を生成する。この場合、二次信号Sum(n)のデータ周期は、1サイクル、すなわち4msである。
なお、第二時間間隔aは、37サイクル分に相当する第一時間長N以下にする。特に、高精度に脈波及び呼吸成分を計測するためには、第二時間間隔aは、37サイクル分に相当する時間長Nに比べて十分に短くするとよい。
二次信号Sum(n)の挙動は、図7Bに示すように、一次信号Dに比べると、ノイズが大幅に低減している。特に、圧力センサセル10aによる計測周期Tbの2分の1以下のノイズ成分は、加算処理によって、除去される。二次信号Sum(n)において、短周期に変化する成分が脈波による変化を表しており、長周期の振動が呼吸による変化を表している。つまり、二次信号Sum(n)は、脈波及び呼吸成分を表す信号となっている。具体的には、図7Bにおいて、丸印の箇所(下方に凹んでいる箇所)が脈波に対応する箇所であり、矢印の箇所(大きく上に突出している箇所)が呼吸成分に対応する箇所である。つまり、所定時間において、丸印の箇所の数が脈拍数となり、矢印の箇所の数が呼吸数となる。
第一差分演算部613及び第二差分演算部614は、第三時間間隔M1、M2(M1、M2≧a)だけずれた二個の二次信号Sum(n)の差分である三次信号Sig3_1、Sig3_2を演算する。
第一差分演算部613は、直前の二個の二次信号Sum(m)、Sum(m−1)の差分である三次信号Sig3_1を演算する。つまり、三次信号Sig3_1は、第三時間間隔M1として、1サイクル分の時間間隔(第二時間間隔aに等しい)、すなわち4msだけずれた二個の二次信号Sum(m)、Sum(m−1)の差分である。
第一差分演算部613により演算される三次信号Sig3_1の挙動は、図7Cに示すように、二次信号Sum(n)の挙動とそれほど変わりない。つまり、三次信号Sig3_1には、脈波及び呼吸成分を表す信号である。具体的には、図7Cにおいて、丸印の箇所(下方に凹んでいる箇所)が脈波に対応する箇所であり、矢印の箇所(大きく上に突出している箇所)が呼吸成分に対応する箇所である。ここで、三次信号Sig3_1は、二次信号Sum(n)に比べて、絶対値が小さくなっている。そのため、三次信号Sig3_1は、二次信号Sum(n)よりも、データ容量が小さく、高速処理を行うことができる。なお、三次信号Sig3_1のデータ周期は、1サイクル、すなわち4msである。
第二差分演算部614は、相対的に長い時間だけずれた二個の二次信号Sum(m)とSum(1)の差分である三次信号Sig3_2を演算する。ここで、第三時間間隔M2は、第一差分演算部613にて用いる第三時間間隔M1よりも十分に長い。第一実施形態では、第三時間間隔M2は、250msである。250msとは、一般的な脈波が1Hz程度であるとすると、その脈波の4分の1である。つまり、この第三時間間隔M2は、脈波の最大振幅の半分程度に相当する周期に等しい。従って、第三時間間隔M2だけずれた二個の二次信号Sum(m)、Sum(1)の差分は、比較的、脈波の影響分を大きく抽出作用がある。
第二差分演算部614により演算される三次信号Sig3_2の挙動は、図7Dに示すような挙動となる。図7Dにおいて、丸印の箇所(下方に凹んでいる箇所)が脈波に対応する箇所であり、矢印の箇所(大きく上に突出している箇所)が呼吸成分に対応する箇所である。三次信号Sig3_2は、二次信号Sum(n)と比べると、呼吸の変化に対して脈波による変化が際立っている。さらに、三次信号Sig3_2は、二次信号Sum(n)に比べて、絶対値が小さくなっている。そのため、三次信号Sig3_1は、二次信号Sum(n)よりも、データ容量が小さく、高速処理を行うことができる。なお、三次信号Sig3_2のデータ周期は、1サイクル、すなわち4msである。
(1−2−3.呼吸数計測部62及び心拍数計測部63の基本構成)
次に、呼吸数計測部62及び心拍数計測部63の基本構成について、図4、図7C及び図7Dを参照して説明する。呼吸数計測部62は、第一差分演算部613により演算された三次信号Sig3_1を、呼吸成分を表すデータとして取得する(図7C)。さらに、呼吸数計測部62は、三次信号Sig3_1に基づいて、呼吸に起因する振動周波数を演算し、1分間当たりの呼吸数を得る。
心拍数計測部63は、第二差分演算部614により演算された三次信号Sig3_2を、脈波成分を表すデータとして取得する(図7D)。さらに、心拍数計測部63は、三次信号Sig3_2に基づいて、心拍に起因する振動周波数を演算し、1分間当たりの心拍数を得る。
呼吸数計測部62及び心拍数計測部63は、対象の周波数が異なるのみで、実質的に同様の処理を行う。呼吸数計測部62及び心拍数計測部63は、単位関連値演算部621,631、総合関連値演算部622,632、周波数特定部623,633を備える。
(1−2−4.単位関連値及び総合関連値の基本的考え方)
次に、単位関連値Q1(a)及び総合関連値R1(a)の基本説明、及び、これらと振動周波数との関係について、図8を参照して説明する。ただし、説明を分かりやすくするために、正弦波信号を用いて説明する。また、図8の横軸は、時間であり、縦軸は、振幅を示す。
単位関連値Q1(a)は、基準時刻θ0の第一データ及び基準時刻θ0から半周期(π)前の第二データに基づいて演算される。第一実施形態では、単位関連値Q1(a)は、式(1)で表される。つまり、単位関連値Q1(a)は、基準時刻θ0における振幅sin(θ0)と、基準時刻θ0から半周期前の時刻(θ0−π)における振幅sin(θ0−π)との差分の絶対値である。正弦波において、半周期ずれた2つの値の中で、位相がπ/4の値と位相が3π/4の値との差分が最大となる。一方、位相が0の値と位相がπ/2の値との差分は、0となる。
Figure 2018023618
そして、総合関連値R1(a)は、基準時刻θ0を所定時間分とする複数の単位関連値Q1(a)に基づいて演算される。第一実施形態では、総合関連値R1(a)は、式(2)で表される。つまり、総合関連値R1(a)は、基準時刻θ0を当該周波数の1周期分とする複数の単位関連値Q1(a)を加算した値である。ここでは、総合関連値R1(a)は、基準時刻θ0を当該周波数の1周期分としたが、半周期分としてもよいし、1.5周期分、2周期分などとしてもよい。
Figure 2018023618
ここで、単位関連値Q1(a)が基準時刻θ0の振幅(第一データ)と基準時刻θ0から半周期(π)前の振幅(第二データ)とに基づいて演算した。仮に、単位関連値Q1として、基準時刻θ0の振幅と、基準時刻θ0から半周期(π)とは異なる時刻前の振幅とに基づいて演算した場合に、総合関連値R1の大きさについて比較検討する。
第一比較例としての単位関連値Q1(b)は、式(3)で表される。つまり、第一比較例としての単位関連値Q1(b)は、基準時刻θ0の振幅sin(θ0)と、基準時刻θ0から(π−α)だけ前の振幅sin(θ0−π+α)との差分の絶対値である。
Figure 2018023618
そして、第一比較例としての総合関連値R1(b)は、式(4)で表される。つまり、第一比較例としての総合関連値R1(b)は、基準時刻θ0を当該周波数の1周期分とする複数の単位関連値Q1(b)を加算した値である。
Figure 2018023618
また、第二比較例としての単位関連値Q1(c)は、式(5)で表される。つまり、第二比較例としての単位関連値Q1(c)は、基準時刻θ0の振幅sin(θ0)と、基準時刻θ0から(π+β)だけ前の振幅sin(θ0−π−β)との差分の絶対値である。
Figure 2018023618
そして、第二比較例としての総合関連値R1(c)は、式(6)で表される。つまり、第二比較例としての総合関連値R1(c)は、基準時刻θ0を当該周波数の1周期分とする複数の単位関連値Q1(c)を加算した値である。
Figure 2018023618
第一実施形態における総合関連値R1(a)と、比較例としての総合関連値R1(b),R1(c)とは、式(7)の関係を有する。つまり、単位関連値Q1(a)が、基準時刻θ0における振幅sin(θ0)と、基準時刻θ0から半周期前の時刻(θ0−π)における振幅sin(θ0−π)との差分の絶対値とすることで、総合関連値R1(a)が最大値となる。
Figure 2018023618
(1−2−5.単位関連値演算部621,631)
次に、単位関連値演算部621,631について、図4、図9及び図10を参照して説明する。呼吸数計測部62の単位関連値演算部621は、三次信号Sig3_1を用いて単位関連値Q1(n,T)を演算する。心拍数計測部63の単位関連値演算部631は、三次信号Sig3_2を用いて単位関連値Q1(n,T)を演算する。
そして、単位関連値演算部621,631は、複数の対象周波数fのそれぞれについて、単位関連値Q1(n,T)を演算する。呼吸数計測部62における単位関連値演算部621は、0.10〜0.45Hzを対象周波数fとする。心拍数計測部63における単位関連値演算部631は、0.8〜3.0Hzを対象周波数fとする。
ここで、単位関連値演算部621,631は、三次信号Sig3_1,Sig3_2が演算される都度、単位関連値Q1(n,T)を演算する。つまり、対象周波数fに関わりなく、単位関連値Q1(n,T)の演算周波数は一定(第一実施形態では、250Hz、周期4ms)である。ただし、対象周波数fに応じて、単位関連値Q1(n,T)の演算周波数を変更してもよい。例えば、対象周波数fが高周波であるほど、単位関連値Q1(n,T)の演算周波数を大きくし、対象周波数fが低周波であるほど、単位関連値Q1(n,T)の演算周波数を小さくする。
そして、単位関連値Q1(n,T)は、式(8)で表される。単位関連値Q1(n,T)は、基準時刻T0の第一データP(T0)、基準時刻T0から半周期前の第二データP(T0−0.5T)、基準時刻T0から1周期前の第三データP(T0−T)、基準時刻T0から1.5周期前の第四データP(T0−1.5T)、基準時刻T0から2周期前の第五データP(T0−2T)、基準時刻T0から2.5周期前の第六データP(T0−2.5T)に基づいて演算される。
換言すると、単位関連値Q1(n,T)は、第一データP(T0)から、第二データP(T0−0.5T)を減算し、第三データP(T0−T)を加算し、第四データP(T0−1.5T)を減算し、第五データP(T0−2T)を加算し、第六データP(T0−2.5T)を減算し、その結果の絶対値をとる。つまり、単位関連値Q1(n,T)は、データの半周期毎に減算と加算を交互に繰り返して得られる。
Figure 2018023618
そして、第一実施形態においては、基準時刻T0が、三次信号Sig3_1,Sig3_2が演算される周期で変化する。第一実施形態においては、カウント値nは、4ms毎に1増加する。つまり、基準時刻T0は、4ms経過の都度、次の演算時刻となる。
式(8)について、詳細には、単位関連値Q1(n,T)は、第一データと第二データの差分、第三データと第四データの差分、及び、第五データと第六データの差分を加算した値の絶対値である。なお、単位関連値Q1(n,T)は、第一データと第二データの差分の絶対値としてもよいし、第一データと第二データの差分と第三データと第四データの差分とを加算した値の絶対値としてもよい。
対象周波数fの周期Tが正弦波の周期Tf0に一致する場合に、第一データから第六データは、図9に示す各位置Pとなる。この場合、第一データ、第三データ及び第五データと、第二データ、第四データ及び第六データとは、正負反対の値となる。
一方、対象周波数fの周期Tが正弦波の周期Tf0に一致しない場合に、第一データから第六データは、図10に示す各位置Pとなる。この場合、第一データ、第三データ及び第五データと、第二データ、第四データ及び第六データとは、常に正負反対の値になるわけではない。
(1−2−6.総合関連値演算部622,632)
次に、総合関連値演算部622,632について説明する。総合関連値演算部622,632は、単位関連値演算部621,631により演算された単位関連値Q1(n,T)を用いて、複数の対象周波数fのそれぞれについて、総合関連値R1(T)を演算する。総合関連値R1(T)は、式(9)で表される。
Figure 2018023618
総合関連値R1(T)は、基準時刻T0を所定時間分(第一実施形態では1周期分)とする複数の単位関連値Q1(n,T)に基づいて演算される。具体例を挙げて、総合関連値R1(T)について説明する。第一実施形態においては、カウント値nは、4ms毎に1増加するため、例えば、対象周波数fが1.0Hzの場合には、1周期分のカウント値nは、250となる。そして、開始時刻Txにおけるカウント値n(Tx)が250とする場合に、1周期前の時刻(Tx−T)は1000msとなり、1周期前のカウント値n(Tx−Tf)は、0となる。そして、開始時刻Txにおける基準時刻T0250における単位関連値Q1(250,1.0)から、基準時刻T00における単位関連値Q1(0、1.0)までを加算した値が、総合関連値R1(1.0)となる。上記の演算を他の対象周波数fについても行い、複数の対象周波数fの全てについての総合関連値R1(T)が得られる。
このようにして、呼吸数計測部62の総合関連値演算部622は、0.10〜0.45Hzの周波数範囲について、例えば、0.01Hz間隔で、総合関連値R1(T)を取得する。また、心拍数計測部63の総合関連値演算部632は、0.8〜3.0Hzの周波数範囲について、例えば、0.1Hz間隔で、総合関連値R1(T)を取得する。
ここで、図9に示すように、対象周波数fの周期Tが正弦波の周期Tf0に一致する場合における総合関連値R1(Tf0)と、図10に示すように、対象周波数fの周期Tが正弦波の周期Tf0に一致しない場合における総合関連値R1(Tf1)とを比較する。この場合、上述した式(7)から分かるように、一致する場合の総合関連値R1(Tf0)と、一致しない場合の総合関連値R1(Tf1)とは、式(10)の関係を有する。
Figure 2018023618
従って、対象周波数fの周期Tが正弦波の周期Tf0に一致する場合における総合関連値R1(Tf0)は、一致しない場合における全ての総合関連値R1(Tf1)よりも必ず大きな値となる。
(1−2−7.周波数特定部623,633による周波数特定処理)
次に、周波数特定部623,633による周波数特定処理について説明する。周波数特定部623,633は、複数の対象周波数fの総合関連値R1(T)に基づいて、検出信号の振動周波数、具体的には、検出信号に含まれる呼吸の周波数及び心拍の周波数を演算する。
呼吸の周波数に60を乗算することで、1分間当たりの呼吸数になる。また、心拍の周波数に60を乗算することで、1分間当たりの心拍数になる。つまり、呼吸数計測部62の周波数特定部623は、呼吸の周波数を演算することで、1分間当たりの呼吸数を得る。心拍数計測部63の周波数特定部633は、心拍の周波数を演算することで、1分間当たりの心拍数を得る。
ここで、周波数特定部623,633による処理は、単位関連値演算部621,631及び総合関連値演算部622,632の処理と相互に関連し合いながら、実行される。単位関連値演算部621,631、総合関連値演算部622,632、周波数特定部623,633による周波数計測処理は、第一周波数計測処理から第四周波数計測処理の何れか一つを採用する。以下に、各周波数計測処理について詳細に説明する。
(1−2−8−1.第一周波数計測処理)
第一周波数計測処理について、図11及び図12を参照して説明する。ただし、図11及び図12は、心拍数計測部63による第一周波数計測処理を例にあげる。この場合、心拍の周波数は、0.8〜3.0Hzの範囲内に含まれるものとする。なお、呼吸数計測部62による第一周波数計測処理も実質的に同様である。
第一周波数計測処理は、一部の所定数の対象周波数fについての総合関連値R1(T)を比較することにより、心拍の周波数を特定する。そして、第一周波数計測処理は、予め設定された開始周波数についての総合関連値R1(T)を演算し、演算した総合関連値R1(T)の傾向に基づいて、対象周波数fを順次変更していく。特に、第一周波数計測処理は、対象周波数fの中で最も高い周波数を開始周波数として、対象周波数fを順次低くしていくことにより、心拍の周波数を特定する。
さらに、第一周波数計測処理は、心拍の周波数が一旦特定された後において、特定された周波数を含む所定数の対象周波数fについての総合関連値R1(T)に基づいて心拍の周波数が変化する傾向を決定する。そして、第一周波数計測処理は、当該傾向に基づいて連続する所定数の対象周波数fについての最新の総合関連値R1(T)を演算し、演算した最新の総合関連値R1(T)を比較することで、変化後の心拍の周波数を特定するようにしている。
まずは、図11を参照して第一周波数計測処理を説明する。図11に示すように、単位関連値演算部631は、最初に、基準となる対象周波数fcを決定する(図11のS1)。ここでは、基準となる対象周波数fcは、(3.0−0.1)Hzすなわち、2.9Hzとする。
続いて、単位関連値演算部631は、対象周波数fを(fc−0.1)、fc、(fc+0.1)とする単位関連値Q1(n,T(fc−0.1))、Q1(n,Tfc)、Q1(n,T(fc+0.1))を演算する(図11のS2)。続いて、総合関連値演算部632は、対象周波数fを(fc−0.1)、fc、(fc+0.1)とする3つの総合関連値R1(T(fc−0.1))、R1(Tfc)、R1(T(fc+0.1))を演算する(図11のS3)。
続いて、周波数特定部633は、3つの総合関連値R1(T(fc−0.1))、R1(Tfc)、R1(T(fc+0.1))を比較して、大小関係を取得する(図11のS4)。総合関連値R1(Tfc)が最も大きい場合には(図11のS5:Y)、得られた総合関連値R1(T)の中で、最大値R1(Tfmax)の周波数fmaxを、検出信号に含まれる振動周波数、すなわち心拍の周波数とする(図11のS9)。
続いて、単位関連値演算部631は、前回のデータPを取得してから4ms経過して、次のデータPを取得したか否かを判定する(図11のS10)。次の取得データPが無い場合には(図11のS10:N)、次のデータPを取得する時間に達していない、すなわち4ms経過していないことになるため、単位関連値演算部631は、次のデータPを取得するまで待機する。一方、次の取得データPが有る場合には(図11のS10:Y)、単位関連値演算部631は、再び、S2に戻り処理を繰り返す。
図11のS5において、取得した3つの中で総合関連値R1(Tfc)が最大ではない場合には(図11のS5:N)、周波数特定部633は、3つの総合関連値R1の比較結果に基づいて振動周波数の傾向を決定し、決定された傾向に基づいて連続する3つの対象周波数を順次変化させる。
すなわち、図11のS5において、取得した3つの中で総合関連値R1(Tfc)が最大ではない場合には(図11のS5:N)、周波数特定部633は、高周波側の総合関連値R1(T(fc+0.1))が低周波側の総合関連値R1(T(fc−0.1))より大きいか否かを判定する(図11のS6)。高周波側の総合関連値R1(T(fc+0.1))が低周波側の総合関連値R1(T(fc−0.1))より大きい場合には(図11のS6:Y)、周波数特定部633は、周波数(fc+0.1)を基準となる対象周波数fcに変更する(図11のS6)。つまり、周波数特定部633は、3つの総合関連値R1に基づいて心拍の周波数がもっと高い周波数であるとの傾向があると判断し、基準となる対象周波数fcを前回から0.1Hz大きな値に変更する。
基準となる対象周波数fcが変更された後には、単位関連値演算部631がS10から処理を継続する。つまり、変更された基準となる対象周波数fcを含む3つの対象周波数(fc−0.1)、fc、(fc+0.1)について、単位関連値Q1(n,T)及び総合関連値R1(T)が演算される(S2,S3)。
一方、図11のS6において、低周波側の総合関連値R1(T(fc−0.1))が高周波側の総合関連値R1(T(fc+0.1))より大きい場合には(図11のS6:N)、周波数特定部633は、基準となる対象周波数(fc−0.1)を基準となる対象周波数fcに変更する(図11のS8)。つまり、周波数特定部633は、3つの総合関連値R1に基づいて心拍の周波数がもっと低い周波数であるとの傾向があると判断し、基準となる対象周波数fcを前回から0.1Hz小さな値に変更する。
基準となる対象周波数fcが変更された後には、単位関連値演算部631がS10から処理を継続する。つまり、変更された基準となる対象周波数fcを含む3つの対象周波数(fc−0.1)、fc、(fc+0.1)について、単位関連値Q1(n,T)及び総合関連値R1(T)が演算される(S2,S3)。
第一周波数計測処理について、図12を参照して、具体例をあげて説明する。被検者の心拍の周波数は、1.2〜1.4Hz程度とする。ここでは、対象周波数fのうち最も高い周波数3.0Hzが開始周波数に予め設定される(図11のS1,S2)。そして、最も高い周波数3.0Hzを含む連続する3つ(一部の所定数)の対象周波数fについて、すなわち2.8Hz、2.9Hz、3.0Hzの周波数について、単位関連値Q1(n,T2.8)、Q1(n,T2.9)、Q1(n,T3.0)及び総合関連値R1(T2.8)、R1(T2.9)、R1(T3.0)が演算される(図11のS2,S3)。
ここで、被検者の心拍の周波数は、2.8〜3.0Hzに対して小さいため、これらの中で最小の対象周波数2.8Hzの総合関連値R1(T2.8)が最大となる。つまり、周波数特定部633は、心拍の周波数がもっと低い周波数帯に存在するとの傾向を取得する。従って、図11のS5:N→S6:Nとなり、基準となる対象周波数fcは、0.1減算されて2.8Hzとなる(図11のS8)。
次のデータPが4ms後に取得されると(図11のS10)、連続する3つの対象周波数fについて、すなわち2.7Hz、2.8Hz、2.9Hzの周波数について、単位関連値Q1(n,T2.7)、Q1(n,T2.8)、Q1(n,T2.9)及び総合関連値R1(T2.7)、R1(T2.8)、R1(T2.9)が演算される(図11のS2,S3)。被検者の心拍の周波数は、2.7〜2.9Hzに対して小さいため、これらの中で最小の対象周波数2.7Hzの総合関連値R1(T2.7)が最大となる。従って、再び、図11のS5:N→S6:Nとなり、基準となる対象周波数fcは、0.1減算されて2.7Hzとなる(図11のS8)。
この処理を数回繰り返して開始から64ms経過後には、連続する3つの対象周波数fについて、すなわち1.2Hz、1.3Hz、1.4Hzの周波数について、単位関連値Q1(n,T1.2)、Q1(n,T1.3)、Q1(n,T1.4)及び総合関連値R1(T1.2)、R1(T1.3)、R1(T1.4)が演算される(図11のS2,S3)。これらの中で最小の対象周波数1.2Hzの総合関連値R1(T1.2)が最大となっている。従って、再び、図11のS5:N→S6:Nとなり、基準となる対象周波数fcは、0.1減算されて1.2Hzとなる(図11のS8)。
そして、処理開始から68ms経過後に、連続する3つの対象周波数fについて、すなわち1.1Hz、1.2Hz、1.3Hzの周波数について、単位関連値Q1(n,T1.1)、Q1(n,T1.2)、Q1(n,T1.3)及び総合関連値R1(T1.1)、R1(T1.2)、R1(T1.3)が演算される(図11のS2,S3)。これらの中で中央の対象周波数1.2Hzの総合関連値R1(T1.2)が最大となっている。そうすると、図11のS5:Yとなり、周波数特定部633は、最大値となる総合関連値R1(T1.2)の周波数fmaxを心拍の周波数に特定する。心拍の周波数が1.2Hzであるので、1分間当たりの心拍数は、72回となる。なお、図12において、横軸の時間の数字に枠で囲んだ標記としているものは、心拍の周波数が特定されたタイミングを意味する。
心拍の周波数が一旦特定された後にも、心拍数計測部63は、心拍の周波数の計測を継続して、心拍の周波数の変化を計測し続ける。すなわち、次のデータPを取得した場合、すなわち処理開始から72ms経過後において、前回と同様の対象周波数1.1Hz、1.2Hz、1.3Hzの周波数について、単位関連値Q1(n,T1.1)、Q1(n,T1.2)、Q1(n,T1.3)が演算される(図11のS2)。ここで、現時点(開始時刻Tx)における第一データP(T0n(Tx))が新たに取得される。従って、単位関連値Q1(n,T1.1)、Q1(n,T1.2)、Q1(n,T1.3)のそれぞれにおいて、カウント値nが現時点(開始時刻Tx)のカウント値n(Tx)の場合のみについての演算が実施される。つまり、カウント値nが、現時点(開始時刻Tx)以外のカウント値の場合には、既に演算が実施されている。従って、演算量は非常に少ない。
続いて、総合関連値R1(T1.1)、R1(T1.2)、R1(T1.3)が演算される(図11のS3)。これらの中で最大の対象周波数1.3Hzの総合関連値R1(T1.3)が最大となっている。つまり、周波数特定部633は、心拍の周波数がもっと高い周波数帯に存在するとの傾向を取得する。従って、再び、図11のS5:N→S6:Yとなり、基準となる対象周波数fcは、0.1加算されて1.3Hzとなる(図11のS7)。
そして、処理開始から76ms経過後に、連続する3つの対象周波数fについて、すなわち1.2Hz、1.3Hz、1.4Hzの周波数について、単位関連値Q1(n,T1.2)、Q1(n,T1.3)、Q1(n,T1.4)及び総合関連値R1(T1.2)、R1(T1.3)、R1(T1.4)が演算される(図11のS2,S3)。これらの中で中央の対象周波数1.3Hzの総合関連値R1(T1.3)が最大となっている。そうすると、図11のS5:Yとなり、周波数特定部633は、最大値となる総合関連値R1(T1.3)の周波数fmaxを心拍の周波数に特定する。この場合、心拍の周波数が1.3Hzであるので、1分間当たりの心拍数は、78回となる。このように、データPを取得するたびに、上記処理を繰り返し、心拍の周波数が変化した場合に、現在の心拍の周波数を特定することができる。
なお、上記処理では、一旦、心拍の周波数が特定された後において、3つの中で中央の対象周波数fcの総合関連値R1(Tfc)が最大となる場合に、心拍の周波数が特定されることとした。これに限られず、一旦、心拍の周波数が特定された後においては、3つの中で中央の対象周波数fcの総合関連値R1(Tfc)が最大とならない場合であっても、3つの中で最大となる総合関連値R1に対応する対象周波数fが心拍の周波数であると特定してもよい。
(1−2−8−2.第二周波数計測処理)
第二周波数計測処理について、図13及び図14を参照して説明する。ただし、図13及び図14は、心拍数計測部63による第二周波数計測処理を例にあげる。なお、呼吸数計測部62による第二周波数計測処理も実質的に同様である。
第二周波数計測処理は、一部の所定数の対象周波数fについての総合関連値R1(T)を比較することにより、心拍の周波数を特定する。そして、第二周波数計測処理は、予め設定された開始周波数についての総合関連値R1(T)を演算し、演算した総合関連値R1(T)の傾向に基づいて、対象周波数fを順次変更していく。特に、第二周波数計測処理は、対象周波数fの中で最も低い周波数を開始周波数として、対象周波数fを順次高くしていくことにより、心拍の周波数を特定する。
さらに、第二周波数計測処理は、心拍の周波数が一旦特定された後において、特定された周波数を含む所定数の対象周波数についての総合関連値R1に基づいて心拍の周波数が変化する傾向を決定する。そして、第二周波数計測処理は、当該傾向に基づいて連続する所定数の対象周波数fについての最新の総合関連値R1(T)を演算し、演算した最新の総合関連値R1を比較することで、変化後の心拍の周波数を特定するようにしている。
まずは、図13を参照して第二周波数計測処理を説明する。図13に示すように、単位関連値演算部631は、最初に、基準となる対象周波数fcを決定する(図13のS11)。ここでは、基準となる対象周波数fcは、(0.8+0.1)Hzすなわち、0.9Hzとする。
続いて、単位関連値演算部631は、対象周波数fを(fc−0.1)、fc、(fc+0.1)とする単位関連値Q1(n,T(fc−0.1))、Q1(n,Tfc)、Q1(n,T(fc+0.1))を演算する(図11のS12)。続いて、総合関連値演算部632は、対象周波数fを(fc−0.1)、fc、(fc+0.1)とする3つの総合関連値R1(T(fc−0.1))、R1(Tfc)、R1(T(fc+0.1))を演算する(図13のS13)。
続いて、周波数特定部633は、3つの総合関連値R1(T(fc−0.1))、R1(Tfc)、R1(T(fc+0.1))を比較して、大小関係を取得する(図13のS14)。そして、総合関連値R1(Tfc)が最も大きい場合には(図13のS15:Y)、得られた総合関連値R1(T)の中で、最大値R1(Tfmax)の周波数fmaxを、検出信号に含まれる振動周波数、すなわち心拍の周波数とする(図13のS19)。
続いて、単位関連値演算部631は、前回のデータPを取得してから4ms経過して、次のデータPを取得したか否かを判定する(図13のS20)。次の取得データPが無い場合には(図13のS20:N)、次のデータPを取得する時間に達していない、すなわち4ms経過していないことになるため、単位関連値演算部631は、次のデータPを取得するまで待機する。一方、次の取得データPが有る場合には(図13のS20:Y)、単位関連値演算部631は、再び、S12に戻り処理を繰り返す。
図13のS14において、取得した3つの中で総合関連値R1(Tfc)が最大ではない場合には(図13のS14:N)、3つの総合関連値R1の比較結果に基づいて振動周波数の傾向を取得し、傾向に基づいて連続する3つの対象周波数を順次変化させる。
すなわち、図13のS15において、取得した3つの中で総合関連値R1(Tfc)が最大ではない場合には(図13のS15:N)、周波数特定部633は、高周波側の総合関連値R1(T(fc+0.1))が低周波側の総合関連値R1(T(fc−0.1))より大きいか否かを判定する(図13のS16)。高周波側の総合関連値R1(T(fc+0.1))が低周波側の総合関連値R1(T(fc−0.1))より大きい場合には(図13のS16:Y)、周波数特定部633は、周波数(fc+0.1)を基準となる対象周波数fcに変更する(図13のS17)。つまり、周波数特定部633は、3つの総合関連値R1(n,T)に基づいて心拍の周波数がもっと高い周波数であるとの傾向があると判断し、基準となる対象周波数fcを前回から0.1Hz大きな値に変更する。
基準となる対象周波数fcが変更された後には、単位関連値演算部631がS20から処理を継続する。つまり、変更された基準となる対象周波数fcを含む3つの対象周波数(fc−0.1)、fc、(fc+0.1)について、単位関連値Q1(n,T)及び総合関連値R1(T)が演算される(S12,S13)。
一方、図13のS16において、低周波側の総合関連値R1(T(fc−0.1))が高周波側の総合関連値R1(T(fc+0.1))より大きい場合には(図13のS16:N)、周波数特定部633は、周波数(fc−0.1)を基準となる対象周波数fcに変更する(図11のS18)。つまり、周波数特定部633は、3つの総合関連値R1に基づいて心拍の周波数がもっと低い周波数であるとの傾向があると判断し、基準となる対象周波数fcを前回から0.1Hz小さな値に変更する。
基準となる対象周波数fcが変更された後には、単位関連値演算部631がS20から処理を継続する。つまり、変更された基準となる対象周波数fcを含む3つの対象周波数(fc−0.1)、fc、(fc+0.1)について、単位関連値Q1(n,T)及び総合関連値R1(T)が演算される(S12,S13)。
第二周波数計測処理について、図14を参照して、具体例をあげて説明する。被検者の心拍の周波数は、1.2〜1.4Hz程度とする。ここでは、対象周波数fのうち最も低い周波数0.8Hzが開始周波数に予め設定される(図13のS11,S12)。そして、最も低い周波数0.8Hzを含む連続する3つ(一部の所定数)の対象周波数fについて、すなわち0.8Hz、0.9Hz、1.0Hzの周波数について、単位関連値Q1(n,T0.8)、Q1(n,T0.9)、Q1(n,T1.0)及び総合関連値R1(T0.8)、R1(T0.9)、R1(T1.0)が演算される(図13のS12,S13)。
ここで、被検者の心拍の周波数は、0.8〜1.0Hzに対して大きいため、これらの中で最大の対象周波数1.0Hzの総合関連値R1(T1.0)が最大となる。つまり、周波数特定部633は、心拍の周波数がもっと高い周波数帯に存在するとの傾向を取得する。従って、図13のS15:N→S16:Yとなり、基準となる対象周波数fcは、0.1加算されて1.0Hzとなる(図13のS17)。
処理を繰り返して、処理開始から12ms経過後に、連続する3つの対象周波数fについて、すなわち1.1Hz、1.2Hz、1.3Hzの周波数について、単位関連値Q1(n,T1.1)、Q1(n,T1.2)、Q1(n,T1.3)及び総合関連値R1(T1.1)、R1(T1.2)、R1(T1.3)が演算される(図13のS12,S13)。これらの中で中央の対象周波数1.2Hzの総合関連値R1(T1.2)が最大となっている。
そうすると、図13のS15:Yとなり、最大値となる総合関連値R1(T1.2)の周波数1.2Hzが、心拍の周波数に特定される。心拍の周波数が1.2Hzであるので、1分間当たりの心拍数は、72回となる。なお、図14において、横軸の時間の数字に枠で囲んだ標記としているものは、心拍の周波数が特定されたタイミングを意味する。
心拍の周波数が一旦特定された後にも、心拍数計測部63は、心拍の周波数の計測を継続して、心拍の周波数の変化を計測し続ける。すなわち、データPを取得するたびに、上記処理を繰り返し、心拍の周波数が変化した場合に、現在の心拍の周波数を特定することができる。
なお、上記処理では、一旦、心拍の周波数が特定された後において、3つの中で中央の対象周波数fcの総合関連値R1(Tfc)が最大となる場合に、心拍の周波数が特定されることとした。これに限られず、一旦、心拍の周波数が特定された後においては、3つの中で中央の対象周波数fcの総合関連値R1(Tfc)が最大とならない場合であっても、3つの中で最大となる総合関連値R1に対応する対象周波数fが心拍の周波数であると特定してもよい。
(1−2−8−3.第三周波数計測処理)
第三周波数計測処理について、図15及び図16を参照して説明する。ただし、図15及び図16は、心拍数計測部63による第三周波数計測処理を例にあげる。なお、呼吸数計測部62による第三周波数計測処理も実質的に同様である。
第三周波数計測処理は、まず、全ての対象周波数fについての単位関連値Q1(n,T)及び総合関連値R1(T)を演算し、演算した総合関連値R1(T)を比較することにより、心拍の周波数を特定する。さらに、第三周波数計測処理は、心拍の周波数が一旦特定された後において、特定された周波数を含む一部の所定数の周波数についての総合関連値R1を演算し、演算した総合関連値R1(T)を比較することにより、変化後の心拍の周波数を特定するようにしている。
まずは、図15を参照して第三周波数計測処理を説明する。図15に示すように、単位関連値演算部631が、0.8〜3.0Hz全ての対象周波数fについて、単位関連値Q1(n,T)を演算する(図15のS21)。続いて、総合関連値演算部632が、0.8〜3.0Hz全ての対象周波数fについて、総合関連値R1(T)を演算する(図15のS22)。
続いて、周波数特定部633は、取得した全ての総合関連値R1(T)を比較して、大小関係を取得する(図15のS23)。そして、周波数特定部633は、得られた総合関連値R1(T)の中で最大値R1(Tfmax)の周波数fmaxを心拍の周波数とする(図15のS24)。このようにして、心拍の周波数は一旦特定される。続いて、周波数特定部633は、心拍の周波数とされた周波数fmaxを、基準となる対象周波数fcに設定する(図15のS25)。
続いて、単位関連値演算部631は、前回のデータPを取得してから4ms経過して、次のデータPを取得したか否かを判定する(図15のS26)。次の取得データPが無い場合には(図15のS26:N)、次のデータPを取得する時間に達していない、すなわち4ms経過していないことになるため、単位関連値演算部631は、次のデータPを取得するまで待機する。
一方、次の取得データPが有る場合には(図15のS26:Y)、単位関連値演算部631は、対象周波数fを(fc−0.1)、fc、(fc+0.1)とする単位関連値Q1(n,T(fc−0.1))、Q1(n,Tfc)、Q1(n,T(fc+0.1))を演算する(図15のS27)。つまり、第三周波数計測処理において、最初は、全ての対象周波数fについて単位関連値Q1(n,T)を演算したが、一旦心拍の周波数が特定された後は、一部の対象周波数fについて単位関連値Q1(n,T)を演算する。
続いて、総合関連値演算部632は、対象周波数fを(fc−0.1)、fc、(fc+0.1)とする3つの総合関連値R1(T(fc−0.1))、R1(Tfc)、R1(T(fc+0.1))を演算する(図15のS28)。続いて、S23に戻り処理が繰り返される。つまり、周波数特定部633は、取得した3つの総合関連値R1(T(fc−0.1))、R1(Tfc)、R1(T(fc+0.1))を比較して、大小関係を取得する(図15のS23)。そして、周波数特定部633は、3つの総合関連値R1(T(fc−0.1))、R1(Tfc)、R1(T(fc+0.1))の中で最大値R1(Tfmax)の周波数fmaxを心拍の周波数とする(図15のS24)。
第三周波数計測処理について、図16を参照して、具体例をあげて説明する。0.8〜3.0Hzの全ての周波数について、単位関連値Q1(n,T)及び総合関連値R1(T)が演算される(図15のS21,S22)。これらの中で、1.2Hzの周波数についての総合関連値R1(T1.2)が最大となっている。従って、最大値となる総合関連値R1(T1.2)の周波数1.2Hzが、心拍の周波数に特定される(図15のS24)。
このとき、特定された心拍の周波数fmaxが基準となる対象周波数fcに設定される。従って、次は、1.1Hz、1.2Hz、1.3Hzの周波数について、単位関連値Q1(n,T1.1)、Q1(n,T1.2)、Q1(n,T1.3)及び総合関連値R1(T1.1)、R1(T1.2)、R1(T1.3)が演算される(図15のS27,S28)。このとき、総合関連値R1(T1.3)が最大値となると、この時点では、心拍の周波数は、1.3Hzとされる(図15のS24)。上記処理が順次繰り返されて、現在の心拍の周波数、及び、これまでの心拍の周波数の変化が得られる。
(1−2−8−4.第四周波数計測処理)
第四周波数計測処理について、図17及び図18を参照して説明する。ただし、図17及び図18は、心拍数計測部63による第四周波数計測処理を例にあげる。なお、呼吸数計測部62による第四周波数計測処理も実質的に同様である。
第四周波数計測処理は、まず、全ての対象周波数fについての単位関連値Q1(n,T)及び総合関連値R1(T)を演算し、演算した総合関連値R1(T)を比較することにより、心拍の周波数を特定する。さらに、第三周波数計測処理は、心拍の周波数が一旦特定された後において、特定された周波数を含む一部の所定数の周波数についての総合関連値R1に基づいて、心拍の周波数が変化する傾向を決定する。そして、第四周波数計測処理は、当該傾向に基づいて連続する所定数の対象周波数fについての最新の総合関連値R1(T)を演算し、演算した最新の総合関連値R1(T)を比較することで、変化後の心拍の周波数を特定するようにしている。
まずは、図17を参照して第四周波数計測処理を説明する。図17に示すように、単位関連値演算部631が、0.8〜3.0Hz全ての対象周波数fについて、単位関連値Q1(n,T)を演算する(図17のS31)。続いて、総合関連値演算部632が、0.8〜3.0Hz全ての対象周波数fについて、総合関連値R1(T)を演算する(図17のS32)。
続いて、周波数特定部633は、取得した全ての総合関連値R1(T)を比較して、大小関係を取得する(図17のS33)。そして、周波数特定部633は、得られた総合関連値R1(T)の中で最大値R1(Tfmax)の周波数fmaxを心拍の周波数とする(図17のS34)。このようにして、心拍の周波数は一旦特定される。
続いて、周波数特定部633は、心拍の周波数とされた周波数fmaxを、一時的に、基準となる対象周波数fcに設定する。そして、基準となる対象周波数fc及びその前後の周波数fc−0.1、fc+0.1についての3つの総合関連値R1(Tfc)を用いて、振動周波数の傾向を取得し、傾向に基づいて演算するための連続する3つの対象周波数を順次変化させる。
すなわち、図17のS36において、周波数特定部633は、基準となる対象周波数fcより1つ高周波側の総合関連値R1(T(fc+0.1))が1つ低周波側の総合関連値R1(T(fc−0.1))より大きいか否かを判定する(図17のS36)。高周波側の総合関連値R1(T(fc+0.1))が低周波側の総合関連値R1(T(fc−0.1))より大きい場合には(図17のS36:Y)、周波数特定部633は、周波数(fc+0.1)を基準となる対象周波数fcに変更する(図17のS37)。つまり、周波数特定部633は、3つの総合関連値R1(n,T)に基づいて、心拍の周波数が高周波側に変化する傾向(可能性)があると判断し、基準となる対象周波数fcを前回から0.1Hz大きな値に変更する。
一方、低周波側の総合関連値R1(T(fc−0.1))が高周波側の総合関連値R1(T(fc+0.1))より大きい場合には(図17のS36:N)、周波数特定部633は、周波数(fc−0.1)を基準となる対象周波数fcに変更する(図17のS38)。つまり、周波数特定部633は、3つの総合関連値R1(n,T)に基づいて、心拍の周波数が低周波側に変化する傾向(可能性)があると判断し、基準となる対象周波数fcを前回から0.1Hz小さな値に変更する。
続いて、単位関連値演算部631は、前回のデータPを取得してから4ms経過して、次のデータPを取得したか否かを判定する(図17のS39)。次の取得データPが無い場合には(図17のS39:N)、次のデータPを取得する時間に達していない、すなわち4ms経過していないことになるため、単位関連値演算部631は、次のデータPを取得するまで待機する。
一方、次の取得データPが有る場合には(図17のS39:Y)、単位関連値演算部631は、対象周波数fを(fc−0.1)、fc、(fc+0.1)とする単位関連値Q1(n,T(fc−0.1))、Q1(n,Tfc)、Q1(n,T(fc+0.1))を演算する(図17のS40)。つまり、第四周波数計測処理において、最初は、全ての対象周波数fについて単位関連値Q1(n,T)を演算したが、一旦心拍の周波数が特定された後は、一部の対象周波数fについて単位関連値Q1(n,T)を演算する。
このとき、基準となる対象周波数fcは、心拍の周波数として特定された周波数ではなく、その前後の何れか一方である。つまり、心拍の変化に確実に追従するために、心拍の周波数として特定された周波数を、次に演算する3つのうちの中央の周波数にせずに、その前後の何れか一方を、次に演算する3つのうちの中央の周波数としている。
続いて、総合関連値演算部632は、対象周波数fを(fc−0.1)、fc、(fc+0.1)とする3つの総合関連値R1(T(fc−0.1))、R1(Tfc)、R1(T(fc+0.1))を演算する(図17のS41)。
続いて、周波数特定部633は、取得した3つの総合関連値R1(T(fc−0.1))、R1(Tfc)、R1(T(fc+0.1))を比較して、大小関係を取得する(図17のS42)。そして、周波数特定部633は、3つの総合関連値R1(T(fc−0.1))、R1(Tfc)、R1(T(fc+0.1))の中で最大値R1(Tfmax)の周波数fmaxを心拍の周波数とする(図17のS43)。続いて、S36に戻り処理が繰り返される。
第四周波数計測処理について、図18を参照して、具体例をあげて説明する。0.8〜3.0Hzの全ての周波数について、単位関連値Q1(n,T)及び総合関連値R1(T)が演算される(図17のS31,S32)。これらの中で、1.2Hzの周波数についての総合関連値R1(T1.2)が最大となっている。従って、最大値となる総合関連値R1(T1.2)の周波数1.2Hzが、心拍の周波数に特定される(図17のS34)。
このとき、特定された心拍の周波数fmaxが基準となる対象周波数fcに設定される。そして、基準となる対象周波数fc(1.2Hz)の前後の周波数についての総合関連値R1(T1.1)、R1(T1.3)を比較する(図17のS36)。1.3Hzの総合関連値R1(T1.3)の方が大きいため、周波数特定部633は、心拍の周波数が1.2Hzよりも高周波側に変化する傾向があると判断する。
そこで、処理開始から4ms経過後に、1.3Hzを中央とする連続する3つの対象周波数1.2Hz、1.3Hz、1.4Hzの周波数について、単位関連値Q1(n,T1.2)、Q1(n,T1.3)、Q1(n,T1.4)及び総合関連値R1(T1.2)、R1(T1.3)、R1(T1.4)が演算される(図17のS40,S41)。これらの中で1.2Hzの総合関連値R1(T1.2)が最大となっている。そうすると、現時点においても、心拍の周波数は1.2Hzとされる。
このとき、3つの周波数のうちで、高周波側の周波数1.4Hzの総合関連値R1(T1.4)と低周波側の周波数1.2Hzの総合関連値R1(T1.2)とを比較すると、1.2Hzの総合関連値R1(T1.2)の方が大きい(図17のS36)。
そこで、処理開始から8ms経過後には、1.2Hzを中央とする連続する3つの対象周波数1.1Hz、1.2Hz、1.3Hzの周波数について、単位関連値Q1(n,T1.1)、Q1(n,T1.2)、Q1(n,T1.3)及び総合関連値R1(T1.1)、R1(T1.2)、R1(T1.3)が演算される(図17のS40,S41)。上記処理が順次繰り返されて、現在の心拍の周波数、及び、これまでの心拍の周波数の変化が得られる。
<2.第二実施形態>
第二実施形態の生体情報計測システム1について説明する。第二実施形態の生体情報計測システム1においては、第一実施形態の生体情報計測システム1に対して、単位関連値Q2(n,T)及び総合関連値R2(T)が相違する。以下に相違点について説明する。
(2−1.単位関連値及び総合関連値の基本的考え方)
単位関連値Q2(a)及び総合関連値R2(a)の基本説明、及び、これらと振動周波数との関係について、図8を参照して説明する。
単位関連値Q2(a)は、基準時刻θ0の第一データ及び基準時刻θ0から半周期(π)前の第二データに基づいて演算される。第二実施形態では、単位関連値Q2(a)は、式(11)で表される。つまり、単位関連値Q2(a)は、基準時刻θ0における振幅sin(θ0)の絶対値と、基準時刻θ0から半周期前の時刻(θ0−π)における振幅sin(θ0−π)の絶対値との加算値である。正弦波において、半周期ずれた2つの値の中で、位相がπ/4の絶対値と位相が3π/4の絶対値とが最大となる。一方、位相が0の絶対値と位相がπ/2の絶対値とは、0となる。
Figure 2018023618
そして、総合関連値R2(a)は、基準時刻θ0を所定時間分とする複数の単位関連値Q2(a)に基づいて演算される。第二実施形態では、総合関連値R2(a)は、式(12)で表される。つまり、総合関連値R2(a)は、基準時刻θ0を当該周波数の1周期分とする複数の単位関連値Q2(a)を加算した値である。ここでは、総合関連値R2(a)は、基準時刻θ0を当該周波数の1周期分としたが、半周期分としてもよいし、1.5周期分、2周期分などとしてもよい。
Figure 2018023618
ここで、単位関連値Q2(a)が基準時刻θ0の振幅(第一データ)と基準時刻θ0から半周期(π)前の振幅(第二データ)とに基づいて演算した。仮に、単位関連値Q2として、基準時刻θ0の振幅と、基準時刻θ0から半周期(π)とは異なる時刻前の振幅とに基づいて演算した場合に、総合関連値R2の大きさについて比較検討する。
第一比較例としての単位関連値Q2(b)は、式(13)で表される。つまり、第一比較例としての単位関連値Q2(b)は、基準時刻θ0の振幅sin(θ0)の絶対値と、基準時刻θ0から(π−α)だけ前の振幅sin(θ0−π+α)の絶対値との加算値である。
Figure 2018023618
そして、第一比較例としての総合関連値R2(b)は、式(14)で表される。つまり、第一比較例としての総合関連値R2(b)は、基準時刻θ0を当該周波数の1周期分とする複数の単位関連値Q2(b)を加算した値である。
Figure 2018023618
また、第二比較例としての単位関連値Q2(c)は、式(15)で表される。つまり、第二比較例としての単位関連値Q2(c)は、基準時刻θ0の振幅sin(θ0)の絶対値と、基準時刻θ0から(π+β)だけ前の振幅sin(θ0−π−β)の絶対値との加算値である。
Figure 2018023618
そして、第二比較例としての総合関連値R2(c)は、式(16)で表される。つまり、第二比較例としての総合関連値R2(c)は、基準時刻θ0を当該周波数の1周期分とする複数の単位関連値Q2(c)を加算した値である。
Figure 2018023618
第二実施形態における総合関連値R2(a)と、比較例としての総合関連値R2(b),R1(c)とは、式(17)の関係を有する。つまり、単位関連値Q2(a)が、基準時刻θ0における振幅sin(θ0)の絶対値と、基準時刻θ0から半周期前の時刻(θ0−π)における振幅sin(θ0−π)の絶対値との加算値とすることで、総合関連値R2(a)が最大値となる。
Figure 2018023618
(2−2.単位関連値演算部621,631)
次に、単位関連値演算部621,631について説明する。呼吸数計測部62の単位関連値演算部621は、三次信号Sig3_1を用いて単位関連値Q2(n,T)を演算する。心拍数計測部63の単位関連値演算部631は、三次信号Sig3_2を用いて単位関連値Q2(n,T)を演算する。
単位関連値Q2(n,T)は、式(18)で表される。単位関連値Q2(n,T)は、基準時刻T0の第一データP(T0)、基準時刻T0から半周期前の第二データP(T0−0.5T)、基準時刻T0から1周期前の第三データP(T0−T)、基準時刻T0から1.5周期前の第四データP(T0−1.5T)、基準時刻T0から2周期前の第五データP(T0−2T)、基準時刻T0から2.5周期前の第六データP(T0−2.5T)に基づいて演算される。
換言すると、単位関連値Q2(n,T)は、第一データP(T0)の絶対値、第二データP(T0−0.5T)の絶対値、第三データP(T0−T)の絶対値、第四データP(T0−1.5T)の絶対値、第五データP(T0−2T)の絶対値、及び、第六データP(T0−2.5T)の絶対値を、加算した値をとる。つまり、単位関連値Q2(n,T)は、データの半周期毎のデータの絶対値を加算して得られる。
Figure 2018023618
そして、第二実施形態においては、基準時刻T0が、三次信号Sig3_1,Sig3_2が演算される周期で変化する。第二実施形態においては、カウント値nは、4ms毎に1増加する。つまり、基準時刻T0は、4ms経過の都度、次の演算時刻となる。
式(18)について、詳細には、単位関連値Q2(n,T)は、第一データ、第二データ、第三データ、第四データ、第五データ及び第六データのそれぞれの絶対値を加算した値である。なお、単位関連値Q2(n,T)は、第一データと第二データのそれぞれの絶対値を加算した値としてもよいし、第一データ、第二データ、第三データ及び第四データのそれぞれの絶対値を加算した値としてもよい。
(2−3.総合関連値演算部622,632)
次に、総合関連値演算部622,632について説明する。総合関連値演算部622,632は、単位関連値演算部621,631により演算された単位関連値Q2(n,T)を用いて、複数の対象周波数fのそれぞれについて、総合関連値R2(T)を演算する。総合関連値R2(T)は、式(19)で表される。
Figure 2018023618
総合関連値R2(T)は、基準時刻T0を所定時間分(第二実施形態では1周期分)とする複数の単位関連値Q2(n,T)に基づいて演算される。総合関連値R2(T)は、第一実施形態の総合関連値R1(T)と基本的には同様である。
図9に示すように、対象周波数fの周期Tが正弦波の周期Tf0に一致する場合における総合関連値R2(Tf0)と、図10に示すように、対象周波数fの周期Tが正弦波の周期Tf0に一致しない場合における総合関連値R2(Tf1)とを比較する。この場合、上述した式(17)から分かるように、一致する場合の総合関連値R2(Tf0)と、一致しない場合の総合関連値R2(Tf1)とは、式(20)の関係を有する。
Figure 2018023618
従って、対象周波数fの周期Tが正弦波の周期Tf0に一致する場合における総合関連値R2(Tf0)は、一致しない場合における全ての総合関連値R2(Tf1)よりも必ず大きな値となる。式(20)の関係から、第一実施形態にて説明したように、第一周波数計測処理から第四周波数計測処理を同様に実行することにより、1分間当たりの呼吸数及び心拍数を計測することができる。
(3.その他)
第一実施形態及び第二実施形態において、単位関連値演算部621,631は、三次信号Sig3_1、Sig3_2そのものである各データPそのものを用いて、単位関連値Q1(n,T)、Q2(n,T)を演算した。この他に、単位関連値演算部621,631は、各データPの相関値を用いて、単位関連値Q1(n,T)、Q2(n,T)を演算するようにしてもよい。
例えば、第一データP(T0)の相関値とは、第一データP(T0)の累乗値、第一データP(T0)の所定値を積算した値などである。第二データP(T0−0.5T)の相関値とは、第二データP(T0−0.5T)の累乗値、第二データP(T0−0.5T)の所定値を積算した値などである。他のデータについても同様である。そして、単位関連値演算部621,631は、式(8)(18)において、各データPを各データの相関値に置換して、単位関連値Q1(n,T)、Q2(n,T)を演算するようにしてもよい。
(4.実施形態の効果)
第一実施形態及び第二実施形態の振動周波数計測装置60は、振動センサ(10a)の検出信号に含まれる振動の周波数を計測する。計測装置60は、検出信号から所定周波数範囲を抽出するフィルタ処理部61と、フィルタ処理部61で処理された処理信号に対して、複数の対象周波数fについて、基準時刻T0の第一データP(T0)及び基準時刻T0から半周期前の第二データP(T0−0.5T)に基づいて単位関連値Q1(n,T),Q2(n,T)を演算する単位関連値演算部621,631と、複数の対象周波数fについて、基準時刻T0を所定時間分とする複数の単位関連値Q1(n,T),Q2(n,T)に基づいて、総合関連値R1(T),R2(T)を演算する総合関連値演算部622,632と、複数の対象周波数fの総合関連値R1(T),R2(T)に基づいて、検出信号の振動周波数を特定する周波数特定部623,633とを備える。
上記計測装置60によれば、単位関連値Q1(n,T),Q2(n,T)は、基準時刻T0の第一データP(T0)と基準時刻T0から半周期前の第二データP(T0−0.5T)とに基づいて演算される。ここで、第一データP(T0)と第二データP(T0−0.5T)との関係について説明する。
図8に示すように、特定周波数の正弦波において、位相(π/4)における振幅が最大値を示し、位相(3π/4)における振幅が最小値を示す。従って、位相(π/4)における振幅と位相(3π/4)における振幅との差が最大となる。また、位相(π/4)における振幅の絶対値及び位相(3π/4)における振幅の絶対値が、最大値を示す。位相(π/4)と位相(3π/4)とは、半周期(π/2)ずれている。
従って、正弦波の特定周波数に一致する周波数f1においては、半周期(π/2)ずれた第一データP(T0)及び第二データP(T0−0.5T)を用いることにより、単位関連値Q1(n,T),Q2(n,T)の一つは、大きな値となる。一方、正弦波の特定周波数に一致しない周波数f2においては、半周期(π/2)ずれた第一データP(T0)及び第二データP(T0−0.5T)を用いた場合に、単位関連値Q1(n,T),Q2(n,T)の一つは、上記の場合に比べて小さな値となる。
ただし、半周期(π/2)ずれた位相であっても、位相0°、π/2における振幅はゼロとなるため、両者の差及び両者の絶対値は、最小となる。そこで、基準時刻T0を所定時間分(例えば、1周期分)とする複数の単位関連値Q1(n,T),Q2(n,T)に基づいた総合関連値R1(T),R2(T)が、演算されている。総合関連値R1(T),R2(T)は、単位関連値Q1(n,T),Q2(n,T)と同様に、周波数f1,f2に対して同様の関係を有する。つまり、正弦波の特定周波数に一致する周波数f1においては、総合関連値R1(T),R2(T)は、大きな値となる。一方、正弦波の特定周波数に一致しない周波数f2においては、総合関連値R1(T),R2(T)は、小さな値となる。従って、周波数特定部623,633は、総合関連値R1(T),R2(T)に基づいて振動周波数を特定することができる。
つまり、相互に半周期ずれた第一データP(T0)と第二データP(T0−0.5T)との関係を利用することにより、振動周波数を特定している。従って、単位関連値Q1(n,T),Q2(n,T)の演算及び総合関連値R1(T),R2(T)の演算は、フーリエ変換処理や自己相関関数の演算処理に比べると、非常に簡易な処理となる。そのため、高速な演算処理装置を用いずに、振動周波数を計測することができる。
さらに、計測装置60は、上述したように、相互に半周期ずれた第一データP(T0)と第二データP(T0−0.5T)との関係を利用することにより、高精度に振動周波数を特定することができる。ただし、特定周波数に一致しない周波数f1における第一データP(T0)と第二データP(T0−0.5T)とが上記関係を有するためには、ある程度ノイズが除去された信号を用いる必要がある。そこで、単位関連値演算部621,631は、振動センサ(10a)の検出信号そのものを用いるのではなく、フィルタ処理部61で処理された処理信号を用いることとしている。このように、計測装置60は、検出信号をフィルタ処理することにより、確実に振動周波数を特定することができる。
また、第一実施形態においては、単位関連値演算部621,631は、例えば、式(8)に示すように、第一データP(T0)と第二データP(T0−0.5T)との差分値に基づいて、単位関連値Q1(n,T)を演算する。又は、単位関連値演算部621,631は、第一データP(T0)の相関値と第二データP(T0−0.5T)の相関値との差分値に基づいて、前記単位関連値Q1(n,T)を演算してもよい。つまり、差分値から得られる単位関連値Q1(n,T)を用いることで、総合関連値R1(T)は、式(7)の関係を有することが明らかである。これにより、計測装置60は、確実に振動周波数を特定することができる。
また、第一実施形態において、単位関連値演算部621,631は、式(8)に示すように、第一データP(T0)に対して、第二データP(T0−0.5T)を減算し、基準時刻T0から1周期前の第三データP(T0−T)を加算し、基準時刻T0から1.5周期前の第四データP(T0−1.5T)を減算することに基づいて、単位関連値Q1(n,T)を演算する。これにより、単位関連値Q1(n,T)が、複数周期分の情報を考慮した値となり、より高精度に振動周波数を特定することができる。
また、単位関連値演算部621,631は、第一データP(T0)の相関値に対して、第二データP(T0−0.5T)の相関値を減算し、基準時刻T0から1周期前の第三データP(T0−T)の相関値を加算し、基準時刻T0から1.5周期前の第四データP(T0−1.5T)の相関値を減算することに基づいて、単位関連値Q1(n,T)を演算するようにしてもよい。この場合も同様の効果を奏する。
また、第二実施形態において、単位関連値演算部621,631は、式(18)に示すように、第一データP(T0)の絶対値と第二データP(T0−0.5T)の絶対値との加算値に基づいて、単位関連値Q1(n,T)を演算する。又は、単位関連値演算部621,631は、第一データP(T0)の絶対値の相関値と第二データP(T0−0.5T)の絶対値の相関値との加算値に基づいて、単位関連値Q1(n,T)を演算してもよい。つまり、差分値から得られる単位関連値Q2(n,T)を用いることで、総合関連値R2(T)は、式(17)の関係を有することが明らかである。これにより、計測装置60は、確実に振動周波数を特定することができる。
また、第一周波数計測処理(図11,12)又は第二周波数計測処理(図13,14)においては、周波数特定部623,633は、一部の所定数の対象周波数fについての総合関連値R1(T),R2(T)の比較を行った後に、異なる一部の所定数の対象周波数fについての総合関連値R1(T),R2(T)の比較を1回以上行うことにより、振動周波数を特定する。これにより、1回当たりの演算量が大幅に少なくできる。
また、第一周波数計測処理(図11,12)又は第二周波数計測処理(図13,14)においては、周波数特定部623,633は、連続する所定数(例えば3つ)の対象周波数fについての総合関連値R1(T),R2(T)を比較することにより振動周波数の傾向を取得し、傾向に基づいて連続する所定数の対象周波数を順次変化させ、得られた総合関連値R1(T),R2(T)の比較により振動周波数を特定する。これにより、早期に振動周波数に到達する。
また、第一周波数計測処理(図11,12)又は第二周波数計測処理(図13,14)においては、周波数特定部623,633は、予め設定された開始周波数を含む一部の所定数の対象周波数fについての総合関連値R1(T),R2(T)の比較を最初に行う。特に、開始周波数は、対象周波数fのうち最も高い周波数又は最も低い周波数とする。これにより、振動周波数に到達するまでの処理が容易となる。なお、開始周波数は、対象周波数fのうち最も高い周波数又は最も低い周波数に限定されるものではなく、対象周波数fの中央値や、心拍や呼吸の一般的な正常範囲の上限値、下限値又は中央値など、対象周波数fの範囲から任意の周波数を設定することもできる。
また、第三周波数計測処理(図15,16)又は第二周波数計測処理(図17,18)においては、周波数特定部623,633は、全ての対象周波数fについての総合関連値R1(T),R2(T)を比較することにより振動周波数を特定するようにしてもよい。これにより、早期に振動周波数を取得することができる。この場合であっても、フーリエ変換処理や自己相関関数の演算処理に比べると、非常に簡易な処理となる。
また、第一周波数計測処理(図11,12)、第二周波数計測処理(図13,14)又は第三周波数計測処理(図15,16)においては、単位関連値演算部621,631及び総合関連値演算部622,632は、周波数特定部623,633により振動周波数fmaxが一旦特定された後において、特定された振動周波数fmax及び振動周波数fmaxの前後の周波数について単位関連値Q1(n,T),Q2(n,T)及び総合関連値R1(T),R2(T)の演算を繰り返している。そして、周波数特定部623,633は、演算された最新の総合関連値R1(T),R2(T)に基づいて、最新の振動周波数fmaxを繰り返して特定している。
これにより、振動周波数fmaxが変化した場合であっても、最新の振動周波数fmaxを特定することができる。さらに、振動周波数fmaxの変化の推移を得ることもできる。これにより、振動周波数fmaxがどのように変化しているかを把握することができる。
また、第四周波数計測処理(図17,18)においては、周波数特定部623,633は、振動周波数が一旦特定された後において、特定された振動周波数を含む所定数の対象周波数についての総合関連値R1(T),R2(T)に基づいて変化する傾向を決定している。そして、単位関連値演算部621,631及び総合関連値演算部622,632は、決定された傾向に基づいて連続する所定数の対象周波数についての単位関連値Q1(n,T),Q2(n,T)及び総合関連値R1(T),R2(T)を演算し、周波数特定部623,633は、演算された最新の総合関連値R1(T),R2(T)に基づいて、最新の振動周波数を繰り返して特定する。
これにより、演算する所定数の対象周波数を、振動周波数の変化に追従することができる。従って、振動周波数の変化が早い場合であっても、変化後の振動周波数を早期に特定することができる。
また、第一周波数計測処理から第四周波数計測処理において、総合関連値演算部622,632は、基準時刻T0を所定時間分とする複数の単位関連値Q1(n,T),Q2(n,T)を取得する際に、対象周波数fに関わらず演算周波数を同一の周波数(250Hz、1サイクルの4ms)として、複数の単位関連値Q1(n,T),Q2(n,T)の演算を行っている。これにより、演算が簡易となる。
一方、総合関連値演算部622,632は、基準時刻T0を所定時間分とする複数の単位関連値Q1(n,T),Q2(n,T)を取得する際に、対象周波数fに応じて演算周波数を異なる周波数として、複数の単位関連値Q1(n,T),Q2(n,T)の演算を行うようにしてもよい。この場合、対象周波数に応じたデータ数にすることができるため、単位関連値Q1(n,T),Q2(n,T)の対象周波数毎の精度のばらつきを抑制できる。
また、第一実施形態及び第二実施形態において、振動センサ(10a)は、身体の表面に配置され、面状に形成される圧力センサセル10aであり、計測装置60は、圧力センサセル10aによる検出信号に基づいて、圧力センサセル10aに面する身体の生体情報としての心拍又は呼吸を計測する。つまり、圧力センサセル10aを用い、且つ、上記処理を行うことにより、非常に微小な振動である心拍又は呼吸を確実に計測できる。
1:生体情報計測システム、 10:センサユニット、 10a:圧力センサセル(振動センサ)、 11:第一電極、 12:第二電極、 13:誘電層、 20:電源装置、 40:スイッチ回路、 50:切替制御装置、 60:振動周波数計測装置、 61:フィルタ処理部、 62:呼吸数計測部、 63:心拍数計測部、 621,631:単位関連値演算部、 622,632:総合関連値演算部、 623,633:周波数特定部、 Q1(n,T),Q2(n,T):単位関連値、 R1(T),R1(T):総合関連値、 T0:基準時刻

Claims (13)

  1. 振動センサの検出信号に含まれる振動の周波数を計測する振動周波数計測装置であって、
    前記検出信号から所定周波数範囲を抽出するフィルタ処理部と、
    前記フィルタ処理部で処理された処理信号に対して、複数の対象周波数について、基準時刻の第一データ及び前記基準時刻から半周期前の第二データに基づいて単位関連値を演算する単位関連値演算部と、
    前記複数の対象周波数について、前記基準時刻を所定時間分とする複数の前記単位関連値に基づいて、総合関連値を演算する総合関連値演算部と、
    前記複数の対象周波数の前記総合関連値に基づいて、前記検出信号の振動周波数を特定する周波数特定部と、
    を備える、振動周波数計測装置。
  2. 前記単位関連値演算部は、前記第一データと前記第二データとの差分値、又は、前記第一データの相関値と前記第二データの相関値との差分値に基づいて、前記単位関連値を演算する、請求項1に記載の振動周波数計測装置。
  3. 前記単位関連値演算部は、
    前記第一データに対して、前記第二データを減算し、前記基準時刻から1周期前の第三データを加算し、前記基準時刻から1.5周期前の第四データを減算することに基づいて、前記単位関連値を演算する、
    又は、前記第一データの相関値に対して、前記第二データの相関値を減算し、前記基準時刻から1周期前の第三データの相関値を加算し、前記基準時刻から1.5周期前の第四データの相関値を減算することに基づいて、前記単位関連値を演算する、請求項2に記載の振動周波数計測装置。
  4. 前記単位関連値演算部は、前記第一データの絶対値と前記第二データの絶対値との加算値、又は、前記第一データの絶対値の相関値と前記第二データの絶対値の相関値との加算値に基づいて、前記単位関連値を演算する、請求項1に記載の振動周波数計測装置。
  5. 前記周波数特定部は、一部の所定数の対象周波数についての前記総合関連値の比較を行った後に、異なる一部の所定数の対象周波数についての前記総合関連値の比較を1回以上行うことにより、前記振動周波数を特定する、請求項1−4の何れか一項に記載の振動周波数計測装置。
  6. 前記周波数特定部は、連続する所定数の対象周波数についての前記総合関連値を比較することにより前記振動周波数の傾向を取得し、前記傾向に基づいて前記連続する所定数の対象周波数を順次変化させ、得られた前記総合関連値の比較により前記振動周波数を特定する、請求項5に記載の振動周波数計測装置。
  7. 前記周波数特定部は、予め設定された開始周波数を含む一部の所定数の対象周波数についての前記総合関連値の比較を最初に行う、請求項5又は6に記載の振動周波数計測装置。
  8. 前記開始周波数は、前記対象周波数のうち最も高い周波数又は最も低い周波数である、請求項7に記載の振動周波数計測装置。
  9. 前記周波数特定部は、全ての対象周波数についての前記総合関連値を比較することにより前記振動周波数を特定する、請求項1−4の何れか一項に記載の振動周波数計測装置。
  10. 前記単位関連値演算部及び前記総合関連値演算部は、前記周波数特定部により前記振動周波数が一旦特定された後において、特定された前記振動周波数及び前記振動周波数の前後の周波数について前記単位関連値及び前記総合関連値の演算を繰り返し、
    前記周波数特定部は、演算された最新の前記総合関連値に基づいて、最新の前記振動周波数を繰り返して特定する、請求項1−9の何れか一項に記載の振動周波数計測装置。
  11. 前記周波数特定部は、前記振動周波数が一旦特定された後において、特定された前記振動周波数を含む所定数の対象周波数についての前記総合関連値に基づいて変化する傾向を決定し、
    前記単位関連値演算部及び前記総合関連値演算部は、決定された前記傾向に基づいて連続する所定数の対象周波数についての前記単位関連値及び前記総合関連値を演算し、
    前記周波数特定部は、演算された最新の前記総合関連値に基づいて、最新の前記振動周波数を繰り返して特定する、請求項1−9の何れか一項に記載の振動周波数計測装置。
  12. 前記総合関連値演算部は、前記基準時刻を所定時間分とする前記複数の単位関連値を取得する際に、対象周波数に応じて演算周波数を異なる周波数として、前記複数の単位関連値の演算を行う、請求項1−11の何れか一項に記載の振動周波数計測装置。
  13. 前記振動センサは、身体の表面に配置され、面状に形成される圧力センサセルであり、
    前記振動周波数計測装置は、前記圧力センサセルによる検出信号に基づいて、前記圧力センサセルに面する身体の生体情報としての心拍又は呼吸を計測する、請求項1−12の何れか一項に記載の振動周波数計測装置。
JP2016157671A 2016-08-10 2016-08-10 振動周波数計測装置 Active JP6783091B2 (ja)

Priority Applications (1)

Application Number Priority Date Filing Date Title
JP2016157671A JP6783091B2 (ja) 2016-08-10 2016-08-10 振動周波数計測装置

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
JP2016157671A JP6783091B2 (ja) 2016-08-10 2016-08-10 振動周波数計測装置

Publications (2)

Publication Number Publication Date
JP2018023618A true JP2018023618A (ja) 2018-02-15
JP6783091B2 JP6783091B2 (ja) 2020-11-11

Family

ID=61193749

Family Applications (1)

Application Number Title Priority Date Filing Date
JP2016157671A Active JP6783091B2 (ja) 2016-08-10 2016-08-10 振動周波数計測装置

Country Status (1)

Country Link
JP (1) JP6783091B2 (ja)

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20070293781A1 (en) * 2003-11-04 2007-12-20 Nathaniel Sims Respiration Motion Detection and Health State Assesment System
JP2011015887A (ja) * 2009-07-10 2011-01-27 Mitsubishi Electric Corp 生体状態取得装置、生体状態取得プログラム、生体状態取得装置を備えた機器及び空気調和機
JP2011193886A (ja) * 2010-03-17 2011-10-06 Seiko Epson Corp 生体情報測定装置、生体情報測定方法、および生体情報測定プログラム
JP2012105762A (ja) * 2010-11-16 2012-06-07 Seiko Epson Corp バイタルサイン計測装置及び体動検出装置
JP2014226272A (ja) * 2013-05-21 2014-12-08 トヨタ自動車株式会社 物体変位検知装置、及び物体変位検知方法

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20070293781A1 (en) * 2003-11-04 2007-12-20 Nathaniel Sims Respiration Motion Detection and Health State Assesment System
JP2011015887A (ja) * 2009-07-10 2011-01-27 Mitsubishi Electric Corp 生体状態取得装置、生体状態取得プログラム、生体状態取得装置を備えた機器及び空気調和機
JP2011193886A (ja) * 2010-03-17 2011-10-06 Seiko Epson Corp 生体情報測定装置、生体情報測定方法、および生体情報測定プログラム
JP2012105762A (ja) * 2010-11-16 2012-06-07 Seiko Epson Corp バイタルサイン計測装置及び体動検出装置
JP2014226272A (ja) * 2013-05-21 2014-12-08 トヨタ自動車株式会社 物体変位検知装置、及び物体変位検知方法

Also Published As

Publication number Publication date
JP6783091B2 (ja) 2020-11-11

Similar Documents

Publication Publication Date Title
JP7054710B2 (ja) 電子機器及びシステム
KR102291934B1 (ko) 심장 모니터링 시스템
US20160242672A1 (en) Vital signal measuring apparatus and method for estimating contact condition
FI126631B (en) Method and apparatus for monitoring stress
CN106073731A (zh) 用于确定哺乳动物心血管分量的方法和系统
JP6579890B2 (ja) 疲労度計
JP2012210236A (ja) 測定装置、測定方法、情報処理装置、情報処理方法、およびプログラム
JP7265022B2 (ja) 心臓血管信号の取得、融合、およびノイズの軽減
Azimi et al. Breathing signal combining for respiration rate estimation in smart beds
JPWO2018105447A1 (ja) 接触状態推定装置及び生体信号測定装置
JP2017537710A (ja) 睡眠段階分類のスペクトル境界を決定するシステム及び方法
EP3200691B1 (en) System and method for detecting slow waves in a subject during a sleep session
JP2017221647A (ja) 筋肉疲労出力装置、筋肉疲労出力方法及びプログラム
JP6608320B2 (ja) 生体情報計測装置
CN106999086B (zh) 活体内信号源位置检测方法及活体内信号源位置检测装置
JP6783091B2 (ja) 振動周波数計測装置
Massaroni et al. Influence of motion artifacts on a smart garment for monitoring respiratory rate
CN109745020A (zh) 基于血流磁电效应的非侵入式血管弹性评估方法
WO2017169993A1 (ja) 体圧分布及び生体情報の計測装置
Banerjee et al. HeartSense: smart phones to estimate blood pressure from photoplethysmography
JP7307185B2 (ja) 心臓血管信号の取得、融合、およびノイズの軽減のための方法
JP6501643B2 (ja) 信号処理方法、生体信号処理方法、信号処理装置及び生体信号処理装置
Gomez-Clapers et al. A novel method to obtain proximal plethysmographic information from distal measurements using the impedance plethysmogram
Rhee et al. Analysis of both hands' two pulse waveforms using a clip-type pulsimeter equipped with magnetic sensing hall device
EP3501392A1 (en) Apparatus for measurement of the galvanic skin response of a person working with a computer and the method of measuring and performing diagnostics with this device

Legal Events

Date Code Title Description
A621 Written request for application examination

Free format text: JAPANESE INTERMEDIATE CODE: A621

Effective date: 20190722

A977 Report on retrieval

Free format text: JAPANESE INTERMEDIATE CODE: A971007

Effective date: 20200710

A131 Notification of reasons for refusal

Free format text: JAPANESE INTERMEDIATE CODE: A131

Effective date: 20200714

A521 Request for written amendment filed

Free format text: JAPANESE INTERMEDIATE CODE: A523

Effective date: 20200910

TRDD Decision of grant or rejection written
A01 Written decision to grant a patent or to grant a registration (utility model)

Free format text: JAPANESE INTERMEDIATE CODE: A01

Effective date: 20201006

A61 First payment of annual fees (during grant procedure)

Free format text: JAPANESE INTERMEDIATE CODE: A61

Effective date: 20201021

R150 Certificate of patent or registration of utility model

Ref document number: 6783091

Country of ref document: JP

Free format text: JAPANESE INTERMEDIATE CODE: R150