JP2016145763A - 信号解析システム、方法およびプログラム - Google Patents

Abstract

【課題】物体の破壊や崩壊の予兆を捉えることができる信号解析システムを提供する。
【解決手段】本発明によれば、解析対象信号のピーク時刻を取得するピーク時刻取得部と、前記ピーク時刻に基づいてＤＦＡ対象データを生成するＤＦＡ対象データ生成部と、前記ＤＦＡ対象データに基づいてＤＦＡを実行してスケーリング指数を取得するＤＦＡ実行部と、取得した前記スケーリング指数を積算するスケーリング指数積算部とを含み、前記スケーリング指数積算部は、取得順にロードした前記スケーリング指数が所定の閾値を下回っていない場合には該スケーリング指数を最新の積算値に加算し、ロードした該スケーリング指数が該所定の閾値を下回っている場合には最新の積算値をリセットする処理を繰り返し実行する信号解析システムが提供される。
【選択図】図２

Description

本発明は、信号解析システムに関し、より詳細には、周期信号のゆらぎを解析するシステムに関する。

従来、転がり軸受の損傷状態をアコースティック・エミッション法（Acoustic Emission, AE）を用いて診断する方法が知られている。例えば、特開平８−１５９１５１号公報（特許文献１）は、軸受から発生する超音波をマイクロフォンで計測し、その音響信号の注目周波数の成分強度を予め用意したマスターカーブと比較することによって、軸受の損傷状態を予測し余寿命を推定する軸受診断方法を開示する。

特開平８−１５９１５１号公報

本発明は、上記従来技術に鑑みてなされたものであり、本発明は、物体の破壊や崩壊の予兆を捉えることができる新規な信号解析システムを提供することを目的とする。

本発明者は、物体の破壊や崩壊の予兆を捉えることができる新規な信号解析システムにつき検討する過程で、振動する物体が破壊や崩壊の途上にある場合、物体の内部に生じたひずみに起因して物体から発する振動波形のピーク間隔にゆらぎが生じることを発見した。この発見から、本発明者は、物体から発する振動波形のゆらぎをＤＦＡ（Detrended Fluctuation Analysis）で連続的に解析する仕組みを作れば、スケーリング指数の時間変動から物体の破壊や崩壊の予兆を捉えることができるのではないかという着想を得た。そして、この着想の具体化につき鋭意検討した結果、以下の構成に想到し、本発明に至ったのである。

すなわち、本発明によれば、解析対象信号のピーク時刻を取得するピーク時刻取得部と、前記ピーク時刻取得部が取得した前記ピーク時刻に基づいてＤＦＡ対象データを生成するＤＦＡ対象データ生成部と、前記ＤＦＡ対象データに基づいてＤＦＡを実行してスケーリング指数を取得するＤＦＡ実行部と、取得した前記スケーリング指数を積算するスケーリング指数積算部と、を含み、前記ＤＦＡ対象データ生成部は、前記ピーク時刻取得部が取得した前記ピーク時刻を時系列順に（Ａ＋１）個ずつロードし、ロードした該（Ａ＋１）個のピーク時刻に基づいてＡ個のピーク間隔を算出し、該Ａ個のピーク間隔を要素とする時系列データに基づいてＭ個の要素を含む擬似的な時系列データを前記ＤＦＡ対象データとして生成する処理を繰り返し実行し（Ａ、Ｍはともに自然数であり、ＭはＡより大きい。以下同様。）、前記ＤＦＡ実行部は、前記ＤＦＡ対象データ生成部が生成した前記ＤＦＡ対象データを生成順にロードし、ロードした該ＤＦＡ対象データに基づいて前記スケーリング指数を取得する処理を繰り返し実行し、前記スケーリング指数積算部は、前記ＤＦＡ実行部が取得した前記スケーリング指数を取得順にロードし、ロードした該スケーリング指数が所定の閾値を下回っていない場合には該スケーリング指数を最新の積算値に加算し、ロードした該スケーリング指数が該所定の閾値を下回っている場合には最新の積算値をリセットする処理を繰り返し実行する、信号解析システムが提供される。

上述したように、本発明によれば、物体の破壊や崩壊の予兆を捉えることができる新規な信号解析システムが提供される。

本実施形態の信号解析システムの概要構成図。 本実施形態の信号解析システム機能構成図。 第１の測定方法における解析対象信号の波形を示す図。 第２の測定方法における解析対象信号の波形を示す図。 第３の測定方法における解析対象信号の波形を示す図。 本実施形態のピーク時刻取得処理を表すフローチャート。 本実施形態のピーク時刻取得処理を説明するための概念図。 本実施形態のＤＦＡ対象データ生成処理を表すフローチャート。 本実施形態のＤＦＡ対象データ生成処理を説明するための概念図。 本実施形態のＤＦＡ対象データ生成処理を説明するための概念図。 本実施形態のＤＦＡ対象データ生成処理を説明するための概念図。 本実施形態のＤＦＡを表すフローチャート。 本実施形態のＤＦＡを表すフローチャート。 本実施形態のボックスサイズ・データを示す図。 実施形態のＤＦＡを説明するための概念図。 本実施形態のスケーリング指数の積算処理を表すフローチャート。 本実施形態のスケーリング指数の積算値の時系列グラフを示す図。

以下、本発明を図面に示した実施の形態をもって説明するが、本発明は、図面に示した実施の形態に限定されるものではない。なお、以下に参照する各図においては、共通する要素について同じ符号を用い、適宜、その説明を省略するものとする。

図１は、本発明の実施形態である信号解析システム１００の概要構成を示す。図１に示すように、本実施形態の信号解析システム１００は、振動センサ１０と、加振器１２と、コンピュータ２０とを含んで構成される。振動センサ１０は、測定対象となる物体が発する振動を電圧信号に変換する振動検出手段である。振動センサ１０は、ピエゾ素子などを用いた接触式のものであってもよく、レーザなどを用いた非接触式のものであってもよい。また、加振器１２は、測定対象となる物体に周期振動を入力するための手段である。加振器１２は、測定対象の規模に応じた適切な振動子を用いて構成することができ、接触式および非接触式のいずれであってもよい。

一方、コンピュータ２０は、振動センサ１０から入力されるセンサ信号のゆらぎを解析するための情報処理装置であり、好ましくは、加振器１２に対して加振信号（正弦波信号）を出力するための構成を併せて備える。なお、図１では、コンピュータ２０としてノート型ＰＣを例示的に示しているが、コンピュータ２０は、タブレット型やデスクトップ型のＰＣであってもよいし、静止物体の破壊や崩壊の予兆を捉える用途に特化したワークステーションやサーバであってもよいし、当該用途に特化した端末装置に搭載される組み込みコンピュータであってもよい。また、コンピュータ２０と振動センサ１０および加振器１２との間の信号の送受信は、有線通信および無線通信のいずれで行ってもよい。

以上、本実施形態の信号解析システム１００概要構成について説明してきたが、続いて、コンピュータ２０の機能構成を図２に示す機能ブロック図に基づいて説明する。

図２に示すように、コンピュータ２０は、ピーク時刻取得部２１と、測定モード設定部２２と、ＤＦＡ対象データ生成部２３と、記憶部２４と、ＤＦＡ実行部２５と、スケーリング指数積算部２６と、出力情報生成部２７とを含む。上述した各機能部は、コンピュータ２０に専用のプログラムをインストールすることによって実現される。

ピーク時刻取得部２１は、入力される周期信号のピーク時刻を取得するための機能部であり、測定モード設定部２２は、信号解析システム１００が有する測定モード（後述する）を選択可能にユーザに提示し、ユーザが選択した測定モードに応じて各種パラメータの設定等を行うための機能部である。

ＤＦＡ対象データ生成部２３は、ピーク時刻取得部２１が取得したピーク時刻に基づいてＤＦＡ（Detrended Fluctuation Analysis）の対象データを生成するための機能部であり、ＤＦＡ実行部２５は、ＤＦＡ対象データ生成部２３が生成した対象データに基づいてＤＦＡを実行するための機能部である。

スケーリング指数積算部２６は、ＤＦＡで取得したスケーリング指数を積算するための機能部であり、出力情報生成部２７は、スケーリング指数の積算値に基づいて出力情報を生成するための機能部である。

コンピュータ２０は、さらに、正弦波発振回路３４と、信号重畳回路３６とを含む。正弦波発振回路３４は、指定された周波数の正弦波信号を生成・発振するための回路であり、信号重畳回路３６は、入力される２つの信号を重畳して出力するための回路である。

ここで、信号解析システム１００を使用した測定方法を図１に基づいて説明する。本実施形態においては、測定対象に応じて、以下の３種類の測定方法のいずれかを適用することができる。

（第１の測定方法）
正常時に所定の周波数近傍で周期振動する運動物体を測定する場合、第１の測定方法を適用することができる。第１の測定方法においては、運動物体の運動に付随して発生する周期振動を振動センサ１０で検出し、振動センサ１０から送信されるセンサ信号をコンピュータ２０で解析する。第１の測定方法が適用される運動物体としては、モーターや内燃機関などの回転駆動装置や当該回転駆動装置に接続されて回転する転がり軸受などの回転体を例示することができる。

第１の測定方法において、コンピュータ２０は、振動センサ１０のセンサ信号がピーク時刻取得部２１に直接入力されるように、図１に示すスイッチＳＷをＯＮにする。

図３は、モーターの周期振動を検出する振動センサ１０からピーク時刻取得部２１に入力されるセンサ信号の振動波形を例示する。この場合、ピーク時刻取得部２１に入力される振動波形は、図３（ａ）、（ｂ）に示すように、モーターの状態変化に応じて逐次変化する。

（第２の測定方法）
それ自体は振動しないが外部から振動を入力することが可能な静止物体を測定する場合、第２の測定方法を適用することができる。第２の測定方法においては、加振器１２から静止物体に周期振動を入力した状態で、加振された静止物体の周期振動を振動センサ１０で検出し、振動センサ１０から送信されるセンサ信号をコンピュータ２０で解析する。第２の測定方法が適用される静止物体としては、躯体（建物の柱、梁、壁など）を例示することができる。

第２の測定方法において、コンピュータ２０は、加振信号生成部３２が生成する加振信号が加振器１２に送信されるように制御するとともに、振動センサ１０のセンサ信号がピーク時刻取得部２１に直接入力されるように、図１に示すスイッチＳＷをＯＮにする。

図４は、アルミ鋼材に周期振動を入力した場合に、振動センサ１０からピーク時刻取得部２１に入力されるセンサ信号の振動波形を例示する。この場合、ピーク時刻取得部２１に入力される振動波形は、図４（ａ）、（ｂ）に示すように、アルミ鋼材の状態変化に応じて逐次変化する。

（第３の測定方法）
微弱な振動を発生する静止物体であって、外部から振動を入力することができない静止物体を測定する場合、第３の測定方法を適用することができる。第３の測定方法においては、静止物体の微弱な振動を振動センサ１０で検出し、振動センサ１０から送信されるセンサ信号と所定の正弦波信号を重畳した重畳信号をコンピュータ２０で解析する。第３の測定方法が適用される静止物体としては、地盤を例示することができる。

第３の方法においては、コンピュータ２０は、正弦波発振回路３４が発振する所定周波数の正弦波信号が信号重畳回路３６に入力されるように制御するとともに、振動センサ１０のセンサ信号が信号重畳回路３６に入力されるように、図１に示すスイッチＳＷをＯＦＦにする。その結果、正弦波発振回路３４から入力される正弦波信号と振動センサ１０から入力されるセンサ信号が信号重畳回路３６において重畳され、その重畳信号がピーク時刻取得部２１に入力されることになる。

図５（ａ）は、地盤に設置された地震計（振動センサ１０に相当）から信号重畳回路３６に入力されるセンサ信号の振動波形を例示し、図５（ｂ）は、正弦波発振回路３４から信号重畳回路３６に入力される正弦波信号を例示する。この場合、信号重畳回路３６において、図５（ａ）、（ｂ）に示す２つの信号が重畳される結果、図５（ｃ）に示す重畳信号がピーク時刻取得部２１に入力される。

上述したように、本実施形態においては、適用する測定方法によって、振動センサ１０のセンサ信号または振動センサ１０のセンサ信号と正弦波信号を重畳した重畳信号のいずれかがピーク時刻取得部２１に入力されることになる。ピーク時刻取得部２１に入力される信号はいずれも周期信号であり、コンピュータ２０は、この周期信号のゆらぎを解析する。なお、以下においては、ピーク時刻取得部２１に入力されるセンサ信号および重畳信号を総称して“解析対象信号”という。

以上、コンピュータ２０の機能構成について概説してきたが、続いて、上述した各機能部が実行する処理の具体的な内容を説明する。なお、以下の説明においては、適宜、図２を参照するものとする。

まず、ピーク時刻取得部２１が実行する処理を図６に示すフローチャートおよび図７に示す解析対象信号の模式的な波形図に基づいて説明する。

ステップ１０１では、入力される解析対象信号の電圧監視を開始する。電圧の監視は、解析対象信号の信号電圧が所定の閾値電圧Ｖ_Ｔｈを超えるまで継続され（Ｓ１０１、Ｎｏ）、解析対象信号の信号電圧が閾値電圧Ｖ_Ｔｈを超えた時点で（Ｓ１０１、Ｙｅｓ）、処理はステップ１０２に進む。

ステップ１０２では、検出期間Ｔ１にわたって解析対象信号の電圧値を所定のサンプリングレートでサンプリングし、検出期間Ｔ１が経過した時点でサンプリングを停止する（Ｓ１０３）。続くステップ１０４では、サンプリングした電圧値の中の最大値Ｖ_ｍａｘをサンプリングした時刻をピーク時刻ｔとして取得し、ピーク時刻ｔの時系列データ｛ｔ_ｉ｝を記憶部２４に保存する。

ここで、解析対象信号のゆらぎを精度よく解析するためには、解析対象信号のピークを正確に検知することが肝要となる。そのため、ステップ１０２におけるサンプリングレートは、解析対象信号に生じるゆらぎのタイムスケールに応じた適切な値に設定する必要があり、解析対象信号のピーク時刻ｔがβミリ秒単位でゆらぐ場合、サンプリングレートは、概ね、（１／β）ｋＨｚ以上に設定することが好ましい。例えば、転がり軸受に代表されるような精密な人工の回転体を測定対象とする場合、サンプリングレートを２００ｋＨｚ以上に設定することが好ましく、回転体を含む複数の部品が非線形結合してなる複合体（例えば、自動車）を測定対象とする場合は、サンプリングレートを２０ｋＨｚ以上に設定することが好ましく、さらに、生体の心拍信号を測定対象とする場合は、サンプリングレートを１ｋＨｚ以上に設定することが好ましい。

続くステップ１０５では、サンプリング停止後、待機期間Ｔ２が経過したか否を判断する。以降、待機期間Ｔ２が経過するまでこの判断処理をループする（Ｓ１０５、Ｎｏ）。

待機期間Ｔ２が経過すると（Ｓ１０５、Ｙｅｓ）、処理は再びステップ１０１に戻って、解析対象信号の電圧監視を再開する。以降、上述した一連の処理（Ｓ１０１〜Ｓ１０５）を繰り返し実行する。

上述したように、本実施形態においては、検出期間Ｔ１に連続する待機期間Ｔ２において解析対象信号の電圧監視を中断することによって、待機期間Ｔ２に発生する予期しないノイズがピークとしてカウントされることを防止する。なお、本実施形態においては、待機期間Ｔ２として、正常時の測定対象から得られる解析対象信号のピーク間隔の略半分の時間長と検出期間Ｔ１の時間長の差分に相当する時間長を設定することが好ましい。

以上、ピーク時刻取得部２１が実行する処理について説明してきたが、続いて、ＤＦＡ対象データ生成部２３が実行する処理を図８に示すフローチャートに基づいて説明する。

ステップ２０１では、ピーク時刻取得部２１が記憶部２４に保存した時系列データ｛ｔ_ｉ｝の中から、時系列順に（Ａ＋１）個のピーク時刻ｔをロードする。ここで、「Ａ」は、後述する時系列データ｛ｐ_ｉ｝の要素数を決定するパラメータである。本実施形態において、「Ａ」の値は、２００〜３００の範囲の自然数に設定することが好ましく、「２５０」に設定することが望ましい。

続くステップ２０２では、記憶部２４からロードした（Ａ＋１）個のピーク時刻ｔからＡ個のピーク間隔ｘを算出し、算出したＡ個のピーク間隔ｘを要素とする時系列データ｛ｐ_ｉ｝を生成する。具体的には、記憶部２４からロードした（Ａ＋１）個のピーク時刻ｔを時系列に並べ、（ｉ＋１）番目の時刻と（ｉ）番目の時刻の時間差分をピーク間隔ｘとして算出する。その結果、｛ｐ_ｉ｝は、Ａ個のピーク間隔ｘからなる時系列データとして生成される。

続くステップ２０３では、ステップ２０２で生成した時系列データ｛ｐ_ｉ｝をＢ個のセグメントに分割する。このとき、各セグメントに連番を付与しておく。ここで、「Ｂ」は、時系列データ｛ｐ_ｉ｝の分割数を決定するパラメータである。本実施形態において、「Ｂ」の値は、４〜６の範囲の自然数に設定することが好ましく、「５」に設定することが望ましい。例えば、Ａ＝２５０、Ｂ＝５の設定の場合、図９（ａ）に示すように、時系列データ｛ｐ_ｉ｝をそれぞれが５０個の要素からなる５種類のセグメントに分割し、各セグメントに１〜５の番号を付する。

続くステップ２０４では、Ｂ種類のセグメントの中からランダムに１つのセグメントを選出し、選出したセグメントを時系列データ｛ｐ_ｉ｝の後方に結合する。例えば、図９（ａ）に示す時系列データ｛ｘ_ｉ｝の場合、乱数発生器（例えば、メルセンヌ・ツイスター擬似乱数発生器）がランダムに出力する１から５まで数字に応じて、当該数字を付したセグメントを選出し、選出したセグメントを時系列データ｛ｐ_ｉ｝の最後尾に結合する。以下、セグメントの結合により延長されるこの擬似的な時系列データを時系列データ｛ｘ_ｉ｝という。

続くステップ２０５では、時系列データ｛ｘ_ｉ｝の要素数が目標数Ｍに達しているか否かを判断する。その結果、時系列データ｛ｘ_ｉ｝の要素数が目標数Ｍに達していない場合は（ステップ２０５、Ｎｏ）、ステップ２０４に戻って、ランダムに選出したセグメントを時系列データ｛ｘ_ｉ｝の最後尾に結合する。以降、｛ｘ_ｉ｝の要素数が目標数Ｍに達するまで、図９（ｂ）に示すように、セグメントの結合が繰り返され（ステップ２０５、Ｎｏ→ステップ２０４）、｛ｘ_ｉ｝の要素数が目標数Ｍに達した時点で（ステップ２０５、Ｙｅｓ）、処理はステップ２０６に進む。ここで、「Ｍ」は、ＤＦＡの対象となる連続データの要素数を決定するパラメータである。本実施形態において、「Ｍ」の値は、２５０００〜３５０００の範囲の自然数に設定することが好ましく、「３００００」に設定することが望ましい。

続くステップ２０６では、時系列データ｛ｘ_ｉ｝を構成するＭ個のピーク間隔ｘの平均値ｘ_ａｖｅを算出した後、時系列データ｛ｘ_ｉ｝を構成する各要素から平均値ｘ_ａｖｅを差し引くことによって、時系列データ｛(ｘ_ｉ-ｘ_ａｖｅ)_ｉ｝を生成する。図１０（ａ）および（ｂ）は、それぞれ、時系列データ｛ｘ_ｉ｝のグラフおよび時系列データ｛(ｘ_ｉ-ｘ_ａｖｅ)_ｉ｝のグラフを例示的に示す。

続くステップ２０７では、時系列データ｛(ｘ_ｉ-ｘ_ａｖｅ)_ｉ｝を積分して時系列データ｛ｙ_ｉ｝を生成する。下記式（１）は、時系列データ｛ｙ_ｉ｝の算出式を示す。

具体的には、時系列データ｛(ｘ_ｉ-ｘ_ａｖｅ)_ｉ｝の各要素を時系列順に足し合わせることによって時系列データ｛ｙ_ｉ｝を生成する。図１０（ｃ）は、ステップ２０７で生成される時系列データ｛ｙ_ｉ｝のグラフを例示的に示す。

続くステップ２０８では、上述した一連の手順で生成した擬似的な時系列データ｛ｙ_ｉ｝をＤＦＡ対象データとして生成順にナンバリングして記憶部２４に保存する。その後、処理は、再びステップ２０１に戻り、以降、ステップ２０１〜２０８を繰り返し実行する。このとき、二巡目以降のステップ２０１では、記憶部２４に保存されている時系列データ｛ｔ_ｉ｝の中から、直前のステップ２０１でロードした（Ａ＋１）個のピーク時刻ｔにおける（Ａ＋１）番目のピーク時刻tを重複する形で（Ａ＋１）個の要素をロードすることが好ましい。

また、二巡目以降のステップ２０１において、直前のステップ２０１でロードした（Ａ＋１）個の要素の中の２以上の要素を重複する形で（Ａ＋１）個の要素をロードするようにしてもよい。例えば、図１１に示すように、Ａ＝２５０に設定した場合、最初のステップ２０１で、２５１個の要素（第１データ）をロードした後、二巡目のステップ２０１で、第１データの後半の１２６個の要素に未ロードの１２５個の要素を加えた２５１個の要素（第２データ）をロードし、三巡目のステップ２０１で、第２データの後半の１２６個の要素に未ロードの１２５個の要素を加えた２５１個の要素（第３データ）をロードするといった手順を繰り返してもよい。

以上、ＤＦＡ対象データ生成部２３が実行する処理について説明してきたが、続いて、ＤＦＡ実行部２５が実行する処理を図１２および図１３に示すフローチャートに基づいて説明する。

まず、処理開始時において、時系列データ｛ｙ_ｉ｝の生成順を示すカウンタ[ｃ]を初期値「１」にセットする（ステップ３０１）。

続くステップ３０２では、ｃ番目（開始時においては１番目）の時系列データ｛ｙ_ｉ｝を記憶部２４からロードする。

続くステップ３０３では、記憶部２４に格納されるボックスサイズ・データから設定されたボックスサイズ範囲のボックスサイズ・データをロードする。ここで、ボックスサイズ・データとは、ＤＦＡにおいて使用する複数のボックスサイズ（整数）のセットを意味し、ボックスサイズとは、データの要素数を意味し、ボックスサイズ範囲とは、ＤＦＡにおいて使用するボックスサイズ（整数）がとる範囲を意味する。

図１４は、記憶部２４に格納されるボックスサイズ・データを例示的に示す。図１４に示すボックスサイズ・データは、１０〜１０００の範囲に存在する１３６個の整数をボックスサイズとして有しており、１３６個の整数は、１０〜１００の間は１ずつ増分し、１００〜５００の間は１０ずつ増分し、５００〜１０００の間は１００ずつ増分している。

ステップ３０３では、図１４に示す整数群の中から設定されたボックスサイズ範囲の整数群をロードする。ここで、ボックスサイズ範囲は、測定対象に応じた適切な範囲を予め設定しておく。例えば、転がり軸受に代表されるような精密な人工の回転体を測定対象とする場合は、ボックスサイズ範囲を１０〜３０を中心とした範囲に設定することが好ましい。仮に、当該回転体に係る解析対象信号の周波数を１００Ｈｚとすれば、上記ボックスサイズ範囲は、０．１秒〜０．３秒を中心とした時間長の範囲に置き換えて把握することができる。また、回転体を含む複数の部品が非線形結合してなる複合体（例えば、自動車）を測定対象とする場合は、ボックスサイズ範囲を３０〜６０を中心とした範囲に設定することが好ましい。仮に、当該複合体に係る解析対象信号の周波数を１００Ｈｚとすれば、上記ボックスサイズ範囲は、０．３秒〜０．６秒を中心とした時間長の範囲に置き換えて把握することができる。さらに、心拍信号に代表されるような周期的な生体信号を測定対象とする場合は、ボックスサイズ範囲を３０〜２７０を中心とした範囲に設定することが好ましい。仮に、当該生体信号の周波数を１Ｈｚとすれば、上記ボックスサイズ範囲は、３０秒〜２７０秒を中心とした時間長の範囲に置き換えて把握することができる。

続くステップ３０４では、ロードしたボックスサイズ・データの中から最初のボックスサイズ[Ｎ]（例えば[１０]）をセットする。

続くステップ３０５では、ステップ３０２でロードした時系列データ｛ｙ_ｉ｝をその時点でセットされているボックスサイズ[Ｎ]で分割する。例えば、その時点でセットされているボックスサイズが[１０]であった場合、時系列データ｛ｙ_ｉ｝を１０個の要素を含む小区間（以下、ボックスという）に分割する。この場合、時系列データ｛ｙ_ｉ｝はＭ個の要素からなるので、ステップ３０５において、時系列データ｛ｙ_ｉ｝は、Ｍ／Ｎ個のボックスに分割されることになる。図１５（ａ）は、ボックスサイズ[１０]で分割した時系列データ｛ｙ_ｉ｝を示す。この場合、各ボックス（ＢＯＸ（１）、ＢＯＸ（２）、ＢＯＸ（３）…）は、１０個の要素を含む。

続くステップ３０６では、分割後のボックス（ＢＯＸ（１）、ＢＯＸ（２）、ＢＯＸ（３）…）のそれぞれにつき、当該ボックス内に存在するＮ個のデータに対して近似曲線をフィッティングし、当該近似曲線上の値を各ボックスの局所トレンドｙ_ｖとして決定する。ここで、近似曲線のフィッティングは１次関数から４次関数までを用いた最小二乗法によって行うことができる。なお、ここでいう近似曲線は、直線を含む概念である。図１５（ｂ）は、各ボックス（ＢＯＸ（１）、ＢＯＸ（２）、ＢＯＸ（３）…）に対して近似曲線ｙ_ｖ（ｙ_ｖ（１）、ｙ_ｖ（２）、ｙ_ｖ（３）…）をフィッティングした状態を示している。なお、図１５（ｂ）においては、説明の便宜上、近似曲線を１次関数で示している。

続くステップ３０８では、ボックス（ＢＯＸ（１）、ＢＯＸ（２）、ＢＯＸ（３）…）のそれぞれにつき、ボックスを構成する各要素から当該ボックスについて決定した局所トレンドｙ_ｖを差し引くことにより、時系列データ｛ｚ_ｉ｝を生成する。下記式（２）は、時系列データ｛ｚ_ｉ｝の式を示し、図１５（ｃ）は、ステップ３０６で生成される時系列データ｛ｚ_ｉ｝のグラフを示す。

続くステップ３０８では、時系列データ｛ｚ_ｉ｝のボックス（ＢＯＸ（１）、ＢＯＸ（２）、ＢＯＸ（３）…）のそれぞれにつき、ボックスを構成する先頭の要素（点線の丸で囲んで示す）の値と末尾の要素（点線の四角で囲んで示す）の値の差分を求める。

なお、上述したステップ３０７および３０８の処理に代えて、以下の方法を採用して計算量を軽減してもよい。すなわち、各ボックス（ＢＯＸ（１）、ＢＯＸ（２）、ＢＯＸ（３）…）の先頭の要素から対応する局所トレンドｙ_ｖを差し引いた値と末尾の要素から対応する局所トレンドｙ_ｖを差し引いた値の差分を求める方法である。

続くステップ３０９では、全てのボックスについて求めた差分（先頭と末尾の差分）の二乗平均平方根[Ｓ]を算出する。

続くステップ３１０では、ステップ３０９で算出した二乗平均平方根[Ｓ]とその時点でセットされているボックスサイズ[Ｎ]からなる数値の組（Ｎ，Ｓ）を記録する。

続くステップ３１１では、ステップ３０３でロードしたボックスサイズ・データに含まれる全てのボックスサイズ[Ｎ]について、組（Ｎ，Ｓ）を記録したか否かを判断する。その結果、全てのボックスサイズ[Ｎ]について組（Ｎ，Ｓ）の記録が完了していない場合には（ステップ３１１、Ｎｏ）、処理はステップ３１２に進む。

ステップ３１２では、ステップ３０３でロードしたボックスサイズ・データに含まれる値の中から次のボックスサイズ[Ｎ]（例えば[１１]）を新たにセットする。その後、処理は、再びステップ３０５に戻り、以降、ステップ３１１において、全てのボックスサイズ[Ｎ]について組（Ｎ，Ｓ）を記録したと判断されるまでの間、上述した一連の処理（Ｓ３０５〜Ｓ３１２）を繰り返す。

つまり、本実施形態では、ステップ３０３でロードしたボックスサイズ・データに含まれる全てのボックスサイズについて、上述した一連の処理（Ｓ３０４〜Ｓ３１２）を実行することによって、ボックスサイズの数だけ数値の組（Ｎ，Ｓ）が記録される。

ここで、ボックスサイズ[Ｎ]と二乗平均平方根[Ｓ]の関係は、下記式（３）に示す関数Ｓ(Ｎ)として定義される。なお、下記式（３）において、“Ｍ”は時系列データ｛ｙ_ｉ｝の要素数を示し、“Ｎ”はボックスサイズを示し、“ｚ_ｊＮ＋Ｎ−ｚ_ｊＮ＋１“は、ｊ番目のボックスの先頭の要素（ｚ_ｊＮ＋１）と末尾の要素（ｚ_ｊＮ＋Ｎ）の差分（変位）を示す。

一方、ステップ３１１において、全てのボックスサイズ[Ｎ]について組（Ｎ，Ｓ）を記録したと判断した場合（ステップ３１１、Ｙｅｓ）、処理は、図１３に示すステップ３１３に進む。

ステップ３１３では、記録した組（Ｎ，Ｓ）を両対数スケールでプロットし、これに１次関数をフィッティングする。

続くステップ３１４では、フィッティングした１次関数の傾きをスケーリング指数αとして取得する。

続くステップ３１５では、このスケーリング指数αを取得順にナンバリングして記憶部２４に保存する。

続くステップ３１６では、[ｃ]を１増分し、その後、処理は、再び、図１２に示すステップ３０２に戻る。以降、上述した一連の処理（Ｓ３０２〜Ｓ３１６）を繰り返し実行する。

以上、ＤＦＡ実行部２５が時系列データ｛ｙ_ｉ｝に基づいてスケーリング指数αを求める処理について説明してきた。上述したＤＦＡの手法は、解析対象信号のゆらぎの分布関数がガウス分布である場合だけでなく、レヴィ分布である場合にも有効である点で好ましい手法である。ただし、解析対象信号のゆらぎの分布関数がガウス分布である場合には、他の既知のＤＦＡの手法でスケーリング指数を求めるようにしてもよい。例えば、下記式（４）、（５）に基づいて時系列データ｛ｙ_ｉ｝の各ボックスにおける分散の平均Ｆ（Ｎ）を算出し、ボックスサイズＮおよび分散Ｆを両対数スケールでプロットし、これにフィッティングした１次関数の傾きをスケーリング指数αとして求めてもよい。なお、下記式（４）、（５）において、“Ｎ”はボックスサイズ示し、“Ｍ”は時系列データ｛ｙ_ｉ｝の要素数を示し、“ｙ_ｖ（ｊＮ+ｋ）”は居所トレンドのフィッティング関数を示す。

以上、ＤＦＡ実行部２５が実行する処理について説明してきたが、続いて、スケーリング指数積算部２６が実行する処理を図１６に示すフローチャートに基づいて説明する。

まず、処理開始時において、閾値α_ｔｈをセットする（ステップ４０１）。ここで、閾値α_ｔｈは、測定対象となる物体の正常時の解析対象信号のＤＦＡにより取得されるスケーリング指数を基準として設定する。具体的には、正常時のスケーリング指数に適切なマージンを加えた値を閾値α_ｔｈとする。

例えば、精密軸受に代表される精密な人工の回転体の正常時のスケーリング指数は、概ね「０」に近い値になるので、このような回転体を測定対象とする場合は、０．２〜０．５の範囲の適切な値を閾値α_ｔｈとして設定することができる。また、複数の部品が非線形結合してなる複合体（モーターや内燃機関、これらを搭載した自動車など）のスケーリング指数は、概ね「０．５」を中心とした値になるので、このような複合体を測定対象とする場合は、０．６〜１．０の範囲の適切な値を閾値α_ｔｈとして設定することができる。

続くステップ４０２では、スケーリング指数αの取得順を示すカウンタ[ｄ]を初期値「１」にセットするとともに、スケーリング指数αの値を積算するための積算カウンタ[Ｉ]を初期値「０」にセットする（ステップ４０２）。

続くステップ４０３では、ｄ番目のスケーリング指数αを記憶部２４からロードする。

続くステップ４０４では、ロードしたスケーリング指数αが閾値α_ｔｈを下回っているか否かを判断する。その結果、閾値α_ｔｈを下回っていない場合には（ステップ４０４、Ｎｏ）、最新の[Ｉ]の値にロードしたスケーリング指数αの値を加算した後（ステップ４０５）、ステップ４０７に進む。一方、閾値α_ｔｈを下回っている場合には（ステップ４０４、Ｙｅｓ）、最新の[Ｉ]を初期値[０]にリセットした後（ステップ４０６）、ステップ４０７に進む。

続くステップ４０７では、[ｄ]を１増分し、その後、処理は、再びステップ４０３に戻る。以降、上述した一連の処理（Ｓ４０３〜Ｓ４０７）を繰り返し実行する。

ここで、ロードしたスケーリング指数αが閾値α_ｔｈを下回っていた場合に最新の[Ｉ]を初期値[０]にリセットする理由について説明する。

図１７（ａ）は、ＤＦＡ実行部２５が取得するスケーリング指数αの時系列グラフ２００を示す。図１７（ａ）に示すように、ＤＦＡ実行部２５が取得するスケーリング指数αは、時間経過の中で大きく増減するため、スケーリング指数αの値そのものに着目するやり方では、測定対象の物体内部に生じたひずみを評価することは難しい。この点につき、本発明は、時間経過の中でスケーリング指数αの値が閾値α_ｔｈを超えて高止まりする状態を測定対象の物体内部にひずみが蓄積している状態として捉えること特徴とする。つまり、本実施形態では、上述した積算値をリセットするステップ（ステップ４０４→ステップ４０６）を採用することにより、スケーリング指数αが高止まりする期間の長さを積算値の大きさで表すようになっている。

上述した点を理解するために、図１７（ｂ）にスケーリング指数積算部２６が導出するスケーリング指数αの積算値の時系列グラフ３００を示す。なお、時系列グラフ３００は、図１７（ａ）に示すスケーリング指数αについて閾値α_ｔｈを[１．０]にセットした場合の積算値の時系列グラフである。ここで、時系列グラフ２００と時系列グラフ３００を対比すれば、時間経過の中で、スケーリング指数αが高止まりする期間（時系列グラフ２００おいて破線で囲む）の継続に従って、積算グラフ２００において積算値が増大していることが理解されるであろう。

以上、スケーリング指数積算部２６が実行する処理について説明してきたが、続いて、出力情報生成部２７が実行する処理を説明する。

出力情報生成部２７は、スケーリング指数積算部２６が導出するスケーリング指数αの最新の積算値に基づいて出力情報を生成する。ここでいう出力情報とは、ディスプレイ４０に出力する画像データ、スピーカー５０に出力する音声データ、他の装置６０を制御する制御データを含む概念である。

好ましい実施形態では、出力情報生成部２７は、カウンタ[Ｉ]の値に基づいて、図１７（ｂ）に示すようなスケーリング指数αの積算値の時系列グラフを表示するための画像データを生成する。出力情報生成部２７は、スケーリング指数αの積算値が更新される度に最新の積算値を反映した時系列グラフを表示する画像データを生成してディスプレイ４０に出力する。この場合、ユーザは、ディスプレイ４０に表示される時系列グラフにおいて、折れ線グラフの三角形状の波形が大きく成長する様子から測定対象となる物体の内部でひずみが蓄積していることを直観的に把握することができ、これを物体の破壊や崩壊の予兆として捉えることができる。なお、別法では、スケーリング指数αの最新の積算値を棒グラフで表現することもできる。この場合、ユーザは、ディスプレイ４０に表示される棒グラフが伸長する様子を物体の破壊や崩壊の予兆として捉えることができる。

また、出力情報生成部２７は、カウンタ[Ｉ]の値に基づいて、スケーリング指数αの最新の積算値の大きさに応じた警報音を発音するための音声データを生成する。具体的には、積算値について１以上の適切な閾値を定めておき、積算値が所定の閾値に達した時点でその閾値の大きさに応じた音色や音量の警報音を発音するための音声データを生成してスピーカー５０に出力する。この場合、ユーザは、スピーカー５０が発音する警報音を測定対象となる物体の破壊や崩壊の予兆として捉えることができる。

さらに、出力情報生成部２７は、カウンタ[Ｉ]の値に基づいて、他の装置６０を制御するための制御データを生成する。ここで、他の装置６０とは、例えば、測定対象の運動物体の制御装置６０である。この場合、出力情報生成部２７は、積算値が所定の閾値に達した時点でその閾値の大きさに応じて、測定対象の運動物体の制御データを生成して制御装置６０に出力する。ここでいう制御データは、例えば、運動物体の速度を制御する速度制御データであり、好ましくは、緊急停止信号である。この場合、運動物体は、当該信号を受けて緊急停止することで破壊や崩壊を免れる。

以上、出力情報生成部２７が実行する処理について説明してきたが、続いて、測定モード設定部２２が実行する処理を説明する。

本実施形態においては、予定される測定対象の類型に応じた複数の測定モードを用意し、測定モード毎に上述した各種パラメータのプリセットデータを記憶部２４に格納しておく。本実施形態におけるプリセットデータとしては、ピーク時刻取得部２１が使用するパラメータ（信号サンプリングレート、検出期間Ｔ１、待機期間Ｔ２、閾値電圧Ｖ_Ｔｈ）、ＤＦＡ対象データ生成部２３が使用するパラメータ（Ａ、Ｂ、Ｍ）、ＤＦＡ実行部２５が使用するパラメータ（ボックスサイズ範囲）、スケーリング指数積算部２６が使用するパラメータ（閾値α_ｔｈ）、加振信号生成部３２に生成させる加振信号の周波数、正弦波発振回路３４に発振させる正弦波信号の周波数、などを例示することができる。本実施形態において、測定モード設定部２２は、ユーザが選択した測定モードに対応するプリセットデータを記憶部２４からロードし、ロードしたプリセットデータに含まれる各パラメータを対応する各機能部にセットする。

以上、説明したように、本実施形態の信号解析システムによれば、物体が発する振動波形のゆらぎをＤＦＡで連続的に解析することによって、物体内部に発生したひずみの蓄積状態をスケーリング指数の積分値によって定量的に評価することができるようになる。例えば、周波数３００Ｈｚの解析対象信号に対して上述したパラメータ「Ａ」を３００に設定した場合、本実施形態の信号解析システムは、約１秒の時間分解能をもって物体内部に発生したひずみの蓄積状態を評価することができる。

なお、本実施形態の信号解析システムによれば、測定対象となる物体に振動センサを常設することによって、当該物体の状態の常時監視が可能になる。この場合、振動センサのセンサ信号を常時連続的に取得するようにしてもよいし、所定期間に渡ってセンサ信号を取得することを間欠的に繰り返すようにしてもよい。

以上、本発明について実施形態をもって説明してきたが、本発明は上述した実施形態に限定されるものではなく、当業者が推考しうる実施態様の範囲内において、本発明の作用・効果を奏する限り、本発明の範囲に含まれるものである。

また、上述したコンピュータ２０の各機能は、Ｃ、Ｃ＋＋、Ｃ＃、Ｊａｖａ（登録商標）などで記述された装置実行可能なプログラムにより実現することができ、本実施形態のプログラムは、ハードディスク装置、ＣＤ−ＲＯＭ、ＭＯ、ＤＶＤ、フレキシブルディスク、ＥＥＰＲＯＭ、ＥＰＲＯＭなどの装置可読な記録媒体に格納して頒布することができ、また他装置が可能な形式でネットワークを介して伝送することができる。

以下、本発明の信号解析システムについて、実施例を用いてより具体的に説明を行なうが、本発明は、後述する実施例に限定されるものではない。

（実験１：周期振動する運動物体の測定）
固定台の上にモーターを設置した上で、モーターを約３００Ｈｚで駆動させた。なお、本実験においては、モーターの破壊を促すためにモーターをガラスウールで被覆してあえて内部に熱がこもるようにした。この状態で、固定台に取り付けた振動センサ（ピエゾ素子）のセンサ信号を本発明の信号解析システムを使用して解析した。下記表１は、実験１で信号解析システムに設定したパラメータの値をまとめて示す。

（実験１の結果）
モーターの駆動を開始してしばらくすると、信号解析システムが表示する折れ線グラフにおいて三角形状の波形が大きく成長し始めたが、この時点では、特に異変は感じ取れなかった。その後、しばらくしてモーターから異音が発生し、最終的にモーターは停止した。

（実験２：静止物体の測定）
固定台の上にアルミ製の棒材を固定した上で、この棒材に周波数３００Ｈｚの音を発信するスピーカーを固定することで振動を入力した。その上で、棒材の一方の端につるした水桶に徐々に水を注入することによって棒材に応力が生じるようにした。この状態で、棒材に取り付けた振動センサ（ピエゾ素子）のセンサ信号を本発明の信号解析システムを使用して解析した。下記表２は、実験２で信号解析システムに設定したパラメータの値をまとめて示す。

（実験２の結果）
水桶への水の注入を開始してしばらくすると、本発明の信号解析システムが表示する折れ線グラフにおいて三角形状の波形が大きく成長し始めた。その後、しばらくして棒材が水桶の重みに負けて折れ曲がった。

以上の実験結果から、本発明の信号解析システムによれば、スケーリング指数の積分値に基づいて物体の破壊や崩壊の予兆を捉えることができることが示された。

１０…振動センサ
１２…加振器
２０…コンピュータ
２１…ピーク時刻取得部
２２…測定モード設定部
２３…ＤＦＡ対象データ生成部
２４…記憶部
２５…ＤＦＡ実行部
２６…スケーリング指数積算部
２７…出力情報生成部
３２…加振信号生成部
３４…正弦波発振回路
３６…信号重畳回路
４０…ディスプレイ
５０…スピーカー
６０…他の装置
１００…信号解析システム
２００…スケーリング指数の時系列データ
３００…スケーリング指数の積算値の時系列データ

Claims (22)

1. 解析対象信号のピーク時刻を取得するピーク時刻取得部と、
   前記ピーク時刻取得部が取得した前記ピーク時刻に基づいてＤＦＡ対象データを生成するＤＦＡ対象データ生成部と、
   前記ＤＦＡ対象データに基づいてＤＦＡを実行してスケーリング指数を取得するＤＦＡ実行部と、
   取得した前記スケーリング指数を積算するスケーリング指数積算部と、
   を含み、
   前記ＤＦＡ対象データ生成部は、
   前記ピーク時刻取得部が取得した前記ピーク時刻を時系列順に（Ａ＋１）個ずつロードし、ロードした該（Ａ＋１）個のピーク時刻に基づいてＡ個のピーク間隔を算出し、該Ａ個のピーク間隔を要素とする時系列データに基づいてＭ個の要素を含む擬似的な時系列データを前記ＤＦＡ対象データとして生成する処理を繰り返し実行し（Ａ、Ｍはともに自然数であり、ＭはＡより大きい。以下同様。）、
   前記ＤＦＡ実行部は、
   前記ＤＦＡ対象データ生成部が生成した前記ＤＦＡ対象データを生成順にロードし、ロードした該ＤＦＡ対象データに基づいて前記スケーリング指数を取得する処理を繰り返し実行し、
   前記スケーリング指数積算部は、
   前記ＤＦＡ実行部が取得した前記スケーリング指数を取得順にロードし、ロードした該スケーリング指数が所定の閾値を下回っていない場合には該スケーリング指数を最新の積算値に加算し、ロードした該スケーリング指数が該所定の閾値を下回っている場合には最新の積算値をリセットする処理を繰り返し実行する、
   信号解析システム。
2. 前記ＤＦＡ対象データ生成部は、
   前記Ａ個のピーク間隔を要素とする時系列データを複数のセグメントに分割し、該時系列データの後方にランダムに選出した１つのセグメントを足し合わせることを繰り返してＭ個のピーク間隔を要素とする擬似的な時系列データを生成し、該擬似的な時系列データの各要素から平均値を引いた値を積分してなる擬似的な時系列データを前記ＤＦＡ対象データとして生成する、
   請求項１に記載の信号解析システム。
3. 前記積算値に基づいて出力情報を生成する出力情報生成部をさらに含む、請求項１または２に記載の信号解析システム。
4. 前記出力情報は、最新の前記積算値を反映したグラフを表示するための画像データである、請求項３に記載の信号解析システム。
5. 前記グラフは、前記積算値の時系列グラフである、請求項４に記載の信号解析システム。
6. 前記出力情報は、警報音を発音するための音声データである、請求項３〜５のいずれか一項に記載の信号解析システム。
7. 前記出力情報は、他の装置を制御するための制御データである、請求項３〜６のいずれか一項に記載の信号解析システム。
8. 前記Ａは、２００〜３００の範囲の自然数であり、前記Ｍは、２５０００〜３５０００の範囲の自然数である、請求項１〜６のいずれか一項に記載の信号解析システム。
9. 前記解析対象信号は、運動する物体の振動を検出するセンサ信号である、請求項１〜７のいずれか一項に記載の信号解析システム。
10. 前記運動する物体は、軸受であり、
    前記ＤＦＡ実行部がスケーリング指数を取得するためのボックスサイズの範囲は、１０〜３０を中心とする範囲である、請求項９に記載の信号解析システム。
11. 前記解析対象信号は、外部から周期振動が入力された物体の振動を検出するセンサ信号である、請求項１〜７のいずれか一項に記載の信号解析システム。
12. 前記解析対象信号は、物体の振動を検出するセンサ信号と所定の正弦波信号の重畳信号である、請求項１〜７のいずれか一項に記載の信号解析システム。
13. 前記所定の閾値は、正常時の前記物体に係る前記解析対象信号のＤＦＡにより取得されるスケーリング指数を基準として設定される、請求項９〜１２のいずれか一項に記載の信号解析システム。
14. 選択された測定モードに応じて、前記Ａ、前記Ｍ、前記ＤＦＡ実行部がスケーリング指数を取得するためのボックスサイズの範囲および前記所定の閾値のそれぞれに係るパラメータのうち１以上のパラメータを自動的に設定する測定モード設定部を含む、請求項１〜１３のいずれか一項に記載の信号解析システム。
15. 信号を解析する方法であって、
    解析対象信号のピーク時刻を取得するステップと、
    取得した前記ピーク時刻に基づいてＤＦＡ対象データを生成するステップと、
    前記ＤＦＡ対象データに基づいてＤＦＡを実行してスケーリング指数を取得するステップと、
    取得した前記スケーリング指数を積算するステップと、
    を含み、
    前記ＤＦＡ対象データを生成するステップは、
    前記ピーク時刻取得部が取得した前記ピーク時刻を時系列順に（Ａ＋１）個ずつロードし、ロードした該（Ａ＋１）個のピーク時刻に基づいてＡ個のピーク間隔を算出し、該Ａ個のピーク間隔を要素とする時系列データに基づいてＭ個の要素を含む擬似的な時系列データを前記ＤＦＡ対象データとして生成する処理を繰り返し実行するステップであり（Ａ、Ｍは自然数であり、ＭはＡより大きい。以下同様。）、
    前記スケーリング指数を取得するステップは、
    生成した前記ＤＦＡ対象データを生成順にロードし、ロードした該前記ＤＦＡ対象データに基づいて前記スケーリング指数を取得する処理を繰り返し実行するステップであり、
    前記スケーリング指数を積算するステップは、
    取得した前記スケーリング指数を取得順にロードし、ロードした該スケーリング指数が所定の閾値を下回っていない場合には該スケーリング指数を最新の積算値に加算し、ロードした該スケーリング指数が該所定の閾値を下回っている場合には最新の積算値をリセットする処理を繰り返し実行するステップである、
    方法。
16. 前記解析対象信号は、運動する物体の振動を検出するセンサ信号である、請求項１５に記載の方法。
17. 前記運動する物体は、軸受であり、
    前記ＤＦＡにおいてスケーリング指数を取得するためのボックスサイズの範囲は、１０〜３０を中心とする範囲である、請求項１６に記載の方法。
18. 前記運動する物体は、複数の部品が非線形結合してなる複合体であり、
    前記ＤＦＡにおいてスケーリング指数を取得するためのボックスサイズの範囲は、３０〜６０を中心とする範囲である、請求項１６に記載の方法。
19. 前記解析対象信号は、外部から周期振動が入力された物体の振動を検出するセンサ信号である、請求項１５に記載の方法。
20. 前記解析対象信号は、物体の振動を検出するセンサ信号と所定の正弦波信号の重畳信号である、請求項１５に記載の方法。
21. 前記解析対象信号は、周期的な生体信号であり、
    前記ＤＦＡにおいてスケーリング指数を取得するためのボックスサイズの範囲は、３０〜２７０を中心とする範囲である、請求項１５に記載の方法。
22. コンピュータに、請求項１５〜２１のいずれか一項に記載の方法の各ステップを実行させるためのコンピュータ実行可能なプログラム。