(57)【要約】      （修正有）【課題】波浪モデルの波浪周波数と波浪減衰係数を同定することができ、該同定された波浪周波数と減衰係数を利用して推定方位と推定角速度を出力する推定器を有する船舶用自動操舵装置を提供する。【解決手段】設定針路と検出方位に基づいて命令舵角を出力するために、推定方位2006188095.tif 000110と推定角速度2006188095.tif 000111とを出力する推定器と、フィードバックループと、検出方位から波浪周波数と波浪減衰係数とを同定する波浪同定器とを有し、波浪同定器は、スペクトル演算部と、スペクトルデータから包絡点群である同定データを抽出する同定データ抽出部２８Ｂと、中心周波数から波浪周波数ｆwnを求める中心周波数算出部２８Ｃと、その波浪周波数に対応する波浪モデルのスペクトル値（Ｓi）との偏差に基づく評価値を最小にする波浪モデルのスペクトルＳの減衰係数を、波浪減衰係数ζwnとする減衰係数算出部２８Ｄと、を備える。【選択図】図４
2006188095.tif 000002
【特許請求の範囲】
【請求項１】
設定針路と検出方位に基づいて命令舵角を出力するために、検出方位
2006188095.tif 000094
と命令舵角δｃとを入力とし推定方位
2006188095.tif 000095
と推定角速度
2006188095.tif 000096
とを出力する推定器と、推定器の出力を設定針路にフィードバックするフィードバックループと、検出方位から波浪周波数と波浪減衰係数とを同定する波浪同定器と、を有する船舶用自動操舵装置において、前記波浪同定器は、
検出方位信号の時系列データをスペクトルデータに変換するスペクトル演算部と、
スペクトルデータから包絡点群である同定データ
2006188095.tif 000097
（ｆiは周波数、Ｓi−はスペクトル値を表す）を抽出する同定データ抽出部と、
同定データに基づき中心周波数から波浪周波数ｆwnを求める中心周波数算出部と、
同定データのスペクトル値（Ｓi−）と、前記中心周波数算出部によって求められた波浪周波数に基づきその波浪周波数に対応する波浪モデルのスペクトル値（Ｓi）との偏差に基づく評価値を最小にする波浪モデルのスペクトルＳの減衰係数を、波浪減衰係数ζwnとする減衰係数算出部と、
を備えることを特徴とする船舶用自動操舵装置。
【請求項２】
前記波浪同定器は、さらに、
前記スペクトル演算部で変換されたスペクトルデータを正規化する正規化部を備えることを特徴とする請求項１記載の船舶用自動操舵装置。
【請求項３】
前記同定データ抽出部は、スペクトルデータの所定の領域において、任意の開始周波数から、その領域において最大スペクトル値をとる周波数に向って、スペクトル値が大きくなるまたは等しくなる包絡点の集まりからなる包絡点群を求め、該最大スペクトル値をとる周波数よりも低域周波数側及び高域周波数側のそれぞれにおいて求めた複数の包絡点群のうちで、包絡点の個数が最大となる包絡点群の中で、開始周波数が最大スペクトル値をとる周波数に最も近い包絡点群を、同定データとして抽出するものであることを特徴とする請求項１または２記載の船舶用自動操舵装置。
【請求項４】
前記中心周波数算出部は、同定データの中でスペクトル値が上部しきい値ＳH以上の値となる同定データ（ｆi，Ｓi−）に対して、
【数１】
2006188095.tif 000098
で波浪周波数ｆwnを求めることを特徴とする請求項１ないし３のいずれか１項に記載の船舶用自動操舵装置。
【請求項５】
前記減衰係数算出部は、同定データの中でスペクトル値が下部しきい値ＳL以上の値となる同定データ（ｆi，Ｓi−）に対して、そのスペクトル値とその周波数に対応する波浪モデルのスペクトルのスペクトル値Ｓとの偏差の２乗和Ｊζを最小値とする、
【数２】
2006188095.tif 000099
を満足する波浪モデルのスペクトルの減衰係数ζを、波浪減衰係数ζwnとすることを特徴とする請求項１ないし４のいずれか１項に記載の船舶用自動操舵装置。
【請求項６】
設定針路と検出方位に基づいて命令舵角を出力するために、検出方位
2006188095.tif 000100
と命令舵角δcとを入力とし推定方位
2006188095.tif 000101
と推定角速度
2006188095.tif 000102
とを出力する推定器と、推定器の出力を設定針路にフィードバックするフィードバックループと、検出方位から波浪周波数と波浪減衰係数とを同定する波浪同定器と、を有する船舶用自動操舵装置において、波浪同定器は、
検出方位信号の時系列データのスペクトルデータから波浪周波数と波浪減衰係数とを同定する波浪同定部と、
時系列データと前記波浪同定部で同定した波浪周波数ｆwn＝ωwn／２πと波浪減衰係数ζwnとを用いて波浪状態を判定し、該波浪状態に適した波浪減衰係数に修正する波浪判定部と、
を備えることを特徴とする船舶用自動操舵装置。
【請求項７】
波浪判定部は、波浪から命令舵角までの伝達関数Ｇeである
【数３】
2006188095.tif 000103
に対して、ζwnに前記波浪同定部で同定した波浪減衰係数を代入したときのＧeである
2006188095.tif 000104
と、ζwnに１/√２を代入したときのＧeである
2006188095.tif 000105
のそれぞれの定数を求める定数設定部と、
前記定数設定部で求めた
2006188095.tif 000106
に対して、
【数４】
2006188095.tif 000107
で表されるそれぞれの評価値Ｊon，Ｊoffを求める評価量算出部と、
前記評価値Ｊon，Ｊoffのうちで小さい方に用いたζwnの値を修正した波浪減衰係数
2006188095.tif 000108
として採用し、即ち、
【数５】
2006188095.tif 000109
とする減衰係数決定部と、
を備えることを特徴とする請求項６記載の船舶用自動操舵装置。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、推定器を有する船舶用自動操舵装置において、その推定器で利用する波浪モデルの波浪パラメータの同定を行なうために波浪同定器を有する船舶用自動操舵装置に関する。
【背景技術】
【０００２】
船舶用自動操舵装置は、設定針路に、ジャイロコンパスからの方位を追従させるために舵を制御する装置であり、制御系は、設定針路と船首方位との入力から偏差と旋回角速度とを求めフィードバックゲインを乗じて制御量である命令舵角を操舵機に出力する。操舵機は命令舵角に比例した舵角を動かすため、船体は舵角によって角速度を生じ方位が変化する。
【０００３】
船舶用自動操舵装置の閉ループ系は図９に示すように船体モデルと波浪モデルとからなる制御対象と、推定器２０とフィードバックループ１５とからなる制御系とから構成される。簡単化のため操舵機を省くので、命令舵角は舵角と等価になる。同図においてｓはラプラス演算子を、
2006188095.tif 000003
はそれぞれ設定針路，検出方位，方位，波浪を、ｒは角速度を、δｃは命令舵角を、νは白色ノイズＮ（０，１）を、それぞれ示す。方位，舵角及びノイズの単位は[deg]とし、角速度のそれは[deg/s]とする。Ｋsn，Ｔsn，Ｔs3nは、船体パラメータのノミナル値で、それぞれ旋回力ゲイン[1/s]，船体時定数[s]，舵感度時定数[s]を、ζwn，ωwn，Ｋwnは波浪パラメータのノミナル値でそれぞれ減衰係数，固有周波数[rad/s]，ゲイン[1/s]を、ＫＰ,ＫDはフィードバックゲインでそれぞれ比例ゲインと微分ゲイン[s]を、それぞれ示す。なお添字(・)ｓ，(・)ｗはそれぞれ船体と波浪とを、(・)ｎ  ，(・)ａはそれぞれノミナル値と実際値とを、アクセント
2006188095.tif 000004
は推定値を、それぞれ意味する。
【０００４】
荒天時、方位に含まれる外乱成分が制御系を経て舵角に入力するために、操舵機は無効な操舵を行い機器の疲労や磨耗、燃料の損失を招くことになる。かかる無駄舵を防止するために各種のフィルタや天候調節機能が提案されている（非特許文献１、非特許文献２）。その中で外乱モデルとして波浪モデルを組み込んだ推定器の外乱除去効果の有効性が報告されている（非特許文献３、非特許文献４）。
【０００５】
また、特許文献１では、定常カルマンフィルタを用いて、旋回角速度及び方位誤差の最適値を推定して天候調整を行うようにしたカルマンフィルタを含む自動操舵装置が開示されている。
【０００６】
しかしながら、制御対象の船体モデルは、載荷の変化によるノミナル値のパラメータ不確かさを持ち、図９中のＫsn，Ｔsn，Ｔs3nはＫsa，Ｔsa，Ｔs3aに変化する。そのためノミナル値のモデルベースで構成する推定器は載荷の変化によって推定値に誤差を生じる。このため、推定値の状態フィードバックによる閉ループ安定性は、推定誤差によって劣化し、船体を蛇行航行させる（ヨーイング）おそれがある。
【０００７】
しかしながら、従来の推定器の設計にあたっては、かかるパラメータ不確かさに対する検討については、何らなされていなかった。
【０００８】
本発明者らは、このような状況に鑑みて、船体パラメータのノミナル値のパラメータ不確かさを積極的に考慮に入れて、該パラメータ不確かさを起因とする推定誤差を小さくすることができるようにした船舶用自動操舵装置及びその推定器の設計方法を特願２００４−２２６６７０で提案している。
【０００９】
この特願２００４−２２６６７０で提案する推定器の設計方法の１つは、まず波浪モデル無しの推定器として船体モデルの状態量を推定するための固有周波数
2006188095.tif 000005
と減衰係数ζｅとを設計し、次に波浪モデル有りの推定器として船体パラメータによる推定誤差が波浪モデル無しの場合とほぼ同一になるような船体モデルの状態量を推定するための固有周波数
2006188095.tif 000006
を設計し、波浪モデルの状態量を推定するための固有周波数ωewと減衰係数ζewとを設計するようにしている。
【００１０】
本発明者らはこのように推定器の設計に波浪モデルを取り込んで、推定器で推定する推定方位と推定角速度に含まれる波浪成分を除去することにより、波浪モデル無しで設計した推定器に対して外乱除去性能の向上を図ることができることを見出した。
【００１１】
【非特許文献１】大津、長谷川、IX. オートパイロットの評価と展望、第3回操縦性シンホジウムテキスト、日本造船学会試験水槽委員会、p.243/279（1981）
【非特許文献２】Fossen, T. I., "Guidance and Control of Ocean Vehicles", John Wiley & Sons Ltd., p.222/245(1994)
【非特許文献３】Fossen, T. I., "High Performance Ship Autopilot With Wave Filter", Proceedings of the 10th International Ship Control Systems Symposium (SCSS’93), p.2.271/2.285 (1993)
【非特許文献４】Grimble, M. J., Patton, R. J., and Wise, D. A., "The Design of Dynamic Positioning Control System Using Stochastic Optimal Control Thory", Optimal Control Applications and Methods, Vol.1, p.167/202 (1980)
【特許文献１】特公平４−５５９８号公報
【発明の開示】
【発明が解決しようとする課題】
【００１２】
上述のように、推定器に波浪モデルを取り込むためには、波浪モデルの波浪周波数と減衰係数といった波浪パラメータを特定する必要がある。波浪モデルは、海象と船体との相対運動に基づくものであるから、上記波浪パラメータの値は船体によって異なるはずである。
【００１３】
本発明は、精度良く波浪モデルの波浪周波数と波浪減衰係数を同定することができ、該同定された波浪周波数と減衰係数を利用して推定方位と推定角速度を出力する推定器を有する船舶用自動操舵装置を提供することをその目的とする。
【課題を解決するための手段】
【００１４】
かかる目的を達成するために、本発明の請求項１に記載の発明は、設定針路と検出方位に基づいて命令舵角を出力するために、検出方位
2006188095.tif 000007
と命令舵角δcとを入力とし推定方位
2006188095.tif 000008
と推定角速度
2006188095.tif 000009
とを出力する推定器と、推定器の出力を設定針路にフィードバックするフィードバックループと、検出方位から波浪周波数と波浪減衰係数とを同定する波浪同定器と、を有する船舶用自動操舵装置において、前記波浪同定器は、
検出方位信号の時系列データをスペクトルデータに変換するスペクトル演算部と、
スペクトルデータから包絡点群である同定データ
2006188095.tif 000010
（ｆiは周波数、Ｓi−はスペクトル値を表す）を抽出する同定データ抽出部と、
同定データに基づき中心周波数から波浪周波数ｆwnを求める中心周波数算出部と、
同定データのスペクトル値（Ｓi−）と、前記中心周波数算出部によって求められた波浪周波数に基づきその波浪周波数に対応する波浪モデルのスペクトル値（Ｓi）との偏差に基づく評価値を最小にする波浪モデルのスペクトルＳの減衰係数を、波浪減衰係数ζwnとする減衰係数算出部と、を備えることを特徴とする。
【００１５】
請求項２記載の発明は、請求項１記載の前記波浪同定器が、さらに、前記スペクトル演算部で変換されたスペクトルデータを正規化する正規化部を備えることを特徴とする。
【００１６】
請求項３記載の発明は、請求項１または２記載のものにおいて、前記同定データ抽出部が、スペクトルデータの所定の領域において、任意の開始周波数から、その領域において最大スペクトル値をとる周波数に向って、スペクトル値が大きくなるまたは等しくなる包絡点の集まりからなる包絡点群を求め、該最大スペクトル値をとる周波数よりも低域周波数側及び高域周波数側のそれぞれにおいて求めた複数の包絡点群のうちで、包絡点の個数が最大となる包絡点群の中で、開始周波数が最大スペクトル値をとる周波数に最も近い包絡点群を、同定データとして抽出するものであることを特徴とする。
【００１７】
請求項４記載の発明は、請求項１ないし３のいずれか１項に記載のものにおいて、前記中心周波数算出部は、同定データの中でスペクトル値が上部しきい値ＳH以上の値となる同定データ（ｆi，Ｓi−）に対して、
【００１８】
【数１】
2006188095.tif 000011
で波浪周波数ｆwnを求めることを特徴とする。
【００１９】
請求項５記載の発明は、請求項１ないし４のいずれか１項に記載のものにおいて、前記減衰係数算出部は、同定データの中でスペクトル値が下部しきい値ＳL以上の値となる同定データ（ｆi，Ｓi−）に対して、そのスペクトル値とその周波数に対応する波浪モデルのスペクトルのスペクトル値Ｓとの偏差の２乗和Ｊζを最小値とする、
【００２０】
【数２】
2006188095.tif 000012
を満足する波浪モデルのスペクトルの減衰係数ζを、波浪減衰係数ζwnとすることを特徴とする。
【００２１】
請求項６記載の発明は、設定針路と検出方位に基づいて命令舵角を出力するために、検出方位
2006188095.tif 000013
と命令舵角δｃとを入力とし推定方位
2006188095.tif 000014
と推定角速度
2006188095.tif 000015
とを出力する推定器と、推定器の出力を設定針路にフィードバックするフィードバックループと、検出方位から波浪周波数と波浪減衰係数とを同定する波浪同定器と、を有する船舶用自動操舵装置において、前記波浪同定器は、
検出方位信号の時系列データのスペクトルデータから波浪周波数と波浪減衰係数とを同定する波浪同定部と、
時系列データと前記波浪同定部で同定した波浪周波数ｆwn＝ωwn／２πと波浪減衰係数ζwnとを用いて波浪状態を判定し、該波浪状態に適した波浪減衰係数に修正する波浪判定部と、を備えることを特徴とする。
【００２２】
請求項７記載の発明は、請求項６記載のものにおいて、波浪判定部は、波浪から命令舵角までの伝達関数Ｇeである
【００２３】
【数３】
2006188095.tif 000016
に対して、ζwnに前記波浪同定部で同定した波浪減衰係数を代入したときのＧeである
2006188095.tif 000017
と、ζwnに１/√２を代入したときのＧeである
2006188095.tif 000018
のそれぞれの定数を求める定数設定部と、
前記定数設定部で求めた
2006188095.tif 000019
に対して、
【００２４】
【数４】
2006188095.tif 000020
で表されるそれぞれの評価値Ｊon，Ｊoffを求める評価量算出部と、
前記評価値Ｊon，Ｊoffのうちで小さい方に用いたζwnの値を修正した波浪減衰係数
2006188095.tif 000021
として採用し、即ち、
【００２５】
【数５】
2006188095.tif 000022
とする減衰係数決定部と、
を備えることを特徴とする。
【発明の効果】
【００２６】
本発明によれば、波浪モデルが反映されているはずである検出方位信号の時系列データをスペクトルデータに変換する。検出方位信号のスペクトルにおいて、低域周波数は操舵成分であり、高域周波数は波浪成分となる。従って、所定領域のスペクトルデータを求めて、その波浪成分が反映しているスペクトルデータのスペクトル解析を行なうことで、波浪パラメータである波浪周波数及び波浪減衰係数を求めることができる。
【００２７】
請求項２記載の発明によれば、検出方位信号の時系列データから得られたスペクトルデータを正規化することにより、波浪モデルの理論的なスペクトルとの一致度合いを比較することが可能になり、波浪モデルのスペクトルとの一致度合いからそのスペクトルに対応する波浪減衰係数値を決定することができるようになる。
【００２８】
請求項３記載の発明によれば、任意の開始周波数から、所定の領域において最大スペクトル値をとる周波数に向って、スペクトル値が大きくなるまたは等しくなる包絡点の集まりからなる包絡点群を求めることで、凹凸のあるスペクトルデータから包絡的なデータを抽出することができる。また、複数の包絡点群のうちで、包絡点の個数が最大となる包絡点群の中で、開始周波数が最大スペクトル値をとる周波数に最も近い包絡点群を、同定データとすることで、データ数を確保し、且つデータ数／データ周波数範囲を大きくとることができて、同定データの精度を高めることができる。
【００２９】
請求項４記載の発明によれば、波浪モデルの理論上のスペクトルは、中心周波数に対して左右対称にはならないが、中心周波数に近いスペクトル値が高い領域では、その非対称性が小さいので、スペクトル値が上部しきい値以上となる同定データのみを用いて中心周波数を求め、それを波浪周波数とすることにより、非対称性の影響を低減させることができる。
【００３０】
請求項５記載の発明によれば、時系列データを変換したスペクトルデータでは、スペクトルの裾野においてノイズ成分が主要となるため、スペクトル値が下部しきい値ＳL以上の値となる同定データ（ｆi，Ｓi−）に対してのみ処理を行なうことで、ノイズ成分の影響を低減させることができる。また、同定データのスペクトル値と波浪モデルのスペクトル値との偏差の２乗和を最小とする波浪モデルのスペクトルの減衰係数から波浪減衰係数を求めることができる。
【００３１】
請求項６記載の発明によれば、海象が凪状態といった波浪状態においては、波浪同定部における検出方位信号の時系列データのスペクトルのＳ／Ｎ比が低下するので、波浪周波数と波浪減衰係数の同定誤差が増加するが、このような凪状態といった波浪状態を判定し、波浪状態に適した波浪減衰係数に修正することにより、フィードバックループ及び推定器による外乱除去性能の悪化を防ぎ、外乱が増幅されるといった状況を回避することができるようになる。
【００３２】
請求項７記載の発明によれば、波浪から命令舵角までの伝達関数Ｇeにノッチフィルタ特性を持たせた場合と、ノッチフィルタ特性を持たせない場合（ζwn＝１/√２）のそれぞれの伝達関数
2006188095.tif 000023
を求め、それぞれの伝達関数
2006188095.tif 000024
に対して、それぞれ検出方位信号との２乗和を評価値とし、その評価量が小さい方に用いたζwnの値を波浪減衰係数
2006188095.tif 000025
として採用することにより、波浪から命令舵角までの伝達関数Ｇeが大きくなり外乱が増幅されることを回避することができる。
【発明を実施するための最良の形態】
【００３３】
以下、図面を用いて本発明の実施の形態を説明する。
１  全体概略構成
本発明で対象とする船舶用自動操舵装置１０の制御系は、図１に示したように、設定針路と検出方位とに基づいて命令舵角δｃを出力する自動操舵装置の制御系において、閉ループを安定化するフィードバックループ１５と、検出方位に含まれる波浪成分を除去する推定器２０とから構成され、推定器２０は、検出方位
2006188095.tif 000026
と命令舵角δcとを入力とし、検出方位の推定値である推定方位
2006188095.tif 000027
と旋回角速度の推定値である推定角速度
2006188095.tif 000028
とを出力し、フィードバックループ１５は推定器２０の出力
2006188095.tif 000029
を設定針路にフィードバックする。このとき、推定器２０は、波浪同定器２２によって同定される波浪周波数ωwnと波浪減衰係数ζwnを利用して、推定方位と推定角速度を求める。
【００３４】
推定器２０の構成の理解を図るために、まず特願２００４−２２６６７０で本発明者らが提案する波浪モデル無しの推定器（即ち、２次推定器）と、波浪モデル有りの推定器（即ち、４次推定器）について概略説明する。
【００３５】
２次推定器は、図２に示すように、ノミナル値の船体モデル（野本モデルを基に考案したもの）と、検出方位と推定方位との差
【００３６】
【数６】
2006188095.tif 000030
を用いると、ここで
【００３７】
【数７】
2006188095.tif 000031
とし、
【００３８】
【数８】
2006188095.tif 000032
を満足する。ここでｋ1，ｋ2は推定ゲインを示す。上式と実際値の船体モデル関係式
【００３９】
【数９】
2006188095.tif 000033
とを整理することによって、２次推定器の伝達関数は、
【００４０】
【数１０】
2006188095.tif 000034
として得られる。ここで、
2006188095.tif 000035
は、それぞれ方位と波浪とから推定方位までの伝達関数、
2006188095.tif 000036
はそれぞれ方位と波浪とから推定角速度までの伝達関数で、
2006188095.tif 000037
は２次推定器を意味し、
【００４１】
【数１１】
2006188095.tif 000038
をそれぞれ示す。λe2は２次推定器の特性多項式で、
【００４２】
【数１２】
2006188095.tif 000039
を示し、Δ2，Δ1は船体パラメータのノミナル値のパラメータ不確かさで
【００４３】
【数１３】
2006188095.tif 000040
をそれぞれ示す。
【００４４】
上記パラメータ不確かさを表すΔ2，Δ1による伝達関数への影響を調べ、閉ループ安定性が劣化しないように上記(１）式の２次推定器の特性多項式における設計パラメータである減衰係数ζeと固有周波数ωe[rad/s]を決定する。
【００４５】
検討の結果、固有周波数ωeについては、ωeを高くすると、
2006188095.tif 000041
の伝達誤差が低減することが分かるので、操舵系の固有周波数ωnと推定係数ρ＞１とを用いて、固有周波数ωeは、
【００４６】
【数１４】
2006188095.tif 000042
と設定する。ここで、ρの設定は閉ループ安定性を満足する適切な値に選択される。また、減衰係数ζeは１／√２とすると、Δ１による伝達関数
2006188095.tif 000043
のパラメータ不確かさ感度を最小にすることができる。
【００４７】
次に、波浪モデル有りの推定器を考える。波浪モデル有りの推定器（即ち、４次推定器）は、２次推定器に波浪モデルを組み込むことで、検出方位に含まれる波浪成分を推定波浪として分離抽出することができる。そのため推定方位と推定角速度とに含まれる波浪成分が除去され無駄舵が防止できる。しかし波浪モデルの導入が推定誤差を増加させる状況は避けなければならない。よって２次推定器の推定誤差を保持しつつ、かつ効果的な外乱除去性能を実現する４次推定器を設計する必要がある。４次推定器は、例えば、図３に示すように、ノミナル値の船体モデル（野本モデルを基に考案したもの）と、検出方位
2006188095.tif 000044
から推定方位
2006188095.tif 000045
と推定波浪
2006188095.tif 000046
とを差し引いた
【００４８】
【数１５】
2006188095.tif 000047
とを用いて、
【００４９】
【数１６】
2006188095.tif 000048
と表すことができる。ここで、ｋi，i＝１〜４は推定ゲインを、
2006188095.tif 000049
は状態量を、それぞれ示す。波浪モデルの固有周波数ωwn、減衰係数ζwnは、波浪同定器２２によって、検出方位
2006188095.tif 000050
から同定される。波浪モデルの未知入力νは白色ノイズを仮定しているので平均値ゼロになる。したがってψwおよびψwに対応する
2006188095.tif 000051
もまた平均値ゼロになるから、上式でνの入力の必要はない。上式と前記δcについての船体モデル関係式とを整理すると、４次推定器の伝達関数は、
【００５０】
【数１７】
2006188095.tif 000052
として得られる。ここで、
2006188095.tif 000053
はそれぞれ方位と波浪とから推定波浪までの伝達関数で、
2006188095.tif 000054
は４次推定器を意味し、
【００５１】
【数１８】
2006188095.tif 000055
をそれぞれ示す。λe4は４次推定器の特性多項式で
【００５２】
【数１９】
2006188095.tif 000056
を示す。λw，ｃewは波浪モデルのそれぞれ多項式とゲイン多項式とで
【００５３】
【数２０】
2006188095.tif 000057
を示す。
【００５４】
前記４次推定器の特性多項式λe4を２つの２次標準系を用いて
【００５５】
【数２１】
2006188095.tif 000058
に定める。ここでλeは船体モデルの状態推定に、λewは波浪モデルのそれにそれぞれ対応すると仮定する。ζe，ωeは２次推定器で設計したものである。ζew，ωewは４次推定器の設計パラメータでそれぞれ減衰係数と固有周波数[rad/s]とを示す。よって推定ゲインは
【００５６】
【数２２】
2006188095.tif 000059
になる。
【００５７】
検討の結果、波浪除去に対応する伝達関数
2006188095.tif 000060
に波浪周波数ωwnでのノッチフィルタ特性を持たせると共に、推定誤差の悪化を防ぐために、ωew≧ωeを満足するように、波浪モデルが無い（λew＝１）としたときの推定器の船体パラメータに対する感度が、波浪モデルが有るときの推定器の船体パラメータに対する感度とほぼ等しくなるときの、ωe＝ωwn（ωwn：波浪モデルの固有周波数である波浪周波数を表す）を満足するときのωeの値をωx、ωe＝ωwnを満足しないときのωeの値をωyとしたときに、波浪周波数ωwnがωxよりも大きいときには、固有周波数ωeをωyに設定し、固有周波数ωewを波浪周波数ωwnと等しく設定し、波浪周波数ωwnがωx以下のときには、固有周波数ωeをωxに設定し、波浪周波数ωwnを修正の波浪周波数ωwn＊＝ωxに修正し、固有周波数ωewを修正の波浪周波数ωwn＊と等しく設定するとよいことが分かっている。
【００５８】
推定器の構成は、以上に概略説明したものに限るものではなく、任意の構成とすることができるが、波浪モデルを導入する推定器とする場合には、波浪モデルの波浪周波数ωwn及び波浪減衰係数ζwnを特定する必要がある。
【００５９】
上記波浪モデルの波浪周波数ωwn及び波浪減衰係数ζwnを同定するための波浪同定器２２の構成を、図４に示す。この実施形態では、バッチ方式で波浪パラメータを同定するもので、波浪同定器２２は、データストア部２４、スペクトル演算部２６、波浪同定部２８及び波浪判別部３０とからなる。
【００６０】
データストア部２４は、Ａ／Ｄ変換された検出方位信号を時系列データとして蓄積するものであり、スペクトル演算部２６は、検出方位信号の時系列データをスペクトルデータに変換するものであり、波浪同定部２８はスペクトル演算部２６からのスペクトルデータから波浪パラメータである中心周波数と減衰係数とを演算するものであり、波浪判別部３０は、検出方位信号の時系列データと波浪パラメータとから波浪状態に適する波浪パラメータの最適値に修正するものである。
【００６１】
波浪同定部２８は、さらに正規化部２８Ａ、同定データ抽出部２８Ｂ、中心周波数算出部２８Ｃ及び減衰係数算出部２８Ｄからなり、波浪判別部３０は、さらに定数設定部３０Ａ、評価量算出部３０Ｂ、減衰係数決定部３０Ｃ及び減衰係数算出部３０Ｄからなる。これら波浪同定部２８及び波浪判別部３０は、コンピュータソフトウェアで構成することができる。この波浪同定部２８と波浪判別部３０との詳細構成及び作用を以下説明する。
【００６２】
２  波浪同定部
２．１  正規化部
波浪モデルの伝達関数Ｇwは
【００６３】
【数２３】
2006188095.tif 000061
を用いる。波浪モデルのスペクトルＳは伝達関数のゲインの２乗であるから、
【００６４】
【数２４】
2006188095.tif 000062
になる。ここでωは周波数を示す。（４）式より、ω＝ωwnのときにＳ（ω）は最大値をとり、その大きさは１になる。
【００６５】
一方、検出方位
2006188095.tif 000063
に含まれる波浪成分
2006188095.tif 000064
の時系列データは
【００６６】
【数２５】
2006188095.tif 000065
で与えられる。ここでσwは波浪の強さを表す定数を、νは平均値ゼロのランダムノイズを、ｔは時間をそれぞれ示す。検出方位の時系列データから変換したスペクトルＳhの最大値は
2006188095.tif 000066
であるから、Ｓhを最大値
2006188095.tif 000067
で割った
2006188095.tif 000068
は波浪モデルのスペクトルＳと等価になる。即ち、
【００６７】
【数２６】
2006188095.tif 000069
になる。以下、時系列データから求めた最大値１に正規化されたスペクトル
2006188095.tif 000070
を実スペクトルと呼ぶことにする。
【００６８】
正規化部２８Ａは、スペクトル演算部２６で検出方位
2006188095.tif 000071
の時系列データを高速フーリエ変換（ＦＦＴ）を用いて変換して得られたスペクトルを正規化して実スペクトルに関するデータを得る。尚、ＦＦＴ変換した実スペクトルは帯域周波数を網羅していないため、滑らかな山形ではなく凹凸状となっている。凹凸状によるパラメータ同定誤差を低減するためには、適宜、平滑処理を行なうと好ましい。エイリアシング対策として、時系列データの移動平均をとることも可能である。
【００６９】
２．２  同定データ抽出部
同定データ抽出部２８Ｂは、正規化部２８Ａによって正規化された実スペクトルから同定に用いる同定データ
2006188095.tif 000072
を抽出する。
【００７０】
図５の実スペクトル例を参照して説明すると、実スペクトルにおいて低域周波数帯域は操舵成分、高域周波数帯域は波浪成分となるので、高域の波浪周波数領域において、最大スペクトル値をとる周波数を境にして低域側と高域側とに別け、それぞれから包絡的なデータを抽出する。
【００７１】
即ち、任意の開始周波数から最大スペクトル値をとる周波数に向かって、スペクトル値がそれよりも大きくなるまたは等しくなる包絡点
2006188095.tif 000073
の集まり(包絡点群)を求める。図において×と○とがそれぞれ包絡点群を示す。
【００７２】
次に、複数の包絡点群のうちで、包絡点の個数が最大となる包絡点群の中で開始周波数が最大値の周波数に最も近接した包絡点群を同定データにする。包絡点の個数を最大とする包絡点群を選ぶ理由はデータ数を多くするためで、開始周波数が最大スペクトル値の周波数に近い包絡点群を選ぶ理由は（データ数÷データ周波数範囲）を大きくするためである。図において低域側と高域側においてそれぞれ×と○との個数は同じであるが、○の開始周波数は×のそれより最大値の周波数に近い。よって○の包絡点群が同定データになる。
【００７３】
２．３  中心周波数算出部
中心周波数算出部２８Ｃは、同定データから中心周波数を算出する。
【００７４】
式（４）の波浪モデルのスペクトルを図６に示す。ある減衰係数を持つスペクトルにおいて、あるスペクトル値をとる周波数範囲を［ω1、ω2］（ω2＞ω1）とする。同図よりスペクトル形状は減衰係数により決定され、ω／ωwn＝１に対してスペクトルは左右対称にならないことが分かる。即ち、
【００７５】
【数２７】
2006188095.tif 000074
である。この非対称性は減衰係数に比例する。
【００７６】
図７はスペクトルに対する非対称度合を求めたものである。同図より減衰係数が０．３のとき非対称度合を３％以内にするためには、上部しきい値ＳH＝０．８とするとよいことが分かる。実スペクトルが林立する場合どのしきい値以上の同定データを対象にすれば良いかの基準として、この上部しきい値ＳHを用いる。最大スペクトル値をとる周波数がその分解能の関係より必ずしも中心周波数にならない場合が生じるので、中心周波数ｆwn＝ωwn／２πは上部しきい値ＳH以上のスペクトル値
2006188095.tif 000075
とその周波数ｆi [Hz]とから加重平均により
【００７７】
【数２８】
2006188095.tif 000076
で求める。ここでΣiは上部しきい値ＳH 以上の同定データの個数で総和を取ることを意味する。
【００７８】
２．４  減衰係数算出部
減衰係数算出部２８Ｄは、同定データから減衰係数を算出する。
【００７９】
図６に示したように、スペクトルは、そのゼロ付近の裾野が急激に広がっている。実スペクトルの裾野領域はノイズ成分が主要になるために同定範囲から除くことにする。その有効領域と無効領域を適当に分ける必要があるが、例えば、スペクトルの有効領域を上部の９５％にして、無効領域を下部の５％にする。その下部しきい値ＳL＝０．０５における周波数範囲［ω1、ω2］は式（４）を用いて
【００８０】
【数２９】
2006188095.tif 000077
の２根
2006188095.tif 000078
として得られる。ＳL＝０．０５による周波数範囲［ω1、ω2］にある同定データを用いて減衰係数の計算を行う。
【００８１】
減衰係数ζwnは同定データのスペクトル値とその周波数に対応する波浪モデルのスペクトル値との偏差の２乗和を用いた評価関数を最小にする波浪モデルのスペクトルの減衰係数ζによって与える。即ち、
【００８２】
【数３０】
2006188095.tif 000079
で求める。ここでＳは波浪モデルのスペクトル値を、Σiは同定範囲にある同定データの個数で総和を取ることを意味する。上式は、既にωwn＝２πｆwnが得られているからζwnの関数になり、１変数による最小値問題に帰着する。その解法は黄金分割を用いると良い。
【００８３】
以上による波浪同定部２８によって検出方位信号の時系列データのスペクトルデータから波浪パラメータを同定することができる。
【００８４】
３  波浪判別部
３．１  伝達関数の導出
次に、波浪から命令舵角までの伝達関数を検討する。図１及び図３に示した操舵系のブロック図を元にし、簡単化のため操舵機は省き、設定針路（参照針路）、パラメータ不確かさ（船体パラメータのノミナル値と実際値との誤差）及び舵角オフセットなど等はゼロとする。同図において、推定器２０は波浪モデルが組み込まれた４次推定器であり、Ｐは船体モデルを、ＫP，ＫDはフィードバックループ１５のフィードバックゲインでそれぞれ比例ゲインと微分ゲインとを表している。
【００８５】
図３に示した４次推定器の伝達関数は、１項で既述の通りであるが、
【００８６】
【数３１】
2006188095.tif 000080
として得られる。（２）式と上式とを整理すると、
2006188095.tif 000081
は、
【００８７】
【数３２】
2006188095.tif 000082
として表すことができる。ここで
【００８８】
【数３３】
2006188095.tif 000083
を表し、Ｋsn、Ｔsn，Ｔs3nは船体パラメータのノミナル値でそれぞれ旋回力ゲインと２つの船体時定数を、ωwn，ζwnは波浪パラメータのノミナル値でそれぞれ中心周波数と減衰周波数を、それぞれ示す。
また、命令舵角の関連式は
【００８９】
【数３４】
2006188095.tif 000084
になることから、波浪から命令舵角までの伝達関数は
【００９０】
【数３５】
2006188095.tif 000085
として得られる。ここで、Ｇp、Ｇeはそれぞれ操舵系関連と推定器関連との伝達関数を示す。
【００９１】
伝達関数Ｇｐは２次操舵系の特性多項式をもち
【数３６】
2006188095.tif 000086
になる。ここでζn、ωnは操舵系のそれぞれ減衰係数と固有周波数とを示す。Ｇpのカットオフ周波数ωcは
【００９２】
【数３７】
2006188095.tif 000087
で与えられる。波浪周波数はωcより高域側であると仮定できるから、Ｇpのωcより高域周波数のゲインは、
【００９３】
【数３８】
2006188095.tif 000088
の一定値になる。ただしＫD＞０とする。
一方、伝達関数Ｇeは４次推定器の特性多項式とフィードバックゲインとからなり、
【００９４】
【数３９】
2006188095.tif 000089
になる。
【００９５】
３．２  伝達関数の特性
伝達関数Ｇp、Ｇeのゲイン特性の数値例を図８に示す。同図はＫsn＝０．０４７，Ｔsn＝６１，Ｔs3n＝８．２，ＫP＝１，ζn＝０．９，ωn＝０．０２４６，ＫD＝４３．８，ωe＝ρωn，ρ＝３，ωew＝ωwn，ζe＝ζew＝１/√２とｋi，ｉ＝１〜４を用いて、減衰係数ζwnと波浪周波数ωwnとを変化させる。
【００９６】
同図よりＧPはω＞ωeでほぼ一定値になるので、波浪パラメータによる影響はＧeを調べればよいことが分かる。波浪周波数ωwnが小さくなると、Ｇeのゲインがω＞ωeで１以上になり外乱が増幅される領域を生じる。即ち、波浪状態によっては、外乱成分を増幅させるために、外乱除去性能の効果が期待できない。この状況は操舵周波数ωnと波浪周波数ωwnとが近くなる小型船で発生しやすい。
【００９７】
３．３  波浪状態
海象が凪状態のとき、波浪同定部において波浪スペクトルのＳＮ比が低下するために、波浪パラメータの同定誤差が増加する。その結果凪状態にも拘わらず、上記外乱除去性能の効果が期待できない状況が生じる場合がある。この状況を避けるために、凪またはそれと同等の状態のときは、波浪モデルによるノッチフィルタ特性が働かないようにする。即ち、ζwn＝ζew＝１/√２とする。このとき式（２８）で波浪モデルに関係する項は
【００９８】
【数４０】
2006188095.tif 000090
になり、ゲインはωwnに無関係にω＞ωeで−２０[dB/dec]になる（図８のζwn＝１/√２）。このため、ζwn＝１/√２のＧeのゲインはζwn＝０．１、ωwn＝０．１５の場合のような外乱増幅領域を持たないので、ノイズ成分を増幅することはない。
【００９９】
この性質を利用し、波浪同定部２８で同定された波浪パラメータと上記伝達関数Ｇeを用いて、波浪状態の判別を行なうことができる。
【０１００】
以下、波浪判定部３０の作用を図４を参照して説明する。
【０１０１】
定数設定部３０Ａは、式（２８）のζwnに同定値（ノッチフィルタ有り）と１/√２（ノッチフィルタ無し）とを代入し、Ｇeのパラメータ、即ち定数をそれぞれ設定する。具体的には、
【０１０２】
【数４１】
2006188095.tif 000091
のそれぞれのパラメータ（定数）を決定する。
【０１０３】
次いで、評価量算出部３０Ｂは、波浪状態の評価を行なうための評価量の算出を行なう。波浪状態の評価量は、ノッチフィルタの有りの場合と無しの場合のそれぞれのＧeと、検出方位時系列データを用いた命令舵角相当の２乗和を用いる。初期値はゼロとする。
【０１０４】
【数４２】
2006188095.tif 000092
ここでＪon，Ｊoffはノッチフィルタのそれぞれ有りと無しの場合の評価量を示す。尚、総和をとるデータ数ｎは波浪の周期の１０倍程度にわたる期間中の時系列データ数とするとよい。
【０１０５】
次いで、減衰係数決定部３０Ｃは、２つの評価量の内で小さい方がその時の波浪状態に対する外乱除去効果のより高い設定であるとして、小さい方に用いた減衰係数を波浪パラメータとして決定する。中心周波数ωwnは同定値のままである。即ち波浪係数は
【０１０６】
【数４３】
2006188095.tif 000093
によって定める。
【０１０７】
このように、波浪判定部３０は、波浪同定部２８によって同定された波浪減衰係数に対して修正を行なうことにより、波浪状態によっては、波浪から命令舵角までの伝達関数Ｇeが大きくなり、外乱が増幅されることを回避することができる。
【図面の簡単な説明】
【０１０８】
【図１】本発明による船舶用自動操舵装置の全体構成を表すブロック図である。
【図２】２次推定器の構成を表すブロック図である。
【図３】４次推定器の構成を表すブロック図である。
【図４】波浪同定器の構成を表すブロック図である。
【図５】同定データの抽出を表す実スペクトル例を表すグラフである。
【図６】波浪モデルのスペクトルを表すグラフである。
【図７】図６のスペクトルに対して非対称度合を求めたグラフである。
【図８】伝達関数Ｇn、Ｇeのゲイン特性を表すグラフである。
【図９】船舶用自動操舵装置の全体構成を表すブロック図である。
【符号の説明】
【０１０９】
１０  船舶用自動操舵装置
１５  フィードバックループ
２０  推定器
２２  波浪同定器
２６  スペクトル演算部
２８  波浪同定部
２８Ａ  正規化部
２８Ｂ  同定データ抽出部
２８Ｃ  中心周波数算出部
２８Ｄ  減衰係数算出部
３０  波浪判別部
３０Ａ  定数設定部
３０Ｂ  評価量算出部
３０Ｃ  減衰係数決定部
【図１】
2006188095.tif 000112
【図２】
2006188095.tif 000113
【図３】
2006188095.tif 000114
【図４】
2006188095.tif 000115
【図５】
2006188095.tif 000116
【図６】
2006188095.tif 000117
【図７】
2006188095.tif 000118
【図８】
2006188095.tif 000119
【図９】
2006188095.tif 000120