誤差伝播による出力データの精度評価方法
【課題】出力データの精度評価を、理論的で、簡単な計算で行え、計算量が少なく、収束計算を行う必要が無く、データ処理を行う多くのソフトウェアに適用可能であり、誤差の要因の分析、要因毎の大きさの比較等の多くの解析が可能となるようにする。
【解決手段】データ処理を入力データによって数値微分することで、データ処理のヤコビ行列Ｊを求める手順（ステップ１１０）と、入力データの持つ誤差の分散共分散を推定して、入力データの誤差行列Ｄを求める手順（ステップ１２０）と、前記ヤコビ行列Ｊ及び入力データの誤差行列Ｄより、出力データの持つ誤差の分散共分散を表す誤差行列Ｒを計算する手順（ステップ１３０）と、を含む。
【解決手段】データ処理を入力データによって数値微分することで、データ処理のヤコビ行列Ｊを求める手順（ステップ１１０）と、入力データの持つ誤差の分散共分散を推定して、入力データの誤差行列Ｄを求める手順（ステップ１２０）と、前記ヤコビ行列Ｊ及び入力データの誤差行列Ｄより、出力データの持つ誤差の分散共分散を表す誤差行列Ｒを計算する手順（ステップ１３０）と、を含む。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、誤差伝播による出力データの精度評価方法に係り、特に、三次元座標測定機（以下、単に三次元測定機又はＣＭＭとも称する）における座標測定値の不確かさを推定する際に用いるのに好適な、誤差伝播による出力データの精度評価方法に関する。
【背景技術】
【０００２】
品質保証の観点から、全ての測定機器の測定結果には不確かさが必要になっている。
【０００３】
その際、複雑なデータ処理、シミュレーションなど（以後、データ処理）の計算結果がどのくらいの精度を持っているかを評価することが重要である。データ処理では、計算の入力となる入力データ、計算処理に必要な入力パラメータなど（以後、パラメータも入力データに含める）により、計算結果である出力データが得られる。出力データがどのくらいの誤差を持つかは、データ処理の前提となる計算手法やモデルがどのくらい正しいかによる場合が多い。しかし、計算手法やモデルが十分正しく、正しい入力データに対して、正しい出力データが得られるような場合でも、入力データの持つ誤差により、出力データが誤差を持つ。以後、正しい入力データに対しては、正しい出力データが得られるようなデータ処理を前提とする。
【０００４】
このような条件では、従来、入力データ（パラメータを含む）を単に少しだけ変化させ、出力データの影響を見たり、入力データをモンテカルロシミュレーションや実験計画法によって変化させ、出力誤差を推定したり、経験により出力データの誤差を評価することが行われている。
【０００５】
特に、近年、精密計測における形状計測の中で、三次元測定機を用いる座標計測の重要性が増加している。この中では、三次元測定機によって測定された形体上の複数の測定点から、最小二乗法により実際形体が計算される。しかし、三次元測定機本体の不確かさ以外にも、測定対象、測定環境、測定戦略等の影響により、不確かさが変わるため、不確かさを正確に推定することが難しい。
【０００６】
従来は、三次元測定機の測定結果の不確かさを、経験的な方法、一連の実験による方法、バーチャルＣＭＭ手法（モンテカルロシミュレーションによる方法）で推定していた（特許文献１、非特許文献１乃至８参照）。
【先行技術文献】
【特許文献】
【０００７】
【特許文献１】特開２００６−２１５０１１号公報
【非特許文献】
【０００８】
【非特許文献１】高増潔：三次元座標測定概論−精密測定の条件と三次元座標測定の不確かさ−、精密工学会講習会「差のつく三次元座標測定機活用法」、２００７年１２月
【非特許文献２】高増潔：ＧＰＳ（製品の幾何特性仕様）に基づいた検証方法−幾何公差の解釈と三次元測定機の不確かさ評価、設計工学、４０（２），２００５，７９−８５
【非特許文献３】高増潔：三次元測定機の不確かさ推定方法とその国際規格化、機械と工具、４５（３），２００１，６６−７０
【非特許文献４】高増潔，佐藤理，下嶋賢，古谷涼秋：座標測定機のアーティファクト校正（第３報）−校正後の測定の不確かさの推定−，精密工学会誌，７１（７），２００５，８９０−８９４
【非特許文献５】高増潔，阿部誠，古谷涼秋，大園成夫：形体計測における不確かさの見積り（第１報）−校正作業で生じる系統誤差の寄与−，精密工学会誌，６７（１），２００１，ｐｐ． ９１−９５
【非特許文献６】高増潔，阿部誠，古谷涼秋，大園成夫：形体計測における不確かさの見積り（第３報）−区分多項式による直線形体の形状偏差モデル−，精密工学会誌，７１（１１），２００５，１４５４−１４５８
【非特許文献７】高増潔，野坂健一郎，阿部誠，古谷涼秋，大園成夫：形体計測における不確かさの見積り（第２報）−円形体の計測における形状偏差の寄与−，精密工学会誌，６９（５），２００３，６９３−６９７
【非特許文献８】高増潔，古谷涼秋，大園成夫：最小領域法の座標計測における統計的評価，精密工学会誌，６４（４），１９９８，ｐｐ．５５７−５６１
【発明の概要】
【発明が解決しようとする課題】
【０００９】
しかし、これらの手法のうち、モンテカルロシミュレーションでは、誤差を含んだ入力データを効率よく発生することが難しく、計算量も多く、結果の信頼性も低い。他の方法では、経験的な側面が多く、一般的に評価する手法となっていない。
【００１０】
特に三次元測定機のような複雑な測定機では、不確かさの推定が難しく、実験のコスト、時間がかかり、推定精度も低いことから、まだ、実際には有効な不確かさ推定手法は確立されていなかった。
【００１１】
本発明は、前記従来の問題点を解消するべくなされたもので、理論的で、簡単な計算で行え、計算量が少なく、収束計算を行う必要が無く、データ処理を行う多くのソフトウェアに適用可能であり、誤差の要因の分析、要因毎の大きさの比較等の多くの解析が可能な出力データの精度評価方法を提供することを課題とする。
【課題を解決するための手段】
【００１２】
本発明は、図１に示す如く、データ処理を入力データによって数値微分することで、データ処理のヤコビ行列Ｊを求める手順（ステップ１１０）と、入力データの持つ誤差の分散共分散を推定して、入力データの誤差行列Ｄを求める手順（ステップ１２０）と、前記ヤコビ行列Ｊ及び入力データの誤差行列Ｄより、出力データの持つ誤差の分散共分散を表す誤差行列Ｒを計算する手順（ステップ１３０）と、を含むことにより、前記課題を解決したものである。
【００１３】
以下、各手順について説明する。
【００１４】
１．データ処理のヤコビ行列Ｊを求める手法（ステップ１１０）
データ処理では、入力データ列（入力パラメータを含める）に対して、出力データ列が得られる。まず、計算したい条件に合った入力データ列Ｉを用意し、その出力データ列Ｏを得る。ここで、入力データは、ベクトル量についてはスカラー量として分割し、全部でｎ個ある場合は、入力データ列Ｉはｎ個のスカラー量の並びとして表現する。同様に出力データも、ｍ個のスカラー量の並びとして表現する。次式参照。
【００１５】
【数１】
【００１６】
データ処理のヤコビ行列Ｊは、ｎ個の入力データとｍ個の出力データの関係を示すため、ｎ行ｍ列の行列として表現できる。ヤコビ行列Ｊのｎ行ｍ列の要素jn,mは、以下の計算で求めることができる。
【００１７】
１）入力データ列Ｉのｋ番目の入力だけを微小量ｄｋだけ変化させた入力データ列Ｉｋを作る。
【００１８】
２）入力データ列Ｉｋに対してデータ処理を行い、出力データ列Ｏ＊を求める。
【００１９】
３）ヤコビ行列Ｊのｋ行の要素jk,1〜jk,mは、出力データ列Ｏ＊とＯの各要素の差とｄｋを利用して、次式で求めることができる。
【００２０】
【数２】
【００２１】
４）同様にヤコビ行列の全ての行（１行からｎ行）を求めることができる。
【００２２】
上記の計算は、機械的に行うことができる。上記処理を行うプログラムからデータ処理をｎ＋１回行うことで、自動的にヤコビ行列Ｊを計算できる。
【００２３】
ここで、微小量ｄｋは、以下のことを考慮して決定する。
【００２４】
１）データ処理の非線形性が小さい場合は、ｄｋをどのように決めても、結果に大きな影響はない。
【００２５】
２）次の２．項で推定する入力データの誤差行列Ｄを利用すれば、その入力データの分散の平方根に対して、何分の一かを利用するなど、小さい量を設定することで、自動的に決定できる。
【００２６】
３）入力データの性質がある程度分かっている場合には、その誤差などに対して小さい量を設定する。
【００２７】
２．入力データの誤差行列Ｄの推定（ステップ１２０）
入力データが持つ誤差の性質を表す誤差行列Ｄを推定するためには、まず、入力データに影響する独立した誤差要因に対して、それぞれの誤差行列を求めて単純に和を求めればよいので、独立した1つの誤差要因に対する誤差行列の求め方を示す。
【００２８】
まず、入力データの誤差がそれぞれ独立で、その標準偏差ｔが分かっている場合は、誤差行列Ｄの対角要素がそれぞれの誤差の標準偏差の二乗（分散）となり、非対角要素は０となる。次式参照。
【００２９】
【数３】
【００３０】
入力データの誤差が独立でない場合は、２つの入力データ間の相関係数ｓ及び誤差の標準偏差ｔが推定できれば、相関係数ｓと標準偏差ｔから共分散を計算することができる。誤差行列Ｄの対角要素はそれぞれの誤差の標準偏差の二乗（分散）となり、非対角要素は対応する２つの入力データ間の共分散となる。次式参照。
【００３１】
【数４】
【００３２】
より複雑なデータの関係に関しては、ｉ）理論的な推定、ii）この部分だけをモンテカルロシミュレーションで行う、iii）経験によるなどの方法で推定する必要がある。又、次の３．項に示すように、本手法をより単純な入力データへ適用することで、誤差行列を推定することも可能である。
【００３３】
３．出力データの誤差行列Ｒの計算（ステップ１３０）
上記手法により求めたヤコビ行列Ｊ及び入力データの誤差行列Ｄより、次式の計算方法で出力データの誤差行列Ｒを計算できる。
【００３４】
【数５】
ここで、Ｊｔは行列Ｊの転置行列を示す。
【発明の効果】
【００３５】
本発明によれば、次のような効果を奏する。
【００３６】
１．出力データの誤差行列Ｒの利用
出力データの誤差行列Ｒは、このデータ処理において計算された結果の誤差を推定している。誤差行列Ｒのｊ行ｊ列の要素rj,jは、出力データ列Ｏのｊ番目の出力値ｏｊの誤差の分散を示している。そこで、出力データ列Ｏのｊ番目の出力値誤差は、標準偏差で誤差行列Ｒのｊ行ｊ列の要素rj,jの平方根tjとなる。次式参照。
【００３７】
【数６】
【００３８】
入力データの持つ誤差が正規分布であることが仮定でき、データ処理の非線形性が少ない場合には、出力データの誤差も正規分布であることが仮定できる。このような場合には、ｏｊ±ｔｊの範囲に出力データは、６０％の確率で入ることになり、ｏｊ±２ｔｊの範囲に９５％の確率で入ることになる。このように、出力データの誤差を評価することができる。
【００３９】
出力データの誤差行列Ｒのｊ行ｌ列の要素rj,lは、ｊ番目の出力とｌ番目の出力の共分散を示している。rj,lとそれぞれの出力の誤差の分散rj,j及びrl,lにより、式（３）を使って、２つの出力の誤差の相関係数sj,lを計算することができる。この相関係数は、±１の範囲で、絶対値が０に近い場合は２つの出力の誤差に相関がないこと、絶対値が１に近い場合は２つの出力の誤差に正又は負の相関が大きいことを示す。相関が大きい場合には、２つの出力値の演算は単に出力値を使うことはせずに、演算するデータ処理に対して、本手法を再び適用することで正しい評価が可能となる。
【００４０】
上記のように、出力データが別のデータ処理の入力データとなる場合には、出力データの誤差行列を次のデータ処理の入力データの誤差行列として使用することができる。
【００４１】
２．ヤコビ行列Ｊの利用
ヤコビ行列Ｊは、入力データと出力データの関係を示している。入力データのｉｋだけがｄｋだけ変化した場合の出力データolの変化はjk,l × dkで評価できる。これが十分小さい場合には、入力データｉｋは出力データolにほとんど関係ないことが分かる。この解析により、入力データと出力データが無関係な部分を抽出することができ、データ処理を独立したいくつかの部分に分離することが可能となる。
【００４２】
３．入力データの誤差行列Ｄと出力データの誤差行列Ｒの関係
入力データの誤差行列がｐ個の複数の独立した誤差要因から構成されているとすると、それぞれの要因毎の誤差行列Ｄ1〜Ｄpが推定され、それぞれに対応した出力データの誤差行列Ｒ1〜Ｒpが次式のように計算される。
【００４３】
【数７】
【００４４】
それぞれのＲに対して１．項に示した解析を行うことができる。又、Ｒの大きさを比較することで、誤差要因のどれが重要かどうかを判断できる。
【００４５】
以上をまとめると、従来の複雑なデータ処理やシミュレーションでは、計算結果を得るまでで終わっている場合が多く、計算結果の評価ができていない。これに対して、本発明では、モンテカルロシミュレーションではなく、誤差伝播を利用することによって、どのようなデータ処理に対しても汎用的に適用できる。本発明によれば、自動的に誤差評価が可能で、その結果に対して種々の分析が行なえ、データ解析に対する強力なツールとなる。
【図面の簡単な説明】
【００４６】
【図１】本発明による基本的な処理手順を示す流れ図
【図２】本発明の第１実施形態の適用対象である三次元測定機の一例を示す斜視図
【図３】第１実施形態における測定点と測定値の関係（２次元での例）を示す図
【図４】同じくソフトウェアの誤差伝播による不確かさ推定の基本的な流れを示す図
【図５】同じく三次元測定機の運動学による不確かさの概念を示す図
【図６】同じくプロービングによる誤差の推定を示す図
【図７】同じく円形体の形状偏差による不確かさを示す図
【図８】フーリエ変換に適用した本発明の第２実施形態におけるフーリエ変換プログラムを示す図
【図９】同じくフーリエ変換プログラムの入出力データ及び行列の例を示す図
【図１０】同じく出力の例を示す図
【図１１】条件探索による形体計測に適用した本発明の第３実施形態の処理手順を示す流れ図
【図１２】光学シミュレーションに適用した本発明の第４実施形態の処理手順を示す流れ図
【図１３】生体信号計測に適用した本発明の第５実施形態の処理手順を示す流れ図
【発明を実施するための形態】
【００４７】
以下図面を参照して、本発明の実施形態を詳細に説明する。
【００４８】
本発明の第１実施形態は、本発明を、図２に例示するような三次元測定機１の不確かさ推定に適用したものである。
【００４９】
図において、２は、三次元測定機１を駆動制御すると共に三次元測定機１から必要な測定座標値を取込むための駆動制御装置、３は、駆動制御装置２を介して三次元測定機１を手動操作するためのジョイスティック操作盤、４は、駆動制御装置２での測定手順を指示するパートプログラムを編集・実行すると共に、駆動制御装置２を介して取込まれた測定座標値に幾何形状を当てはめるための計算を行ったり、パートプログラムを記録、送信する機能を備えたホストシステムである。
【００５０】
前記三次元測定機１は、次のように構成されている。即ち、除振台１０の上には、定盤１１がその上面をベース面として水平面と一致するように載置され、この定盤１１の両側端から立設されたアーム支持体１２ａ、１２ｂの上端でＸ軸ガイド１３を支持している。アーム支持体１２ａは、その下端がＹ軸駆動機構１４によってＹ軸方向に駆動され、アーム支持体１２ｂは、その下端がエアーベアリングによって定盤１１上にＹ軸方向に移動可能に支持されている。Ｘ軸ガイド１３は、垂直方向に延びるＺ軸ガイド１５をＸ軸方向に駆動する。Ｚ軸ガイド１５には、Ｚ軸アーム１６がＺ軸ガイド１５に沿って駆動されるように設けられ、Ｚ軸アーム１６の下端に接触式のプローブ１７が装着されている。このプローブ１７が定盤１１上に載置されたワーク５に接触したときに、プローブ１７から駆動制御装置２にタッチ信号が出力され、そのときのＸＹＺ軸座標値を駆動制御装置２が取り込むようになっている。
【００５１】
１．基本的な手法の概要
以下、２次元での例として図３を利用して説明するが、手法そのものは、簡単に三次元へ拡張可能である。図３でＰ１からＰ１７までは、測定点であり、Ｄ１からＤ１３までが測定値である。この例では、Ｐ１からＰ４によりワーク座標系を設定し、Ｐ５からＰ８で２つの点の位置と辺の角度、Ｐ９及びＰ１０で辺の位置、Ｐ１１からＰ１４及びＰ１５からＰ１７で２つの円の中心と直径を測定している。測定値は例えば、Ｄ１とＤ２は点のＸ座標及びＹ座標、Ｄ３は２つの辺の角度、Ｄ８、Ｄ９及びＤ１０は円の中心のＸ座標、Ｙ座標及び直径を表す。
【００５２】
本実施形態は、図４に示す如く、基本的に誤差伝播に基づいて測定の不確かさを推定する。次式のように、測定値Ｄは、測定処理をするソフトウェア（測定用ソフトウェアと称する）Ｆの出力として測定点座標Ｐを入力することで得られる。
【００５３】
【数８】
【００５４】
このソフトウェアを次式のようにＰで偏微分することでヤコビ行列Ａを計算できる。
【００５５】
【数９】
【００５６】
この計算には、測定用ソフトウェアをブラックボックスとして利用する。このため、測定用ソフトウェアがどんなものであれ、オフラインで入力として測定点座標、出力として測定値が得られればよい。測定点座標の不確かさ行列Ｓが推定できれば、測定値の不確かさは次式の誤差伝播によって計算できる。
【００５７】
【数１０】
【００５８】
ヤコビ行列Ａは、式（９）によりソフトウェアをオフラインで利用して、機械的に計算できる。ヤコビ行列の計算方法は後述する。
【００５９】
そこで、この手法のキーとなる技術は、誤差行列Ｓの推定である。測定点座標の持つ不確かさが互いに独立で同じであるような場合は単純であるが、低レベルの推定にしかならない。３．項で後述するように、多くの不確かさ要因を解析することで、できるだけ実際に近い誤差行列の推定を行う。
【００６０】
ヤコビ行列Ａが計算でき、不確かさ行列Ｓが推定できれば、測定値Ｄの不確かさ行列Ｔ（それぞれの測定値の分散、共分散）は、式（１０）で計算できる。式（１０）の形は、Ｆが観測方程式の場合とは違うことに注意する必要がある。Ｆが観測方程式の場合は、Ｆの入力は測定値だが、出力は測定結果ではなく求めたい形体からの偏差になっている。そして、偏差の二乗和を最小にするように最小二乗法を解いた結果が、ここでのＦと等しくなる。そこでＦが観測方程式の場合のヤコビ行列をＢとするとＡｔ＝（ＢｔＳ−１Ｂ）−１ＢｔＳ−１となる。しかし、一般の三次元測定機のプログラムは最小二乗法を解くときにＳを考慮しないので、実際はＡ＝（ＢｔＢ)−１Ｂとなる。
【００６１】
２．ヤコビ行列の推定
ヤコビ行列Ａは、次式に示すように、各測定点ＰiのＸ座標Ｐix，Ｙ座標Ｐiyと測定値Ｄjの関係について、数値的な差分を取ることで求める。
【００６２】
【数１１】
【００６３】
これにより、測定手順及びソフトウェア部分をブラックボックスＦとして扱うことができる。従来は、理論的な計算又はモンテカルロ法で計算を行なっていたが、測定手順及びソフトウェア部分をオフラインモード（三次元測定機を使わずに、仮想の測定を行なうモード）で実施し、ブラックボックスとして扱うことができ、以下の計算を、機械的に行えるため理論的な難しさはない。
【００６４】
ヤコビ行列Ａは、測定点の座標の不確かさと測定値の不確かさの伝播を示すため、ｎ個の測定値（２次元なので２ｎ個の座標値）とｍ個の測定値に対応して、２ｎ行ｍ列の行列となる。図３の例では、ｎ＝１７、ｍ＝１３となるので、ヤコビ行列Ａは３４行１３列の行列である。ｉ番目の測定点Ｐiの実際の測定座標を（Ｐix,Ｐiy）とし、ｊ番目の計算された測定値をＤjとする。オフラインにおいて、ｉ番目の測定点のＸ座標だけをｄだけ増やしてＰiの座標を（Ｐix＋ｄ，Ｐiy）とし、同じ計算を行ない、そのときの計算された測定値をＤj,ixとする。この計算から、ヤコビ行列Ａの２ｉ−ｌ行の要素が（Ｄj,ix−Ｄj）／ｄとして計算できる。
【００６５】
実際の測定座標は、三次元測定機で測定しなくても、名目的な値（例えば図面上に記載されている値）を使っても構わない。又、ｄの値は数値誤差が問題にならない範囲で、十分小さい値を使うことが重要である。一般的には、例えば0.000000001mmといった値を使うことができる。
【００６６】
以上の計算を各測定点のＸ座標、Ｙ座標について計算すれば、ヤコビ行列Ａを計算できる。この場合、ソフトウェア部分は２ｉ回呼び出されることになる。実際には、ソフトウェアは、部分的に独立なため、全ての計算を繰り返す必要はないが、ソフトウェアの独立性についての判断方法は今後の課題とする。ヤコビ行列を計算した後なら、ヤコビ行列の共分散を見ることによって、自動的にソフトウェアの独立部分を抜き出すことができる。
【００６７】
３．測定点座標の不確かさ行列の推定
各測定点のＸ座標、Ｙ座標のそれぞれが持つ不確かさ（分散）と各測定点間の座標の関係（共分散）を推定することが、提案手法のキーとなる技術である。まず、各測定点の座標が持つ不確かさとして、以下の要因を考えることとする。このうち、いくつかはすでに定式化されている。
【００６８】
（１）ＣＭＭの運動学補正による不確かさ
（２）プロービング（プロービング、プローブ補正及びプローブ交換）による不確かさ
（３）測定物の形状偏差による不確かさ
（４）測定環境（温度ドリフトなど）による不確かさ
【００６９】
基本的な考え方は、それぞれの要因によって与えられる不確かさを求め、その要因がそれぞれの座標値に与える相関係数から共分散を推定する。以下、ひとつひとつの要因に対して、分散と共分散の推定方法の概要を示す。
【００７０】
（１）ＣＭＭの運動学補正による不確かさ
ＣＭＭは運動学に基づいた補正を行っているが、補正しきれていない高周波のランダムな不確かさと低周波の座標間の相関を持つ不確かさが存在する。この不確かさでは、Ｘ座標とＹ座標との相関は考える必要はないと思われる。座標の不確かさは全ての場所で一定の値を持つ。又、２つの座標間の相関は、２つの座標の差の絶対値の関数となる（非特許文献４参照）。
【００７１】
差し当たって、それぞれの測定点のＸ座標、Ｙ座標に不確かさscmmを与え、２つの測定点ｈとｉの座標の差の絶対値によって相関rcmm(|Phx - Pix|)を与える。Ｘ座標とＹ座標との間の相関は０とする。scmmとrcmmは、ＣＭＭの長さ測定における最大許容誤差のＭＰＥ線図、もしくはステップゲージ等の測定結果から推定する。
【００７２】
ＣＭＭの運動学による不確かさの概念を図５に示す
【００７３】
（２）プロービング（プロービング、プローブ補正およびプローブ交換）による不確かさ
プロービングによる不確かさは、ランダムな不確かさと測定ベクトルによる測定方向に依存する不確かさがある。又、２つの点を同じプローブで測定した場合、プローブの補正による不確かさの相関は２つの測定ベクトルの方向の差の関数によって記述できる（非特許文献５参照）。
【００７４】
まず、プロービングについては不確かさspを与え、プローブ径の補正を含めて、２つの測定点ｈとｉの測定ベクトルのなす角との関係で相関rp(Vh・Vi)を与える。これらの値は、理論的な考察及び（１）のscmmとrcmmから推定できる。
【００７５】
プロービングによる誤差の推定を図６に示す。
【００７６】
（３）被測定物の形状偏差による不確かさ
被測定物に形状偏差がある場合、形状偏差の持つランダムな不確かさと形状の自己相関による不確かさがある。２つの点が同じ被測定物を測定した場合、測定ベクトルの方向及び測定位置の関数として、不確かさの相関が記述できる（非特許文献６、７参照）。
【００７７】
又、形状偏差の持つランダムな誤差及び相関により、測定点の不確かさsf及び相関rfを推定できる（非特許文献６、７参照）。相関rfは、円の場合は測定ベクトルから推定される測定点の円上での角度位置の差で決まり、直線の場合は、直線上での測定点の位置の差で決まる。
【００７８】
円形体の形状偏差による不確かさを図７に示す。
【００７９】
（４）測定環境（温度ドリフトなど）による不確かさ
環境の影響による不確かさは評価が難しいが、温度ドリフトなどの環境変化による不確かさを時間の関数として評価することが考えられる。２つの測定点の測定した時刻の差が大きいと、環境が変化し相関が低くなる。一方、２つの測定点を連続して測定した場合には、環境の相関が高いため、環境による不確かさが打ち消しあうことになる。
【００８０】
本実施形態によれば、ヤコビ行列及び不確かさ行列から誤差伝播により、測定結果の不確かさを推定する手法で、理論的に高い推定精度が期待できる。又、ヤコビ行列の計算手法が、ソフトウェアをブラックボックスとして数値差分で行なえるため、高速で精度が高く、経験や実験が不要である。更に、測定座標の不確かさ行列を推定する手法が、理論的な解析に基づき、多くの要因をカバーできるため、最終的な測定結果の不確かさの推定精度が高い。又、測定座標の不確かさ行列を推定する方法が、理論的な方法のため、実験等がほとんど不要で、コストが低く、推定時間が速い。
【００８１】
なお、測定座標の不確かさ行列Ｓの推定要因として、次の６つを考えることもできる。
【００８２】
（１）三次元測定機の直角誤差による不確かさＳ１
（２）三次元測定機のスケール誤差による不確かさＳ２
（３）三次元測定機の直角誤差及びスケール誤差以外の運動学誤差による不確かさＳ３
（４）プロービングによる不確かさＳ４
（５）被測定物の形状偏差による不確かさＳ５
（６）測定環境（運動ドリフト等）による不確かさＳ６
【００８３】
次に、フーリエ変換に適用した本発明の第２実施形態を説明する。
【００８４】
この第２実施形態では、８個のデータをフーリエ変換して、０次から４次のパワースペクトルを求める。ここで、データに誤差がある場合、求めたパワースペクトルの誤差を推定する。
【００８５】
図８のプログラムで、ｆｆｔ８（５４行〜６１行）は、８個の入力データに対して、５個の出力データを求めるデータ処理で、中身の処理を知る必要は無い。入力データｉｉが与えられ、入力データの分散共分散行列ｄｄが推定できれば、後は機械的に計算できる。ｄｄとしては、分散が０．１２で共分散が２つのデータの距離に関係した相関を与えている（１４行〜２２行）。
【００８６】
２８行〜５２行までは、完全に機械的に計算できる部分で、対象となるデータ処理が何でも同じプログラムで計算が可能である。このプログラムの入出力データ及び行列の例を図９に、出力の例を図１０に示す。結果は、出力データが大きい次数に大きな誤差があること、相対的に低次の誤差が大きいこと等が分かる。出力データの誤差の相関係数がそれほど大きくないのはＦＦＴが直交的な処理だからと思われる。
【００８７】
図１１に、条件検索による形態計測に適用した本発明の第３実施形態の処理手順を示す。
【００８８】
例えば三次元測定機や形状測定機を用いて形体計測を行う場合、最小二乗法による評価の他に、最小領域法による評価を行うことがある。最小領域法では、条件探索を通して採択された数少ない代表点によってのみ結果が決まり、他の棄却された測定点は結果に寄与しない。そのため計算プロセスの統計的な評価が難しく、測定の不確かさを理論的に検討した研究は少ない。わずかに非特許文献８が挙げられる。この研究によると、最小領域法の不確かさを統計的に推定する手法が提案され、理論的な方法論は確立されている。
【００８９】
しかしながら最小領域法を実装した市販のソフトウェアにおいては、内部のアルゴリズムや数値演算の詳細についての情報は開示されていない。このような場合においても、本願発明による不確かさ推定法を適用することによって、最小領域法の測定の不確かさを容易に推定することができる。又最小領域法に限らず、複数の計算ステップを組み合わせた複雑なデータ処理プログラムに対しても適用することができる。
【００９０】
図１２に、光学シミュレーションに適用した本発明の第４実施形態の処理手順を示す。
【００９１】
例えば、レンズ等の光学素子の設計、及びこれらを組み合わせたシステム設計は、装置の小型化や高機能化の流れを受けて複雑化している。又、非線型光学効果を有する種々の素子の実用化や複雑な形状と光学特性を有する光共振器等についても多様なものが研究され、実用に供されている。これらの光学素子の設計や、それらを組み込んだシステム設計を適切に行うためには、光源から最終的な光学的出力にいたるまでの光電磁場の挙動を物理光学的に記述し得るシミュレーションツールが不可欠である。
【００９２】
一方、市販されている光学シミュレーションツールの内部では、解析対象に特有な光学場の条件設定に依存し、あるいはオペレータが指示する入力パラメータに依存し、より高速で、且つ数値的な精度が十分な、異なるアルゴリズムや解法を自動や手動で使いわけている。こうしたシミュレーションプログラムの内部情報は、全てがユーザに開示されているわけではなく、ユーザは自身が実施しているシミュレーションの内部処理を把握せずに、いわばブラックボックスとして使っていることになる。
【００９３】
そこで、「誤差伝播により計算結果に対する精度評価を行う」本願発明の特徴の一つである、ブラックボックス処理についての入力／出力間のヤコビ行列を数値差分により求め、誤差行列と併せて出力値の不確かさの推定を適用することができる。
【００９４】
一例として結像レンズ系の組み立てにおける許容差の設定を挙げることができる。
【００９５】
図１３に、生体信号計測に適用した本発明の第５実施形態の処理手順を示す。
【００９６】
例えば、自動車運転者の運転中のモニタリングは、居眠りや疲労等に起因する事故を未然に防止する観点において重要な技術的課題であると考えられている。しかしながら人間の生体信号の検出を非侵襲で実現することが強い要求として存在する一方、安心安全の実現の手段としては、高い信頼性が一方では求められている。外乱に対する複雑な感受性を示し、本質的に揺らぎを示す生体信号を元に、信頼性高く運転者のモニタリングを行い、例えば運転者の覚醒度を定量化するためのモデリングと検証が行われている。
【００９７】
非侵襲の条件下では、運転者のモニタリングに供され得る生体信号の種類は非常に限られている。例えば、心拍、呼吸、皮膚温度、筋電位、瞳や瞼、あるいは首から頭部にかけての挙動などが挙げられる。これらを単独、あるいは組み合わせて生体モデルの入力信号とし、例えば覚醒度を出力とするモデリング、及びそれに必要な信号処理アルゴリズムを構築することが考えられる。
【００９８】
そこで、「誤差伝播により計算結果に対する精度評価を行う」本願発明の特徴の一つである、ブラックボックス処理についての入力／出力間のヤコビ行列を数値差分により求め、誤差行列と併せて出力値の不確かさの推定を適用することができる。
【符号の説明】
【００９９】
１…三次元（座標）測定機
２…駆動制御装置
５…ワーク
１１…定盤
１３…Ｘ軸ガイド
１４…Ｙ軸駆動機構
１５…Ｚ軸ガイド
１６…Ｚ軸アーム
１７…プローブ
【技術分野】
【０００１】
本発明は、誤差伝播による出力データの精度評価方法に係り、特に、三次元座標測定機（以下、単に三次元測定機又はＣＭＭとも称する）における座標測定値の不確かさを推定する際に用いるのに好適な、誤差伝播による出力データの精度評価方法に関する。
【背景技術】
【０００２】
品質保証の観点から、全ての測定機器の測定結果には不確かさが必要になっている。
【０００３】
その際、複雑なデータ処理、シミュレーションなど（以後、データ処理）の計算結果がどのくらいの精度を持っているかを評価することが重要である。データ処理では、計算の入力となる入力データ、計算処理に必要な入力パラメータなど（以後、パラメータも入力データに含める）により、計算結果である出力データが得られる。出力データがどのくらいの誤差を持つかは、データ処理の前提となる計算手法やモデルがどのくらい正しいかによる場合が多い。しかし、計算手法やモデルが十分正しく、正しい入力データに対して、正しい出力データが得られるような場合でも、入力データの持つ誤差により、出力データが誤差を持つ。以後、正しい入力データに対しては、正しい出力データが得られるようなデータ処理を前提とする。
【０００４】
このような条件では、従来、入力データ（パラメータを含む）を単に少しだけ変化させ、出力データの影響を見たり、入力データをモンテカルロシミュレーションや実験計画法によって変化させ、出力誤差を推定したり、経験により出力データの誤差を評価することが行われている。
【０００５】
特に、近年、精密計測における形状計測の中で、三次元測定機を用いる座標計測の重要性が増加している。この中では、三次元測定機によって測定された形体上の複数の測定点から、最小二乗法により実際形体が計算される。しかし、三次元測定機本体の不確かさ以外にも、測定対象、測定環境、測定戦略等の影響により、不確かさが変わるため、不確かさを正確に推定することが難しい。
【０００６】
従来は、三次元測定機の測定結果の不確かさを、経験的な方法、一連の実験による方法、バーチャルＣＭＭ手法（モンテカルロシミュレーションによる方法）で推定していた（特許文献１、非特許文献１乃至８参照）。
【先行技術文献】
【特許文献】
【０００７】
【特許文献１】特開２００６−２１５０１１号公報
【非特許文献】
【０００８】
【非特許文献１】高増潔：三次元座標測定概論−精密測定の条件と三次元座標測定の不確かさ−、精密工学会講習会「差のつく三次元座標測定機活用法」、２００７年１２月
【非特許文献２】高増潔：ＧＰＳ（製品の幾何特性仕様）に基づいた検証方法−幾何公差の解釈と三次元測定機の不確かさ評価、設計工学、４０（２），２００５，７９−８５
【非特許文献３】高増潔：三次元測定機の不確かさ推定方法とその国際規格化、機械と工具、４５（３），２００１，６６−７０
【非特許文献４】高増潔，佐藤理，下嶋賢，古谷涼秋：座標測定機のアーティファクト校正（第３報）−校正後の測定の不確かさの推定−，精密工学会誌，７１（７），２００５，８９０−８９４
【非特許文献５】高増潔，阿部誠，古谷涼秋，大園成夫：形体計測における不確かさの見積り（第１報）−校正作業で生じる系統誤差の寄与−，精密工学会誌，６７（１），２００１，ｐｐ． ９１−９５
【非特許文献６】高増潔，阿部誠，古谷涼秋，大園成夫：形体計測における不確かさの見積り（第３報）−区分多項式による直線形体の形状偏差モデル−，精密工学会誌，７１（１１），２００５，１４５４−１４５８
【非特許文献７】高増潔，野坂健一郎，阿部誠，古谷涼秋，大園成夫：形体計測における不確かさの見積り（第２報）−円形体の計測における形状偏差の寄与−，精密工学会誌，６９（５），２００３，６９３−６９７
【非特許文献８】高増潔，古谷涼秋，大園成夫：最小領域法の座標計測における統計的評価，精密工学会誌，６４（４），１９９８，ｐｐ．５５７−５６１
【発明の概要】
【発明が解決しようとする課題】
【０００９】
しかし、これらの手法のうち、モンテカルロシミュレーションでは、誤差を含んだ入力データを効率よく発生することが難しく、計算量も多く、結果の信頼性も低い。他の方法では、経験的な側面が多く、一般的に評価する手法となっていない。
【００１０】
特に三次元測定機のような複雑な測定機では、不確かさの推定が難しく、実験のコスト、時間がかかり、推定精度も低いことから、まだ、実際には有効な不確かさ推定手法は確立されていなかった。
【００１１】
本発明は、前記従来の問題点を解消するべくなされたもので、理論的で、簡単な計算で行え、計算量が少なく、収束計算を行う必要が無く、データ処理を行う多くのソフトウェアに適用可能であり、誤差の要因の分析、要因毎の大きさの比較等の多くの解析が可能な出力データの精度評価方法を提供することを課題とする。
【課題を解決するための手段】
【００１２】
本発明は、図１に示す如く、データ処理を入力データによって数値微分することで、データ処理のヤコビ行列Ｊを求める手順（ステップ１１０）と、入力データの持つ誤差の分散共分散を推定して、入力データの誤差行列Ｄを求める手順（ステップ１２０）と、前記ヤコビ行列Ｊ及び入力データの誤差行列Ｄより、出力データの持つ誤差の分散共分散を表す誤差行列Ｒを計算する手順（ステップ１３０）と、を含むことにより、前記課題を解決したものである。
【００１３】
以下、各手順について説明する。
【００１４】
１．データ処理のヤコビ行列Ｊを求める手法（ステップ１１０）
データ処理では、入力データ列（入力パラメータを含める）に対して、出力データ列が得られる。まず、計算したい条件に合った入力データ列Ｉを用意し、その出力データ列Ｏを得る。ここで、入力データは、ベクトル量についてはスカラー量として分割し、全部でｎ個ある場合は、入力データ列Ｉはｎ個のスカラー量の並びとして表現する。同様に出力データも、ｍ個のスカラー量の並びとして表現する。次式参照。
【００１５】
【数１】
【００１６】
データ処理のヤコビ行列Ｊは、ｎ個の入力データとｍ個の出力データの関係を示すため、ｎ行ｍ列の行列として表現できる。ヤコビ行列Ｊのｎ行ｍ列の要素jn,mは、以下の計算で求めることができる。
【００１７】
１）入力データ列Ｉのｋ番目の入力だけを微小量ｄｋだけ変化させた入力データ列Ｉｋを作る。
【００１８】
２）入力データ列Ｉｋに対してデータ処理を行い、出力データ列Ｏ＊を求める。
【００１９】
３）ヤコビ行列Ｊのｋ行の要素jk,1〜jk,mは、出力データ列Ｏ＊とＯの各要素の差とｄｋを利用して、次式で求めることができる。
【００２０】
【数２】
【００２１】
４）同様にヤコビ行列の全ての行（１行からｎ行）を求めることができる。
【００２２】
上記の計算は、機械的に行うことができる。上記処理を行うプログラムからデータ処理をｎ＋１回行うことで、自動的にヤコビ行列Ｊを計算できる。
【００２３】
ここで、微小量ｄｋは、以下のことを考慮して決定する。
【００２４】
１）データ処理の非線形性が小さい場合は、ｄｋをどのように決めても、結果に大きな影響はない。
【００２５】
２）次の２．項で推定する入力データの誤差行列Ｄを利用すれば、その入力データの分散の平方根に対して、何分の一かを利用するなど、小さい量を設定することで、自動的に決定できる。
【００２６】
３）入力データの性質がある程度分かっている場合には、その誤差などに対して小さい量を設定する。
【００２７】
２．入力データの誤差行列Ｄの推定（ステップ１２０）
入力データが持つ誤差の性質を表す誤差行列Ｄを推定するためには、まず、入力データに影響する独立した誤差要因に対して、それぞれの誤差行列を求めて単純に和を求めればよいので、独立した1つの誤差要因に対する誤差行列の求め方を示す。
【００２８】
まず、入力データの誤差がそれぞれ独立で、その標準偏差ｔが分かっている場合は、誤差行列Ｄの対角要素がそれぞれの誤差の標準偏差の二乗（分散）となり、非対角要素は０となる。次式参照。
【００２９】
【数３】
【００３０】
入力データの誤差が独立でない場合は、２つの入力データ間の相関係数ｓ及び誤差の標準偏差ｔが推定できれば、相関係数ｓと標準偏差ｔから共分散を計算することができる。誤差行列Ｄの対角要素はそれぞれの誤差の標準偏差の二乗（分散）となり、非対角要素は対応する２つの入力データ間の共分散となる。次式参照。
【００３１】
【数４】
【００３２】
より複雑なデータの関係に関しては、ｉ）理論的な推定、ii）この部分だけをモンテカルロシミュレーションで行う、iii）経験によるなどの方法で推定する必要がある。又、次の３．項に示すように、本手法をより単純な入力データへ適用することで、誤差行列を推定することも可能である。
【００３３】
３．出力データの誤差行列Ｒの計算（ステップ１３０）
上記手法により求めたヤコビ行列Ｊ及び入力データの誤差行列Ｄより、次式の計算方法で出力データの誤差行列Ｒを計算できる。
【００３４】
【数５】
ここで、Ｊｔは行列Ｊの転置行列を示す。
【発明の効果】
【００３５】
本発明によれば、次のような効果を奏する。
【００３６】
１．出力データの誤差行列Ｒの利用
出力データの誤差行列Ｒは、このデータ処理において計算された結果の誤差を推定している。誤差行列Ｒのｊ行ｊ列の要素rj,jは、出力データ列Ｏのｊ番目の出力値ｏｊの誤差の分散を示している。そこで、出力データ列Ｏのｊ番目の出力値誤差は、標準偏差で誤差行列Ｒのｊ行ｊ列の要素rj,jの平方根tjとなる。次式参照。
【００３７】
【数６】
【００３８】
入力データの持つ誤差が正規分布であることが仮定でき、データ処理の非線形性が少ない場合には、出力データの誤差も正規分布であることが仮定できる。このような場合には、ｏｊ±ｔｊの範囲に出力データは、６０％の確率で入ることになり、ｏｊ±２ｔｊの範囲に９５％の確率で入ることになる。このように、出力データの誤差を評価することができる。
【００３９】
出力データの誤差行列Ｒのｊ行ｌ列の要素rj,lは、ｊ番目の出力とｌ番目の出力の共分散を示している。rj,lとそれぞれの出力の誤差の分散rj,j及びrl,lにより、式（３）を使って、２つの出力の誤差の相関係数sj,lを計算することができる。この相関係数は、±１の範囲で、絶対値が０に近い場合は２つの出力の誤差に相関がないこと、絶対値が１に近い場合は２つの出力の誤差に正又は負の相関が大きいことを示す。相関が大きい場合には、２つの出力値の演算は単に出力値を使うことはせずに、演算するデータ処理に対して、本手法を再び適用することで正しい評価が可能となる。
【００４０】
上記のように、出力データが別のデータ処理の入力データとなる場合には、出力データの誤差行列を次のデータ処理の入力データの誤差行列として使用することができる。
【００４１】
２．ヤコビ行列Ｊの利用
ヤコビ行列Ｊは、入力データと出力データの関係を示している。入力データのｉｋだけがｄｋだけ変化した場合の出力データolの変化はjk,l × dkで評価できる。これが十分小さい場合には、入力データｉｋは出力データolにほとんど関係ないことが分かる。この解析により、入力データと出力データが無関係な部分を抽出することができ、データ処理を独立したいくつかの部分に分離することが可能となる。
【００４２】
３．入力データの誤差行列Ｄと出力データの誤差行列Ｒの関係
入力データの誤差行列がｐ個の複数の独立した誤差要因から構成されているとすると、それぞれの要因毎の誤差行列Ｄ1〜Ｄpが推定され、それぞれに対応した出力データの誤差行列Ｒ1〜Ｒpが次式のように計算される。
【００４３】
【数７】
【００４４】
それぞれのＲに対して１．項に示した解析を行うことができる。又、Ｒの大きさを比較することで、誤差要因のどれが重要かどうかを判断できる。
【００４５】
以上をまとめると、従来の複雑なデータ処理やシミュレーションでは、計算結果を得るまでで終わっている場合が多く、計算結果の評価ができていない。これに対して、本発明では、モンテカルロシミュレーションではなく、誤差伝播を利用することによって、どのようなデータ処理に対しても汎用的に適用できる。本発明によれば、自動的に誤差評価が可能で、その結果に対して種々の分析が行なえ、データ解析に対する強力なツールとなる。
【図面の簡単な説明】
【００４６】
【図１】本発明による基本的な処理手順を示す流れ図
【図２】本発明の第１実施形態の適用対象である三次元測定機の一例を示す斜視図
【図３】第１実施形態における測定点と測定値の関係（２次元での例）を示す図
【図４】同じくソフトウェアの誤差伝播による不確かさ推定の基本的な流れを示す図
【図５】同じく三次元測定機の運動学による不確かさの概念を示す図
【図６】同じくプロービングによる誤差の推定を示す図
【図７】同じく円形体の形状偏差による不確かさを示す図
【図８】フーリエ変換に適用した本発明の第２実施形態におけるフーリエ変換プログラムを示す図
【図９】同じくフーリエ変換プログラムの入出力データ及び行列の例を示す図
【図１０】同じく出力の例を示す図
【図１１】条件探索による形体計測に適用した本発明の第３実施形態の処理手順を示す流れ図
【図１２】光学シミュレーションに適用した本発明の第４実施形態の処理手順を示す流れ図
【図１３】生体信号計測に適用した本発明の第５実施形態の処理手順を示す流れ図
【発明を実施するための形態】
【００４７】
以下図面を参照して、本発明の実施形態を詳細に説明する。
【００４８】
本発明の第１実施形態は、本発明を、図２に例示するような三次元測定機１の不確かさ推定に適用したものである。
【００４９】
図において、２は、三次元測定機１を駆動制御すると共に三次元測定機１から必要な測定座標値を取込むための駆動制御装置、３は、駆動制御装置２を介して三次元測定機１を手動操作するためのジョイスティック操作盤、４は、駆動制御装置２での測定手順を指示するパートプログラムを編集・実行すると共に、駆動制御装置２を介して取込まれた測定座標値に幾何形状を当てはめるための計算を行ったり、パートプログラムを記録、送信する機能を備えたホストシステムである。
【００５０】
前記三次元測定機１は、次のように構成されている。即ち、除振台１０の上には、定盤１１がその上面をベース面として水平面と一致するように載置され、この定盤１１の両側端から立設されたアーム支持体１２ａ、１２ｂの上端でＸ軸ガイド１３を支持している。アーム支持体１２ａは、その下端がＹ軸駆動機構１４によってＹ軸方向に駆動され、アーム支持体１２ｂは、その下端がエアーベアリングによって定盤１１上にＹ軸方向に移動可能に支持されている。Ｘ軸ガイド１３は、垂直方向に延びるＺ軸ガイド１５をＸ軸方向に駆動する。Ｚ軸ガイド１５には、Ｚ軸アーム１６がＺ軸ガイド１５に沿って駆動されるように設けられ、Ｚ軸アーム１６の下端に接触式のプローブ１７が装着されている。このプローブ１７が定盤１１上に載置されたワーク５に接触したときに、プローブ１７から駆動制御装置２にタッチ信号が出力され、そのときのＸＹＺ軸座標値を駆動制御装置２が取り込むようになっている。
【００５１】
１．基本的な手法の概要
以下、２次元での例として図３を利用して説明するが、手法そのものは、簡単に三次元へ拡張可能である。図３でＰ１からＰ１７までは、測定点であり、Ｄ１からＤ１３までが測定値である。この例では、Ｐ１からＰ４によりワーク座標系を設定し、Ｐ５からＰ８で２つの点の位置と辺の角度、Ｐ９及びＰ１０で辺の位置、Ｐ１１からＰ１４及びＰ１５からＰ１７で２つの円の中心と直径を測定している。測定値は例えば、Ｄ１とＤ２は点のＸ座標及びＹ座標、Ｄ３は２つの辺の角度、Ｄ８、Ｄ９及びＤ１０は円の中心のＸ座標、Ｙ座標及び直径を表す。
【００５２】
本実施形態は、図４に示す如く、基本的に誤差伝播に基づいて測定の不確かさを推定する。次式のように、測定値Ｄは、測定処理をするソフトウェア（測定用ソフトウェアと称する）Ｆの出力として測定点座標Ｐを入力することで得られる。
【００５３】
【数８】
【００５４】
このソフトウェアを次式のようにＰで偏微分することでヤコビ行列Ａを計算できる。
【００５５】
【数９】
【００５６】
この計算には、測定用ソフトウェアをブラックボックスとして利用する。このため、測定用ソフトウェアがどんなものであれ、オフラインで入力として測定点座標、出力として測定値が得られればよい。測定点座標の不確かさ行列Ｓが推定できれば、測定値の不確かさは次式の誤差伝播によって計算できる。
【００５７】
【数１０】
【００５８】
ヤコビ行列Ａは、式（９）によりソフトウェアをオフラインで利用して、機械的に計算できる。ヤコビ行列の計算方法は後述する。
【００５９】
そこで、この手法のキーとなる技術は、誤差行列Ｓの推定である。測定点座標の持つ不確かさが互いに独立で同じであるような場合は単純であるが、低レベルの推定にしかならない。３．項で後述するように、多くの不確かさ要因を解析することで、できるだけ実際に近い誤差行列の推定を行う。
【００６０】
ヤコビ行列Ａが計算でき、不確かさ行列Ｓが推定できれば、測定値Ｄの不確かさ行列Ｔ（それぞれの測定値の分散、共分散）は、式（１０）で計算できる。式（１０）の形は、Ｆが観測方程式の場合とは違うことに注意する必要がある。Ｆが観測方程式の場合は、Ｆの入力は測定値だが、出力は測定結果ではなく求めたい形体からの偏差になっている。そして、偏差の二乗和を最小にするように最小二乗法を解いた結果が、ここでのＦと等しくなる。そこでＦが観測方程式の場合のヤコビ行列をＢとするとＡｔ＝（ＢｔＳ−１Ｂ）−１ＢｔＳ−１となる。しかし、一般の三次元測定機のプログラムは最小二乗法を解くときにＳを考慮しないので、実際はＡ＝（ＢｔＢ)−１Ｂとなる。
【００６１】
２．ヤコビ行列の推定
ヤコビ行列Ａは、次式に示すように、各測定点ＰiのＸ座標Ｐix，Ｙ座標Ｐiyと測定値Ｄjの関係について、数値的な差分を取ることで求める。
【００６２】
【数１１】
【００６３】
これにより、測定手順及びソフトウェア部分をブラックボックスＦとして扱うことができる。従来は、理論的な計算又はモンテカルロ法で計算を行なっていたが、測定手順及びソフトウェア部分をオフラインモード（三次元測定機を使わずに、仮想の測定を行なうモード）で実施し、ブラックボックスとして扱うことができ、以下の計算を、機械的に行えるため理論的な難しさはない。
【００６４】
ヤコビ行列Ａは、測定点の座標の不確かさと測定値の不確かさの伝播を示すため、ｎ個の測定値（２次元なので２ｎ個の座標値）とｍ個の測定値に対応して、２ｎ行ｍ列の行列となる。図３の例では、ｎ＝１７、ｍ＝１３となるので、ヤコビ行列Ａは３４行１３列の行列である。ｉ番目の測定点Ｐiの実際の測定座標を（Ｐix,Ｐiy）とし、ｊ番目の計算された測定値をＤjとする。オフラインにおいて、ｉ番目の測定点のＸ座標だけをｄだけ増やしてＰiの座標を（Ｐix＋ｄ，Ｐiy）とし、同じ計算を行ない、そのときの計算された測定値をＤj,ixとする。この計算から、ヤコビ行列Ａの２ｉ−ｌ行の要素が（Ｄj,ix−Ｄj）／ｄとして計算できる。
【００６５】
実際の測定座標は、三次元測定機で測定しなくても、名目的な値（例えば図面上に記載されている値）を使っても構わない。又、ｄの値は数値誤差が問題にならない範囲で、十分小さい値を使うことが重要である。一般的には、例えば0.000000001mmといった値を使うことができる。
【００６６】
以上の計算を各測定点のＸ座標、Ｙ座標について計算すれば、ヤコビ行列Ａを計算できる。この場合、ソフトウェア部分は２ｉ回呼び出されることになる。実際には、ソフトウェアは、部分的に独立なため、全ての計算を繰り返す必要はないが、ソフトウェアの独立性についての判断方法は今後の課題とする。ヤコビ行列を計算した後なら、ヤコビ行列の共分散を見ることによって、自動的にソフトウェアの独立部分を抜き出すことができる。
【００６７】
３．測定点座標の不確かさ行列の推定
各測定点のＸ座標、Ｙ座標のそれぞれが持つ不確かさ（分散）と各測定点間の座標の関係（共分散）を推定することが、提案手法のキーとなる技術である。まず、各測定点の座標が持つ不確かさとして、以下の要因を考えることとする。このうち、いくつかはすでに定式化されている。
【００６８】
（１）ＣＭＭの運動学補正による不確かさ
（２）プロービング（プロービング、プローブ補正及びプローブ交換）による不確かさ
（３）測定物の形状偏差による不確かさ
（４）測定環境（温度ドリフトなど）による不確かさ
【００６９】
基本的な考え方は、それぞれの要因によって与えられる不確かさを求め、その要因がそれぞれの座標値に与える相関係数から共分散を推定する。以下、ひとつひとつの要因に対して、分散と共分散の推定方法の概要を示す。
【００７０】
（１）ＣＭＭの運動学補正による不確かさ
ＣＭＭは運動学に基づいた補正を行っているが、補正しきれていない高周波のランダムな不確かさと低周波の座標間の相関を持つ不確かさが存在する。この不確かさでは、Ｘ座標とＹ座標との相関は考える必要はないと思われる。座標の不確かさは全ての場所で一定の値を持つ。又、２つの座標間の相関は、２つの座標の差の絶対値の関数となる（非特許文献４参照）。
【００７１】
差し当たって、それぞれの測定点のＸ座標、Ｙ座標に不確かさscmmを与え、２つの測定点ｈとｉの座標の差の絶対値によって相関rcmm(|Phx - Pix|)を与える。Ｘ座標とＹ座標との間の相関は０とする。scmmとrcmmは、ＣＭＭの長さ測定における最大許容誤差のＭＰＥ線図、もしくはステップゲージ等の測定結果から推定する。
【００７２】
ＣＭＭの運動学による不確かさの概念を図５に示す
【００７３】
（２）プロービング（プロービング、プローブ補正およびプローブ交換）による不確かさ
プロービングによる不確かさは、ランダムな不確かさと測定ベクトルによる測定方向に依存する不確かさがある。又、２つの点を同じプローブで測定した場合、プローブの補正による不確かさの相関は２つの測定ベクトルの方向の差の関数によって記述できる（非特許文献５参照）。
【００７４】
まず、プロービングについては不確かさspを与え、プローブ径の補正を含めて、２つの測定点ｈとｉの測定ベクトルのなす角との関係で相関rp(Vh・Vi)を与える。これらの値は、理論的な考察及び（１）のscmmとrcmmから推定できる。
【００７５】
プロービングによる誤差の推定を図６に示す。
【００７６】
（３）被測定物の形状偏差による不確かさ
被測定物に形状偏差がある場合、形状偏差の持つランダムな不確かさと形状の自己相関による不確かさがある。２つの点が同じ被測定物を測定した場合、測定ベクトルの方向及び測定位置の関数として、不確かさの相関が記述できる（非特許文献６、７参照）。
【００７７】
又、形状偏差の持つランダムな誤差及び相関により、測定点の不確かさsf及び相関rfを推定できる（非特許文献６、７参照）。相関rfは、円の場合は測定ベクトルから推定される測定点の円上での角度位置の差で決まり、直線の場合は、直線上での測定点の位置の差で決まる。
【００７８】
円形体の形状偏差による不確かさを図７に示す。
【００７９】
（４）測定環境（温度ドリフトなど）による不確かさ
環境の影響による不確かさは評価が難しいが、温度ドリフトなどの環境変化による不確かさを時間の関数として評価することが考えられる。２つの測定点の測定した時刻の差が大きいと、環境が変化し相関が低くなる。一方、２つの測定点を連続して測定した場合には、環境の相関が高いため、環境による不確かさが打ち消しあうことになる。
【００８０】
本実施形態によれば、ヤコビ行列及び不確かさ行列から誤差伝播により、測定結果の不確かさを推定する手法で、理論的に高い推定精度が期待できる。又、ヤコビ行列の計算手法が、ソフトウェアをブラックボックスとして数値差分で行なえるため、高速で精度が高く、経験や実験が不要である。更に、測定座標の不確かさ行列を推定する手法が、理論的な解析に基づき、多くの要因をカバーできるため、最終的な測定結果の不確かさの推定精度が高い。又、測定座標の不確かさ行列を推定する方法が、理論的な方法のため、実験等がほとんど不要で、コストが低く、推定時間が速い。
【００８１】
なお、測定座標の不確かさ行列Ｓの推定要因として、次の６つを考えることもできる。
【００８２】
（１）三次元測定機の直角誤差による不確かさＳ１
（２）三次元測定機のスケール誤差による不確かさＳ２
（３）三次元測定機の直角誤差及びスケール誤差以外の運動学誤差による不確かさＳ３
（４）プロービングによる不確かさＳ４
（５）被測定物の形状偏差による不確かさＳ５
（６）測定環境（運動ドリフト等）による不確かさＳ６
【００８３】
次に、フーリエ変換に適用した本発明の第２実施形態を説明する。
【００８４】
この第２実施形態では、８個のデータをフーリエ変換して、０次から４次のパワースペクトルを求める。ここで、データに誤差がある場合、求めたパワースペクトルの誤差を推定する。
【００８５】
図８のプログラムで、ｆｆｔ８（５４行〜６１行）は、８個の入力データに対して、５個の出力データを求めるデータ処理で、中身の処理を知る必要は無い。入力データｉｉが与えられ、入力データの分散共分散行列ｄｄが推定できれば、後は機械的に計算できる。ｄｄとしては、分散が０．１２で共分散が２つのデータの距離に関係した相関を与えている（１４行〜２２行）。
【００８６】
２８行〜５２行までは、完全に機械的に計算できる部分で、対象となるデータ処理が何でも同じプログラムで計算が可能である。このプログラムの入出力データ及び行列の例を図９に、出力の例を図１０に示す。結果は、出力データが大きい次数に大きな誤差があること、相対的に低次の誤差が大きいこと等が分かる。出力データの誤差の相関係数がそれほど大きくないのはＦＦＴが直交的な処理だからと思われる。
【００８７】
図１１に、条件検索による形態計測に適用した本発明の第３実施形態の処理手順を示す。
【００８８】
例えば三次元測定機や形状測定機を用いて形体計測を行う場合、最小二乗法による評価の他に、最小領域法による評価を行うことがある。最小領域法では、条件探索を通して採択された数少ない代表点によってのみ結果が決まり、他の棄却された測定点は結果に寄与しない。そのため計算プロセスの統計的な評価が難しく、測定の不確かさを理論的に検討した研究は少ない。わずかに非特許文献８が挙げられる。この研究によると、最小領域法の不確かさを統計的に推定する手法が提案され、理論的な方法論は確立されている。
【００８９】
しかしながら最小領域法を実装した市販のソフトウェアにおいては、内部のアルゴリズムや数値演算の詳細についての情報は開示されていない。このような場合においても、本願発明による不確かさ推定法を適用することによって、最小領域法の測定の不確かさを容易に推定することができる。又最小領域法に限らず、複数の計算ステップを組み合わせた複雑なデータ処理プログラムに対しても適用することができる。
【００９０】
図１２に、光学シミュレーションに適用した本発明の第４実施形態の処理手順を示す。
【００９１】
例えば、レンズ等の光学素子の設計、及びこれらを組み合わせたシステム設計は、装置の小型化や高機能化の流れを受けて複雑化している。又、非線型光学効果を有する種々の素子の実用化や複雑な形状と光学特性を有する光共振器等についても多様なものが研究され、実用に供されている。これらの光学素子の設計や、それらを組み込んだシステム設計を適切に行うためには、光源から最終的な光学的出力にいたるまでの光電磁場の挙動を物理光学的に記述し得るシミュレーションツールが不可欠である。
【００９２】
一方、市販されている光学シミュレーションツールの内部では、解析対象に特有な光学場の条件設定に依存し、あるいはオペレータが指示する入力パラメータに依存し、より高速で、且つ数値的な精度が十分な、異なるアルゴリズムや解法を自動や手動で使いわけている。こうしたシミュレーションプログラムの内部情報は、全てがユーザに開示されているわけではなく、ユーザは自身が実施しているシミュレーションの内部処理を把握せずに、いわばブラックボックスとして使っていることになる。
【００９３】
そこで、「誤差伝播により計算結果に対する精度評価を行う」本願発明の特徴の一つである、ブラックボックス処理についての入力／出力間のヤコビ行列を数値差分により求め、誤差行列と併せて出力値の不確かさの推定を適用することができる。
【００９４】
一例として結像レンズ系の組み立てにおける許容差の設定を挙げることができる。
【００９５】
図１３に、生体信号計測に適用した本発明の第５実施形態の処理手順を示す。
【００９６】
例えば、自動車運転者の運転中のモニタリングは、居眠りや疲労等に起因する事故を未然に防止する観点において重要な技術的課題であると考えられている。しかしながら人間の生体信号の検出を非侵襲で実現することが強い要求として存在する一方、安心安全の実現の手段としては、高い信頼性が一方では求められている。外乱に対する複雑な感受性を示し、本質的に揺らぎを示す生体信号を元に、信頼性高く運転者のモニタリングを行い、例えば運転者の覚醒度を定量化するためのモデリングと検証が行われている。
【００９７】
非侵襲の条件下では、運転者のモニタリングに供され得る生体信号の種類は非常に限られている。例えば、心拍、呼吸、皮膚温度、筋電位、瞳や瞼、あるいは首から頭部にかけての挙動などが挙げられる。これらを単独、あるいは組み合わせて生体モデルの入力信号とし、例えば覚醒度を出力とするモデリング、及びそれに必要な信号処理アルゴリズムを構築することが考えられる。
【００９８】
そこで、「誤差伝播により計算結果に対する精度評価を行う」本願発明の特徴の一つである、ブラックボックス処理についての入力／出力間のヤコビ行列を数値差分により求め、誤差行列と併せて出力値の不確かさの推定を適用することができる。
【符号の説明】
【００９９】
１…三次元（座標）測定機
２…駆動制御装置
５…ワーク
１１…定盤
１３…Ｘ軸ガイド
１４…Ｙ軸駆動機構
１５…Ｚ軸ガイド
１６…Ｚ軸アーム
１７…プローブ
【特許請求の範囲】
【請求項１】
データ処理を入力データによって数値微分することで、データ処理のヤコビ行列を求める手順と、
入力データの持つ誤差の分散共分散を推定して、入力データの誤差行列を求める手順と、
前記ヤコビ行列及び入力データの誤差行列より、出力データの持つ誤差の分散共分散を表す誤差行列を計算する手順と、
を含むことを特徴とする誤差伝播による出力データの精度評価方法。
【請求項２】
前記データ処理のヤコビ行列を、
入力データ列Ｉのｋ番目の入力だけを微小量ｄｋだけ変化させた入力データＩｋを作り、
入力データ列Ｉｋに対して、データ処理を行い出力データ列Ｏ＊を求め、
出力データ列Ｏ＊とＯの各要素の差とｄｋを利用してヤコビ行列Ｊのｋ行の要素ｊｋ，ｌ〜ｊｋ，ｍを求めることにより求めることを特徴とする請求項１に記載の出力データの精度評価方法。
【請求項３】
前記入力データの誤差行列を、入力データの誤差がそれぞれ独立で、その標準偏差ｔが分かっている場合は、誤差行列Ｄの対角要素がそれぞれの誤差の標準偏差の二乗（分散）となり、非対角要素は０となるように求めることを特徴とする請求項１に記載の出力データの精度評価方法。
【請求項４】
前記入力データの誤差行列を、入力データの誤差が独立でないが、２つの入力データ間の相関係数及び誤差の標準偏差が推定できる場合は、相関係数ｓと標準偏差ｔから共分散を計算し、誤差行列Ｄの対角要素をそれぞれの誤差の標準偏差の二乗（分散）とし、非対角要素は対応する２つの入力データ間の共分散とするように求めることを特徴とする請求項１に記載の出力データの精度評価方法。
【請求項５】
前記ヤコビ行列Ｊ及び入力データの誤差行列Ｄより、次式で出力データの誤差行列Ｒ
Ｒ＝ＪｔＤＪ（ここで、Ｊｔは行列Ｊの対角行列を示す）
を計算することを特徴とする請求項１に記載の出力データの精度評価方法。
【請求項１】
データ処理を入力データによって数値微分することで、データ処理のヤコビ行列を求める手順と、
入力データの持つ誤差の分散共分散を推定して、入力データの誤差行列を求める手順と、
前記ヤコビ行列及び入力データの誤差行列より、出力データの持つ誤差の分散共分散を表す誤差行列を計算する手順と、
を含むことを特徴とする誤差伝播による出力データの精度評価方法。
【請求項２】
前記データ処理のヤコビ行列を、
入力データ列Ｉのｋ番目の入力だけを微小量ｄｋだけ変化させた入力データＩｋを作り、
入力データ列Ｉｋに対して、データ処理を行い出力データ列Ｏ＊を求め、
出力データ列Ｏ＊とＯの各要素の差とｄｋを利用してヤコビ行列Ｊのｋ行の要素ｊｋ，ｌ〜ｊｋ，ｍを求めることにより求めることを特徴とする請求項１に記載の出力データの精度評価方法。
【請求項３】
前記入力データの誤差行列を、入力データの誤差がそれぞれ独立で、その標準偏差ｔが分かっている場合は、誤差行列Ｄの対角要素がそれぞれの誤差の標準偏差の二乗（分散）となり、非対角要素は０となるように求めることを特徴とする請求項１に記載の出力データの精度評価方法。
【請求項４】
前記入力データの誤差行列を、入力データの誤差が独立でないが、２つの入力データ間の相関係数及び誤差の標準偏差が推定できる場合は、相関係数ｓと標準偏差ｔから共分散を計算し、誤差行列Ｄの対角要素をそれぞれの誤差の標準偏差の二乗（分散）とし、非対角要素は対応する２つの入力データ間の共分散とするように求めることを特徴とする請求項１に記載の出力データの精度評価方法。
【請求項５】
前記ヤコビ行列Ｊ及び入力データの誤差行列Ｄより、次式で出力データの誤差行列Ｒ
Ｒ＝ＪｔＤＪ（ここで、Ｊｔは行列Ｊの対角行列を示す）
を計算することを特徴とする請求項１に記載の出力データの精度評価方法。
【図１】
【図２】
【図５】
【図３】
【図４】
【図６】
【図７】
【図８】
【図９】
【図１０】
【図１１】
【図１２】
【図１３】
【図２】
【図５】
【図３】
【図４】
【図６】
【図７】
【図８】
【図９】
【図１０】
【図１１】
【図１２】
【図１３】
【公開番号】特開２０１１−４７７０３（Ｐ２０１１−４７７０３Ａ）
【公開日】平成２３年３月１０日（２０１１．３．１０）
【国際特許分類】
【出願番号】特願２００９−１９４５３８（Ｐ２００９−１９４５３８）
【出願日】平成２１年８月２５日（２００９．８．２５）
【新規性喪失の例外の表示】特許法第３０条第１項適用申請有り 発行者名：社団法人精密工学会 刊行物名：２００９年度精密工学会春季大会学術講演会講演論文集 発行年月日：２００９年２月２５日
【出願人】（０００１３７６９４）株式会社ミツトヨ (979)
【Ｆターム（参考）】
【公開日】平成２３年３月１０日（２０１１．３．１０）
【国際特許分類】
【出願日】平成２１年８月２５日（２００９．８．２５）
【新規性喪失の例外の表示】特許法第３０条第１項適用申請有り 発行者名：社団法人精密工学会 刊行物名：２００９年度精密工学会春季大会学術講演会講演論文集 発行年月日：２００９年２月２５日
【出願人】（０００１３７６９４）株式会社ミツトヨ (979)
【Ｆターム（参考）】
[ Back to top ]