形状パラメータ決定方法、形状パラメータ決定装置、及びコンピュータプログラム
【課題】 ３次元形状の計測データに対応する理論形状を記述するパラメータを高速に且つ安定して決定することができるようにする。
【解決手段】 パラメータα(k)、β(k)を既知量として、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0(k+1)，ｙ0(k+1)，α(k)，β(k)）を最小とするパラメータｘ0(k+1)、ｙ0(k+1)を決定し、決定したパラメータｘ0(k+1)、ｙ0(k+1)を既知量として、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0(k+1)，ｙ0(k+1)，α(k+1)，β(k+1)）を最小とするパラメータα(k+1)、β(k+1)を決定する。
【解決手段】 パラメータα(k)、β(k)を既知量として、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0(k+1)，ｙ0(k+1)，α(k)，β(k)）を最小とするパラメータｘ0(k+1)、ｙ0(k+1)を決定し、決定したパラメータｘ0(k+1)、ｙ0(k+1)を既知量として、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0(k+1)，ｙ0(k+1)，α(k+1)，β(k+1)）を最小とするパラメータα(k+1)、β(k+1)を決定する。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、形状パラメータ決定方法、形状パラメータ決定装置、及びコンピュータプログラムに関し、特に、計測された３次元形状に対応する理論形状を規定する複数のパラメータを求めるために用いて好適なものである。
【背景技術】
【０００２】
近年、３次元形状の座標を、接触式又は非接触式などで測定する３次元形状座標計測器が発達し、広く普及されるようになった。この３次元形状座標計測器を用いて計測されたデータと、理論形状（設計値）とを比較することにより、部品の加工精度等を直接評価できるようになっている。特許文献１では、３次元座標計測器で計測されたデータに対応する理論形状を記述する理論形状パラメータを決定する方法が提示されている。
【０００３】
【特許文献１】特開平６−２８１４４０号公報
【発明の開示】
【発明が解決しようとする課題】
【０００４】
しかしながら、前述した従来の技術において、理論形状パラメータを決定する問題は、最適値探索問題である。このため、理論形状を記述するためのパラメータの数が多くなる。すなわち、解探索空間の次元が大きくなる。このように、解探索空間の次元が大きくなると、解が局所的な適合値に収束してしまい、真の最適値へ到達できなくなる可能性が大きくなる。前述した特許文献１に記載されている技術においても、局所的な適合値に解が落ち込むことへの対策は明示されていない。
【０００５】
このように、従来の技術では、３次元形状の計測データに対応する理論形状を記述するパラメータを高速に且つ安定して決定することが困難であるという問題点があった。
【０００６】
本発明は、このような問題点に鑑みてなされたものであり、３次元形状の計測データに対応する理論形状を記述するパラメータを高速に且つ安定して決定することができるようにすることを目的とする。
【課題を解決するための手段】
【０００７】
本発明の形状パラメータ決定方法は、計測された３次元形状に対応する理論形状を規定する複数のパラメータを決定する形状パラメータ決定方法であって、前記複数のパラメータを少なくとも２組に分け、前記計測された３次元形状の計測データと、その３次元形状の計測データに対応する理論形状との差を表す量が最小となるようなパラメータを、前記分けた組毎に順番に繰り返し求めることによって、前記複数のパラメータを決定することを特徴とする。
また、本発明の他の特徴とするところは、計測された３次元形状に対応する理論形状を規定する複数のパラメータを第１〜第ｎ（ｎは２以上の自然数）の組に分け、分けた組毎に前記複数のパラメータを決定する形状パラメータ決定方法であって、前記計測された３次元形状の計測データを入力する入力ステップと、前記入力ステップにより入力された計測データと、前記計測データに対応する理論形状との差を表す量が最小となるようなパラメータを、前記分けた組毎に求める第１〜第ｎのパラメータ取得ステップとを有することを特徴とする。
【０００８】
本発明の形状パラメータ決定装置は、計測された３次元形状に対応する理論形状を規定する複数のパラメータを決定する形状パラメータ決定装置であって、前記複数のパラメータを少なくとも２組に分け、前記計測された３次元形状の計測データと、その３次元形状の計測データに対応する理論形状との差を表す量が最小となるようなパラメータを、前記分けた組毎に順番に繰り返し求める手段を有することを特徴とする。
また、本発明の他の特徴とするところは、計測された３次元形状に対応する理論形状を規定する複数のパラメータを第１〜第ｎ（ｎは２以上の自然数）の組に分け、分けた組毎に前記複数のパラメータを決定する形状パラメータ決定装置であって、前記計測された３次元形状の計測データを入力する入力手段と、前記入力手段により入力された計測データと、前記計測データに対応する理論形状との差を表す量が最小となるようなパラメータを、前記分けた組毎に求める第１〜第ｎのパラメータ取得手段とを有することを特徴とする。
【０００９】
本発明のコンピュータプログラムは、計測された３次元形状に対応する理論形状を規定する複数のパラメータを決定することをコンピュータに実行させるためのコンピュータプログラムであって、前記複数のパラメータを少なくとも２組に分け、前記計測された３次元形状の計測データと、その３次元形状の計測データに対応する理論形状との差を表す量が最小となるようなパラメータを、前記分けた組毎に順番に繰り返し求めることをコンピュータに実行させることを特徴とする。
また、本発明の他の特徴とするところは、計測された３次元形状に対応する理論形状を規定する複数のパラメータを第１〜第ｎ（ｎは２以上の自然数）の組に分け、分けた組毎に前記複数のパラメータを決定することをコンピュータに実行させるためのコンピュータプログラムであって、前記計測された３次元形状の計測データを入力する入力ステップと、前記入力ステップにより入力された計測データと、前記計測データに対応する理論形状との差を表す量が最小となるようなパラメータを、前記分けた組毎に求める第１〜第ｎのパラメータ取得ステップとをコンピュータに実行させることを特徴とする。
【発明の効果】
【００１０】
本発明によれば、３次元形状の計測データに対応する理論形状を規定するパラメータが多数の場合であっても、高速に且つ安定してパラメータを決定することができるので、三次元座標計測器を用いて計測した形状の評価が効率的に行えるようになる。
【発明を実施するための最良の形態】
【００１１】
（第１の実施形態）
次に、図面を参照しながら、本発明の第１の実施形態について説明する。
なお、本実施形態では、３次元形状として回転対称体を例に挙げて説明する。また、回転対称体の中心軸（対称軸）の方向は、ｚ軸と略同じ方向を向いているとしても任意性は失われないので、以下、中心軸は、ｚ軸と略同じ方向を向いているとする。
【００１２】
図１は、３次元空間内での形状計測データ（座標値）と、その形状計測データに対応する理論形状を記述するパラメータとの関係の一例を示す図である。なお、形状計測データは、３次元形状座標計測器で測定されたデータである。
【００１３】
３次元形状座標計測器の一例としては、プローブと、測定台にセットされた被測定物の表面に前記プローブをトレースさせるトレース機構と、前記トレース機構によりトレースされたプローブの位置情報を処理して被測定物の形状を測定するコンピュータとを有するものが挙げられる。
【００１４】
図１において、形状計測データ（測定点）はＰi（ｘi,ｙi,ｚi）である。元の理論形状（例えば設計値）１００は、中心軸（ここではｚ軸）に沿って回転半径Ｒ（ｚ）で与えられているとする。この元の理論形状１００の基準点Ｏ（０，０，０）は、測定空間では、基準点Ｃ（ｘ0,ｙ0,ｚ0）へ移動（シフト）している。また、測定空間では、元の理論形状１００の中心軸（ｚ軸）に対し、方向ベクトルｎ（＝（α,β,１））の方向に中心軸１０１が傾いている。
【００１５】
すなわち、図１に示した例では、形状計測データＰiに対応する理論形状を記述するパラメータは、ｘ0、ｙ0、ｚ0、α、βの５つである。これらの５つのパラメータは、理論値（設計値）との差異、つまり、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0,ｙ0,ｚ0,α,β）を最小にするものとして算出する。ここで、形状計測データ（測定点）Ｐiから中心軸１０１へおろした垂線の足をＨiとすると、以下の（１式）が成り立つ。
【００１６】
【数１】
【００１７】
ただし、Ｎは、形状計測データＰiの数である。また、以下の（２式）〜（４式）が成り立つとする。
【００１８】
【数２】
【００１９】
３次元形状の一例である前記回転対称体が円筒の場合は、回転半径Ｒ（ｚ）はｚ軸方向の値によらず、基準点Ｃ（ｘ0,ｙ0,ｚ0）は、中心軸１０１のどこにあっても等価となる。よって、基準点Ｃのｚ軸方向の値ｚ0は定数と考えてよい。そうすると、形状計測データＰiに対応する理論形状を記述するパラメータは、ｘ0、ｙ0、α、βの４つとなり、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0,ｙ0,α,β）を考えればよいことになる。
【００２０】
また、回転半径Ｒ（ｚ）はｚ軸方向の値によらないが、設計値からは変化している場合は、回転半径Ｒも、形状計測データＰiに対応する理論形状を記述するパラメータに加えてもよい。この場合、形状計測データＰiに対応する理論形状を記述するパラメータはｘ0、ｙ0、Ｒ、α、βの５つとなり、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0,ｙ0,Ｒ,α,β）を考えればよいことになる。
【００２１】
いずれの場合（回転対称体が円筒の場合、及び回転対称体が、設計値と異なる回転半径を有する円筒である場合）についても、中心軸１０１までの距離の誤差の２乗和Ｓを最小にするパラメータを決定することは変わらないが、本実施形態では、回転対称体が、前述した円筒であり、回転半径Ｒ（ｚ）が一定である場合を説明する。
【００２２】
以下、図２のフローチャートを参照しながら、パラメータを決定する方法の一例を説明する。前述したように、本実施形態では、回転対称体が、前述した円筒であり、回転半径Ｒ（ｚ）が一定の場合、すなわち、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0,ｙ0,α,β）を最小にする場合を説明する。
【００２３】
なお、ここでは、形状計測データＰiに対応する理論形状を記述するパラメータｘ0、ｙ0、α、βを、２つの組に分けるようにする。例えば、パラメータｘ0、ｙ0の組と、パラメータα、βの組とに分けるようにする。このようにパラメータｘ0、ｙ0、α、βが、２つの組に分けられる設定がなされた上で、以下の処理が実行される。
【００２４】
（１） 形状計測データＰiを入力する（ステップＳ２０１）。
（２） 各パラメータｘ0、ｙ0、α、βの初期値ｘ0(0)、ｙ0(0)、α(0)、β(0)を設定する（ステップＳ２０２）。
ここで、右肩の添え字は反復回数を示すこととする。０回目の値だということで初期値の右肩には（０）が付いている。
【００２５】
（３） パラメータを反復計算する（ｋ（ｋは自然数）回目の値をもとに、（ｋ＋１）回目の値を決定する）。
（３−１） パラメータα(k)、β(k)を既知量として、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0(k+1)，ｙ0(k+1)，α(k)，β(k)）を最小とするパラメータｘ0(k+1)、ｙ0(k+1)を決定する。すなわち、以下の（５式）の連立方程式からパラメータｘ0(k+1)、ｙ0(k+1)を求める（ステップＳ２０３）。
【００２６】
【数３】
【００２７】
（３−２） 求めたパラメータｘ0(k+1)、ｙ0(k+1)を既知量として、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0(k+1)，ｙ0(k+1)，α(k+1)，β(k+1)）を最小とするパラメータα(k+1)、β(k+1)を決定する。すなわち、以下の（６式）の連立方程式から、パラメータα(k+1)、β(k+1)を求める（ステップＳ２０４）。
【００２８】
後述するように、本実施形態では、被測定物の形状を測定するときに工夫をして、中心軸１０１がほぼｚ軸に平行となる位置で測定されたとする。そうすると、中心軸１０１の方向を記述するパラメータα、βは、１に比べて十分小さな量となり、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0，ｙ0，α，β）は、α、βの２次形式で表される。したがって、このステップＳ２０４では、線形連立方程式の問題を解いて、パラメータα(k+1)、β(k+1)を算出することになる。
【００２９】
【数４】
【００３０】
（３−３） ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下であるか否かを判定する（ステップＳ２０５）。この判定の結果、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下となり、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値とが十分に近づいたならば、パラメータｘ0(k+1)、ｙ0(k+1)、α(k+1)、β(k+1)を求める処理を終了して、ステップＳ２０６（（４）結果を出力する）へ進む。
【００３１】
一方、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値とが十分に近づいていない場合には、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値とが十分に近づくまで、ステップＳ２０３〜Ｓ２０５を繰り返し行う。
【００３２】
（４）（ｋ＋１）回目に求めたパラメータｘ0(k+1)、ｙ0(k+1)、α(k+1)、β(k+1)の値を、形状計測データＰiに対応する理論形状を記述するパラメータとして出力する（ステップＳ２０６）。
【００３３】
なお、前述した説明では、４つのパラメータｘ0、ｙ0、α、βを、パラメータｘ0、ｙ0の組と、パラメータα、βの組に分けたが、組の分け方はこれに限ったものではない。複数の組に分けることで各ステップ（ステップＳ２０３、Ｓ２０４）における最適解の探索空間（最適値探索問題の解空間）を狭めることができれば、他の分け方でもよい。
【００３４】
また、前述した（３−１）と、（３−２）のステップＳ２０３、Ｓ２０４におけるパラメータの最適値の探索は、一般には非線形の最小自乗法である。したがって、各ステップＳ２０３、Ｓ２０４内でさらに反復計算が必要となる。
ところが、前述したように、本実施形態では、被測定物の形状を測定するときに工夫をして、中心軸１０１がｚ軸に略平行となる位置で測定されるので、中心軸１０１の方向を記述するパラメータα、βは、１に比べて十分小さな量であることが期待される。その場合、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0，ｙ0，α，β）は、α、βの２次形式で表される。そうすると、（３−２）のステップＳ２０４では、線形連立方程式の問題を解けばよいことなり、極めて簡単にパラメータα、βを算出することができる。
【００３５】
つまり、パラメータの組の分け方としては、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓを、２次形式で表現できるパラメータの組とそうでない組とに分けることが望ましい。このようにしてパラメータを分けることによって、各ステップＳ２０３、Ｓ２０４における最適解の探索空間を狭めるだけでなく、線形探索となり、計算コストと、安定性とが著しく向上する。
【００３６】
また、ここでは、計算順序として、（３−１）のステップＳ２０３を行った後に、（３−２）のステップＳ２０４を行う場合を例に挙げて説明したが、初期値の与え方によっては、（３−２）のステップＳ２０４を行った後に、（３−１）のステップＳ２０３を行ってもよい場合がある。
【００３７】
（実施例）
次に、本実施形態の実施例として、円筒の側面上に配置されるべき形状計測データＰiから中心軸を決定する方法の一例を示す。本実施例では、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0，ｙ0，α，β）が最小となるパラメータｘ0、ｙ0、α、βを決定した。ここで、３次元形状計測器では、中心軸１０１がｚ軸と略平行になるようにして測定されている。よって、パラメータα、βは１よりも十分小さい（α、β≪１）とみなしてよい。そうすると、形状計測データＰiから中心軸までの距離の誤差の２乗和Ｓ（ｘ0，ｙ0，α，β）を、以下の（７式）のように書くことができる。（７式）に示すように、形状計測データＰiから中心軸までの距離の誤差の２乗和Ｓは、４つのパラメータｘ0、ｙ0、α、βのうち、パラメータα、βについての２次形式と見ることができる。
【００３８】
【数５】
【００３９】
そこで、パラメータｘ0、ｙ0、α、βの組の分け方としては、少なくとも１つの組は、パラメータα、βの中から選択されたパラメータのみを含む組にしたい。そうすることによって、その組のパラメータを算出するための計算は線形計算となる。本実施例においては４つのパラメータｘ0、ｙ0、α、βを、パラメータα、βの組と、パラメータｘ0、ｙ0の組とに分けるようにした。そして、パラメータα、βに関しては線形方程式を解くことで求め、パラメータｘ0、ｙ0に関しては非線形方程式を解くことで解を求めた。ここで、非線形方程式の反復解法としては、ニュートン・ラフソン法を用いた。
【００４０】
計算精度を確認するために、パラメータα、β、ｘ0、ｙ0が、それぞれ０．１０、０．０５、１．０、２．０（α＝０．１０、β＝０．０５、ｘ0＝１．０、ｙ0＝２．０）となる半径Ｒが３（Ｒ＝３）の円筒の側面上に配置された６０点の擬似形状計測データを生成し、生成した擬似形状計測データから、パラメータα、β、ｘ0、ｙ0の値を逆算した。図３に、３次元座標上にプロットされた擬似形状計測データの一例を、図４に、入力した擬似形状計測データの値の一例をそれぞれ示す。
【００４１】
各パラメータα、β、ｘ0、ｙ0の初期値α(0)、β(0)、ｘ0(0)、ｙ0(0)を、それぞれ０．０、０．０、１．１９８７６、２．０９９３８（α(0)＝０．０、β(0)＝０．０、ｘ0(0)＝１．１９８７６、ｙ0(0)＝２．０９９３８）として計算した結果、パラメータα、β、ｘ0、ｙ0の値は、それぞれ０．０９９３８７２、０．０４９６９３５、１．００１２２、２．０００６１（α＝０．０９９３８７２、β＝０．０４９６９３５、ｘ0＝１．００１２２、ｙ0＝２．０００６１）となり、本実施形態の手法が有効であることが確認できた。
【００４２】
すなわち、３次元形状の計測データに対応する理論形状を規定するパラメータが多数の場合であっても、高速に且つロバスト（robust）にパラメータを決定することができるので、接触式、非接触式などの３次元座標計測器を用いて計測された３次元形状の評価を効率的に行うことができる。具体的には、回転対称体の形状計測データＰiの中心軸１０１を正しく算出し、傾いて測定された形状計測データＰiに関して高速に且つロバスト（robust）にその傾きを補正することができる。
【００４３】
以上のようにしてパラメータｘ0、ｙ0、α、βを解く処理を実現するハードウェア資源としては、例えば図５に示す情報処理装置が挙げられる。
図５は、前記情報処理装置の構成の一例を示したブロック図である。
図５において、情報処理装置１２０は、ＣＰＵ１２１と、ＲＯＭ１２２と、ＲＡＭ１２３と、キーボード（ＫＢ）１２４のキーボードコントローラ（ＫＢＣ）１２５と、表示部としてのＣＲＴディスプレイ（ＣＲＴ）１２６のＣＲＴコントローラ（ＣＲＴＣ）１２７と、ハードディスク（ＨＤ）１２８及びフレキシブルディスク（ＦＤ）１２９のディスクコントローラ（ＤＫＣ）１３０と、ネットワーク１３１との接続のためのネットワークインターフェースコントローラ（ＮＩＣ）１３２とが、システムバス１１３を介して互いに通信可能に接続された構成としている。
【００４４】
ＣＰＵ１２１は、ＲＯＭ１２２或いはＨＤ１２８に記憶されたソフトウェア、或いはＦＤ１２９より供給されるソフトウェアを実行することで、システムバス１２３に接続された各構成部を総括的に制御する。
すなわち、ＣＰＵ１２１は、所定の処理シーケンスに従った処理プログラムを、ＲＯＭ１２２、或いはＨＤ１２８、或いはＦＤ１２９から読み出して実行することで、後述する動作を実現するための制御を行う。
【００４５】
ＲＡＭ１２３は、ＣＰＵ１２１の主メモリ或いはワークエリア等として機能する。
ＫＢＣ１２５は、ＫＢ１２４や図示していないポインティングデバイス等からの指示入力を制御する。
【００４６】
ＣＲＴＣ１２７は、ＣＲＴ１２６の表示を制御する。
ＤＫＣ１３０は、ブートプログラム、種々のアプリケーション、編集ファイル、ユーザファイル、ネットワーク管理プログラム、及び本実施の形態における所定の処理プログラム等を記憶するＨＤ１２８及びＦＤ１２９とのアクセスを制御する。
ＮＩＣ１３２は、ネットワーク１３１上の装置或いはシステムと双方向にデータをやりとりする。
【００４７】
ＣＰＵ１２１は、ＲＯＭ１２２或いはＨＤ１２８に記憶されたソフトウェア、或いはＦＤ１２９より供給されるソフトウェアを実行することで、システムバス１２３に接続された各構成部を総括的に制御する。
すなわち、ＣＰＵ１２１は、所定の処理シーケンスに従った処理プログラムを、ＲＯＭ１２２、或いはＨＤ１２８、或いはＦＤ１２９から読み出して実行することで、図２に示したフローチャートの処理を実現するための制御を行う。なお、ステップＳ２０６における出力は、例えば、ＣＲＴ１２６に表示することによって実現される。
【００４８】
（第２の実施形態）
次に、本発明の第２の実施形態について説明する。前述した第１の実施形態では、回転対称体が円筒であり、その円筒の回転半径Ｒ（ｚ）が一定である場合について説明したが、本実施形態では、一般的な回転対称体の円筒の場合について説明する。このように、本実施形態と前述した第１の実施形態とは、形状計測データＰiに対応する理論形状を記述するパラメータを求める対象となる３次元形状が異なるだけである。したがって、第１の実施形態と同一の部分については、図１〜図５に付した符号と同一の符号を付すなどして、詳細な説明を省略する。
【００４９】
一般的な回転対称体の円筒の場合、すなわち、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0,ｙ0,ｚ0,α,β）を最小にする場合には、形状計測データＰiに対応する理論形状を記述するパラメータｘ0、ｙ0、ｚ0、α、βを、２つの組に分けるようにする。例えば、パラメータｘ0、ｙ0の組と、パラメータｚ0、α、βの組にパラメータが分けられる設定がなされた上で、図６に示すフローチャートの処理を行う。
（１） 形状計測データＰiを入力する（ステップＳ６０１）。
（２） 各パラメータの初期値ｘ0(0)、ｙ0(0)、ｚ0(0)、α(0)、β(0)を設定する（ステップＳ６０２）。
【００５０】
（３） パラメータを反復計算する（ｋ回目の値をもとに、（ｋ＋１）回目の値を決定する）。
（３−１） パラメータｚ0(k)、α(k)、β(k)を既知量として、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0(k+1)，ｙ0(k+1)，ｚ0(k+1)，α(k)，β(k)）を最小とするパラメータｘ0(k+1)、ｙ0(k+1)を決定する。すなわち、前述した（５式）の連立方程式からパラメータｘ0(k+1)、ｙ0(k+1)を求める（ステップＳ６０３）
【００５１】
（３−２） 求めたパラメータｘ0(k+1)、ｙ0(k+1)を既知量として、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0(k+1)，ｙ0(k+1)，ｚ0(k+1)，α(k+1)，β(k+1)）を最小とするパラメータｚ0(k+1)、α(k+1)、β(k+1)を決定する。すなわち、以下の（８式）の連立方程式から、パラメータｚ0(k+1)、α(k+1)、β(k+1)を求める（ステップＳ６０４）。
【００５２】
後述するように、本実施形態では、被測定物の形状を測定するときに工夫をして、中心軸１０１がほぼｚ軸に平行となる位置で測定されたとする。そうすると、中心軸１０１の方向を記述するパラメータα、βは、１に比べて十分小さな量となり、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0，ｙ0，ｚ0，α，β）は、α、βの２次形式で表される。したがって、このステップＳ６０４では、線形連立方程式の問題を解いて、パラメータα(k+1)、β(k+1)を算出することになる。
【００５３】
【数６】
【００５４】
（３−３） ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下であるか否かを判定する（ステップＳ６０５）。この判定の結果、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下となり、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値とが十分に近づいたならば、パラメータｘ0(k+1)、ｙ0(k+1)、ｚ0(k+1)、α(k+1)、β(k+1)を求める処理を終了して、ステップＳ６０６（（４）結果を出力する）へ進む。
【００５５】
一方、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値とが十分に近づいていない場合には、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値とが十分に近づくまで、ステップＳ６０３〜Ｓ６０５を繰り返し行う。
（４）（ｋ＋１）回目に求めたパラメータｘ0(k+1)、ｙ0(k+1)、ｚ0(k+1)、α(k+1)、β(k+1)の値を、形状計測データＰiに対応する理論形状を記述するパラメータとして出力する（ステップＳ６０６）。
【００５６】
なお、第１の実施形態と同様に、５つのパラメータｘ0、ｙ0、ｚ0、α、βを、複数の組に分けることで各ステップ（ステップＳ６０３、６０４）における最適解の探索空間を狭めることができれば、組の分け方は、図６のフローチャートで示した方法でなくてもよい。
また、（３−１）、（３−２）のステップＳ６０３、Ｓ６０４における最適値の探索は、一般には非線形の最小自乗法であり、各ステップＳ２０３、Ｓ２０４内でさらに反復計算が必要となる。
【００５７】
ところが、前述したように、本実施形態では、被測定物の形状を測定するときに工夫をして、中心軸１０１がほぼｚ軸に平行となる位置で測定されるので、中心軸１０１の方向を記述するパラメータα、β、ｚ0は、１に比べて十分小さな量であることが期待される。その場合、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0，ｙ0，ｚ0，α，β）は、α、β、ｚ0の２次形式で表され、（３−２）のステップＳ６０４では、線形連立方程式の問題を解けばよいことになり極めて簡単にパラメータα、β、ｚ0を算出することができる。このようにしてパラメータを分けることによって、各ステップＳ６０３、６０４における最適解の探索空間を狭めるだけでなく、線形探索となり、計算コストと、安定性とが著しく向上する。
【００５８】
また、初期値の与え方によっては、（３−２）のステップＳ６０４を行った後に、（３−１）のステップＳ６０３を行ってもよい場合があることも第１の実施形態と同様である。
【００５９】
（実施例）
次に、本実施形態の実施例として、円錐の側面上に配置されるべき形状計測データＰiから中心軸を決定する方法の一例を示す。本実施例では、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0，ｙ0，ｚ0，α，β）が最小となるパラメータｘ0、ｙ0、ｚ0、α、βを決定した。ここで、３次元形状計測器では、中心軸１０１がz軸と略平行になるようにすると共に、基準点のｚ座標も略ゼロ（０）となるように測定されている。よって、パラメータｚ0、α、βは１よりも十分小さい（ｚ0、α、β≪１）とみなしてよい。そうすると、形状計測データＰiから中心軸までの距離の誤差の２乗和Ｓ（ｘ0，ｙ0，ｚ0，α，β）を、以下の（９式）のように書くことができる。（９式）に示すように、形状計測データＰiから中心軸までの距離の誤差の２乗和Ｓは、５つのパラメータｘ0、ｙ0、ｚ0、α、βのうち、パラメータｚ0、α、βについての２次形式と見ることができる。
【００６０】
【数７】
【００６１】
そこで、パラメータｘ0、ｙ0、ｚ0、α、βの組の分け方としては、少なくとも１つの組は、パラメータｚ0、α、βの中から選択されたパラメータのみを含む組にしたい。そうすることによって、その組のパラメータを算出するための計算は線形計算となる。本実施例においては５つのパラメータｘ0、ｙ0、ｚ0、α、βを、パラメータｚ0、α、βの組と、パラメータｘ0，ｙ0の組とに分けるようにした。そして、パラメータｚ0、α、βに関しては線形方程式を解くことで求め、パラメータｘ0、ｙ0に関しては非線形方程式を解くことで解を求めた。ここで、非線形方程式の反復解法としては、ニュートン・ラフソン法を用いた。
【００６２】
計算精度を確認するために、パラメータα、β、ｘ0、ｙ0、ｚ0が、それぞれ０．０２、０．０１、１．０、２．０、０．０２（α＝０．０２、β＝０．０１、ｘ0＝１．０、ｙ0＝２．０、ｚ0＝０．０２）となる円錐の面上にらせん状に配置された５０点の擬似形状計測データを生成し、生成した擬似形状計測データからパラメータα、β、ｘ0、ｙ0、ｚ0の値を逆算した。図７に、３次元座標上にプロットされた擬似形状計測データの一例を、図８に、入力した擬似形状計測データの値の一例をそれぞれ示す。
【００６３】
各パラメータα、β、ｘ0、ｙ0、ｚ0の初期値α(0)、β(0)、ｘ0(0)、ｙ0(0)、ｚ0(0)を、それぞれ０．０（α(0)＝０．０、β(0)＝０．０、ｘ0(0)＝０．０、ｙ0(0)＝０．０、ｚ0(0)＝０．０）として計算した結果、パラメータα、β、ｘ0、ｙ0、ｚ0の値は、それぞれ０．０２０５９５５、０．００７３０５９９、０．９９９６６１、２．０００６８、０．０１８６３５９（α＝０．０２０５９５５、β＝０．００７３０５９９、ｘ0＝０．９９９６６１、ｙ0＝２．０００６８、ｚ0＝０．０１８６３５９）となり、本実施形態の手法が有効であることが確認できた。
【００６４】
すなわち、３次元形状の計測データに対応する理論形状を規定するパラメータが多数の場合であっても、高速に且つロバスト（robust）にパラメータを決定することができるので、接触式、非接触式などの３次元座標計測器を用いて計測された３次元形状の評価を効率的に行うことができる。具体的には、回転対称体の形状計測データＰiの中心軸１０１を正しく算出し、傾いて測定された形状計測データＰiに関して高速に且つロバスト（robust）にその傾きを補正することができる。
なお、第１の実施形態と同様に、例えば図５に示した情報処理装置を用いることにより、以上のようにしてパラメータｘ0、ｙ0、ｚ0、α、βを解く処理を実現することができる。
【００６５】
（第３の実施形態）
次に、本発明の第３の実施形態について説明する。前述した第１の実施形態では、回転対称体が円筒であり、その円筒の回転半径Ｒ（ｚ）が一定である場合について説明したが、本実施形態では、回転対称体は円筒であるが、その円筒の回転半径Ｒ（ｚ）もパラメータとする場合について説明する。このように、本実施形態と前述した第１の実施形態とは、形状計測データＰiに対応する理論形状を記述するパラメータを求める対象となる３次元形状の一部が異なるだけである。したがって、第１の実施形態と同一の部分については、図１〜図５に付した符号と同一の符号を付すなどして、詳細な説明を省略する。
【００６６】
回転半径Ｒ（ｚ）をパラメータとする円筒の場合、すなわち、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0,ｙ0,Ｒ，α,β）を最小にする場合には、形状計測データＰiに対応する理論形状を記述するパラメータｘ0、ｙ0、Ｒ、α、βを２つの組に分けるようにする。例えば、パラメータｘ0、ｙ0の組と、パラメータα、β、Ｒの組にパラメータが分けられる設定がなされた上で、図９に示すフローチャートの処理を行う。
（１）形状計測データＰiを入力する（ステップＳ９０１）。
（２）各パラメータの初期値ｘ0(0)、ｙ0(0)、Ｒ(0)、α(0)、β(0)を設定する（ステップＳ９０２）。
（３）パラメータを反復計算する（ｋ回目の値をもとに、（ｋ＋１）回目の値を決定する）。
（３−１） パラメータα(k)、β(k)を既知量として、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0(k+1)，ｙ0(k+1)，Ｒ(k+1)，α(k)，β(k)）を最小とするパラメータｘ0(k+1)、ｙ0(k+1)、Ｒ(k+1)を決定する。すなわち、以下の（１０式）の連立方程式からパラメータｘ0(k+1)、ｙ0(k+1)、Ｒ(k+1)を求める（ステップＳ９０３）。
【００６７】
【数８】
【００６８】
（３−２） 求めたパラメータｘ0(k+1)、ｙ0(k+1)、Ｒ(k+1)を既知量として、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0(k+1)，ｙ0(k+1)，Ｒ(k+1)，α(k+1)，β(k+1)）を最小とするパラメータα(k+1)、β(k+1)を決定する。すなわち、前述した（６式）の連立方程式から、パラメータα(k+1)、β(k+1)を求める（ステップＳ９０４）。
【００６９】
後述するように、本実施形態では、被測定物の形状を測定するときに工夫をして、中心軸１０１がほぼｚ軸に平行となる位置で測定されたとする。そうすると、中心軸１０１の方向を記述するパラメータα、βは、１に比べて十分小さな量となり、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0，ｙ0，Ｒ，α，β）は、α、βの２次形式で表される。したがって、このステップＳ２０４では、線形連立方程式の問題を解いて、パラメータα(k+1)、β(k+1)を算出することになる。
【００７０】
（３−３） ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下であるか否かを判定する（ステップＳ９０５）。この判定の結果、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下となり、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値とが十分に近づいたならば、パラメータｘ0(k+1)、ｙ0(k+1)、Ｒ(k+1)、α(k+1)、β(k+1)を求める処理を終了して、ステップＳ９０６（（４）結果を出力する）へ進む。
【００７１】
一方、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値とが十分に近づいていない場合には、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値とが十分に近づくまで、ステップＳ９０３〜Ｓ９０５を繰り返し行う。
（４）（ｋ＋１）回目に求めたパラメータｘ0(k+1)、ｙ0(k+1)、Ｒ(k+1)、α(k+1)、β(k+1)の値を、形状計測データＰiに対応する理論形状を記述するパラメータとして出力する（ステップＳ９０６）
【００７２】
なお、第１の実施形態と同様に、５つのパラメータｘ0、ｙ0、Ｒ、α、βを、複数の組に分けることで各ステップ（ステップＳ９０３、９０４）における最適解の探索空間を狭めることができれば、組の分け方は、図９のフローチャートで示した方法でなくてもよい。
また、（３−１）、（３−２）のステップＳ９０３、Ｓ９０４における最適値の探索は、一般には非線形の最小自乗法であり、各ステップ９０３、９０４内でさらに反復計算が必要となる。
【００７３】
ところが、前述したように、本実施形態では、被測定物の形状を測定するときに工夫をして、中心軸１０１がほぼｚ軸に平行となる位置で測定されるので、中心軸１０１の方向を記述するパラメータα、βは、１に比べて十分小さな量であることが期待される。その場合、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0，ｙ0，Ｒ，α，β）は、α、βの２次形式で表され、（３−２）のステップＳ９０４では、線形連立方程式の問題を解けばよいことになり極めて簡単にパラメータα、βを算出することができる。このようにしてパラメータを分けることによって、各ステップＳ９０３、９０４における最適解の探索空間を狭めるだけでなく、線形探索となり、計算コストと、安定性とが著しく向上する。
また、初期値の与え方によっては、（３−２）のステップＳ９０４を行った後に、（３−１）のステップＳ９０３を行ってもよい場合があることも第１の実施形態と同様である。
【００７４】
（実施例）
次に、本実施形態の実施例として、円筒の側面上に配置されるべき形状計測データＰiから中心軸を決定する方法の一例を示す。本実施例では、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0，ｙ0，Ｒ，α，β）が最小となるパラメータｘ0、ｙ0、Ｒ、α、βを決定した。ここで、３次元形状計測器では、中心軸１０１がz軸と略平行になるように測定されている。よって、パラメータα、βは１よりも十分小さい（α、β≪１）とみなしてよい。そうすると、形状計測データＰiから中心軸までの距離の誤差の２乗和Ｓ（ｘ0，ｙ0，Ｒ，α，β）を、以下の（１１式）のように書くことができる。（１１式）に示すように、形状計測データＰiから中心軸までの距離の誤差の２乗和Ｓは、５つのパラメータｘ0、ｙ0、Ｒ、α、βのうち、パラメータα、β、Ｒについての２次形式と見ることができる。
【００７５】
【数９】
【００７６】
そこで、パラメータｘ0、ｙ0、Ｒ、α、βの組の分け方としては、少なくとも１つの組は、パラメータα、β、Ｒの中から選択されたパラメータのみを含む組にしたい。そうすることによって、その組のパラメータを算出するための計算は線形計算となる。本実施例においては５つのパラメータｘ0、ｙ0、Ｒ、α、βを、パラメータα、βの組と、パラメータｘ0，ｙ0、Ｒの組とに分けるようにした。そして、パラメータα、βに関しては線形方程式を解くことで、ｘ0，ｙ0、Ｒに関しては非線形方程式を解くことで解を求めた。ここで、非線形方程式の反復解法としては、ニュートン・ラフソン法を用いた。
【００７７】
計算精度を確認するために、パラメータα、β、ｘ0、ｙ0、Ｒがそれぞれ０．１０、０．０５、１．０、２．０、３．０（α＝０．１０、β＝０．０５、ｘ0＝１．０、ｙ0＝２．０、Ｒ＝３．０）となる円筒の側面上に配置された６０点の擬似形状計測データを生成し、生成した擬似形状計測データからパラメータα、β、ｘ0、ｙ0、Ｒの値を逆算した。３次元座標上にプロットされた擬似形状計測データと、入力した擬似形状計測データの値は、それぞれ図３及び図４に示した通りである。
【００７８】
各パラメータα、β、ｘ0、ｙ0、Ｒの初期値α(0)、β(0)、ｘ0(0)、ｙ0(0)、Ｒ(0)を、それぞれ０．０（α(0)＝０．０、β(0)＝０．０、ｘ0(0)＝０．０、ｙ0(0)＝０．０、Ｒ(0)＝０．０）として計算した結果、パラメータα、β、Ｒ、ｘ0、ｙ0の値は、それぞれ０．０９９６５４４、０．０４９８２７１、３．００９２９、１．０００６９、２．０００３４（α＝０．０９９６５４４、β＝０．０４９８２７１、Ｒ＝３．００９２９、ｘ0＝１．０００６９、ｙ0＝２．０００３４）となり、本実施形態の手法が有効であることが確認できた。
【００７９】
すなわち、３次元形状の計測データに対応する理論形状を規定するパラメータが多数の場合であっても、高速に且つロバスト（robust）にパラメータを決定することができるので、接触式、非接触式などの３次元座標計測器を用いて計測された３次元形状の評価を効率的に行うことができる。具体的には、回転対称体の形状計測データＰiの中心軸１０１を正しく算出し、傾いて測定された形状計測データＰiに関して高速に且つロバスト（robust）にその傾きを補正することができる。
なお、第１の実施形態と同様に、例えば図５に示した情報処理装置を用いることにより、以上のようにしてパラメータｘ0、ｙ0、Ｒ、α、βを解く処理を実現することができる。
【００８０】
（第４の実施形態）
次に、本発明の第４の実施形態について説明する。前述した第１の実施形態では、ステップＳ２０４において、線形連立方程式の問題を解くことで、極めて簡単にパラメータα、βの最適値を算出することができるようにしたが、本実施形態では、ステップＳ２０３、Ｓ２０４におけるパラメータの最適値の探索を、それぞれ非線形の最小自乗法を用いて行うようにした。このように、本実施形態と前述した第１の実施形態とは、パラメータα、βを求める方法が異なるだけである。したがって、第１の実施形態と同一の部分については、図１〜図５に付した符号と同一の符号を付すなどして、詳細な説明を省略する。
【００８１】
本実施形態においても、第１の実施形態と同様に、回転対称体が、前述した円筒であり、回転半径Ｒ（ｚ）が一定であるので、形状計測データＰiに対応する理論形状を記述するパラメータｘ0、ｙ0、α、βが、パラメータｘ0、ｙ0の組と、パラメータα、βの組とに分けられる設定がなされた上で、図２のフローチャートに従った処理を行う。
（１）形状計測データＰiを入力する（ステップＳ２０１）。
（２）各パラメータｘ0、ｙ0、α、βの初期値ｘ0(0)、ｙ0(0)、α(0)、β(0)を設定する（ステップＳ２０２）。
（３）パラメータを反復計算する（ｋ（ｋは自然数）回目の値をもとに、（ｋ＋１）回目の値を決定する）。
（３−１） 前述した（５式）の連立方程式からパラメータｘ0(k+1)、ｙ0(k+1)を求める（ステップＳ２０３）
【００８２】
（３−２） 前述した（６式）の連立方程式からパラメータα(k+1)、β(k+1)を求める（ステップＳ２０４）。後述するように、本実施形態では、非線形の最小自乗法を用いて、パラメータα(k+1)、β(k+1)を求める。
【００８３】
（３−３） ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下であるか否かを判定する（ステップＳ２０５）。この判定の結果、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下であれば、パラメータｘ0(k+1)、ｙ0(k+1)、α(k+1)、β(k+1)を求める処理を終了して、ステップＳ２０６（（４）結果を出力する）へ進む。
【００８４】
一方、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下でない場合には、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下になるまで、ステップＳ２０３〜Ｓ２０５を繰り返し行う。
【００８５】
（４）（ｋ＋１）回目に求めたパラメータｘ0(k+1)、ｙ0(k+1)、α(k+1)、β(k+1)の値を、形状計測データＰiに対応する理論形状を記述するパラメータとして出力する（ステップＳ２０６）。
【００８６】
なお、第１の実施形態と同様に、４つのパラメータｘ0、ｙ0、α、βを、複数の組に分けることで各ステップ（ステップＳ２０３、Ｓ２０４）における最適解の探索空間（最適値探索問題の解空間）を狭めることができれば、組の分け方は、図２のフローチャートで示した方法でなくてもよい。
【００８７】
また、本実施形態の（３−１）、（３−２）のステップＳ２０３、Ｓ２０４におけるパラメータの最適値の探索は、非線形の最小自乗法を用いて行われ、ステップＳ２０３、Ｓ２０４内でさらに反復計算が必要となる。
【００８８】
以上のようにしてパラメータを複数の組に分けることによって、各ステップＳ２０３、２０４における最適値の探索空間を狭めることで、計算コストと、安定性とが著しく向上する。
また、初期値の与え方によっては、（３−２）のステップＳ２０４を行った後に、（３−１）のステップＳ２０３を行ってもよい場合があることも第１の実施形態と同様である。
【００８９】
（実施例）
次に、本実施形態の実施例として、円筒の側面上に配置されるべき形状計測データＰiから中心軸を決定する方法の一例を示す。本実施例では、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0，ｙ0，α，β）が最小となるパラメータｘ0、ｙ0、α、βを決定した。ここで、４つのパラメータｘ0，ｙ0，α，βを、パラメータα、βの組と、パラメータｘ0、ｙ0の組とに分けて、それぞれ非線形の最小自乗法を用いて解を求めた。ここで、非線形方程式の反復解法としては、ニュートン・ラフソン法を用いた。
【００９０】
計算精度を確認するために、パラメータα、β、ｘ0、ｙ0が、それぞれ０．１０、０．０５、１．０、２．０（α＝０．１０、β＝０．０５、ｘ0＝１．０、ｙ0＝２．０）となる半径Ｒが３（Ｒ＝３）の円筒の側面上に配置された６０点の擬似形状計測データを生成し、生成した擬似形状計測データから、パラメータα、β、ｘ0、ｙ0の値を逆算した。３次元座標上にプロットされた擬似形状計測データと、入力した擬似形状計測データの値は、それぞれ図３、図４に示した通りである。
【００９１】
各パラメータα、β、ｘ0、ｙ0の初期値α(0)、β(0)、ｘ0(0)、ｙ0(0)を、それぞれ０．０、０．０、１．１９８７６、２．０９９３８（α(0)＝０．０、β(0)＝０．０、ｘ0(0)＝１．１９８７６、ｙ0(0)＝２．０９９３８）として計算した結果、パラメータα、β、ｘ0、ｙ0の値は、それぞれ０．０９９９９９８、０．０４９９９９８、１．０００００、２．０００００（α＝０．０９９９９９８、β＝０．０４９９９９８、ｘ0＝１．０００００、ｙ0＝２．０００００）となり、本実施形態の手法がパラメータを高精度に計算するのに極めて有効であることが確認できた。
【００９２】
すなわち、３次元形状の計測データに対応する理論形状を規定するパラメータが多数の場合であっても、高速に且つロバスト（robust）にパラメータを決定することができるので、接触式、非接触式などの３次元座標計測器を用いて計測された３次元形状の評価を効率的に行うことができる。具体的には、回転対称体の形状計測データＰiの中心軸１０１を正しく算出し、傾いて測定された形状計測データＰiに関して高速に且つロバスト（robust）にその傾きを補正することができる。
なお、第１の実施形態と同様に、例えば図５に示した情報処理装置を用いることにより、以上のようにしてパラメータｘ0、ｙ0、α、βを解く処理を実現することができる。
【００９３】
（第５の実施形態）
次に、本発明の第５の実施形態について説明する。前述した第３の実施形態では、ステップＳ９０４において、線形連立方程式の問題を解くことで、極めて簡単にパラメータα、βの最適値を算出することができるようにしたが、本実施形態では、ステップＳ９０３、Ｓ９０４におけるパラメータの最適値の探索を、それぞれ非線形の最小自乗法を用いて行うようにした。このように、本実施形態と前述した第３の実施形態とは、パラメータα、βを求める方法が異なるだけである。したがって、第３の実施形態と同一の部分については、図１〜図５、図９に付した符号と同一の符号を付すなどして、詳細な説明を省略する。
【００９４】
本実施形態においても、第３の実施形態と同様に、回転対称体は、回転半径Ｒ（ｚ）をパラメータとする円筒であるので、形状計測データＰiに対応する理論形状を記述するパラメータｘ0、ｙ0、Ｒ、α、βが、パラメータｘ0、ｙ0の組と、パラメータα、β、Ｒの組とに分けられる設定がなされた上で、図９のフローチャートに従った処理を行う。
（１）形状計測データＰiを入力する（ステップＳ９０１）。
（２）各パラメータの初期値ｘ0(0)、ｙ0(0)、Ｒ(0)、α(0)、β(0)を設定する（ステップＳ９０２）。
（３）パラメータを反復計算する（ｋ回目の値をもとに、（ｋ＋１）回目の値を決定する）。
（３−１） 前述した（１０式）の連立方程式から、パラメータｘ0(k+1)、ｙ0(k+1)、Ｒ(k+1)を求める（ステップＳ９０３）。
【００９５】
（３−２） 前述した（６式）の連立方程式から、パラメータα(k+1)、β(k+1)を求める（ステップＳ９０４）。後述するように、本実施形態では、非線形の最小自乗法を用いて、パラメータα(k+1)、β(k+1)を求める。
【００９６】
（３−３）ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下であるか否かを判定する（ステップＳ９０５）。この判定の結果、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下であれば、パラメータｘ0(k+1)、ｙ0(k+1)、Ｒ(k+1)、α(k+1)、β(k+1)を求める処理を終了して、ステップＳ９０６（（４）結果を出力する）へ進む。
【００９７】
一方、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下でない場合には、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下になるまで、ステップＳ９０３〜Ｓ９０５を繰り返し行う。
【００９８】
（４）（ｋ＋１）回目に求めたパラメータｘ0(k+1)、ｙ0(k+1)、Ｒ(k+1)、α(k+1)、β(k+1)の値を、形状計測データＰiに対応する理論形状を記述するパラメータとして出力する（ステップＳ９０６）。
【００９９】
なお、第３の実施形態と同様に、５つのパラメータｘ0、ｙ0、Ｒ、α、βを、複数の組に分けることで各ステップ（ステップＳ９０３、Ｓ９０４）における最適解の探索空間（最適値探索問題の解空間）を狭めることができれば、組の分け方は、図９のフローチャートで示した方法でなくてもよい。
【０１００】
また、（３−１）、（３−２）のステップＳ９０３、Ｓ９０４における最適値の探索は、非線形の最小自乗法を用いて行われ、ステップＳ９０３、Ｓ９０４内でさらに反復計算が必要となる。
【０１０１】
以上のようにしてパラメータを複数の組に分けることによって、各ステップＳ９０３、Ｓ９０４における最適値の探索空間を狭めることで、計算コスト、安定性が著しく向上する。
また、初期値の与え方によっては、（３−２）のステップＳ９０４を行った後に、（３−１）のステップＳ９０３を行ってもよい場合があることも第３の実施形態と同様である。
【０１０２】
（実施例）
次に、本実施形態の実施例として、円筒の側面上に配置されるべき形状計測データＰiから中心軸を決定する方法の一例を示す。本実施例では、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0，ｙ0，Ｒ，α，β）が最小となるパラメータｘ0、ｙ0、Ｒ、α、βを決定した。ここで、５つのパラメータｘ0、ｙ0、Ｒ、α、βを、パラメータα、βの組と、パラメータｘ0、ｙ0、Ｒの組とに分けて、それぞれ非線形の最小自乗法を用いて解を求めた。ここで、非線形方程式の反復解法としては、ニュートン・ラフソン法を用いた。
【０１０３】
計算精度を確認するために、パラメータα、β、ｘ0、ｙ0、Ｒがそれぞれ０．１０、０．０５、１．０、２．０、３．０（α＝０．１０、β＝０．０５、ｘ0＝１．０、ｙ0＝２．０、Ｒ＝３．０）となる円筒の側面上に配置された６０点の擬似形状計測データを生成し、生成した擬似形状計測データからパラメータα、β、ｘ0、ｙ0、Ｒの値を逆算した。３次元座標上にプロットされた擬似形状形状データと、入力した擬似形状計測データの値は、それぞれ図３及び図４に示した通りである。
【０１０４】
各パラメータα、β、ｘ0、ｙ0、Ｒの初期値α(0)、β(0)、ｘ0(0)、ｙ0(0)、Ｒ(0)を、それぞれ０．０（α(0)＝０．０、β(0)＝０．０、ｘ0(0)＝０．０、ｙ0(0)＝０．０、Ｒ(0)＝０．０）として計算した結果、パラメータα、β、Ｒ、ｘ0、ｙ0の値は、それぞれ０．０９９９９９８、０．０４９９９９８、３．０００００、１．０００００、２．０００００（α＝０．０９９９９９８、β＝０．０４９９９９８、Ｒ＝３．０００００、ｘ0＝１．０００００、ｙ0＝２．０００００）となり、本実施形態の手法がパラメータを高精度に計算するのに極めて有効であることが確認できた。
【０１０５】
すなわち、３次元形状の計測データに対応する理論形状を規定するパラメータが多数の場合であっても、高速に且つロバスト（robust）にパラメータを決定することができるので、接触式、非接触式などの３次元座標計測器を用いて計測された３次元形状の評価を効率的に行うことができる。具体的には、回転対称体の形状計測データＰiの中心軸１０１を正しく算出し、傾いて測定された形状計測データＰiに関して高速に且つロバスト（robust）にその傾きを補正することができる。
なお、第３の実施形態と同様に、例えば図５に示した情報処理装置を用いることにより、以上のようにしてパラメータｘ0、ｙ0、Ｒ、α、βを解く処理を実現することができる。
【０１０６】
（第６の実施形態）
次に、本発明の第６の実施形態について説明する。前述した第２の実施形態では、ステップＳ６０４において、線形連立方程式の問題を解くことで、極めて簡単にパラメータα、βの最適値を算出することができるようにしたが、本実施形態では、ステップＳ６０３、Ｓ６０４におけるパラメータの最適値の探索を、それぞれ非線形の最小自乗法を用いて行うようにした。このように、本実施形態と前述した第２の実施形態とは、パラメータα、βを求める方法が異なるだけである。したがって、第２の実施形態と同一の部分については、図１〜図２、図６〜図８に付した符号と同一の符号を付すなどして、詳細な説明を省略する。
【０１０７】
本実施形態においても、第２の実施形態と同様に、回転対称体は、一般的な回転対称体の円筒であるので、形状計測データＰiに対応する理論形状を記述するパラメータｘ0、ｙ0、ｚ0、α、βが、パラメータｘ0、ｙ0の組と、パラメータｚ0、α、βの組にすなわちＳ（ｘ0,ｙ0,ｚ0,α,β）の場合を例に説明する。
まず、理論形状記述パラメータを２つの組、ｘ0、ｙ0の組とに分けられる設定がなされた上で、図６のフローチャートに従った処理を行う。
（１）形状計測データＰiを入力する（ステップＳ６０１）。
（２）各パラメータの初期値ｘ0(0)、ｙ0(0)、ｚ0(0)、α(0)、β(0)を設定する（ステップＳ６０２）。
（３）パラメータを反復計算する（ｋ回目の値をもとに、（ｋ＋１）回目の値を決定する）。
（３−１）前述した（５式）の連立方程式から、パラメータｘ0(k+1)、ｙ0(k+1)を求める（ステップＳ６０３）。
【０１０８】
（３−２）前述した（６式）の連立方程式から、パラメータα(k+1)、β(k+1)を求める（ステップＳ６０４）。後述するように、本実施形態では、非線形の最小自乗法を用いて、パラメータα(k+1)、β(k+1)を求める。
【０１０９】
（３−３）ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下であるか否かを判定する（ステップＳ９０５）。この判定の結果、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下であれば、パラメータｘ0(k+1)、ｙ0(k+1)、Ｒ(k+1)、α(k+1)、β(k+1)を求める処理を終了して、ステップＳ６０６（（４）結果を出力する）へ進む。
【０１１０】
一方、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下でない場合には、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下になるまで、ステップＳ６０３〜Ｓ６０５を繰り返し行う。
【０１１１】
（４）（ｋ＋１）回目に求めたパラメータｘ0(k+1)、ｙ0(k+1)、Ｒ(k+1)、α(k+1)、β(k+1)の値を、形状計測データＰiに対応する理論形状を記述するパラメータとして出力する（ステップＳ６０６）。
【０１１２】
なお、第２の実施形態と同様に、５つのパラメータｘ0、ｙ0、ｚ0、α、βを、複数の組に分けることで各ステップ（ステップＳ６０３、６０４）における最適解の探索空間を狭めることができれば、組の分け方は、図６のフローチャートで示した方法でなくてもよい。
【０１１３】
また、（３−１）、（３−２）のステップＳ６０３、Ｓ６０４における最適値の探索は、非線形の最小自乗法を用いて行われ、ステップＳ６０３、Ｓ６０４内でさらに反復計算が必要となる。
【０１１４】
以上のようにしてパラメータを複数の組に分けることによって、各ステップＳ６０３、Ｓ６０４における最適値の探索空間を狭めることで、計算コスト、安定性が著しく向上する。
また、初期値の与え方によっては、（３−２）のステップＳ６０４を行った後に、（３−１）のステップＳ６０３を行ってもよい場合があることも第２の実施形態と同様である。
【０１１５】
（実施例）
次に、本実施形態の実施例として、円錐の側面上に配置されるべき形状計測データＰiから中心軸を決定する方法の一例を示す。本実施例では、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0，ｙ0，ｚ0，α，β）が最小となるパラメータｘ0、ｙ0、ｚ0、α、βを決定した。ここで、５つのパラメータｘ0、ｙ0、ｚ0、α、βを、パラメータｚ0、α、βの組と、パラメータｘ0，ｙ0の組とに分けて、それぞれ非線形の最小自乗法を用いて解を求めた。非線形方程式の反復解法としては、ニュートン・ラフソン法を用いた。
【０１１６】
計算精度を確認するために、パラメータα、β、ｘ0、ｙ0、ｚ0がそれぞれ０．０２、０．０１、１．０、２．０、０．０２（α＝０．０２、β＝０．０１、ｘ0＝１．０、ｙ0＝２．０、ｚ0＝０．０２）となる円錐面上にらせん状に配置された５０点の擬似形状計測データを生成し、生成した擬似形状計測データからパラメータα、β、ｘ0、ｙ0、ｚ0の値を逆算した。３次元座標上にプロットされた擬似形状形状データと、入力した擬似形状計測データの値は、それぞれ図７及び図８に示した通りである。
【０１１７】
各パラメータα、β、ｘ0、ｙ0、ｚ0の初期値α(0)、β(0)、ｘ0(0)、ｙ0(0)、ｚ0(0)を、それぞれ０．０（α(0)＝０．０、β(0)＝０．０、ｘ0(0)＝０．０、ｙ0(0)＝０．０、ｚ0(0))＝０．０）として計算した結果、パラメータα(0)、β(0)、ｘ0(0)、ｙ0(0)、ｚ0の値は、それぞれ０．０１９９９９９８、０．００９９９９９８、０．９９９９９９９９、２．０００００００、０．０１９９９９８（α＝０．０１９９９９９８、β＝０．００９９９９９８、ｘ0＝０．９９９９９９９９、ｙ0＝２．０００００００、ｚ0＝０．０１９９９９８）となり、本実施形態の手法がパラメータを高精度に計算するのに極めて有効であることが確認できた。
【０１１８】
すなわち、３次元形状の計測データに対応する理論形状を規定するパラメータが多数の場合であっても、高速に且つロバスト（robust）にパラメータを決定することができるので、接触式、非接触式などの３次元座標計測器を用いて計測された３次元形状の評価を効率的に行うことができる。具体的には、回転対称体の形状計測データＰiの中心軸１０１を正しく算出し、傾いて測定された形状計測データＰiに関して高速に且つロバスト（robust）にその傾きを補正することができる。
なお、第２の実施形態と同様に、例えば図５に示した情報処理装置を用いることにより、以上のようにしてパラメータｘ0、ｙ0、Ｒ、α、βを解く処理を実現することができる。
【０１１９】
（本発明の他の実施形態）
上述した実施形態の機能を実現するべく各種のデバイスを動作させるように、該各種デバイスと接続された装置あるいはシステム内のコンピュータに対し、前記実施形態の機能を実現するためのソフトウェアのプログラムコードを供給し、そのシステムあるいは装置のコンピュータ（ＣＰＵあるいはＭＰＵ）に格納されたプログラムに従って前記各種デバイスを動作させることによって実施したものも、本発明の範疇に含まれる。
【０１２０】
また、この場合、前記ソフトウェアのプログラムコード自体が上述した実施形態の機能を実現することになり、そのプログラムコード自体、及びそのプログラムコードをコンピュータに供給するための手段、例えば、かかるプログラムコードを格納した記録媒体は本発明を構成する。かかるプログラムコードを記憶する記録媒体としては、例えばフレキシブルディスク、ハードディスク、光ディスク、光磁気ディスク、ＣＤ−ＲＯＭ、磁気テープ、不揮発性のメモリカード、ＲＯＭ等を用いることができる。
【０１２１】
また、コンピュータが供給されたプログラムコードを実行することにより、上述の実施形態の機能が実現されるだけでなく、そのプログラムコードがコンピュータにおいて稼働しているＯＳ（オペレーティングシステム）あるいは他のアプリケーションソフト等と共同して上述の実施形態の機能が実現される場合にもかかるプログラムコードは本発明の実施形態に含まれることは言うまでもない。
【０１２２】
さらに、供給されたプログラムコードがコンピュータの機能拡張ボードやコンピュータに接続された機能拡張ユニットに備わるメモリに格納された後、そのプログラムコードの指示に基づいてその機能拡張ボードや機能拡張ユニットに備わるＣＰＵ等が実際の処理の一部または全部を行い、その処理によって上述した実施形態の機能が実現される場合にも本発明に含まれることは言うまでもない。
【図面の簡単な説明】
【０１２３】
【図１】本発明の第１の実施形態を示し、３次元空間内での形状計測データと、その形状計測データに対応する理論形状を記述するパラメータとの関係の一例を示す図である。
【図２】本発明の第１の実施形態を示し、パラメータを決定する方法の一例を説明するフローチャートである。
【図３】本発明の第１の実施形態を示し、３次元座標上にプロットされた擬似形状計測データを示す図である。
【図４】本発明の第１の実施形態を示し、入力した擬似形状計測データの値の一例を示す図である。
【図５】本発明の第１の実施形態を示し、情報処理装置の構成の一例を示したブロック図である。
【図６】本発明の第２の実施形態を示し、パラメータを決定する方法の一例を説明するフローチャートである。
【図７】本発明の第２の実施形態を示し、３次元座標上にプロットされた擬似形状計測データの一例を示す図である。
【図８】本発明の第２の実施形態を示し、入力した擬似形状計測データの値の一例をそれぞれ示す図である。
【図９】本発明の第３の実施形態を示し、パラメータを決定する方法の一例を説明するフローチャートである。
【符号の説明】
【０１２４】
１００ 理論形状
１０１ 中心軸
Ｃ 測定空間での基準点
Ｏ 理論形状の基準点
Ｐi 形状計測データ
ｘ0、ｙ0、ｚ0、α、β、Ｒ パラメータ
【技術分野】
【０００１】
本発明は、形状パラメータ決定方法、形状パラメータ決定装置、及びコンピュータプログラムに関し、特に、計測された３次元形状に対応する理論形状を規定する複数のパラメータを求めるために用いて好適なものである。
【背景技術】
【０００２】
近年、３次元形状の座標を、接触式又は非接触式などで測定する３次元形状座標計測器が発達し、広く普及されるようになった。この３次元形状座標計測器を用いて計測されたデータと、理論形状（設計値）とを比較することにより、部品の加工精度等を直接評価できるようになっている。特許文献１では、３次元座標計測器で計測されたデータに対応する理論形状を記述する理論形状パラメータを決定する方法が提示されている。
【０００３】
【特許文献１】特開平６−２８１４４０号公報
【発明の開示】
【発明が解決しようとする課題】
【０００４】
しかしながら、前述した従来の技術において、理論形状パラメータを決定する問題は、最適値探索問題である。このため、理論形状を記述するためのパラメータの数が多くなる。すなわち、解探索空間の次元が大きくなる。このように、解探索空間の次元が大きくなると、解が局所的な適合値に収束してしまい、真の最適値へ到達できなくなる可能性が大きくなる。前述した特許文献１に記載されている技術においても、局所的な適合値に解が落ち込むことへの対策は明示されていない。
【０００５】
このように、従来の技術では、３次元形状の計測データに対応する理論形状を記述するパラメータを高速に且つ安定して決定することが困難であるという問題点があった。
【０００６】
本発明は、このような問題点に鑑みてなされたものであり、３次元形状の計測データに対応する理論形状を記述するパラメータを高速に且つ安定して決定することができるようにすることを目的とする。
【課題を解決するための手段】
【０００７】
本発明の形状パラメータ決定方法は、計測された３次元形状に対応する理論形状を規定する複数のパラメータを決定する形状パラメータ決定方法であって、前記複数のパラメータを少なくとも２組に分け、前記計測された３次元形状の計測データと、その３次元形状の計測データに対応する理論形状との差を表す量が最小となるようなパラメータを、前記分けた組毎に順番に繰り返し求めることによって、前記複数のパラメータを決定することを特徴とする。
また、本発明の他の特徴とするところは、計測された３次元形状に対応する理論形状を規定する複数のパラメータを第１〜第ｎ（ｎは２以上の自然数）の組に分け、分けた組毎に前記複数のパラメータを決定する形状パラメータ決定方法であって、前記計測された３次元形状の計測データを入力する入力ステップと、前記入力ステップにより入力された計測データと、前記計測データに対応する理論形状との差を表す量が最小となるようなパラメータを、前記分けた組毎に求める第１〜第ｎのパラメータ取得ステップとを有することを特徴とする。
【０００８】
本発明の形状パラメータ決定装置は、計測された３次元形状に対応する理論形状を規定する複数のパラメータを決定する形状パラメータ決定装置であって、前記複数のパラメータを少なくとも２組に分け、前記計測された３次元形状の計測データと、その３次元形状の計測データに対応する理論形状との差を表す量が最小となるようなパラメータを、前記分けた組毎に順番に繰り返し求める手段を有することを特徴とする。
また、本発明の他の特徴とするところは、計測された３次元形状に対応する理論形状を規定する複数のパラメータを第１〜第ｎ（ｎは２以上の自然数）の組に分け、分けた組毎に前記複数のパラメータを決定する形状パラメータ決定装置であって、前記計測された３次元形状の計測データを入力する入力手段と、前記入力手段により入力された計測データと、前記計測データに対応する理論形状との差を表す量が最小となるようなパラメータを、前記分けた組毎に求める第１〜第ｎのパラメータ取得手段とを有することを特徴とする。
【０００９】
本発明のコンピュータプログラムは、計測された３次元形状に対応する理論形状を規定する複数のパラメータを決定することをコンピュータに実行させるためのコンピュータプログラムであって、前記複数のパラメータを少なくとも２組に分け、前記計測された３次元形状の計測データと、その３次元形状の計測データに対応する理論形状との差を表す量が最小となるようなパラメータを、前記分けた組毎に順番に繰り返し求めることをコンピュータに実行させることを特徴とする。
また、本発明の他の特徴とするところは、計測された３次元形状に対応する理論形状を規定する複数のパラメータを第１〜第ｎ（ｎは２以上の自然数）の組に分け、分けた組毎に前記複数のパラメータを決定することをコンピュータに実行させるためのコンピュータプログラムであって、前記計測された３次元形状の計測データを入力する入力ステップと、前記入力ステップにより入力された計測データと、前記計測データに対応する理論形状との差を表す量が最小となるようなパラメータを、前記分けた組毎に求める第１〜第ｎのパラメータ取得ステップとをコンピュータに実行させることを特徴とする。
【発明の効果】
【００１０】
本発明によれば、３次元形状の計測データに対応する理論形状を規定するパラメータが多数の場合であっても、高速に且つ安定してパラメータを決定することができるので、三次元座標計測器を用いて計測した形状の評価が効率的に行えるようになる。
【発明を実施するための最良の形態】
【００１１】
（第１の実施形態）
次に、図面を参照しながら、本発明の第１の実施形態について説明する。
なお、本実施形態では、３次元形状として回転対称体を例に挙げて説明する。また、回転対称体の中心軸（対称軸）の方向は、ｚ軸と略同じ方向を向いているとしても任意性は失われないので、以下、中心軸は、ｚ軸と略同じ方向を向いているとする。
【００１２】
図１は、３次元空間内での形状計測データ（座標値）と、その形状計測データに対応する理論形状を記述するパラメータとの関係の一例を示す図である。なお、形状計測データは、３次元形状座標計測器で測定されたデータである。
【００１３】
３次元形状座標計測器の一例としては、プローブと、測定台にセットされた被測定物の表面に前記プローブをトレースさせるトレース機構と、前記トレース機構によりトレースされたプローブの位置情報を処理して被測定物の形状を測定するコンピュータとを有するものが挙げられる。
【００１４】
図１において、形状計測データ（測定点）はＰi（ｘi,ｙi,ｚi）である。元の理論形状（例えば設計値）１００は、中心軸（ここではｚ軸）に沿って回転半径Ｒ（ｚ）で与えられているとする。この元の理論形状１００の基準点Ｏ（０，０，０）は、測定空間では、基準点Ｃ（ｘ0,ｙ0,ｚ0）へ移動（シフト）している。また、測定空間では、元の理論形状１００の中心軸（ｚ軸）に対し、方向ベクトルｎ（＝（α,β,１））の方向に中心軸１０１が傾いている。
【００１５】
すなわち、図１に示した例では、形状計測データＰiに対応する理論形状を記述するパラメータは、ｘ0、ｙ0、ｚ0、α、βの５つである。これらの５つのパラメータは、理論値（設計値）との差異、つまり、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0,ｙ0,ｚ0,α,β）を最小にするものとして算出する。ここで、形状計測データ（測定点）Ｐiから中心軸１０１へおろした垂線の足をＨiとすると、以下の（１式）が成り立つ。
【００１６】
【数１】
【００１７】
ただし、Ｎは、形状計測データＰiの数である。また、以下の（２式）〜（４式）が成り立つとする。
【００１８】
【数２】
【００１９】
３次元形状の一例である前記回転対称体が円筒の場合は、回転半径Ｒ（ｚ）はｚ軸方向の値によらず、基準点Ｃ（ｘ0,ｙ0,ｚ0）は、中心軸１０１のどこにあっても等価となる。よって、基準点Ｃのｚ軸方向の値ｚ0は定数と考えてよい。そうすると、形状計測データＰiに対応する理論形状を記述するパラメータは、ｘ0、ｙ0、α、βの４つとなり、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0,ｙ0,α,β）を考えればよいことになる。
【００２０】
また、回転半径Ｒ（ｚ）はｚ軸方向の値によらないが、設計値からは変化している場合は、回転半径Ｒも、形状計測データＰiに対応する理論形状を記述するパラメータに加えてもよい。この場合、形状計測データＰiに対応する理論形状を記述するパラメータはｘ0、ｙ0、Ｒ、α、βの５つとなり、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0,ｙ0,Ｒ,α,β）を考えればよいことになる。
【００２１】
いずれの場合（回転対称体が円筒の場合、及び回転対称体が、設計値と異なる回転半径を有する円筒である場合）についても、中心軸１０１までの距離の誤差の２乗和Ｓを最小にするパラメータを決定することは変わらないが、本実施形態では、回転対称体が、前述した円筒であり、回転半径Ｒ（ｚ）が一定である場合を説明する。
【００２２】
以下、図２のフローチャートを参照しながら、パラメータを決定する方法の一例を説明する。前述したように、本実施形態では、回転対称体が、前述した円筒であり、回転半径Ｒ（ｚ）が一定の場合、すなわち、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0,ｙ0,α,β）を最小にする場合を説明する。
【００２３】
なお、ここでは、形状計測データＰiに対応する理論形状を記述するパラメータｘ0、ｙ0、α、βを、２つの組に分けるようにする。例えば、パラメータｘ0、ｙ0の組と、パラメータα、βの組とに分けるようにする。このようにパラメータｘ0、ｙ0、α、βが、２つの組に分けられる設定がなされた上で、以下の処理が実行される。
【００２４】
（１） 形状計測データＰiを入力する（ステップＳ２０１）。
（２） 各パラメータｘ0、ｙ0、α、βの初期値ｘ0(0)、ｙ0(0)、α(0)、β(0)を設定する（ステップＳ２０２）。
ここで、右肩の添え字は反復回数を示すこととする。０回目の値だということで初期値の右肩には（０）が付いている。
【００２５】
（３） パラメータを反復計算する（ｋ（ｋは自然数）回目の値をもとに、（ｋ＋１）回目の値を決定する）。
（３−１） パラメータα(k)、β(k)を既知量として、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0(k+1)，ｙ0(k+1)，α(k)，β(k)）を最小とするパラメータｘ0(k+1)、ｙ0(k+1)を決定する。すなわち、以下の（５式）の連立方程式からパラメータｘ0(k+1)、ｙ0(k+1)を求める（ステップＳ２０３）。
【００２６】
【数３】
【００２７】
（３−２） 求めたパラメータｘ0(k+1)、ｙ0(k+1)を既知量として、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0(k+1)，ｙ0(k+1)，α(k+1)，β(k+1)）を最小とするパラメータα(k+1)、β(k+1)を決定する。すなわち、以下の（６式）の連立方程式から、パラメータα(k+1)、β(k+1)を求める（ステップＳ２０４）。
【００２８】
後述するように、本実施形態では、被測定物の形状を測定するときに工夫をして、中心軸１０１がほぼｚ軸に平行となる位置で測定されたとする。そうすると、中心軸１０１の方向を記述するパラメータα、βは、１に比べて十分小さな量となり、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0，ｙ0，α，β）は、α、βの２次形式で表される。したがって、このステップＳ２０４では、線形連立方程式の問題を解いて、パラメータα(k+1)、β(k+1)を算出することになる。
【００２９】
【数４】
【００３０】
（３−３） ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下であるか否かを判定する（ステップＳ２０５）。この判定の結果、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下となり、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値とが十分に近づいたならば、パラメータｘ0(k+1)、ｙ0(k+1)、α(k+1)、β(k+1)を求める処理を終了して、ステップＳ２０６（（４）結果を出力する）へ進む。
【００３１】
一方、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値とが十分に近づいていない場合には、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値とが十分に近づくまで、ステップＳ２０３〜Ｓ２０５を繰り返し行う。
【００３２】
（４）（ｋ＋１）回目に求めたパラメータｘ0(k+1)、ｙ0(k+1)、α(k+1)、β(k+1)の値を、形状計測データＰiに対応する理論形状を記述するパラメータとして出力する（ステップＳ２０６）。
【００３３】
なお、前述した説明では、４つのパラメータｘ0、ｙ0、α、βを、パラメータｘ0、ｙ0の組と、パラメータα、βの組に分けたが、組の分け方はこれに限ったものではない。複数の組に分けることで各ステップ（ステップＳ２０３、Ｓ２０４）における最適解の探索空間（最適値探索問題の解空間）を狭めることができれば、他の分け方でもよい。
【００３４】
また、前述した（３−１）と、（３−２）のステップＳ２０３、Ｓ２０４におけるパラメータの最適値の探索は、一般には非線形の最小自乗法である。したがって、各ステップＳ２０３、Ｓ２０４内でさらに反復計算が必要となる。
ところが、前述したように、本実施形態では、被測定物の形状を測定するときに工夫をして、中心軸１０１がｚ軸に略平行となる位置で測定されるので、中心軸１０１の方向を記述するパラメータα、βは、１に比べて十分小さな量であることが期待される。その場合、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0，ｙ0，α，β）は、α、βの２次形式で表される。そうすると、（３−２）のステップＳ２０４では、線形連立方程式の問題を解けばよいことなり、極めて簡単にパラメータα、βを算出することができる。
【００３５】
つまり、パラメータの組の分け方としては、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓを、２次形式で表現できるパラメータの組とそうでない組とに分けることが望ましい。このようにしてパラメータを分けることによって、各ステップＳ２０３、Ｓ２０４における最適解の探索空間を狭めるだけでなく、線形探索となり、計算コストと、安定性とが著しく向上する。
【００３６】
また、ここでは、計算順序として、（３−１）のステップＳ２０３を行った後に、（３−２）のステップＳ２０４を行う場合を例に挙げて説明したが、初期値の与え方によっては、（３−２）のステップＳ２０４を行った後に、（３−１）のステップＳ２０３を行ってもよい場合がある。
【００３７】
（実施例）
次に、本実施形態の実施例として、円筒の側面上に配置されるべき形状計測データＰiから中心軸を決定する方法の一例を示す。本実施例では、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0，ｙ0，α，β）が最小となるパラメータｘ0、ｙ0、α、βを決定した。ここで、３次元形状計測器では、中心軸１０１がｚ軸と略平行になるようにして測定されている。よって、パラメータα、βは１よりも十分小さい（α、β≪１）とみなしてよい。そうすると、形状計測データＰiから中心軸までの距離の誤差の２乗和Ｓ（ｘ0，ｙ0，α，β）を、以下の（７式）のように書くことができる。（７式）に示すように、形状計測データＰiから中心軸までの距離の誤差の２乗和Ｓは、４つのパラメータｘ0、ｙ0、α、βのうち、パラメータα、βについての２次形式と見ることができる。
【００３８】
【数５】
【００３９】
そこで、パラメータｘ0、ｙ0、α、βの組の分け方としては、少なくとも１つの組は、パラメータα、βの中から選択されたパラメータのみを含む組にしたい。そうすることによって、その組のパラメータを算出するための計算は線形計算となる。本実施例においては４つのパラメータｘ0、ｙ0、α、βを、パラメータα、βの組と、パラメータｘ0、ｙ0の組とに分けるようにした。そして、パラメータα、βに関しては線形方程式を解くことで求め、パラメータｘ0、ｙ0に関しては非線形方程式を解くことで解を求めた。ここで、非線形方程式の反復解法としては、ニュートン・ラフソン法を用いた。
【００４０】
計算精度を確認するために、パラメータα、β、ｘ0、ｙ0が、それぞれ０．１０、０．０５、１．０、２．０（α＝０．１０、β＝０．０５、ｘ0＝１．０、ｙ0＝２．０）となる半径Ｒが３（Ｒ＝３）の円筒の側面上に配置された６０点の擬似形状計測データを生成し、生成した擬似形状計測データから、パラメータα、β、ｘ0、ｙ0の値を逆算した。図３に、３次元座標上にプロットされた擬似形状計測データの一例を、図４に、入力した擬似形状計測データの値の一例をそれぞれ示す。
【００４１】
各パラメータα、β、ｘ0、ｙ0の初期値α(0)、β(0)、ｘ0(0)、ｙ0(0)を、それぞれ０．０、０．０、１．１９８７６、２．０９９３８（α(0)＝０．０、β(0)＝０．０、ｘ0(0)＝１．１９８７６、ｙ0(0)＝２．０９９３８）として計算した結果、パラメータα、β、ｘ0、ｙ0の値は、それぞれ０．０９９３８７２、０．０４９６９３５、１．００１２２、２．０００６１（α＝０．０９９３８７２、β＝０．０４９６９３５、ｘ0＝１．００１２２、ｙ0＝２．０００６１）となり、本実施形態の手法が有効であることが確認できた。
【００４２】
すなわち、３次元形状の計測データに対応する理論形状を規定するパラメータが多数の場合であっても、高速に且つロバスト（robust）にパラメータを決定することができるので、接触式、非接触式などの３次元座標計測器を用いて計測された３次元形状の評価を効率的に行うことができる。具体的には、回転対称体の形状計測データＰiの中心軸１０１を正しく算出し、傾いて測定された形状計測データＰiに関して高速に且つロバスト（robust）にその傾きを補正することができる。
【００４３】
以上のようにしてパラメータｘ0、ｙ0、α、βを解く処理を実現するハードウェア資源としては、例えば図５に示す情報処理装置が挙げられる。
図５は、前記情報処理装置の構成の一例を示したブロック図である。
図５において、情報処理装置１２０は、ＣＰＵ１２１と、ＲＯＭ１２２と、ＲＡＭ１２３と、キーボード（ＫＢ）１２４のキーボードコントローラ（ＫＢＣ）１２５と、表示部としてのＣＲＴディスプレイ（ＣＲＴ）１２６のＣＲＴコントローラ（ＣＲＴＣ）１２７と、ハードディスク（ＨＤ）１２８及びフレキシブルディスク（ＦＤ）１２９のディスクコントローラ（ＤＫＣ）１３０と、ネットワーク１３１との接続のためのネットワークインターフェースコントローラ（ＮＩＣ）１３２とが、システムバス１１３を介して互いに通信可能に接続された構成としている。
【００４４】
ＣＰＵ１２１は、ＲＯＭ１２２或いはＨＤ１２８に記憶されたソフトウェア、或いはＦＤ１２９より供給されるソフトウェアを実行することで、システムバス１２３に接続された各構成部を総括的に制御する。
すなわち、ＣＰＵ１２１は、所定の処理シーケンスに従った処理プログラムを、ＲＯＭ１２２、或いはＨＤ１２８、或いはＦＤ１２９から読み出して実行することで、後述する動作を実現するための制御を行う。
【００４５】
ＲＡＭ１２３は、ＣＰＵ１２１の主メモリ或いはワークエリア等として機能する。
ＫＢＣ１２５は、ＫＢ１２４や図示していないポインティングデバイス等からの指示入力を制御する。
【００４６】
ＣＲＴＣ１２７は、ＣＲＴ１２６の表示を制御する。
ＤＫＣ１３０は、ブートプログラム、種々のアプリケーション、編集ファイル、ユーザファイル、ネットワーク管理プログラム、及び本実施の形態における所定の処理プログラム等を記憶するＨＤ１２８及びＦＤ１２９とのアクセスを制御する。
ＮＩＣ１３２は、ネットワーク１３１上の装置或いはシステムと双方向にデータをやりとりする。
【００４７】
ＣＰＵ１２１は、ＲＯＭ１２２或いはＨＤ１２８に記憶されたソフトウェア、或いはＦＤ１２９より供給されるソフトウェアを実行することで、システムバス１２３に接続された各構成部を総括的に制御する。
すなわち、ＣＰＵ１２１は、所定の処理シーケンスに従った処理プログラムを、ＲＯＭ１２２、或いはＨＤ１２８、或いはＦＤ１２９から読み出して実行することで、図２に示したフローチャートの処理を実現するための制御を行う。なお、ステップＳ２０６における出力は、例えば、ＣＲＴ１２６に表示することによって実現される。
【００４８】
（第２の実施形態）
次に、本発明の第２の実施形態について説明する。前述した第１の実施形態では、回転対称体が円筒であり、その円筒の回転半径Ｒ（ｚ）が一定である場合について説明したが、本実施形態では、一般的な回転対称体の円筒の場合について説明する。このように、本実施形態と前述した第１の実施形態とは、形状計測データＰiに対応する理論形状を記述するパラメータを求める対象となる３次元形状が異なるだけである。したがって、第１の実施形態と同一の部分については、図１〜図５に付した符号と同一の符号を付すなどして、詳細な説明を省略する。
【００４９】
一般的な回転対称体の円筒の場合、すなわち、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0,ｙ0,ｚ0,α,β）を最小にする場合には、形状計測データＰiに対応する理論形状を記述するパラメータｘ0、ｙ0、ｚ0、α、βを、２つの組に分けるようにする。例えば、パラメータｘ0、ｙ0の組と、パラメータｚ0、α、βの組にパラメータが分けられる設定がなされた上で、図６に示すフローチャートの処理を行う。
（１） 形状計測データＰiを入力する（ステップＳ６０１）。
（２） 各パラメータの初期値ｘ0(0)、ｙ0(0)、ｚ0(0)、α(0)、β(0)を設定する（ステップＳ６０２）。
【００５０】
（３） パラメータを反復計算する（ｋ回目の値をもとに、（ｋ＋１）回目の値を決定する）。
（３−１） パラメータｚ0(k)、α(k)、β(k)を既知量として、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0(k+1)，ｙ0(k+1)，ｚ0(k+1)，α(k)，β(k)）を最小とするパラメータｘ0(k+1)、ｙ0(k+1)を決定する。すなわち、前述した（５式）の連立方程式からパラメータｘ0(k+1)、ｙ0(k+1)を求める（ステップＳ６０３）
【００５１】
（３−２） 求めたパラメータｘ0(k+1)、ｙ0(k+1)を既知量として、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0(k+1)，ｙ0(k+1)，ｚ0(k+1)，α(k+1)，β(k+1)）を最小とするパラメータｚ0(k+1)、α(k+1)、β(k+1)を決定する。すなわち、以下の（８式）の連立方程式から、パラメータｚ0(k+1)、α(k+1)、β(k+1)を求める（ステップＳ６０４）。
【００５２】
後述するように、本実施形態では、被測定物の形状を測定するときに工夫をして、中心軸１０１がほぼｚ軸に平行となる位置で測定されたとする。そうすると、中心軸１０１の方向を記述するパラメータα、βは、１に比べて十分小さな量となり、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0，ｙ0，ｚ0，α，β）は、α、βの２次形式で表される。したがって、このステップＳ６０４では、線形連立方程式の問題を解いて、パラメータα(k+1)、β(k+1)を算出することになる。
【００５３】
【数６】
【００５４】
（３−３） ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下であるか否かを判定する（ステップＳ６０５）。この判定の結果、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下となり、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値とが十分に近づいたならば、パラメータｘ0(k+1)、ｙ0(k+1)、ｚ0(k+1)、α(k+1)、β(k+1)を求める処理を終了して、ステップＳ６０６（（４）結果を出力する）へ進む。
【００５５】
一方、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値とが十分に近づいていない場合には、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値とが十分に近づくまで、ステップＳ６０３〜Ｓ６０５を繰り返し行う。
（４）（ｋ＋１）回目に求めたパラメータｘ0(k+1)、ｙ0(k+1)、ｚ0(k+1)、α(k+1)、β(k+1)の値を、形状計測データＰiに対応する理論形状を記述するパラメータとして出力する（ステップＳ６０６）。
【００５６】
なお、第１の実施形態と同様に、５つのパラメータｘ0、ｙ0、ｚ0、α、βを、複数の組に分けることで各ステップ（ステップＳ６０３、６０４）における最適解の探索空間を狭めることができれば、組の分け方は、図６のフローチャートで示した方法でなくてもよい。
また、（３−１）、（３−２）のステップＳ６０３、Ｓ６０４における最適値の探索は、一般には非線形の最小自乗法であり、各ステップＳ２０３、Ｓ２０４内でさらに反復計算が必要となる。
【００５７】
ところが、前述したように、本実施形態では、被測定物の形状を測定するときに工夫をして、中心軸１０１がほぼｚ軸に平行となる位置で測定されるので、中心軸１０１の方向を記述するパラメータα、β、ｚ0は、１に比べて十分小さな量であることが期待される。その場合、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0，ｙ0，ｚ0，α，β）は、α、β、ｚ0の２次形式で表され、（３−２）のステップＳ６０４では、線形連立方程式の問題を解けばよいことになり極めて簡単にパラメータα、β、ｚ0を算出することができる。このようにしてパラメータを分けることによって、各ステップＳ６０３、６０４における最適解の探索空間を狭めるだけでなく、線形探索となり、計算コストと、安定性とが著しく向上する。
【００５８】
また、初期値の与え方によっては、（３−２）のステップＳ６０４を行った後に、（３−１）のステップＳ６０３を行ってもよい場合があることも第１の実施形態と同様である。
【００５９】
（実施例）
次に、本実施形態の実施例として、円錐の側面上に配置されるべき形状計測データＰiから中心軸を決定する方法の一例を示す。本実施例では、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0，ｙ0，ｚ0，α，β）が最小となるパラメータｘ0、ｙ0、ｚ0、α、βを決定した。ここで、３次元形状計測器では、中心軸１０１がz軸と略平行になるようにすると共に、基準点のｚ座標も略ゼロ（０）となるように測定されている。よって、パラメータｚ0、α、βは１よりも十分小さい（ｚ0、α、β≪１）とみなしてよい。そうすると、形状計測データＰiから中心軸までの距離の誤差の２乗和Ｓ（ｘ0，ｙ0，ｚ0，α，β）を、以下の（９式）のように書くことができる。（９式）に示すように、形状計測データＰiから中心軸までの距離の誤差の２乗和Ｓは、５つのパラメータｘ0、ｙ0、ｚ0、α、βのうち、パラメータｚ0、α、βについての２次形式と見ることができる。
【００６０】
【数７】
【００６１】
そこで、パラメータｘ0、ｙ0、ｚ0、α、βの組の分け方としては、少なくとも１つの組は、パラメータｚ0、α、βの中から選択されたパラメータのみを含む組にしたい。そうすることによって、その組のパラメータを算出するための計算は線形計算となる。本実施例においては５つのパラメータｘ0、ｙ0、ｚ0、α、βを、パラメータｚ0、α、βの組と、パラメータｘ0，ｙ0の組とに分けるようにした。そして、パラメータｚ0、α、βに関しては線形方程式を解くことで求め、パラメータｘ0、ｙ0に関しては非線形方程式を解くことで解を求めた。ここで、非線形方程式の反復解法としては、ニュートン・ラフソン法を用いた。
【００６２】
計算精度を確認するために、パラメータα、β、ｘ0、ｙ0、ｚ0が、それぞれ０．０２、０．０１、１．０、２．０、０．０２（α＝０．０２、β＝０．０１、ｘ0＝１．０、ｙ0＝２．０、ｚ0＝０．０２）となる円錐の面上にらせん状に配置された５０点の擬似形状計測データを生成し、生成した擬似形状計測データからパラメータα、β、ｘ0、ｙ0、ｚ0の値を逆算した。図７に、３次元座標上にプロットされた擬似形状計測データの一例を、図８に、入力した擬似形状計測データの値の一例をそれぞれ示す。
【００６３】
各パラメータα、β、ｘ0、ｙ0、ｚ0の初期値α(0)、β(0)、ｘ0(0)、ｙ0(0)、ｚ0(0)を、それぞれ０．０（α(0)＝０．０、β(0)＝０．０、ｘ0(0)＝０．０、ｙ0(0)＝０．０、ｚ0(0)＝０．０）として計算した結果、パラメータα、β、ｘ0、ｙ0、ｚ0の値は、それぞれ０．０２０５９５５、０．００７３０５９９、０．９９９６６１、２．０００６８、０．０１８６３５９（α＝０．０２０５９５５、β＝０．００７３０５９９、ｘ0＝０．９９９６６１、ｙ0＝２．０００６８、ｚ0＝０．０１８６３５９）となり、本実施形態の手法が有効であることが確認できた。
【００６４】
すなわち、３次元形状の計測データに対応する理論形状を規定するパラメータが多数の場合であっても、高速に且つロバスト（robust）にパラメータを決定することができるので、接触式、非接触式などの３次元座標計測器を用いて計測された３次元形状の評価を効率的に行うことができる。具体的には、回転対称体の形状計測データＰiの中心軸１０１を正しく算出し、傾いて測定された形状計測データＰiに関して高速に且つロバスト（robust）にその傾きを補正することができる。
なお、第１の実施形態と同様に、例えば図５に示した情報処理装置を用いることにより、以上のようにしてパラメータｘ0、ｙ0、ｚ0、α、βを解く処理を実現することができる。
【００６５】
（第３の実施形態）
次に、本発明の第３の実施形態について説明する。前述した第１の実施形態では、回転対称体が円筒であり、その円筒の回転半径Ｒ（ｚ）が一定である場合について説明したが、本実施形態では、回転対称体は円筒であるが、その円筒の回転半径Ｒ（ｚ）もパラメータとする場合について説明する。このように、本実施形態と前述した第１の実施形態とは、形状計測データＰiに対応する理論形状を記述するパラメータを求める対象となる３次元形状の一部が異なるだけである。したがって、第１の実施形態と同一の部分については、図１〜図５に付した符号と同一の符号を付すなどして、詳細な説明を省略する。
【００６６】
回転半径Ｒ（ｚ）をパラメータとする円筒の場合、すなわち、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0,ｙ0,Ｒ，α,β）を最小にする場合には、形状計測データＰiに対応する理論形状を記述するパラメータｘ0、ｙ0、Ｒ、α、βを２つの組に分けるようにする。例えば、パラメータｘ0、ｙ0の組と、パラメータα、β、Ｒの組にパラメータが分けられる設定がなされた上で、図９に示すフローチャートの処理を行う。
（１）形状計測データＰiを入力する（ステップＳ９０１）。
（２）各パラメータの初期値ｘ0(0)、ｙ0(0)、Ｒ(0)、α(0)、β(0)を設定する（ステップＳ９０２）。
（３）パラメータを反復計算する（ｋ回目の値をもとに、（ｋ＋１）回目の値を決定する）。
（３−１） パラメータα(k)、β(k)を既知量として、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0(k+1)，ｙ0(k+1)，Ｒ(k+1)，α(k)，β(k)）を最小とするパラメータｘ0(k+1)、ｙ0(k+1)、Ｒ(k+1)を決定する。すなわち、以下の（１０式）の連立方程式からパラメータｘ0(k+1)、ｙ0(k+1)、Ｒ(k+1)を求める（ステップＳ９０３）。
【００６７】
【数８】
【００６８】
（３−２） 求めたパラメータｘ0(k+1)、ｙ0(k+1)、Ｒ(k+1)を既知量として、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0(k+1)，ｙ0(k+1)，Ｒ(k+1)，α(k+1)，β(k+1)）を最小とするパラメータα(k+1)、β(k+1)を決定する。すなわち、前述した（６式）の連立方程式から、パラメータα(k+1)、β(k+1)を求める（ステップＳ９０４）。
【００６９】
後述するように、本実施形態では、被測定物の形状を測定するときに工夫をして、中心軸１０１がほぼｚ軸に平行となる位置で測定されたとする。そうすると、中心軸１０１の方向を記述するパラメータα、βは、１に比べて十分小さな量となり、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0，ｙ0，Ｒ，α，β）は、α、βの２次形式で表される。したがって、このステップＳ２０４では、線形連立方程式の問題を解いて、パラメータα(k+1)、β(k+1)を算出することになる。
【００７０】
（３−３） ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下であるか否かを判定する（ステップＳ９０５）。この判定の結果、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下となり、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値とが十分に近づいたならば、パラメータｘ0(k+1)、ｙ0(k+1)、Ｒ(k+1)、α(k+1)、β(k+1)を求める処理を終了して、ステップＳ９０６（（４）結果を出力する）へ進む。
【００７１】
一方、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値とが十分に近づいていない場合には、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値とが十分に近づくまで、ステップＳ９０３〜Ｓ９０５を繰り返し行う。
（４）（ｋ＋１）回目に求めたパラメータｘ0(k+1)、ｙ0(k+1)、Ｒ(k+1)、α(k+1)、β(k+1)の値を、形状計測データＰiに対応する理論形状を記述するパラメータとして出力する（ステップＳ９０６）
【００７２】
なお、第１の実施形態と同様に、５つのパラメータｘ0、ｙ0、Ｒ、α、βを、複数の組に分けることで各ステップ（ステップＳ９０３、９０４）における最適解の探索空間を狭めることができれば、組の分け方は、図９のフローチャートで示した方法でなくてもよい。
また、（３−１）、（３−２）のステップＳ９０３、Ｓ９０４における最適値の探索は、一般には非線形の最小自乗法であり、各ステップ９０３、９０４内でさらに反復計算が必要となる。
【００７３】
ところが、前述したように、本実施形態では、被測定物の形状を測定するときに工夫をして、中心軸１０１がほぼｚ軸に平行となる位置で測定されるので、中心軸１０１の方向を記述するパラメータα、βは、１に比べて十分小さな量であることが期待される。その場合、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0，ｙ0，Ｒ，α，β）は、α、βの２次形式で表され、（３−２）のステップＳ９０４では、線形連立方程式の問題を解けばよいことになり極めて簡単にパラメータα、βを算出することができる。このようにしてパラメータを分けることによって、各ステップＳ９０３、９０４における最適解の探索空間を狭めるだけでなく、線形探索となり、計算コストと、安定性とが著しく向上する。
また、初期値の与え方によっては、（３−２）のステップＳ９０４を行った後に、（３−１）のステップＳ９０３を行ってもよい場合があることも第１の実施形態と同様である。
【００７４】
（実施例）
次に、本実施形態の実施例として、円筒の側面上に配置されるべき形状計測データＰiから中心軸を決定する方法の一例を示す。本実施例では、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0，ｙ0，Ｒ，α，β）が最小となるパラメータｘ0、ｙ0、Ｒ、α、βを決定した。ここで、３次元形状計測器では、中心軸１０１がz軸と略平行になるように測定されている。よって、パラメータα、βは１よりも十分小さい（α、β≪１）とみなしてよい。そうすると、形状計測データＰiから中心軸までの距離の誤差の２乗和Ｓ（ｘ0，ｙ0，Ｒ，α，β）を、以下の（１１式）のように書くことができる。（１１式）に示すように、形状計測データＰiから中心軸までの距離の誤差の２乗和Ｓは、５つのパラメータｘ0、ｙ0、Ｒ、α、βのうち、パラメータα、β、Ｒについての２次形式と見ることができる。
【００７５】
【数９】
【００７６】
そこで、パラメータｘ0、ｙ0、Ｒ、α、βの組の分け方としては、少なくとも１つの組は、パラメータα、β、Ｒの中から選択されたパラメータのみを含む組にしたい。そうすることによって、その組のパラメータを算出するための計算は線形計算となる。本実施例においては５つのパラメータｘ0、ｙ0、Ｒ、α、βを、パラメータα、βの組と、パラメータｘ0，ｙ0、Ｒの組とに分けるようにした。そして、パラメータα、βに関しては線形方程式を解くことで、ｘ0，ｙ0、Ｒに関しては非線形方程式を解くことで解を求めた。ここで、非線形方程式の反復解法としては、ニュートン・ラフソン法を用いた。
【００７７】
計算精度を確認するために、パラメータα、β、ｘ0、ｙ0、Ｒがそれぞれ０．１０、０．０５、１．０、２．０、３．０（α＝０．１０、β＝０．０５、ｘ0＝１．０、ｙ0＝２．０、Ｒ＝３．０）となる円筒の側面上に配置された６０点の擬似形状計測データを生成し、生成した擬似形状計測データからパラメータα、β、ｘ0、ｙ0、Ｒの値を逆算した。３次元座標上にプロットされた擬似形状計測データと、入力した擬似形状計測データの値は、それぞれ図３及び図４に示した通りである。
【００７８】
各パラメータα、β、ｘ0、ｙ0、Ｒの初期値α(0)、β(0)、ｘ0(0)、ｙ0(0)、Ｒ(0)を、それぞれ０．０（α(0)＝０．０、β(0)＝０．０、ｘ0(0)＝０．０、ｙ0(0)＝０．０、Ｒ(0)＝０．０）として計算した結果、パラメータα、β、Ｒ、ｘ0、ｙ0の値は、それぞれ０．０９９６５４４、０．０４９８２７１、３．００９２９、１．０００６９、２．０００３４（α＝０．０９９６５４４、β＝０．０４９８２７１、Ｒ＝３．００９２９、ｘ0＝１．０００６９、ｙ0＝２．０００３４）となり、本実施形態の手法が有効であることが確認できた。
【００７９】
すなわち、３次元形状の計測データに対応する理論形状を規定するパラメータが多数の場合であっても、高速に且つロバスト（robust）にパラメータを決定することができるので、接触式、非接触式などの３次元座標計測器を用いて計測された３次元形状の評価を効率的に行うことができる。具体的には、回転対称体の形状計測データＰiの中心軸１０１を正しく算出し、傾いて測定された形状計測データＰiに関して高速に且つロバスト（robust）にその傾きを補正することができる。
なお、第１の実施形態と同様に、例えば図５に示した情報処理装置を用いることにより、以上のようにしてパラメータｘ0、ｙ0、Ｒ、α、βを解く処理を実現することができる。
【００８０】
（第４の実施形態）
次に、本発明の第４の実施形態について説明する。前述した第１の実施形態では、ステップＳ２０４において、線形連立方程式の問題を解くことで、極めて簡単にパラメータα、βの最適値を算出することができるようにしたが、本実施形態では、ステップＳ２０３、Ｓ２０４におけるパラメータの最適値の探索を、それぞれ非線形の最小自乗法を用いて行うようにした。このように、本実施形態と前述した第１の実施形態とは、パラメータα、βを求める方法が異なるだけである。したがって、第１の実施形態と同一の部分については、図１〜図５に付した符号と同一の符号を付すなどして、詳細な説明を省略する。
【００８１】
本実施形態においても、第１の実施形態と同様に、回転対称体が、前述した円筒であり、回転半径Ｒ（ｚ）が一定であるので、形状計測データＰiに対応する理論形状を記述するパラメータｘ0、ｙ0、α、βが、パラメータｘ0、ｙ0の組と、パラメータα、βの組とに分けられる設定がなされた上で、図２のフローチャートに従った処理を行う。
（１）形状計測データＰiを入力する（ステップＳ２０１）。
（２）各パラメータｘ0、ｙ0、α、βの初期値ｘ0(0)、ｙ0(0)、α(0)、β(0)を設定する（ステップＳ２０２）。
（３）パラメータを反復計算する（ｋ（ｋは自然数）回目の値をもとに、（ｋ＋１）回目の値を決定する）。
（３−１） 前述した（５式）の連立方程式からパラメータｘ0(k+1)、ｙ0(k+1)を求める（ステップＳ２０３）
【００８２】
（３−２） 前述した（６式）の連立方程式からパラメータα(k+1)、β(k+1)を求める（ステップＳ２０４）。後述するように、本実施形態では、非線形の最小自乗法を用いて、パラメータα(k+1)、β(k+1)を求める。
【００８３】
（３−３） ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下であるか否かを判定する（ステップＳ２０５）。この判定の結果、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下であれば、パラメータｘ0(k+1)、ｙ0(k+1)、α(k+1)、β(k+1)を求める処理を終了して、ステップＳ２０６（（４）結果を出力する）へ進む。
【００８４】
一方、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下でない場合には、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下になるまで、ステップＳ２０３〜Ｓ２０５を繰り返し行う。
【００８５】
（４）（ｋ＋１）回目に求めたパラメータｘ0(k+1)、ｙ0(k+1)、α(k+1)、β(k+1)の値を、形状計測データＰiに対応する理論形状を記述するパラメータとして出力する（ステップＳ２０６）。
【００８６】
なお、第１の実施形態と同様に、４つのパラメータｘ0、ｙ0、α、βを、複数の組に分けることで各ステップ（ステップＳ２０３、Ｓ２０４）における最適解の探索空間（最適値探索問題の解空間）を狭めることができれば、組の分け方は、図２のフローチャートで示した方法でなくてもよい。
【００８７】
また、本実施形態の（３−１）、（３−２）のステップＳ２０３、Ｓ２０４におけるパラメータの最適値の探索は、非線形の最小自乗法を用いて行われ、ステップＳ２０３、Ｓ２０４内でさらに反復計算が必要となる。
【００８８】
以上のようにしてパラメータを複数の組に分けることによって、各ステップＳ２０３、２０４における最適値の探索空間を狭めることで、計算コストと、安定性とが著しく向上する。
また、初期値の与え方によっては、（３−２）のステップＳ２０４を行った後に、（３−１）のステップＳ２０３を行ってもよい場合があることも第１の実施形態と同様である。
【００８９】
（実施例）
次に、本実施形態の実施例として、円筒の側面上に配置されるべき形状計測データＰiから中心軸を決定する方法の一例を示す。本実施例では、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0，ｙ0，α，β）が最小となるパラメータｘ0、ｙ0、α、βを決定した。ここで、４つのパラメータｘ0，ｙ0，α，βを、パラメータα、βの組と、パラメータｘ0、ｙ0の組とに分けて、それぞれ非線形の最小自乗法を用いて解を求めた。ここで、非線形方程式の反復解法としては、ニュートン・ラフソン法を用いた。
【００９０】
計算精度を確認するために、パラメータα、β、ｘ0、ｙ0が、それぞれ０．１０、０．０５、１．０、２．０（α＝０．１０、β＝０．０５、ｘ0＝１．０、ｙ0＝２．０）となる半径Ｒが３（Ｒ＝３）の円筒の側面上に配置された６０点の擬似形状計測データを生成し、生成した擬似形状計測データから、パラメータα、β、ｘ0、ｙ0の値を逆算した。３次元座標上にプロットされた擬似形状計測データと、入力した擬似形状計測データの値は、それぞれ図３、図４に示した通りである。
【００９１】
各パラメータα、β、ｘ0、ｙ0の初期値α(0)、β(0)、ｘ0(0)、ｙ0(0)を、それぞれ０．０、０．０、１．１９８７６、２．０９９３８（α(0)＝０．０、β(0)＝０．０、ｘ0(0)＝１．１９８７６、ｙ0(0)＝２．０９９３８）として計算した結果、パラメータα、β、ｘ0、ｙ0の値は、それぞれ０．０９９９９９８、０．０４９９９９８、１．０００００、２．０００００（α＝０．０９９９９９８、β＝０．０４９９９９８、ｘ0＝１．０００００、ｙ0＝２．０００００）となり、本実施形態の手法がパラメータを高精度に計算するのに極めて有効であることが確認できた。
【００９２】
すなわち、３次元形状の計測データに対応する理論形状を規定するパラメータが多数の場合であっても、高速に且つロバスト（robust）にパラメータを決定することができるので、接触式、非接触式などの３次元座標計測器を用いて計測された３次元形状の評価を効率的に行うことができる。具体的には、回転対称体の形状計測データＰiの中心軸１０１を正しく算出し、傾いて測定された形状計測データＰiに関して高速に且つロバスト（robust）にその傾きを補正することができる。
なお、第１の実施形態と同様に、例えば図５に示した情報処理装置を用いることにより、以上のようにしてパラメータｘ0、ｙ0、α、βを解く処理を実現することができる。
【００９３】
（第５の実施形態）
次に、本発明の第５の実施形態について説明する。前述した第３の実施形態では、ステップＳ９０４において、線形連立方程式の問題を解くことで、極めて簡単にパラメータα、βの最適値を算出することができるようにしたが、本実施形態では、ステップＳ９０３、Ｓ９０４におけるパラメータの最適値の探索を、それぞれ非線形の最小自乗法を用いて行うようにした。このように、本実施形態と前述した第３の実施形態とは、パラメータα、βを求める方法が異なるだけである。したがって、第３の実施形態と同一の部分については、図１〜図５、図９に付した符号と同一の符号を付すなどして、詳細な説明を省略する。
【００９４】
本実施形態においても、第３の実施形態と同様に、回転対称体は、回転半径Ｒ（ｚ）をパラメータとする円筒であるので、形状計測データＰiに対応する理論形状を記述するパラメータｘ0、ｙ0、Ｒ、α、βが、パラメータｘ0、ｙ0の組と、パラメータα、β、Ｒの組とに分けられる設定がなされた上で、図９のフローチャートに従った処理を行う。
（１）形状計測データＰiを入力する（ステップＳ９０１）。
（２）各パラメータの初期値ｘ0(0)、ｙ0(0)、Ｒ(0)、α(0)、β(0)を設定する（ステップＳ９０２）。
（３）パラメータを反復計算する（ｋ回目の値をもとに、（ｋ＋１）回目の値を決定する）。
（３−１） 前述した（１０式）の連立方程式から、パラメータｘ0(k+1)、ｙ0(k+1)、Ｒ(k+1)を求める（ステップＳ９０３）。
【００９５】
（３−２） 前述した（６式）の連立方程式から、パラメータα(k+1)、β(k+1)を求める（ステップＳ９０４）。後述するように、本実施形態では、非線形の最小自乗法を用いて、パラメータα(k+1)、β(k+1)を求める。
【００９６】
（３−３）ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下であるか否かを判定する（ステップＳ９０５）。この判定の結果、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下であれば、パラメータｘ0(k+1)、ｙ0(k+1)、Ｒ(k+1)、α(k+1)、β(k+1)を求める処理を終了して、ステップＳ９０６（（４）結果を出力する）へ進む。
【００９７】
一方、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下でない場合には、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下になるまで、ステップＳ９０３〜Ｓ９０５を繰り返し行う。
【００９８】
（４）（ｋ＋１）回目に求めたパラメータｘ0(k+1)、ｙ0(k+1)、Ｒ(k+1)、α(k+1)、β(k+1)の値を、形状計測データＰiに対応する理論形状を記述するパラメータとして出力する（ステップＳ９０６）。
【００９９】
なお、第３の実施形態と同様に、５つのパラメータｘ0、ｙ0、Ｒ、α、βを、複数の組に分けることで各ステップ（ステップＳ９０３、Ｓ９０４）における最適解の探索空間（最適値探索問題の解空間）を狭めることができれば、組の分け方は、図９のフローチャートで示した方法でなくてもよい。
【０１００】
また、（３−１）、（３−２）のステップＳ９０３、Ｓ９０４における最適値の探索は、非線形の最小自乗法を用いて行われ、ステップＳ９０３、Ｓ９０４内でさらに反復計算が必要となる。
【０１０１】
以上のようにしてパラメータを複数の組に分けることによって、各ステップＳ９０３、Ｓ９０４における最適値の探索空間を狭めることで、計算コスト、安定性が著しく向上する。
また、初期値の与え方によっては、（３−２）のステップＳ９０４を行った後に、（３−１）のステップＳ９０３を行ってもよい場合があることも第３の実施形態と同様である。
【０１０２】
（実施例）
次に、本実施形態の実施例として、円筒の側面上に配置されるべき形状計測データＰiから中心軸を決定する方法の一例を示す。本実施例では、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0，ｙ0，Ｒ，α，β）が最小となるパラメータｘ0、ｙ0、Ｒ、α、βを決定した。ここで、５つのパラメータｘ0、ｙ0、Ｒ、α、βを、パラメータα、βの組と、パラメータｘ0、ｙ0、Ｒの組とに分けて、それぞれ非線形の最小自乗法を用いて解を求めた。ここで、非線形方程式の反復解法としては、ニュートン・ラフソン法を用いた。
【０１０３】
計算精度を確認するために、パラメータα、β、ｘ0、ｙ0、Ｒがそれぞれ０．１０、０．０５、１．０、２．０、３．０（α＝０．１０、β＝０．０５、ｘ0＝１．０、ｙ0＝２．０、Ｒ＝３．０）となる円筒の側面上に配置された６０点の擬似形状計測データを生成し、生成した擬似形状計測データからパラメータα、β、ｘ0、ｙ0、Ｒの値を逆算した。３次元座標上にプロットされた擬似形状形状データと、入力した擬似形状計測データの値は、それぞれ図３及び図４に示した通りである。
【０１０４】
各パラメータα、β、ｘ0、ｙ0、Ｒの初期値α(0)、β(0)、ｘ0(0)、ｙ0(0)、Ｒ(0)を、それぞれ０．０（α(0)＝０．０、β(0)＝０．０、ｘ0(0)＝０．０、ｙ0(0)＝０．０、Ｒ(0)＝０．０）として計算した結果、パラメータα、β、Ｒ、ｘ0、ｙ0の値は、それぞれ０．０９９９９９８、０．０４９９９９８、３．０００００、１．０００００、２．０００００（α＝０．０９９９９９８、β＝０．０４９９９９８、Ｒ＝３．０００００、ｘ0＝１．０００００、ｙ0＝２．０００００）となり、本実施形態の手法がパラメータを高精度に計算するのに極めて有効であることが確認できた。
【０１０５】
すなわち、３次元形状の計測データに対応する理論形状を規定するパラメータが多数の場合であっても、高速に且つロバスト（robust）にパラメータを決定することができるので、接触式、非接触式などの３次元座標計測器を用いて計測された３次元形状の評価を効率的に行うことができる。具体的には、回転対称体の形状計測データＰiの中心軸１０１を正しく算出し、傾いて測定された形状計測データＰiに関して高速に且つロバスト（robust）にその傾きを補正することができる。
なお、第３の実施形態と同様に、例えば図５に示した情報処理装置を用いることにより、以上のようにしてパラメータｘ0、ｙ0、Ｒ、α、βを解く処理を実現することができる。
【０１０６】
（第６の実施形態）
次に、本発明の第６の実施形態について説明する。前述した第２の実施形態では、ステップＳ６０４において、線形連立方程式の問題を解くことで、極めて簡単にパラメータα、βの最適値を算出することができるようにしたが、本実施形態では、ステップＳ６０３、Ｓ６０４におけるパラメータの最適値の探索を、それぞれ非線形の最小自乗法を用いて行うようにした。このように、本実施形態と前述した第２の実施形態とは、パラメータα、βを求める方法が異なるだけである。したがって、第２の実施形態と同一の部分については、図１〜図２、図６〜図８に付した符号と同一の符号を付すなどして、詳細な説明を省略する。
【０１０７】
本実施形態においても、第２の実施形態と同様に、回転対称体は、一般的な回転対称体の円筒であるので、形状計測データＰiに対応する理論形状を記述するパラメータｘ0、ｙ0、ｚ0、α、βが、パラメータｘ0、ｙ0の組と、パラメータｚ0、α、βの組にすなわちＳ（ｘ0,ｙ0,ｚ0,α,β）の場合を例に説明する。
まず、理論形状記述パラメータを２つの組、ｘ0、ｙ0の組とに分けられる設定がなされた上で、図６のフローチャートに従った処理を行う。
（１）形状計測データＰiを入力する（ステップＳ６０１）。
（２）各パラメータの初期値ｘ0(0)、ｙ0(0)、ｚ0(0)、α(0)、β(0)を設定する（ステップＳ６０２）。
（３）パラメータを反復計算する（ｋ回目の値をもとに、（ｋ＋１）回目の値を決定する）。
（３−１）前述した（５式）の連立方程式から、パラメータｘ0(k+1)、ｙ0(k+1)を求める（ステップＳ６０３）。
【０１０８】
（３−２）前述した（６式）の連立方程式から、パラメータα(k+1)、β(k+1)を求める（ステップＳ６０４）。後述するように、本実施形態では、非線形の最小自乗法を用いて、パラメータα(k+1)、β(k+1)を求める。
【０１０９】
（３−３）ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下であるか否かを判定する（ステップＳ９０５）。この判定の結果、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下であれば、パラメータｘ0(k+1)、ｙ0(k+1)、Ｒ(k+1)、α(k+1)、β(k+1)を求める処理を終了して、ステップＳ６０６（（４）結果を出力する）へ進む。
【０１１０】
一方、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下でない場合には、ｋ回目に求めたパラメータの値と、（ｋ＋１）回目に求めたパラメータの値との差が閾値以下になるまで、ステップＳ６０３〜Ｓ６０５を繰り返し行う。
【０１１１】
（４）（ｋ＋１）回目に求めたパラメータｘ0(k+1)、ｙ0(k+1)、Ｒ(k+1)、α(k+1)、β(k+1)の値を、形状計測データＰiに対応する理論形状を記述するパラメータとして出力する（ステップＳ６０６）。
【０１１２】
なお、第２の実施形態と同様に、５つのパラメータｘ0、ｙ0、ｚ0、α、βを、複数の組に分けることで各ステップ（ステップＳ６０３、６０４）における最適解の探索空間を狭めることができれば、組の分け方は、図６のフローチャートで示した方法でなくてもよい。
【０１１３】
また、（３−１）、（３−２）のステップＳ６０３、Ｓ６０４における最適値の探索は、非線形の最小自乗法を用いて行われ、ステップＳ６０３、Ｓ６０４内でさらに反復計算が必要となる。
【０１１４】
以上のようにしてパラメータを複数の組に分けることによって、各ステップＳ６０３、Ｓ６０４における最適値の探索空間を狭めることで、計算コスト、安定性が著しく向上する。
また、初期値の与え方によっては、（３−２）のステップＳ６０４を行った後に、（３−１）のステップＳ６０３を行ってもよい場合があることも第２の実施形態と同様である。
【０１１５】
（実施例）
次に、本実施形態の実施例として、円錐の側面上に配置されるべき形状計測データＰiから中心軸を決定する方法の一例を示す。本実施例では、形状計測データＰiから中心軸１０１までの距離の誤差の２乗和Ｓ（ｘ0，ｙ0，ｚ0，α，β）が最小となるパラメータｘ0、ｙ0、ｚ0、α、βを決定した。ここで、５つのパラメータｘ0、ｙ0、ｚ0、α、βを、パラメータｚ0、α、βの組と、パラメータｘ0，ｙ0の組とに分けて、それぞれ非線形の最小自乗法を用いて解を求めた。非線形方程式の反復解法としては、ニュートン・ラフソン法を用いた。
【０１１６】
計算精度を確認するために、パラメータα、β、ｘ0、ｙ0、ｚ0がそれぞれ０．０２、０．０１、１．０、２．０、０．０２（α＝０．０２、β＝０．０１、ｘ0＝１．０、ｙ0＝２．０、ｚ0＝０．０２）となる円錐面上にらせん状に配置された５０点の擬似形状計測データを生成し、生成した擬似形状計測データからパラメータα、β、ｘ0、ｙ0、ｚ0の値を逆算した。３次元座標上にプロットされた擬似形状形状データと、入力した擬似形状計測データの値は、それぞれ図７及び図８に示した通りである。
【０１１７】
各パラメータα、β、ｘ0、ｙ0、ｚ0の初期値α(0)、β(0)、ｘ0(0)、ｙ0(0)、ｚ0(0)を、それぞれ０．０（α(0)＝０．０、β(0)＝０．０、ｘ0(0)＝０．０、ｙ0(0)＝０．０、ｚ0(0))＝０．０）として計算した結果、パラメータα(0)、β(0)、ｘ0(0)、ｙ0(0)、ｚ0の値は、それぞれ０．０１９９９９９８、０．００９９９９９８、０．９９９９９９９９、２．０００００００、０．０１９９９９８（α＝０．０１９９９９９８、β＝０．００９９９９９８、ｘ0＝０．９９９９９９９９、ｙ0＝２．０００００００、ｚ0＝０．０１９９９９８）となり、本実施形態の手法がパラメータを高精度に計算するのに極めて有効であることが確認できた。
【０１１８】
すなわち、３次元形状の計測データに対応する理論形状を規定するパラメータが多数の場合であっても、高速に且つロバスト（robust）にパラメータを決定することができるので、接触式、非接触式などの３次元座標計測器を用いて計測された３次元形状の評価を効率的に行うことができる。具体的には、回転対称体の形状計測データＰiの中心軸１０１を正しく算出し、傾いて測定された形状計測データＰiに関して高速に且つロバスト（robust）にその傾きを補正することができる。
なお、第２の実施形態と同様に、例えば図５に示した情報処理装置を用いることにより、以上のようにしてパラメータｘ0、ｙ0、Ｒ、α、βを解く処理を実現することができる。
【０１１９】
（本発明の他の実施形態）
上述した実施形態の機能を実現するべく各種のデバイスを動作させるように、該各種デバイスと接続された装置あるいはシステム内のコンピュータに対し、前記実施形態の機能を実現するためのソフトウェアのプログラムコードを供給し、そのシステムあるいは装置のコンピュータ（ＣＰＵあるいはＭＰＵ）に格納されたプログラムに従って前記各種デバイスを動作させることによって実施したものも、本発明の範疇に含まれる。
【０１２０】
また、この場合、前記ソフトウェアのプログラムコード自体が上述した実施形態の機能を実現することになり、そのプログラムコード自体、及びそのプログラムコードをコンピュータに供給するための手段、例えば、かかるプログラムコードを格納した記録媒体は本発明を構成する。かかるプログラムコードを記憶する記録媒体としては、例えばフレキシブルディスク、ハードディスク、光ディスク、光磁気ディスク、ＣＤ−ＲＯＭ、磁気テープ、不揮発性のメモリカード、ＲＯＭ等を用いることができる。
【０１２１】
また、コンピュータが供給されたプログラムコードを実行することにより、上述の実施形態の機能が実現されるだけでなく、そのプログラムコードがコンピュータにおいて稼働しているＯＳ（オペレーティングシステム）あるいは他のアプリケーションソフト等と共同して上述の実施形態の機能が実現される場合にもかかるプログラムコードは本発明の実施形態に含まれることは言うまでもない。
【０１２２】
さらに、供給されたプログラムコードがコンピュータの機能拡張ボードやコンピュータに接続された機能拡張ユニットに備わるメモリに格納された後、そのプログラムコードの指示に基づいてその機能拡張ボードや機能拡張ユニットに備わるＣＰＵ等が実際の処理の一部または全部を行い、その処理によって上述した実施形態の機能が実現される場合にも本発明に含まれることは言うまでもない。
【図面の簡単な説明】
【０１２３】
【図１】本発明の第１の実施形態を示し、３次元空間内での形状計測データと、その形状計測データに対応する理論形状を記述するパラメータとの関係の一例を示す図である。
【図２】本発明の第１の実施形態を示し、パラメータを決定する方法の一例を説明するフローチャートである。
【図３】本発明の第１の実施形態を示し、３次元座標上にプロットされた擬似形状計測データを示す図である。
【図４】本発明の第１の実施形態を示し、入力した擬似形状計測データの値の一例を示す図である。
【図５】本発明の第１の実施形態を示し、情報処理装置の構成の一例を示したブロック図である。
【図６】本発明の第２の実施形態を示し、パラメータを決定する方法の一例を説明するフローチャートである。
【図７】本発明の第２の実施形態を示し、３次元座標上にプロットされた擬似形状計測データの一例を示す図である。
【図８】本発明の第２の実施形態を示し、入力した擬似形状計測データの値の一例をそれぞれ示す図である。
【図９】本発明の第３の実施形態を示し、パラメータを決定する方法の一例を説明するフローチャートである。
【符号の説明】
【０１２４】
１００ 理論形状
１０１ 中心軸
Ｃ 測定空間での基準点
Ｏ 理論形状の基準点
Ｐi 形状計測データ
ｘ0、ｙ0、ｚ0、α、β、Ｒ パラメータ
【特許請求の範囲】
【請求項１】
計測された３次元形状に対応する理論形状を規定する複数のパラメータを決定する形状パラメータ決定方法であって、
前記複数のパラメータを少なくとも２組に分け、前記計測された３次元形状の計測データと、その３次元形状の計測データに対応する理論形状との差を表す量が最小となるようなパラメータを、前記分けた組毎に順番に繰り返し求めることによって、前記複数のパラメータを決定することを特徴とする形状パラメータ決定方法。
【請求項２】
計測された３次元形状に対応する理論形状を規定する複数のパラメータを第１〜第ｎ（ｎは２以上の自然数）の組に分け、分けた組毎に前記複数のパラメータを決定する形状パラメータ決定方法であって、
前記計測された３次元形状の計測データを入力する入力ステップと、
前記入力ステップにより入力された計測データと、前記計測データに対応する理論形状との差を表す量が最小となるようなパラメータを、前記分けた組毎に求める第１〜第ｎのパラメータ取得ステップとを有することを特徴とする形状パラメータ決定方法。
【請求項３】
前記第１〜第ｎのパラメータ取得ステップは、他のパラメータ取得ステップにより求められたパラメータを既知量として、前記入力ステップにより入力された計測データと、前記計測データに対応する理論形状との差を表す量が最小となるようなパラメータを求めることを特徴とする請求項２に記載の形状パラメータ決定方法。
【請求項４】
前記第１〜第ｎのパラメータ取得ステップを、所定の条件が満たされるまで繰り返し行うことを特徴とする請求項３に記載の形状パラメータ決定方法。
【請求項５】
前記３次元形状に対応する理論形状が、回転対称体であることを特徴とする請求項１〜４の何れか１項に記載の形状パラメータ決定方法。
【請求項６】
前記回転対称体が、円筒であることを特徴とする請求項５に記載の形状パラメータ決定方法。
【請求項７】
前記複数のパラメータを、前記計測データと、前記計測データに対応する理論形状との差を表す量を２次形式で記述するパラメータの組と、そうでないパラメータの組とに分け、
前記計測データに対応する理論形状との差を表す量を２次形式で記述するパラメータを、線形の最小自乗法を用いて求め、
前記計測データに対応する理論形状との差を表す量を２次形式で記述しないパラメータを、非線形の最小自乗法を用いて求めることを特徴とする請求項１〜６の何れか１項に記載の形状パラメータ決定方法。
【請求項８】
前記３次元形状に対応する理論形状が、円筒であり、
前記円筒の中心軸の方向を規定するパラメータを、線形の最小自乗法を用いて求め、
前記中心軸の並進位置を規定するパラメータを、非線形の最小自乗法を用いて求めることを特徴とする請求項７に記載の形状パラメータ決定方法。
【請求項９】
前記３次元形状に対応する理論形状が、円筒であり、
前記複数のパラメータは、前記円筒の半径を含み、
前記円筒の中心軸の方向を規定するパラメータを、線形の最小自乗法を用いて求め、
前記中心軸の並進位置を規定するパラメータと、前記円筒の半径とを、非線形の最小自乗法を用いて求めることを特徴とする請求項７に記載の形状パラメータ決定方法。
【請求項１０】
前記複数のパラメータを、非線形の最小自乗法を用いて求めることを特徴とする請求項１〜６の何れか１項に記載の形状パラメータ決定方法。
【請求項１１】
前記３次元形状に対応する理論形状が、円筒であり、
前記複数のパラメータを、前記円筒の中心軸の方向を規定するパラメータの組と、前記中心軸の並進位置を規定するパラメータの組とに分けることを特徴とする請求項１０に記載の形状パラメータ決定方法。
【請求項１２】
前記３次元形状に対応する理論形状が、円筒であり、
前記複数のパラメータは、前記円筒の半径を含み、
前記複数のパラメータを、前記円筒の中心軸の方向を規定するパラメータの組と、前記中心軸の並進位置を規定するパラメータ及び前記円筒の半径の組とに分けることを特徴とする請求項１０に記載の形状パラメータ決定方法。
【請求項１３】
計測された３次元形状に対応する理論形状を規定する複数のパラメータを決定する形状パラメータ決定装置であって、
前記複数のパラメータを少なくとも２組に分け、前記計測された３次元形状の計測データと、その３次元形状の計測データに対応する理論形状との差を表す量が最小となるようなパラメータを、前記分けた組毎に順番に繰り返し求める手段を有することを特徴とする形状パラメータ決定装置。
【請求項１４】
計測された３次元形状に対応する理論形状を規定する複数のパラメータを第１〜第ｎ（ｎは２以上の自然数）の組に分け、分けた組毎に前記複数のパラメータを決定する形状パラメータ決定装置であって、
前記計測された３次元形状の計測データを入力する入力手段と、
前記入力手段により入力された計測データと、前記計測データに対応する理論形状との差を表す量が最小となるようなパラメータを、前記分けた組毎に求める第１〜第ｎのパラメータ取得手段とを有することを特徴とする形状パラメータ決定装置。
【請求項１５】
計測された３次元形状に対応する理論形状を規定する複数のパラメータを決定することをコンピュータに実行させるためのコンピュータプログラムであって、
前記複数のパラメータを少なくとも２組に分け、前記計測された３次元形状の計測データと、その３次元形状の計測データに対応する理論形状との差を表す量が最小となるようなパラメータを、前記分けた組毎に順番に繰り返し求めることをコンピュータに実行させることを特徴とするコンピュータプログラム。
【請求項１６】
計測された３次元形状に対応する理論形状を規定する複数のパラメータを第１〜第ｎ（ｎは２以上の自然数）の組に分け、分けた組毎に前記複数のパラメータを決定することをコンピュータに実行させるためのコンピュータプログラムであって、
前記計測された３次元形状の計測データを入力する入力ステップと、
前記入力ステップにより入力された計測データと、前記計測データに対応する理論形状との差を表す量が最小となるようなパラメータを、前記分けた組毎に求める第１〜第ｎのパラメータ取得ステップとをコンピュータに実行させることを特徴とするコンピュータプログラム。
【請求項１】
計測された３次元形状に対応する理論形状を規定する複数のパラメータを決定する形状パラメータ決定方法であって、
前記複数のパラメータを少なくとも２組に分け、前記計測された３次元形状の計測データと、その３次元形状の計測データに対応する理論形状との差を表す量が最小となるようなパラメータを、前記分けた組毎に順番に繰り返し求めることによって、前記複数のパラメータを決定することを特徴とする形状パラメータ決定方法。
【請求項２】
計測された３次元形状に対応する理論形状を規定する複数のパラメータを第１〜第ｎ（ｎは２以上の自然数）の組に分け、分けた組毎に前記複数のパラメータを決定する形状パラメータ決定方法であって、
前記計測された３次元形状の計測データを入力する入力ステップと、
前記入力ステップにより入力された計測データと、前記計測データに対応する理論形状との差を表す量が最小となるようなパラメータを、前記分けた組毎に求める第１〜第ｎのパラメータ取得ステップとを有することを特徴とする形状パラメータ決定方法。
【請求項３】
前記第１〜第ｎのパラメータ取得ステップは、他のパラメータ取得ステップにより求められたパラメータを既知量として、前記入力ステップにより入力された計測データと、前記計測データに対応する理論形状との差を表す量が最小となるようなパラメータを求めることを特徴とする請求項２に記載の形状パラメータ決定方法。
【請求項４】
前記第１〜第ｎのパラメータ取得ステップを、所定の条件が満たされるまで繰り返し行うことを特徴とする請求項３に記載の形状パラメータ決定方法。
【請求項５】
前記３次元形状に対応する理論形状が、回転対称体であることを特徴とする請求項１〜４の何れか１項に記載の形状パラメータ決定方法。
【請求項６】
前記回転対称体が、円筒であることを特徴とする請求項５に記載の形状パラメータ決定方法。
【請求項７】
前記複数のパラメータを、前記計測データと、前記計測データに対応する理論形状との差を表す量を２次形式で記述するパラメータの組と、そうでないパラメータの組とに分け、
前記計測データに対応する理論形状との差を表す量を２次形式で記述するパラメータを、線形の最小自乗法を用いて求め、
前記計測データに対応する理論形状との差を表す量を２次形式で記述しないパラメータを、非線形の最小自乗法を用いて求めることを特徴とする請求項１〜６の何れか１項に記載の形状パラメータ決定方法。
【請求項８】
前記３次元形状に対応する理論形状が、円筒であり、
前記円筒の中心軸の方向を規定するパラメータを、線形の最小自乗法を用いて求め、
前記中心軸の並進位置を規定するパラメータを、非線形の最小自乗法を用いて求めることを特徴とする請求項７に記載の形状パラメータ決定方法。
【請求項９】
前記３次元形状に対応する理論形状が、円筒であり、
前記複数のパラメータは、前記円筒の半径を含み、
前記円筒の中心軸の方向を規定するパラメータを、線形の最小自乗法を用いて求め、
前記中心軸の並進位置を規定するパラメータと、前記円筒の半径とを、非線形の最小自乗法を用いて求めることを特徴とする請求項７に記載の形状パラメータ決定方法。
【請求項１０】
前記複数のパラメータを、非線形の最小自乗法を用いて求めることを特徴とする請求項１〜６の何れか１項に記載の形状パラメータ決定方法。
【請求項１１】
前記３次元形状に対応する理論形状が、円筒であり、
前記複数のパラメータを、前記円筒の中心軸の方向を規定するパラメータの組と、前記中心軸の並進位置を規定するパラメータの組とに分けることを特徴とする請求項１０に記載の形状パラメータ決定方法。
【請求項１２】
前記３次元形状に対応する理論形状が、円筒であり、
前記複数のパラメータは、前記円筒の半径を含み、
前記複数のパラメータを、前記円筒の中心軸の方向を規定するパラメータの組と、前記中心軸の並進位置を規定するパラメータ及び前記円筒の半径の組とに分けることを特徴とする請求項１０に記載の形状パラメータ決定方法。
【請求項１３】
計測された３次元形状に対応する理論形状を規定する複数のパラメータを決定する形状パラメータ決定装置であって、
前記複数のパラメータを少なくとも２組に分け、前記計測された３次元形状の計測データと、その３次元形状の計測データに対応する理論形状との差を表す量が最小となるようなパラメータを、前記分けた組毎に順番に繰り返し求める手段を有することを特徴とする形状パラメータ決定装置。
【請求項１４】
計測された３次元形状に対応する理論形状を規定する複数のパラメータを第１〜第ｎ（ｎは２以上の自然数）の組に分け、分けた組毎に前記複数のパラメータを決定する形状パラメータ決定装置であって、
前記計測された３次元形状の計測データを入力する入力手段と、
前記入力手段により入力された計測データと、前記計測データに対応する理論形状との差を表す量が最小となるようなパラメータを、前記分けた組毎に求める第１〜第ｎのパラメータ取得手段とを有することを特徴とする形状パラメータ決定装置。
【請求項１５】
計測された３次元形状に対応する理論形状を規定する複数のパラメータを決定することをコンピュータに実行させるためのコンピュータプログラムであって、
前記複数のパラメータを少なくとも２組に分け、前記計測された３次元形状の計測データと、その３次元形状の計測データに対応する理論形状との差を表す量が最小となるようなパラメータを、前記分けた組毎に順番に繰り返し求めることをコンピュータに実行させることを特徴とするコンピュータプログラム。
【請求項１６】
計測された３次元形状に対応する理論形状を規定する複数のパラメータを第１〜第ｎ（ｎは２以上の自然数）の組に分け、分けた組毎に前記複数のパラメータを決定することをコンピュータに実行させるためのコンピュータプログラムであって、
前記計測された３次元形状の計測データを入力する入力ステップと、
前記入力ステップにより入力された計測データと、前記計測データに対応する理論形状との差を表す量が最小となるようなパラメータを、前記分けた組毎に求める第１〜第ｎのパラメータ取得ステップとをコンピュータに実行させることを特徴とするコンピュータプログラム。
【図１】
【図２】
【図３】
【図４】
【図５】
【図６】
【図７】
【図８】
【図９】
【図２】
【図３】
【図４】
【図５】
【図６】
【図７】
【図８】
【図９】
【公開番号】特開２００６−３４９４８０（Ｐ２００６−３４９４８０Ａ）
【公開日】平成１８年１２月２８日（２００６．１２．２８）
【国際特許分類】
【出願番号】特願２００５−１７５４８２（Ｐ２００５−１７５４８２）
【出願日】平成１７年６月１５日（２００５．６．１５）
【出願人】（０００００１００７）キヤノン株式会社 (59,756)
【Ｆターム（参考）】
【公開日】平成１８年１２月２８日（２００６．１２．２８）
【国際特許分類】
【出願日】平成１７年６月１５日（２００５．６．１５）
【出願人】（０００００１００７）キヤノン株式会社 (59,756)
【Ｆターム（参考）】
[ Back to top ]