JP3628479B2 - 目標運動解析装置および目標運動解析方法 - Google Patents

目標運動解析装置および目標運動解析方法 Download PDF

Info

Publication number
JP3628479B2
JP3628479B2 JP14187097A JP14187097A JP3628479B2 JP 3628479 B2 JP3628479 B2 JP 3628479B2 JP 14187097 A JP14187097 A JP 14187097A JP 14187097 A JP14187097 A JP 14187097A JP 3628479 B2 JP3628479 B2 JP 3628479B2
Authority
JP
Japan
Prior art keywords
analysis
target motion
individual
individuals
target
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.)
Expired - Lifetime
Application number
JP14187097A
Other languages
English (en)
Other versions
JPH10332806A (ja
Inventor
昌浩 松本
和夫 佐藤
凖 原
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.)
Hitachi Ltd
Original Assignee
Hitachi 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 Hitachi Ltd filed Critical Hitachi Ltd
Priority to JP14187097A priority Critical patent/JP3628479B2/ja
Publication of JPH10332806A publication Critical patent/JPH10332806A/ja
Application granted granted Critical
Publication of JP3628479B2 publication Critical patent/JP3628479B2/ja
Anticipated expiration legal-status Critical
Expired - Lifetime legal-status Critical Current

Links

Images

Landscapes

  • Complex Calculations (AREA)
  • Measurement Of Velocity Or Position Using Acoustic Or Ultrasonic Waves (AREA)

Description

【0001】
【発明の属する技術分野】
本発明は、水中および水上などを運動する目標から放射される音波を音響センサにより観測し、そこから得られる到来音波の方向および周波数を用いて目標の運動諸元を解析する目標運動解析の技術に関するものである。
【0002】
【従来の技術】
図11に従来の目標運動解析の原理を示す。
【0003】
図中において、観測船は、音響センサを使って目標102が放射する音波を受信することにより、到来音波の方位および周波数を観測する。そして、観測された方位および周波数、観測時刻および観測装置の時系列の状態情報に基づき目標102の状態ベクトルx(距離r、方位By、速力Mt、針路Ct、放射音波の周波数f0)を解析する。この解析は、状態ベクトルxから推定される各観測時刻の解析方位と観測方位の差分二乗を観測方位誤差の分散で規格化したものと、状態ベクトルxから推定される解析周波数と観測周波数との差分二乗を観測周波数誤差の分散で規格化したものを時系列的に合計した、数1の評価関数Φ(x)が最小となるような状態ベクトルを解析する非線形最適化問題である。
【0004】
【数1】
Figure 0003628479
【0005】
このような非線形最適化問題を解析する方法としては、これまで、カルマンフィルタ、最小二乗法によるものなどが知られている。カルマンフィルタを用いた目標運動解析についてはIEEE会報1985年Vol AC−30 No.10「A Stochastic Analysis of a Modified Gain Extended Kalman Filter with Application to Estimation with Bearing Only Measurements」に記載されており、最小二乗法を用いた目標運動解析についてはIEEE会報1984年VolAC−29No9「Fundamental Properties and Performance of Conventional Bearing− only Target Motion Analysis」に記載されている。
【0006】
【発明が解決しようとする課題】
しかし、これらの技術に係るアルゴリズムは、音響センサのランダム誤差のみを考慮しており、音響センサ固有のバイアス誤差による影響については考慮していない。
【0007】
図12に、音響センサの誤差発生要因を示す。
【0008】
図12aに示すようにランダム誤差204は、目標との相対距離に依存して変化する特性を持ち、音響センサ202の音波受信感度および指向性ビーム幅203をパラメータとして正規分布状に発生するのに対して、図12b、cに示すように、固有バイアス誤差210は、海水と隔壁内媒質との相違や船体の材質や形状による音波の屈折によって発生し、音波入射角をパラメータとして増減する特性211を有する。
【0009】
従来、このような音響センサ205の固有バイアス誤差に対しては次のように対処していた。すなわち、距離および方位が明確な目標をいくつか用意し、音響センサによる観測方位と自明な目標方位との差分を求めてテーブル化しておく。そして、実際の観測時には、このテーブルを用いて観測値などを補正することにより、音響センサ205の固有バイアス誤差に対処していた。
【0010】
しかしながら、固有バイアス誤差は、同型船でも音響センサ205の据え付け状態によって異なるうえ、音響センサ205の経年変化や環境条件の変化によっても特性211が異なってくる。このため、前述した対処法では、補正が不充分となる場合があり、この場合に従来のカルマンフィルタ、最小二乗法による解析を、たとえば図13の状況に対して適用すると、図14に示した解析例に示すように、局所最適解への落ち込みが生じ、目標運動解析精度が劣化することがあった。図14において、a)が距離rの解析例、b)が速力Mtの解析例、c)が進路ctの解析例であり、各図において、実線が真の値、波線が従来のカルマンフィルタによる解析結果、白抜きの四角が最小2乗法による解析結果を示している。各図中において、丸で囲んだ箇所などにおいて、局所最適解への落ち込みが生じ、解析結果が真の値からかけはずれてしまっていることが分かる。
【0011】
そこで、本発明は、固有バイアス誤差に起因する目標運動解析精度の劣化を軽減することを目的とする。
【0012】
【課題を解決するための手段】
前記目的達成のために、本発明は、たとえば、運動する目標から放射される音波を音響センサにより観測した観測値に基づいて、目標の運動諸元を解析する目標運動解析装置であって、
音響センサの固有バイアス誤差を特定する情報を遺伝子とする個体の集合に遺伝的操作を行い、順次、次世代の個体の集合を生成する遺伝的操作部と、
生成された各世代の個体の集合に含まれる各個体のそれぞれについて、当該個体に対応する固有バイアス誤差を考慮した目標の運動諸元の解析を行うと共に、当該解析の結果と観測値との適合度を求める非線形解析処理部とを備え、
前記遺伝的操作部は、次世代の個体の集合の生成を、現世代の個体のうち、前記非線形解析部において、より高い前記適合度が当該個体について求められた個体の特徴が、より高い確率で次世代に遺伝されるように、前記遺伝的操作を行うことを特徴とする目標運動解析装置を提供する。
【0013】
このようにして、固有バイアス誤差特性を遺伝的操作処理部において推定しつつ、各個体について非線形解析処理部において目標の運動諸元の解析を実施することにより、各目標に対する解析値が最適解となるような固有バイアス誤差特性の推定が行われるとともに、固有バイアス誤差を補正した観測情報に基づく目標運動諸元の解析が可能となる。
【0014】
【発明の実施の形態】
以下、本発明に係る目標運動解析装置の実施形態を、艦船に装備され、目標運動諸元の解析を行う目標運動解析システムへの適用を例にとり説明する。
【0015】
まず、第1の実施形態について説明する。
【0016】
図1に、本第1実施形態に係る目標運動解析システムの概略構成を示す。
【0017】
図中において、観測船511は目標510から到来する音波を音響センサ501で受信し、音波の観測周波数と、音波の到来方位であるところの観測方位を、受信時刻および受信位置、受信時の観測船511の針路である観測針路、受信時の観測船の速力である観測速力と共に、目標運動解析装置505のセンサデータベース508に入力時系列データとして蓄積する。このとき、受信位置、観測針路、観測速力は、位置演算器504が、観測船に取り付けられた速力検出器502と針路検出器503から周期的に算出する。また、受信位置は、初めて音波を受信した位置を基準とする相対位置座標として蓄積する。
【0018】
目標運動解析装置505において、非線形解析処理部507は、センサデータベース508に蓄積された時系列データを参照し目標の運動諸元を解析する。また、この際、遺伝的処理部506において音響センサ501の固有バイアス誤差を推定する。
【0019】
目標運動解析装置505で求められた運動諸元の解析値および固有バイアス誤差の推定値は、センサデータベース508の各種情報と共に解析データ表示器509に表示される。
【0020】
なお、端末装置509を介して、オペレータは、端末装置509に表示されるセンサデータベース508の入力時系列データの連続性を監視することにより発見した入力時系列データ中の異常データの除去と、非線形解析処理部507による再計算を目標運動解析装置505に指示することができる。
【0021】
ここで、このような目標運動解析システムが行う目標運動解析の原理を示しておく。
【0022】
本実施形態において行う目標運動解析は、固有バイアス誤差特性403が音響センサ402で観測している複数の目標に対して共通であることに着目したものである。
【0023】
すなわち、図2に示すように、遺伝的処理部404で、係数a,b,cををパラメータとして持つ正弦波曲線405を固有バイアス誤差特性403のモデル関数とし、この係数a,b,cを数値遺伝子として、この係数の組合せを個体と定義して遺伝的操作を実施する。遺伝的操作は、非線形解析処理部406で固有バイアス誤差を考慮して各目標対応に算出される解析残差二乗和の合計値Rの逆数を各個体の適応度とすることにより、解析残差二乗和の合計値が小さな個体ほど適応度が大きくなり、固有バイアス誤差推定値の最適解として選択されるように行う。ここで、解析残差二乗和とは、数2に示すように、各目標毎にセンサデータべースに蓄積された音響センサによる入力時系列データ中の観測方位および観測周波数と、運動諸元の解析の結果求まる時系列の解析方位および解析周波数との差分二乗を各々の観測情報の分散で規格化し、合算したものである。
【0024】
【数2】
Figure 0003628479
【0025】
このようにして、固有バイアス誤差特性403を遺伝的処理部404において推定しつつ、非線形解析処理部406において各目標の運動諸元の解析を実施することにより、各目標に対する解析値が最適解となるような固有バイアス誤差特性403の推定が行われるとともに、固有バイアス誤差を補正した観測情報に基づく目標運動諸元の解析が可能となる。
【0026】
以下、本実施形態に係る目標運動解析システムの詳細について説明する。
【0027】
図3に、目標運動解析装置505の詳細な構成を示す。
【0028】
図中、目標運動解析装置505は固有バイアス誤差の推定を行う遺伝的処理部506、目標運動諸元の解析を行う非線形解析処理部507、これらの両処理部の数値を変換する数値変換処理部607、観測データが蓄積されているセンサデータベース508から構成される。
【0029】
遺伝的処理部506では、固有バイアス誤差のモデル曲線の係数の組合せを個体とする個体の初期集合P0を生成603し、選択604、交叉605、突然変異606の遺伝的操作を繰り返すことにより、最大の適応度を持つ個体を生成する。なお、この遺伝的操作により最大の適応度を持つ個体が推定される基本原理については次の文献に詳しく述べられている。
【0030】
Holland, John (1975). 「Adaptation in Natural and Artificial Systems」, MIT press 1992、Goldberg, David (1987). 「Genetic Algorithm in Search, Optimization,and Machine Learning. Reading, Mass.」: Addison−Wesley , Davis,Lawrence (1990). 「Handbook of Genetic Algorithms」, Van Nostrand (邦訳:遺伝アルゴリズムハンドブック 嘉数(かかず)他 森北出版)
遺伝的処理部506における遺伝的操作の過程において生成される各世代の個体集合Pnは数値変換処理部607により、固有バイアス誤差のモデル曲線の係数の解集合Zに変換する。そして、非線形解析処理部507において、最小二乗法または適応型フィルタによる非線形解析609において、解集合Zの各解について各々、その解を観測方位の補正値として使用し、センサデータベース508に蓄積された入力時系列データを参照しながら、解析対象の各目標についての運動諸元の解析を行う。
【0031】
そして、個別目標毎に非線形解析609の内部で計算される解析残差二乗和を評価610で合計し、その合計値を評価値Rとする。
【0032】
さて、評価値Rは、各個体に対応する固有バイアス誤差のモデル曲線の係数の妥当性を表す。そこで、このようにして算出した評価値Rの評価値集合Rを、前述した遺伝的処理部506の個体の選択604操作に使用し、次世代の個体集合Pnを生成する。
【0033】
以上のような、遺伝的操作を終了判定611の条件を満足するまで繰り返し、それまでに得られた、評価値を最小とする係数の組合せを最適解ziとし、最適解ziを用いて求めた運動諸元の解析値を解析値として出力する。これによって、固有バイアス誤差の曲線係数に関する最適解の推定とともに目標の運動諸元に関する最適解の解析が可能となる。
【0034】
終了判定611では、遺伝処理部506で生成された世代の数(世代交代数)が所定の最大数を越えた場合に処理を終了させる。また世代交代数が最大値に達しない場合でも、今回の世代に対する評価値集合R中の評価値の平均値もしくは評価値の最小値が、前回の世代に対する評価値集合R中の評価値の平均値もしくは評価値の最小値より大きくなった場合に世代交代により適応度が劣化したとして処理を終了させるようにする。
【0035】
また、本実施形態では、今回の目標運動解析で得られた最適解ziを、次回解析時の遺伝的処理部506の初期集合の個体の一つとして継承する処理612を行い次回解析の収束が加速されるようにしている。
【0036】
次に、遺伝的処理部506の行う処理の詳細について説明する。
【0037】
図7は遺伝的処理部506の処理の詳細を示したものである。
【0038】
図示するように、この処理では、まず、初期集合生成603において固有バイアス誤差のモデル曲線の係数に該当する数値遺伝子a,b,cを一様乱数により生成し、その組合せである個体をm種類生成して、個体集合P0を形成する。但し、前回の解析時の最適解Ziを継承しているときは、m−1種類の個体を生成し、これに継承した最適解Ziに対応する個体を加えて個体集合P0とし、これを現世代の個体集合とする。
【0039】
次ぎに、現世代に対する非線形解析処理部507の処理が終了したら、選択604で現世代の各個体の適応度fitの大きさに応じた確立でランダムに、現世代の個体集合中より次世代の個体が抽出されるようルーレットルール705と乱数706により世代交代を実施する。適応度fitは、その個体に対して非線形解析処理部507で求められた評価値Rの逆数を用いる。
【0040】
次に、交叉703では個体の適応度fitを標準偏差とする正規分布に従い数値遺伝子の組み替えを実施する。ここで、正規分布を適用することにより、交叉対象の個体の適応度fitが大きいものほど数値遺伝子a,b,cの変化が少なく、適応度fitが小さいものほど数値遺伝子a,b,cの変化が大きくなり、一様分布を適用した場合に比べて解探索の効率化が図れる。たとえば、正規分布交叉707に示すように、数値遺伝子a1および適応度fit1をもつ個体1と、数値遺伝子a2および適応度fit2をもつ個体2の交叉は、数3および数4に従って行うようにする。
【0041】
【数3】
Figure 0003628479
【0042】
【数4】
Figure 0003628479
【0043】
最後に突然変異704では、個体pの数値遺伝子a,b,cの任意の一つを一様乱数により変化させることにより、解探索の多様性を確保し、局所最適解への落ち込みを防止する。
【0044】
そして、この突然変異704が終了した個体集合を、新たな現世代の個体集合とする。
【0045】
以上一連の、遺伝的操作を繰り返すことにより固有バイアス誤差の曲線係数の最適値を推定する。
【0046】
次に、非線形解析処理部507中の、非線形解析609の詳細について説明する。
【0047】
まず、非線形解析609において、Marquardt法を利用した最小二乗法による非線形解析を行う場合について説明する。
【0048】
この場合の処理の手順を、図5a)に示す。
【0049】
この最小二乗法による非線形解析801では、センサデータベース508に蓄積されている全観測時刻のデータを使用し、観測データ誤差が微少であることを仮定した疑似線形解析805で初期解析値を算出し、この解析値をMarquardt法817により反復修正して最終的な解析値を算出する。
【0050】
疑似線形解析805では目標運動諸元の状態ベクトルを数5のように定義し、線形解析準備803で係数行列Fを数6、定数行列Gを数7により計算しておき、線形解析804の数8により初期解析値を算出する。ただし、条件802のとおり、前回非線形解析時における最小二乗法801による解析値が「解有り」として存在する場合は、この解を今回の初期解析値とし、疑似線形解析805をスキップする。これにより、次のMarquardt法解析817における解の収束を早めることができる。
【0051】
【数5】
Figure 0003628479
【0052】
【数6】
Figure 0003628479
【0053】
【数7】
Figure 0003628479
【0054】
【数8】
Figure 0003628479
【0055】
次に、Marquardt法解析817では目標運動諸元の状態ベクトルを数9のように定義し、疑似線形解析805で算出した初期解析値とこれに基づく目標の放射周波数を数10により算出して設定する。
【0056】
【数9】
Figure 0003628479
【0057】
【数10】
Figure 0003628479
【0058】
そして、縮小因子λ設定806ではMarquardt法解析における反復修正を安定的に行うことを目的として縮小因子λの初期値を設定する。
【0059】
処理807では反復修正回数loop1、loop2の初期化を行う。残差二乗行列算出808では数11により、各時刻における観測情報と解析情報との差分の二乗値を算出する。
【0060】
【数11】
Figure 0003628479
【0061】
ヤコビアン行列形成809では数12および数13により方位情報および周波数情報に関する偏微分係数を算出する。
【0062】
【数12】
Figure 0003628479
【0063】
【数13】
Figure 0003628479
【0064】
次に、状態ベクトル修正量Δxの算出810では、ステップ14に従って、処理808で算出した残差二乗行列と処理809で算出したヤコビアン行列を用いて状態ベクトルの修正量を算出する。
【0065】
【数14】
Figure 0003628479
【0066】
そして、loop1が所定の最大回数を越えたか、もしくは、前回算出のΔxが今回算出Δx以上であるかを判定する条件811を満足していれば、処理813で状態ベクトルを修正し、条件811を満足していなければ、処理812で縮小因子λを2倍にし、loop1を1加算し処理808から810を繰り返す。
【0067】
また、loop1が1以上であることを判定する条件814を満足していれば、縮小因子λを半分にして次回の反復修正処理に備える。
【0068】
そして、条件816で状態ベクトルの修正量が所定の極小値minより小さくなるか、反復修正回数loop2が上限値に達した場合、反復を打ち切り状態ベクトル数9を目標運動諸元(距離r,方位By,速力Mt,針路Ct,放射周波数f0)に変換して、処理を終了する。反復修正回数loop2が上限値に達していない場合には、loop2を1加算し808からの処理を繰り返す。
【0069】
なお、条件816で反復修正回数loop2が上限値に達して終了した場合は今回の非線形解析は「解なし」であったとして扱い、条件816で状態ベクトルの修正量が所定の極小値minより小さくなった場合には、今回の非線形解析は「解あり」であったとして扱う。
【0070】
なお、誤差をともなって観測される情報から内部の状態モデルを推定する最小二乗法およびMarquardt法の詳細については以下の文献に記載されている。
【0071】
「最小二乗法による実験データ解析」中川徹 他 東京大学出版会
次に、非線形解析609において、カルマンフィルタによる非線形解析を行う場合について説明する。
【0072】
この場合の処理の手順を、図5b)に示す。
【0073】
この場合の処理であるカルマンフィルタによる適応型フィルタ819では、目標運動諸元の状態ベクトルを数15のように定義し、センサデータベース508に蓄積されている各入力時系列データ中の各観測時刻毎のデータを一つ一つ取り出して状態ベクトルおよび共分散行列の初期値を逐次更新することにより解析値を収束させていく。
【0074】
【数15】
Figure 0003628479
【0075】
まず、条件820で前回の解析値が存在しないと判断されたとき、すなわち、第1回目の処理の場合には、初期値設定821において、数16により状態ベクトルの初期値を設定し、数17により共分散行列の初期値を設定する。
【0076】
【数16】
Figure 0003628479
【0077】
【数17】
Figure 0003628479
【0078】
次に、状態ベクトル状態状態遷移822では、状態ベクトルの初期値をセンサデータべースから、今回取り出したデータの観測時刻まで数18により時間遷移させ、共分散行列状態遷移823では、共分散行列の初期値をセンサデータべース508から今回取り出したデータの観測時刻まで数19により時間遷移させる。
【0079】
【数18】
Figure 0003628479
【0080】
【数19】
Figure 0003628479
【0081】
さらに、状態ベクトルフィルタリング824では数20によりセンサーデータベース508から取り出したデータに基づき状態ベクトルを修正更新し、共分散ベクトルフィルタリング825では数21により共分散行列を修正更新する。
【0082】
【数20】
Figure 0003628479
【0083】
【数21】
Figure 0003628479
【0084】
そして、最後に処理826で数22により状態ベクトルを目標運動諸元に変換し、数23により共分散行列を目標運動諸元解析値の解析誤差標準偏差に変換する。
【0085】
【数22】
Figure 0003628479
【0086】
【数23】
Figure 0003628479
【0087】
そして、条件827に従ってセンサデータベース508から取り出すデータが無くなるまで822から826までの処理を繰り返す。なお、誤差をともなって観測される情報から内部の状態モデルを推定するカルマンフィルタの詳細については以下の文献に記載されている。
【0088】
「応用カルマンフィルタ」片山徹 朝倉書店
次に、図3の数値変換処理部607の行う処理の詳細を説明する。
【0089】
図6に、数値変換処理部607の行う処理を示す。
【0090】
図示するように、数値変換処理部607では、前回解析時の最適解を継承する処理903、個体の数値変換処理904、評価値を適応度に変換する処理905の3つの処理を行う。
【0091】
前回解析時の最適解を継承する処理903では、固有バイアス誤差の曲線係数の最適解zjを遺伝的処理部の個体pjに変換する。これにより、最適解zjの実数値aj,bj,cjは0〜255の均等な数値遺伝子aj,bj,cjに変換され、乱数による遺伝的操作が可能となる。
【0092】
個体の数値変換処理904では、固有バイアス誤差の曲線係数の組合せである個体pjを、振幅、角周波数、初期位相の3つの実数値の組である解候補zjに変換する。これにより、遺伝的処理部で推定した個体pjの適応度として、各目標の運動諸元解析によって算出される残差二乗和の合計値で評価することが可能となる。
【0093】
評価値を適応度に変換する処理905では、非線形解析処理部で算出した評価値Rの逆数を遺伝的処理部における各個体の適応度Fとする。これにより、遺伝的処理部における個体の選択操作が可能となる。
【0094】
以上、目標運動解析装置505の詳細について説明した。
【0095】
次に、目標解析装置505の解析結果を出力する端末装置509について説明する。
【0096】
図7は、端末装置505に備えたディスプレイの表示のようすを表した図である。
【0097】
この表示画面は、グラフィックユーザインタフェースを形成し、オペレータは、この表示から標運動解析装置505の解析結果を視認することができると共に、表示上にキーボードなどから入力を行うことにより所望の操作を行うことができる。
【0098】
図7において、1001は観測船と目標との相対的な位置関係を表示する領域であり、領域1001中、1005が観測船の航跡を示している。また、1003が目標運動解析装置によって解析された目標の航跡を、1004が観測した方位線を、1002が目標運動解析装置によって解析された目標の解析誤差範囲を、各々示している。また、領域1010には現在時刻が表示される。
【0099】
領域1006から1009は横軸を時間とする表示領域であり、これらの表示において最右端が現在時刻である。これらの表示は、新たな表示データが発生する都度左側へスクロールしてゆく。
【0100】
ここで、領域1006は観測方位と解析方位との差分を観測誤差標準偏差で規格化してプロットした表示であり、1007は観測周波数と解析周波数との差分を観測誤差標準偏差で規格化してプロットした表示である。オペレーターは、これらの表示より、解析の収束状況を判断でき、また、目標の針路および速力の変更状況も判断することができる。
【0101】
次に、1008は観測方位をプロットした表示であり、1009は観測周波数をプロットした表示である。これらの表示より、オペレータは観測情報の連続性および目標の運動状況を把握することができる。
【0102】
ここで、1011から1014は異常データを示しており、1011および1013はオペレータの操作により一時消去が指示されたデータであり、1012および1014はオペレータの操作により削除が指示されたデータを示している。
【0103】
このように、オペレータよりデータの一時消去や削除の指示が成されると、入力時系列データ中の一時消去や削除の指示が成されたデータの観測時刻のデータを無視して、再解析を行うよう、端末装置505は、目標運動解析装置505に指示し、その結果得られた値に従って図7の表示を更新する。一時消去と削除の違いは、削除はセンサデータベース中の対応するデータを削除してしまうのに対して、一時消去は対応するデータ自身は削除せずに解析に対して用いないようにする点である。
【0104】
次、1015および1016は解析対象とする範囲の開始点および終了点を示している。これら1015および1016が操作されると、この範囲についてを再解析を行うよう、端末装置505は、目標運動解析装置505に指示し、その結果得られた値に従って図7の表示を更新する。
【0105】
次に、1017は遺伝的処理部によって推定された固有バイアス誤差特性を表示しており、その曲線の係数が1018に表示される。ここで、オペレータ入力部により1019に曲線係数が入力されると、入力された固有バイアス誤差を固定的に用いて、再解析を行うよう、端末装置505は、目標運動解析装置505に指示し、その結果得られた値に従って図7の表示を更新する。
【0106】
以上、本発明の第1の実施形態について説明した。
【0107】
以上のような処理により、本第1実施形態によれば、図14に示すように、局所最適解への落ち込みを防ぐことができる。図14において、a)が距離rの解析例、b)が速力Mtの解析例、c)が進路ctの解析例であり、各図において、X印が、本第1実施形態の遺伝的処理とカルマンフィルタを用いた解析を適用した場合、黒塗りの四角が本第1実施形態の遺伝的処理と最小2乗法による解析を適用した場合を示している。
【0108】
以上のように、本第1実施形態によれば、遺伝処理部に置いて固有バイアス誤差を推定しつつ非線形解析処理部において目標の運動諸元を解析するので、安定的かつ高精度な解析が可能となる。また、非線形解析処理の結果得られる評価値最小の解を次回解析時における遺伝的処理部の個体の一つとして継承したり、遺伝的処理部に置いて正規分布交叉を行ったり、非線形解析処理部におけるMarquardt法解析の初期値を継承することにより解析収束時間の短縮が実現される。また、さらに、オペレータは、任意に、固有バイアス誤差のモデル曲線の係数や解析対象データを変更し、再解析を行わせることができる。
【0109】
以下、本発明の第2の実施形態について説明する。
【0110】
に、本第2実施形態に係る目標運動解析装置の構成を示す。
【0111】
本第2実施形態では、固有バイアス誤差のモデル曲線の係数を含む目標運動諸元の状態ベクトルを遺伝的処理の個体とし、個別目標毎に固有バイアス誤差の推定と目標運動諸元の解析を実施する。
【0112】
本第2実施形態では、個別目標毎に次のような処理を行う。すなわち、第1実施形態と異なる点のみを説明すると、遺伝的処理部1102では、状態ベクトルを個体として扱う。また、非解析処理部では、各個体を用いて、今解析の対象としている目標の運動諸元を解析する。具体的には、非線形解析処理部1108の非線形解析1109において次のように処理する。すなわち、最小二乗法を用いる場合は、図5a)において疑似線形解析805を行わず、個体に対応する値を初期解析値としてMarquard法による解析817を行う。一方、カルマンフィルタによる場合は、図5b)の処理821で、個体に対応する値を初期値として初期状態ベクトルを求める。また、評価1110では、入力時系列データ中の観測方位および観測周波数と、今対象としている目標について解析結果として得られた時系列の解析方位および解析周波数との差分二乗を各々の観測情報の分散で規格化し、合算したもの評価値とする。なお、本第2実施形態で、遺伝処理部1102に前回の解析時より継承される個体(最適解)は、評価値を最小としたものである。
【0113】
なお、本第2実施形態では固有バイアス誤差特性が全観測目標について共通であるということに着目した推定は行っていない。
【0114】
以下、本発明の第3の実施形態について説明する。
【0115】
図9に、本第3実施形態に係る目標運動解析装置の構成を示す。
【0116】
図示するように、本第3実施形態では遺伝的処理を行わない。
【0117】
本第3実施形態では、固有バイアス誤差のモデル曲線の係数を含む目標運動諸元の解析値を状態ベクトルとし、個別目標毎に固有バイアス誤差の推定と目標運動諸元の解析を実施する。
【0118】
本第3実施形態では、個別目標毎に次のような処理を行う。すなわち、非解析処理部1202における解析を回帰的に繰り返し、回帰回数が所定の最大値を越えたか、非解析処理部1202における前回の解析結果と今回の解析結果の差が所定の最小値以下となった時点で処理を終了し、非解析処理部1202で最後に求まった最適解を解析結果として端末装置1209に出力する。ここで、非解析処理部1202における解析を回帰的に繰り返す際に、前回の解析結果の最適解を、非線形解析1204の初期値として与える。具体的には、非線形解析処理部1108の非線形解析1109において次のように処理する。すなわち、最小二乗法を用いる場合は、図5a)において疑似線形解析805を行わず、与えられた初期値を初期解析値として用いてMarquard法による解析817を行う。一方、カルマンフィルタによる場合は、図5b)の処理821では、与えられた初期値を初期状態ベクトルとして用いる。また、さらに、新たに解析を行う場合において、第1回目の非線形処理部1292の処理の際には、前回同一目標について求め、端末装置1209に出力した最適解を、非線形解析1204の初期値として与える。また、前回同一目標について求め、端末装置1209に出力した最適解が存在しない場合には、適当な初期値を非線形解析1204の初期値として与える。
【0119】
なお、本第3実施形態でも固有バイアス誤差特性が全観測目標について共通であるということに着目した推定は行っていない。
【0120】
次に、本発明の第4の実施形態について説明する。
【0121】
図10に、本第4実施形態に係る目標運動解析装置の構成を示す。
【0122】
本第4実施形態では、固有バイアス誤差のモデル曲線の係数を含む目標運動諸元の状態ベクトル遺伝的処理の個体とし、個別目標毎に固有バイアス誤差の推定と目標運動諸元の解析を実施する。各個体の評価は、各個体に対応する状態ベクトルから推定される時系列の解析方位および解析周波数と入力時系列データ中の観測方位および観測周波数との差分二乗を各々の観測情報の分散で規格化し、合算したもの評価値とする。その他の処理は、第1実施形態と同様である。
【0123】
なお、本第4実施形態でも固有バイアス誤差特性が全観測目標について共通であるということに着目した推定は行っていない。
【0124】
以上、本発明の実施形態について説明した。
【0125】
【発明の効果】
以上のように、本発明によれば、固有バイアス誤差に起因する目標運動解析精度の劣化を軽減することができる。
【図面の簡単な説明】
【図1】目標運動解析システムの概略構成を示したブロック図である。
【図2】目標運動解析の原理を示した図である。
【図3】目標解析装置の構成を示したブロック図である。
【図4】遺伝的処理部の行う処理を示した図である。
【図5】非線形解析処理部の行う処理を示したフローチャートである。
【図6】数値変換処理部の行う処理を示した図である。
【図7】端末装置におけるGUIを示した図である。
【図8】目標解析装置の他の構成例を示したブロック図である。
【図9】目標解析装置の他の構成例を示したブロック図である。
【図10】目標解析装置の他の構成例を示したブロック図である。
【図11】従来の目標運動解析の原理を示した図である。
【図12】音響センサの固有バイアス誤差を示した図である。
【図13】従来の問題点が生じ得る状態を示した図である。
【図14】従来の問題点を示した図である。
【符号の説明】
101・・・観測船、102・・・目標
201・・・観測船、202・・・音響センサ、203・・・指向性ビーム、204・・・ランダム誤差標準偏差、205・・・音響センサ、206・・・音響センサ隔壁、207・・・音響的視野角、208・・・目標方位の真値(船首基準)、209・・・観測方位(船首基準)、210・・・固有バイアス誤差、211・・・音響センサ固有バイアス誤差特性
301・・・目標方位の真値の時系列推移、302・・・船首基準の相対方位、303・・・船首基準の相対方位に対する固有バイアス誤差特性、304・・・固有バイアス誤差の時系列推移、305・・・相対距離の時系列推移、306・・・ランダム誤差特性、307・・・真北基準の観測方位
401・・・音響センサ隔壁、402・・・音響センサ、403・・・船首基準の相対方位に対する固有バイアス誤差特性、404・・・遺伝的処理部、405・・・各種係数による特性曲線、406・・・非線形解析処理部、407・・・センサデータベース
501・・・音響センサ、502・・・速力検出器、503・・・針路検出器、504・・・位置演算器、
505・・・目標運動解析装置、506・・・遺伝的処理部、507・・・非線形解析処理部、508・・・センサデータベース、509・・・端末装置
603・・・初期集合生成操作、604・・・選択操作、605・・・交叉操作、606・・・突然変異操作、607・・・数値変換処理部、609・・・非線形解析(最小二乗法または適応型フィルタ)、610・・・評価処理、611・・・解析終了判定処理、612・・・評価値最小解zjの継承処理
705・・・適応度の大きさに比例した面積を各個体に割り当てたルーレット、706・・・現世代の個体から次世代の個体を選択するためにルーレットに対して放つ矢、707・・・正規分布交叉の概念図
801・・・Marquardt法により実現した最小二乗法の処理例、802・・・前回解析時の「解有り」状況判定処理、803・・・線形解析準備処理、804・・・線形解析処理、805・・・疑似線形解析処理部、806・・・縮小因子λ設定処理、807・・・loop1初期化処理、808・・・残差二乗行列算出処理、809・・・ヤコビアン行列形成処理、810・・・状態ベクトル修正量Δxの算出処理、811・・・loop1終了判定処理、812・・・縮小因子更新処理、813・・・状態ベクトル修正処理、814・・・縮小因子更新判定処理、815・・・縮小因子更新処理、816・・・loop2終了判定処理、817・・・Marquardt法解析処理部、818・・・状態ベクトルを目標運動諸元に変換する処理、819・・・カルマンフィルタにより実現した適応型フィルタの処理例、820・・・第一回目解析の判定処理、821・・・初期値設定処理、822・・・状態ベクトル状態遷移処理、823・・・共分散行列状態遷移処理、824・・・状態ベクトルフィルタリング処理、825・・・共分散行列フィルタリング処理、826・・・状態ベクトルを目標運動諸元に変換する処理、827・・・センサデータベースのデータ有無判定処理
901・・・遺伝的処理部の個体集合、902・・・非線形解析処理部の解集合、903・・・数値変換処理部の機能(1)、904・・・数値変換処理部の機能(2)、905・・・数値変換処理部の機能(3)、906・・・個体が保有する数値遺伝子の種類、907・・・数値遺伝子の範囲、908・・・各個体の適応度の集合、909・・・最適解の実数パラメータの種類、910・・・実数パラメータの範囲、911・・・各解候補の評価値の集合
1001・・・目標と観測船との対勢図、1002・・・解析誤差標準偏差に基づく目標存在範囲、1003・・・解析結果に基づく目標運動線、1004・・・観測方位線、1005・・・観測船運動航跡、1006・・・観測方位と解析方位との残差を観測方位誤差標準偏差で規格化した時系列プロット、1007・・・観測周波数と解析周波数との残差を観測周波数誤差標準偏差で規格化した時系列プロット、1008・・・観測方位の時系列プロット、1009・・・観測周波数の時系列プロット、1010・・・現在時刻、1011・・・一時消去データ、1012・・・削除データ、1013・・・一時消去データ、1014・・・削除データ、1015・・・解析対象データ範囲の開始時刻、1016・・・解析対象データ範囲の終了時刻、1017・・・固有バイアス誤差特性表示部、1018・・・固有バイアス誤差特性関数の係数の解析値表示部、1019・・・固有バイアス誤差特性関数の係数に関するオペレータ入力値表示部、1020・・・目標運動諸元の解析結果表示部、1021・・・解析誤差標準偏差の表示部、1022・・・目標運動諸元に関するオペレータ入力値表示部
1101・・・目標運動解析装置、1102・・・遺伝的処理部、1103・・・初期集合生成操作、1104・・・選択操作、1105・・・交叉操作、1106・・・突然変異操作、1107・・・数値変換処理部、1108・・・非線形解析処理部、1109・・・非線形解析(最小二乗法または適応型フィルタ)、1110・・・評価処理、1111・・・解析終了判定処理、1112・・・評価値最小解zjの継承処理、1113・・・センサデータベース、1114・・・解析データ表示器
1201・・・目標運動解析装置、1202・・・非線形解析処理部、1203・・・疑似線形解析処理、1204・・・非線形解析(最小二乗法または適応型フィルタ)、1205・・・評価処理、1206・・・解析終了判定処理、1207・・・解xの継承処理、1208・・・センサデータベース、1209・・・解析データ表示器
1301・・・目標運動解析装置、1302・・・遺伝的処理部、1303・・・初期集合生成操作、1304・・・選択操作、1305・・・交叉操作、1306・・・突然変異操作、1307・・・数値変換処理部、1308・・・評価処理、1309・・・解析終了判定処理、1310・・・評価値最小解xjの継承処理、1311・・・センサデータベース、1312・・・解析データ表示器

Claims (13)

  1. 運動する目標から放射される音波を音響センサにより観測した観測値に基づいて、目標の運動諸元を解析する目標運動解析装置であって、
    音響センサの固有バイアス誤差を特定する情報を遺伝子とする個体の集合に遺伝的操作を行い、順次、次世代の個体の集合を生成する遺伝的操作部と、
    生成された各世代の個体の集合に含まれる各個体のそれぞれについて、当該個体に対応する固有バイアス誤差を考慮した目標の運動諸元の解析を行うと共に、当該解析の結果と観測値との適合度を求める非線形解析処理部とを備え、
    前記遺伝的操作部は、次世代の個体の集合の生成を、現世代の個体のうち、前記非線形解析部において、より高い前記適合度が当該個体について求められた個体の特徴が、より高い確率で次世代に遺伝されるように、前記遺伝的操作を行うことを特徴とする目標運動解析装置。
  2. 請求項1記載の目標運動解析装置であって、
    現世代の各個体について前記非線形解析部が求めた適合度の平均が、前世代の各個体について前記非線形解析部が求めた適合度の平均よりも低くなった場合、もしくは、現世代の各個体について前記非線形解析部が求めた適合度のうち最も高い適合度が、前世代の各個体について前記非線形解析部が求めた適合度のうち最も高い適合度より低くなった場合に、非線形解析部が、それまでに解析した解析の結果の内、最も適合度が高い解析の結果を、目標の運動諸元として出力する終了判定部を備えたことを特徴とする目標運動解析装置。
  3. 請求項1記載の目標運動解析装置であって、
    前記遺伝的操作部は、第1世代の個体の集合を生成する場合に、前回行われた目標の運動諸元の解析において、当該個体について前記非線形解析部において最も高い適合度が求められた個体を、第1世代の個体の集合に含めることを特徴とする目標運動解析装置。
  4. 請求項1記載の目標運動解析装置であって、
    音響センサの固有バイアス誤差の観測装置に対する目標の相対的な方位角に対する特性をモデル化し、前記個体の遺伝子とする音響センサの固有バイアス誤差を特定する情報として、モデル関数の振幅、周期、初期位相の組み合わせを用いることを特徴とする目標運動解析装置。
  5. 請求項1記載の目標運動解析装置であって、
    前記遺伝的操作処理部は、前記遺伝的操作として、個体の選択操作、個体間の交叉操作、個体の突然変異操作を行うことを特徴とする目標運動解析装置。
  6. 請求項1記載の目標運動解析装置であって、
    前記運動する目標から放射される音波を音響センサにより観測した観測値として、音波の到来方位である観測方位と音波の観測された周波数である観測周波数を用い、
    前記非線形解析処理部は、各個体についての前記運動諸元の解析を、複数の目標について、最小二乗法または適応型フィルタによって行い、時系列に得られた観測方位および観測周波数と、解析の結果得られた時系列の解析方位および解析周波数の重み付き残差二乗和の、各目標についての合計値の逆数を前記適合度とすることを特徴とする目標運動解析装置。
  7. 請求項1記載の目標運動解析装置であって、
    前記遺伝的操作処理部の行う遺伝的操作として、個体間の交叉を、前記個体間の距離と、当該個体について求められた適合度とに比例した最大値をとる正規乱数分、各個体を変化させることにより行うことを特徴とする目標運動解析装置。
  8. 請求項6記載の目標運動解析装置であって、
    前記非線形解析処理部は、前記各個体についての前記運動諸元の解析を、Marquardt法解析を用いる最小二乗法によって行い、
    Marquardt法解析の回帰回数が上限値に達するまでに収束した場合に収束した解析値を、次回の目標の運動諸元の解析においてMarquardt法解析の初期値とすることを特徴とする目標運動解析装置。
  9. 請求項1記載の目標運動解析装置であって
    オペレータの操作を受け入れる端末部を備え、
    当該端末部は、個体の遺伝子の入力を受け入れる手段と、
    入力された遺伝子を持つ個体についての、目標の運動諸元の解析を前記非線形解析部に実行させる手段を有することを特徴とする目標運動解析装置。
  10. 請求項1記載の目標運動解析装置であって
    オペレータの操作を受け入れる端末部を備え、
    当該端末部は、解析の対象とする観測値の集合の指定を、観測値の範囲もしくは解析対象観測値の集合からの特定の観測値の除去によって受け入れる手段と、
    受け入れた観測値を対象とする、目標の運動諸元の解析を前記遺伝的操作処理部および前記非線形解析部に実行させる手段を有することを特徴とする目標運動解析装置。
  11. 運動する目標から放射される音波を音響センサにより観測した観測値に基づいて、目標の運動諸元を解析する目標運動解析装置であって、
    音響センサの固有バイアス誤差を特定する情報と目標の運動諸元を遺伝子とする個体の集合に遺伝的操作を行い、順次、次世代の個体の集合を生成する遺伝的操作部と、
    生成された各世代の個体の集合に含まれる各個体のそれぞれについて、当該個体に対応する固有バイアス誤差と目標の運動諸元を初期値として、音響センサの固有バイアス誤差を特定する情報と目標の運動諸元の解析を、最小二乗法または適応型フィルタによって処理によって行うと共に、当該解析の結果と観測値との適合度を求める非線形解析処理部とを備え、
    前記遺伝的操作部は、次世代の個体の集合の生成を、現世代の個体のうち、前記非線形解析部において、より高い前記適合度が当該個体について求められた個体の特徴が、より高い確率で次世代に遺伝されるように、前記遺伝的操作を行うことを特徴とする目標運動解析装置。
  12. 運動する目標から放射される音波を音響センサにより観測した観測値に基づいて、目標の運動諸元を解析する目標運動解析装置であって、
    音響センサの固有バイアス誤差を特定する情報と目標の運動諸元を遺伝子とする個体の集合に遺伝的操作を行い、順次、次世代の個体の集合を生成する遺伝的操作部と、
    生成された各世代の個体の集合に含まれる各個体のそれぞれについて、当該個体に対応する固有バイアス誤差と目標の運動諸元と観測値との適合度を求める評価部とを備え、
    前記遺伝的操作部は、次世代の個体の集合の生成を、現世代の個体のうち、前記評価部において、より高い前記適合度が当該個体について求められた個体の特徴が、より高い確率で次世代に遺伝されるように、前記遺伝的操作を行い、かつ、第1世代の個体の集合を生成する場合に、前回行われた目標の運動諸元の解析において、当該個体について前記評価部において最も高い適合度が求められた個体を、第1世代の個体の集合に含めることを特徴とする目標運動解析装置。
  13. 運動する目標から放射される音波を音響センサにより観測した観測値に基づいて、目標の運動諸元を解析する目標運動解析方法であって、
    音響センサの固有バイアス誤差を特定する情報を遺伝子とする個体の集合に遺伝的操作を行い、順次、次世代の個体の集合を生成する遺伝的操作ステップと、
    生成された各世代の個体の集合に含まれる各個体のそれぞれについて、当該個体に対応する固有バイアス誤差を考慮した目標の運動諸元の解析を行うと共に、当該解析の結果と観測値との適合度を求める非線形解析ステップとを有し、
    前記遺伝的操作ステップにおいて、次世代の個体の集合の生成を、現世代の個体のうち、前記非線形解析部において、より高い前記適合度が当該個体について求められた個体の特徴が、より高い確率で次世代に遺伝されるように、前記遺伝的操作を行うことを特徴とする目標運動解析方法。
JP14187097A 1997-05-30 1997-05-30 目標運動解析装置および目標運動解析方法 Expired - Lifetime JP3628479B2 (ja)

Priority Applications (1)

Application Number Priority Date Filing Date Title
JP14187097A JP3628479B2 (ja) 1997-05-30 1997-05-30 目標運動解析装置および目標運動解析方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
JP14187097A JP3628479B2 (ja) 1997-05-30 1997-05-30 目標運動解析装置および目標運動解析方法

Publications (2)

Publication Number Publication Date
JPH10332806A JPH10332806A (ja) 1998-12-18
JP3628479B2 true JP3628479B2 (ja) 2005-03-09

Family

ID=15302084

Family Applications (1)

Application Number Title Priority Date Filing Date
JP14187097A Expired - Lifetime JP3628479B2 (ja) 1997-05-30 1997-05-30 目標運動解析装置および目標運動解析方法

Country Status (1)

Country Link
JP (1) JP3628479B2 (ja)

Families Citing this family (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US7190633B2 (en) 2004-08-24 2007-03-13 Bbn Technologies Corp. Self-calibrating shooter estimation
KR100856601B1 (ko) * 2004-08-24 2008-09-03 비비엔 테크놀로지스 코포레이션 슈터 위치의 명확화 시스템 및 방법
JP5134314B2 (ja) * 2007-09-03 2013-01-30 株式会社日立製作所 目標運動解析方法及び装置
JP5730506B2 (ja) * 2010-06-25 2015-06-10 日本無線株式会社 到来方向標定装置および位置標定装置
JP5679856B2 (ja) * 2011-02-14 2015-03-04 三菱電機株式会社 測位装置及び測位方法

Also Published As

Publication number Publication date
JPH10332806A (ja) 1998-12-18

Similar Documents

Publication Publication Date Title
Gruen Adaptive least squares correlation: a powerful image matching technique
CN109597864B (zh) 椭球边界卡尔曼滤波的即时定位与地图构建方法及系统
US6054710A (en) Method and apparatus for obtaining two- or three-dimensional information from scanning electron microscopy
Leung et al. Detection and estimation using an adaptive rational function filter
US7461112B2 (en) Signal processor, signal processing method, signal processing program, recording medium with the signal processing program recorded therein and measuring instrument
JP5546244B2 (ja) 二次元画像の類似性を測定するための方法および電子顕微鏡
Petrov et al. Box particle filtering for extended object tracking
US7974816B2 (en) Sensor registration by global optimization procedures
JP3628479B2 (ja) 目標運動解析装置および目標運動解析方法
Huletski et al. Vinyslam: an indoor slam method for low-cost platforms based on the transferable belief model
US20080319558A1 (en) Global optimization by continuous greedy randomized adaptive search procedure
US20100127933A1 (en) Algorithm of collecting and constructing training location data in a positioning system and the positioning method therefor
JP2005201869A (ja) 信号処理方法、信号処理プログラム、この信号処理プログラムを記録した記録媒体および信号処理装置
Premus Modal scintillation index: A physics-based statistic for acoustic source depth discrimination
Keesman et al. Nonlinear set-membership estimation: A support vector machine approach
Rajani et al. Direction of arrival estimation by using artificial neural networks
Tan et al. Evolutionary system identification in the time domain
WO2021192038A1 (ja) 画像解析装置および画像解析方法
CN114252871A (zh) 一种基于机器学习的雷达测量精度补偿方法
EP3657212A1 (en) Method and system of decomposition of composite targets on elements of a radar target signature with a super-resolution, using the total signal spectrum
JP2007306373A (ja) 信号分離装置、信号分離方法、信号分離プログラム及び記録媒体
WO2022070243A1 (ja) 画像解析装置および画像解析方法
CN117148380B (zh) 基于频率下降优化的振动相位补偿卫星isal成像方法
Nguyen et al. Electronic scan strategy for phased array weather radar using a space–time characterization model
Lindqvist et al. Posterior linearisation smoothing with robust iterations

Legal Events

Date Code Title Description
A977 Report on retrieval

Free format text: JAPANESE INTERMEDIATE CODE: A971007

Effective date: 20040802

A131 Notification of reasons for refusal

Free format text: JAPANESE INTERMEDIATE CODE: A131

Effective date: 20040907

A521 Written amendment

Free format text: JAPANESE INTERMEDIATE CODE: A523

Effective date: 20041027

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: 20041130

A61 First payment of annual fees (during grant procedure)

Free format text: JAPANESE INTERMEDIATE CODE: A61

Effective date: 20041208

FPAY Renewal fee payment (event date is renewal date of database)

Free format text: PAYMENT UNTIL: 20071217

Year of fee payment: 3

FPAY Renewal fee payment (event date is renewal date of database)

Free format text: PAYMENT UNTIL: 20081217

Year of fee payment: 4

FPAY Renewal fee payment (event date is renewal date of database)

Free format text: PAYMENT UNTIL: 20081217

Year of fee payment: 4

FPAY Renewal fee payment (event date is renewal date of database)

Free format text: PAYMENT UNTIL: 20091217

Year of fee payment: 5

FPAY Renewal fee payment (event date is renewal date of database)

Free format text: PAYMENT UNTIL: 20101217

Year of fee payment: 6

FPAY Renewal fee payment (event date is renewal date of database)

Free format text: PAYMENT UNTIL: 20101217

Year of fee payment: 6

FPAY Renewal fee payment (event date is renewal date of database)

Free format text: PAYMENT UNTIL: 20111217

Year of fee payment: 7

FPAY Renewal fee payment (event date is renewal date of database)

Free format text: PAYMENT UNTIL: 20111217

Year of fee payment: 7

FPAY Renewal fee payment (event date is renewal date of database)

Free format text: PAYMENT UNTIL: 20121217

Year of fee payment: 8

FPAY Renewal fee payment (event date is renewal date of database)

Free format text: PAYMENT UNTIL: 20131217

Year of fee payment: 9

EXPY Cancellation because of completion of term