制御出力のサンプリング周期が、制御入力のそれの偶数P倍である入力多重マルチレート系の、制御入力側のサンプリング周期における制御対象の周波数応答を同定する方法であって、擬似白色信号の一種であるM系列信号の周期Mpに基づき、データ取得長Mp×P-1分のM系列信号を生成するM系列信号生成ステップと、前記生成されたM系列信号と、前記生成されたM系列信号を入力して前記制御対象から得られる出力データに基づき、前記制御対象のインパルス応答推定値を計算するインパルス応答推定ステップと、前記インパルス応答推定値を離散フーリエ変換することによって、前記制御対象の周波数応答を同定する周波数応答同定ステップとを備えた方法。前記インパルス応答推定ステップは、ukを、前記M系列信号における時刻kの信号、ykを、ukを入力したときの前記制御対象の出力データ、σu2をukの分散としたとき、により、前記インパルス応答推定値を計算する、ことを特徴とする請求項1に記載の方法。制御出力のサンプリング周期が、制御入力のそれの偶数P倍である入力多重マルチレート系の、制御入力側のサンプリング周期における制御対象の周波数応答を同定する装置であって、擬似白色信号の一種であるM系列信号の周期Mpに基づき、データ取得長Mp×P-1分のM系列信号を生成するM系列信号生成部と、前記生成されたM系列信号と、前記生成されたM系列信号を入力して前記制御対象から得られる出力データ基づき、前記制御対象のインパルス応答推定値を計算するインパルス応答推定部と、前記インパルス応答推定値を離散フーリエ変換することによって、前記制御対象の周波数応答を同定する周波数応答同定部とを備えた装置。前記インパルス応答推定部は、ukを、前記M系列信号における時刻kの信号、ykを、ukを入力したときの前記制御対象の出力データ、σu2をukの分散としたとき、により、前記インパルス応答推定値を計算する、ことを特徴とする請求項3に記載の装置。本発明の一側面は、制御出力のサンプリング周期が制御入力のそれの偶数倍である入力多重マルチレート系において、制御入力側のサンプリング周期における制御対象の周波数応答を同定することを目的とする。本発明の一態様としてのマルチレート系の高速周波数応答同定法は、制御出力のサンプリング周期が、制御入力のそれの偶数P倍である入力多重マルチレート系の、制御入力側のサンプリング周期における制御対象の周波数応答を同定する方法であって、M系列信号生成ステップと、インパルス応答推定ステップと、周波数応答同定ステップとを備える。前記M系列信号生成ステップは、擬似白色信号の一種であるM系列信号の周期Mpに基づき、データ取得長Mp×P-1分のM系列信号を生成する。前記インパルス応答推定ステップは、前記生成されたM系列信号と、前記生成されたM系列信号を入力して前記制御対象から得られる出力データに基づき、前記制御対象のインパルス応答推定値を計算する。前記周波数応答同定ステップは、前記インパルス応答推定値を離散フーリエ変換することによって、前記制御対象の周波数応答を同定する。本発明の実施形態に係る同定アルゴリズムは、同定入力として擬似白色信号の一種であるM系列信号、同定モデルとしてFIRFiniteImpulseResponseモデルを用いる。また、マルチレート比は偶数P倍、すなわち、制御出力のサンプリング周期が制御入力のそれの246...倍であることが条件となる。ここで、入力多重マルチレート系の伝達関数は、以下数式で与えられることが知られている。数式において、Yは出力の伝達関数、Uは入力の伝達関数、G0は入力側のサンプリング周期における制御対象の伝達関数、Pはマルチレート比である。本実施形態で同定すべき制御対象の伝達関数は、数式におけるG0となる。数式は、例えばマルチレート比P2のとき、となる。以下のでは、マルチレート比P2を用いることを仮定して、本実施形態の同定アルゴリズムを説明する。まず、本実施形態ではG0の同定モデルを、未知パラメータを用いて、以下数式のFIRモデルとする。また、入力の時系列をu0u1...umとすれば、数式のUは数式と同様に、以下数式のFIRモデルで表すことが出来る。さらに、となる。ここで、数式のG0をに置き換えた場合の出力、すなわちモデル出力の各時刻k012...における時系列は、数式より以下のように単純な畳み込みによって計算することが出来る。以下同様以上より、モデル出力の時系列は、kが奇数のとき出力が間引きダウンサンプリングされた応答になることが分かる。同定モデルの未知パラメータを推定するには、実際のハードウェア例えば磁気ディスク装置のヘッド位置決め制御系に入力信号u0u1...umを印加した場合の出力の時系列ykと、モデル出力の時系列との誤差ekの二乗和を最小化する問題最小二乗法を解けば良い。マルチレート比P2の入力多重マルチレート系では、kが奇数のときはとも間引きされ0になるのでekは、以下数式のようになる。従って、解くべき最小二乗法は、となる。数式の解は、未知パラメータと入力u0u1u2...から成る以下数式のベクトル、を用いて以下数式で与えられる。次に、数式を実際に計算することを考える。本実施形態ではukに擬似白色信号の一種であるM系列を用いるので、数式の総和回数Nを、M系列の周期Mpより、と選ぶと、数式右辺の1番目の括弧内、は、となる。ここで、M系列の周期性を考慮すると数式は、と変形出来、M系列1周期分の自己相関行列となることが分かる。従って、白色信号の自己相関の性質より数式の非対角項は、以下数式のように近似的に0とすることが出来る。数式のσu2はukの分散である。以上より、数式右辺の逆行列演算は単純な除算に置き換えることが出来、計算負荷を大幅に低減することが出来る。すなわち、計算負荷の大きい逆行列演算を行う必要はない。数式を成立させるために必要最小限の同定実験時の取得データ長Lは、となる。同定実験時の取得データ長の概念を図1に示す。図1では、マルチレート比P2、M系列の周期Mp7の場合のデータ取得長が示される。以上より、マルチレート比P2の場合の未知パラメータ推定アルゴリズムは、となる。一般的な偶数倍マルチレート比P246...の場合は、となる。以上の数式が本実施形態の基本的なアルゴリズムである。また、数式によって計算されるはG0のインパルス応答推定値FIRモデルの分子多項式係数であるため、を離散フーリエ変換すればG0の周波数応答を求めることが出来る。離散フーリエ変換はFFTFastFourierTransformを用いて容易に実行可能である。以上示した本実施形態のアルゴリズムは、同定入力信号をM系列信号と規定しているため、一般的に正弦波掃引法などと比較して同定実験に要する時間を短縮することが出来る。また、数式は単純な加算及び乗除算で構成されているため、本実施形態を適用するハードウェアに搭載されるマイクロプロセッサや外部の計算機に、簡易な記述で実装することが出来る。本実施形態の実行手順を示したフローチャートは図2となる。この実行手順に従った処理は、ハードウェアに搭載されるマイクロプロセッサや、外部の計算機で実行される。ステップS1において、取得したデータ長Lが、Mp×P-1に達したかを判定する。Mp×P-1に達していないときは、次の時刻のM系列信号を生成しS2、同定対象例:磁気ディスク装置にM系列信号を入力するS3。なお、M系列信号の生成アルゴリズムは事前に与えておく。制御対象同定対象の出力を観測しS4、M系列と出力データを、内部メモリまたは外部記憶装置等の記憶部に保存するS5。データ長Lが、Mp×P-1に達するまで、以上のステップS1~S4を繰り返す。データ長LがMp×P1に達したら、すなわち、データ長L分の信号および出力データを取得したら、制御対象のインパルス応答推定値を計算する。具体的に、制御対象の同定モデルであるFIRモデルのインパルス応答推定値を、数式により計算するS6。そして、を離散フーリエ変換するS7。離散フーリエ変換の結果が、求めるべき同定対象の周波数応答となる。図6に、本実施形態に係わる装置の機能ブロック図を示す。M系列信号生成部11は、擬似白色信号の一種であるM系列信号を、データ取得長L=Mp×P-1分、生成する。生成したM系列信号の各時刻の信号がそれぞれ同定対象12に入力され、出力データが制御対象同定対象12から出力される。記憶部15は、生成部11で生成されたM系列信号と、各時刻の信号に対応して同定対象12から得られた出力データを記憶する。インパルス応答推定部13は、記憶部15に記憶されたM系列信号および出力データに基づき、制御対象12のインパルス応答推定値を計算する。この計算は、数式Aにより行う。数式は、前述した数式と同じものである。前述した通り、ukは、時刻kk番目で取得したM系列信号の成分、ykは、ukを同定対象12の入力としたときに得られる出力データ、σu2はukの分散である。フーリエ変換部14は、を離散フーリエ変換することによって、制御対象の周波数応答を同定する。以上、本実施形態によれば、入力多重マルチレート系において、制御入力側のサンプリング周期における制御対象の周波数応答を、高速かつ簡易な演算によって同定することが可能となる。以下、本発明の実施例を計算機シミュレーションによって示す。まず、複数の共振モードを持つ機械システムを想定した、以下数式の離散時間系伝達関数を、計算機シミュレーションにおける同定対象とする。サンプリング周期は50μsとする。この同定対象にサンプリング周期50μsでM系列信号を入力し、マルチレート比P2でダウンサンプリングされた出力サンプリング周期100μsを計算機シミュレーションによって得る。このときの同定対象の入出力関係を示したブロック線図を図3に示す。M系列のシフトレジスタ段数は13段周期Mp8191とした。また、数式で推定するインパルス応答の長さは1024とした。これは、数式においてn1023、すなわち、とすることに相当する。最後に、数式を実行し得られたインパルス応答推定値に対して1024点のFFTを実行し同定対象の周波数応答を得る。ここで、入出力の取得データ長を本実施形態式に従い、とした場合と、数式式に従わず任意に設定した値、とした場合とで数式の計算を行い、取得データ長が同定精度に与える影響を調べる。以上の設定の下で計算機シミュレーションを実行した結果同定対象の周波数応答を図4および図5に示す。図4および図5において、実線が同定された周波数応答、点線が同定対象の真の周波数応答である。また、図中5000Hzに表示されている破線はサンプリング周期100μs時のナイキスト周波数を示している。本実施形態の同定手法と異なり、マルチレート系の特性を陽に考慮しない同定手法では、ナイキスト周波数以上破線右側の周波数応答を同定することは困難である。図4より、本実施形態の数式を満たす取得データ長Lを設定した場合は、同定された周波数応答が真の周波数応答と良く一致しており、ナイキスト周波数以上の帯域においても正確に同定出来ていることが分かる。一方、図5より、数式を満たさない任意の取得データ長Lを設定した場合、周波数応答の大まかな一致は見られるものの、数式を満たす場合に比べて、真の周波数応答との誤差が大きくなっていることが分かる。以上の計算機シミュレーションによる実施例により、本実施形態式の効果が検証された。なお、磁気ディスク装置のヘッド位置決め制御系のように、ナイキスト周波数以上に共振モードを有する制御対象に対して、本実施形態のアルゴリズムによって得られる制御入力側のサンプリング周期における周波数応答は、共振モードを安定化するために用いられるマルチレートノッチフィルターの設計に利用することが出来る。