画像処理ソリューション
これを見れば画像処理の入門から基礎~応用まで全てがわかるのを目指して!
   
翻訳(Translate)

プロフィール

Akira

ニックネーム:Akira
東京都の町田事業所に勤務
画像処理ソフトの開発を行っています。リンクフリーです!
詳細プロフィールは こちら
お問い合わせは、こちら↓

【補助HP】
画像処理ソリューションWeb版 【Newブログ】
イメージングソリューション

スポンサーリンク


カテゴリ

最近のコメント

カレンダー

09 | 2017/10 | 11
S M T W T F S
1 2 3 4 5 6 7
8 9 10 11 12 13 14
15 16 17 18 19 20 21
22 23 24 25 26 27 28
29 30 31 - - - -

趣味のブログ

iPhone萬歳!
iPhoneの情報いろいろ。
ブログ学習帳
ブログ、SEO、アフィリエイト情報など(まだまだこれから)
俺流クラフト日記
ハンドメイド作品の記録(現在、放置中)

スポンサーリンク 最近の記事
(09/18)  計測測定展に光切断のデモを出展しました
(08/17)  ディジタル画像技術事典200に記事が載りました
(06/09)  光切断を画像センシング展で公開
(05/14)  中国(上海)へ行って来ました
(04/12)  韓国へ行って来ました
(03/10)  私の求める新人像
(01/18)  エレクトロテストジャパンにカラー光切断法のデモを出展しました。
(12/23)  ユニークアクセス200万達成!
(12/10)  【カラー光切断法】YouTube動画まとめ
(11/04)  国際画像機器展2014にカラー光切断法を出展します。
(10/05)  第25回コンピュータビジョン勉強会@関東に参加してきました。
(09/08)  フーリエ変換の記事を追加しました。
(08/09)  【画像処理】ランキング低下中
(07/06)  記事の更新が停滞中...
(06/08)  画像センシング展2014でカラー光切断法のデモを行います。
(05/17)  カラー光切断法の動画を公開しました。
(04/30)  ソニーα NEX-5Rで星空撮影
(04/10)  カラー光切断法の取込結果を追加しました
(03/08)  Korea Vision Show 2014へ行ってきました
(02/05)  フーリエ変換シリーズを始めます。
(01/06)  2014年、あけましておめでとうございます。
(12/04)  カラー光切断法を公開(国際画像機器展2013にて)
(11/13)  国際画像機器展2013に出展します
(10/14)  「画像処理のためのC#」はじめます。
(09/16)  【C#,VB.NET】高速描画コントロールをバージョンアップしました。
(09/04)  拡大鏡に輝度値表示、ルーラー機能を追加した個人ツールを公開
(08/05)  7月の拍手Top5
(07/06)  2013年6月人気記事Top5
(05/12)  SONY α NEX-5Rレビュー
(04/24)  SONY α NEX-5RY購入

最小二乗法による楕円近似

楕円についても、「一般式による最小二乗法」で説明したように、最小二乗法を用いて、
近似することができるのですが、若干、計算が複雑なので、簡単に説明をしたいと思います。

この楕円近似の計算をエクセルで解いたファイルは
こちらからダウンロードできます。
楕円近似の図

上図のように
 楕円の中心を(X,Y
 X軸方向の長さをa
 Y軸方向の長さをb
 楕円の傾きをθ
としたときの一般式は

楕円の公式
となります。
(この一般式の作成方法はグラフ(関数)の拡大縮小、平行移動のページを参照下さい。)
この式を展開し、Xi、Yiに関して整理すると、

最小二乗法による楕円近似

となり、右辺の1を左辺に移項して未知数部分を変数で置き換えると
楕円の一般式1
としたくなりますが、これだと全ての項に変数がかかっているため、変数A~Fの定数倍の
組み合わせができてしまい、 A~Fの解が不定になってしまうので、式全体をA~Fの
どれかの値で割ると変数を1つ減らす事ができ、

楕円の一般式2
と置きます。

この式の2乗の和が最小になるよにA~Eの値を決めればいいので、
楕円の総和
の式をA~Eで偏微分し、行列で現すと

楕円の行列式1

となり、A~Eを求めるために両辺に逆行列を掛け合わせると

楕円の逆行列

となり、未知数A~Eを求める事ができます。

ここから先の計算が複雑で自信の無い部分なので、自己責任で参照願います。

A~Eの値を元に本来求めるべき楕円の中心、傾き、長軸、短軸の長さは

中心の座標(X,Y

楕円の中心

楕円の傾きθ
楕円の傾き

X軸方向の長さa、Y軸方向の長さb
楕円の長軸、短軸の長さ

となります。

 傾いた楕円の描画方法については別途記事にまとめました。参照下さい。

傾いた楕円の描画(Graphicsクラスを用いた方法)

傾いた楕円の描画(楕円上の点の座標を求める方法)


Loading...
スポンサーリンク

この記事に対するコメント

いつも大変参考にさせて頂いています。(特に画像処理)
楕円近似の二元二次式から半径a,bを求めるとき、楕円の角度θがジャスト45°の時は、正しい値が得られないようです。
角度は正しいので、0°に座標変換、二元二次の近似式を求め、半径a,bを計算すると正しい値が求められるようです。
異常値が出た時に、試行錯誤の結果ジャスト45°の時だけ、正しい値が計算されないことに気が付きました。
【2017/08/13 22:09】 URL | hamu #- [ 編集]

楕円にならない場合は・・・
X^2 + AXY + BY^2 + CX + DY + E = 0
は一般的な"楕円"というより"2次曲線"なので
係数の取り方によっては楕円にならない(双曲線や放物線になる)
場合がある気がするのですが
いかがでしょう。
【2017/07/06 00:00】 URL | cia_rana #- [ 編集]


ありがとうございます!
【2016/10/20 13:28】 URL | Kai #- [ 編集]

Re: タイトルなし
直接的に”楕円の”最小二乗法を見た事がなかったので、記事にしたぐらいなので、楕円近似に関する参考文献というのはありません。
一般的な最小二乗法であれば
http://imagingsolution.net/math/robust-tukeys-biweight/
のページの下の方に記載している本とはを参考にしています。

近似というか、最小二乗法は残差の2乗が最小になるように解くものなので、それを応用しています。

未知数を求めるという意味では、最近は疑似逆行列を用いることが多いでしょうか。。
http://imagingsolution.net/math/pseudoinversematrix/

以上、ご参考まで。
【2016/10/20 12:28】 URL | Akira #- [ 編集]


こんにちは。
お世話になっております。
こちらの楕円近似で何か参考にされた書籍や文献などがありますでしょうか。
教えていただけると幸いです。
【2016/10/20 11:48】 URL | Kai #- [ 編集]


ご返信ありがとうございます。
中心座標を(0,0)にする…思いつきませんでした!
さっそく取り掛かってみたいと思います!!
【2016/10/19 13:28】 URL | Kai #- [ 編集]

Re: タイトルなし
Kaiさん、こんにちは。
中心座標が既知の場合、X0、Y0に値を代入するというよりも、近似する座標(Xi,Yi)から既知の中心座標(X0,Y0)を引いて、中心座標が原点に来るようにしてから、ここの記事の中心座標(X0,Y0)を(0,0)とすると、X0,Y0に関する項が全て除外できるので、計算的にはかなり、すっきりできると思います。
【2016/10/19 12:22】 URL | Akira #- [ 編集]


こんにちは。
楕円の近似参考になります。
質問があるのですが,中心点が分かっており(自分で指定),中心点周りの測定点の近似を行いたい場合,X0とY0の値は既知であり,Eの値は算出することができるためA,B,C,Dの4つの未知数を求めるようにA,B,C,D,Eの行列を直せばよいのでしょうか。
よろしくお願いいたします。
【2016/10/19 09:47】 URL | Kai #- [ 編集]

ありがとうございます
ロバスト推定ですか・・・ハードルが高くなりますが
トライしてみたいと思います。
初期値を求めるにしても、複素数にならないものから
始めないといけないわけですね。
ありがとうございます。参考になりました。
【2013/01/04 20:46】 URL | ダークナイト #- [ 編集]

Re: 教えてください
ダークナイトさん、コメント頂きありがとうございます。
大きくズレたデータがある場合ですが、そのような場合には、個人的にはロバスト推定法に持ち込む場合が多いです。
(参考)
http://imagingsolution.blog107.fc2.com/blog-entry-32.html
その場合にしても、初期値が必要となるのですが、まずは、データの重心を計算して中心座標の初期値とするとか、円の近似で、中心と半径を初期値とするなど、何かしらの工夫が必要になって来ると思いますが、いきなり楕円で近似するのは、難しいかも?しれません。
【2013/01/04 00:50】 URL | Akira #- [ 編集]

教えてください
参考にさせていただいております。
大きく楕円からずれたような入力に対して計算すると、長径、短径の値が複素数となってしまいました。これはa, bを求める際のルートの中身がマイナスになるからですが、根本的な原因としては、楕円の式よりも、x^2-y^2=r^2系の式にフィットするからだとおもいます。
しかし、点をグラフにプロットしてみた場合には、楕円としてみた場合にもどこかに最も近い答えがあるように思います。
もしなにか工夫があるようでしたらご教授願えないでしょうか?
【2013/01/03 13:48】 URL | ダークナイト #- [ 編集]

有難うございました!
左辺に移行した右辺の"1”を忘れていた為に解けずにいましたが,教えて頂き全ての計算過程を確認することが出来ました.
【2012/04/11 23:51】 URL | ken #mQop/nM. [ 編集]

Re: ご教示ください
kenさん。コメント頂きありがとうございます。
このへんの計算は記事中にもわざわざ赤文字で書いてあるぐらい、複雑(大学ノートで8ページ使いました)で私自身もあまり自信のない部分なのですが...
本文中では
AX^2 + BXY + CY^2 + DX + EY + F = 0

X^2 + AXY + BY^2 + CX + DY + E = 0
と置く。と書いてありますが、この部分でA~Fを同じ文字を使ってA~Eへ置き換えてしまっている説明がちょっと都合が悪く、さらに楕円の一般式の右辺(=1)を本文中のFの中に入れちゃっているのも良くないので、
PX^2 + QXY + RY^2 + SX + TY + U - 1 = 0  ・・・①

X^2 + AXY + BY^2 + CX + DY + E = 0   ・・・②
の2式を使って解きました。
この①、②式と楕円の一般式をそれぞれ代入しながら計算すると解けました。
ただ、なにせ自信がありません。
【2012/04/11 01:36】 URL | Akira #- [ 編集]

ご教示ください
はじめまして.楕円近似参考にさせて頂いております.

順を追って楕円の一般式から計算過程を確認しているのですが,長軸a,短軸bだけ求める事が出来ずにいます.媒介変数A~Eの値を元に求めるかと思いますが,算出するのにどの変数を使用されたか教えて頂けないでしょうか?
【2012/04/11 00:26】 URL | ken #- [ 編集]


tamagoさん、はじめまして。
誤記についてご連絡頂きありがとうございました。
ご指摘頂いた部分を修正しましたので、ご確認願います。
この式の部分ですが、この記事を書いてしばらくしてから書き加えたものなので、その後の部分が
必ずしも間違えているとは限らないのですが、本文中にも赤字で書いてあるように、この先の
計算は結構複雑なので、あまり自信がありません。ただし近似した結果を見る限りでは、
合っていると思われます。
【2009/05/14 22:27】 URL | Akira #- [ 編集]


はじめまして。楕円近似参考にさせてもらっています。
1つ教えてください。2つ目の式で、展開してX、Yでまとめていますが、XYの内容が「-」「+」が違っていますか?その後の計算にも関係しますか?
【2009/05/14 14:08】 URL | tamago #1owpcyVw [ 編集]


sus4さん。こんばんは。
確かに引き継いだソフトをそのまま使うのは不安になるときもありますよね?
ただ、楕円上の座標の計算はここ↓
http://imagingsolution.blog107.fc2.com/blog-entry-87.html
でも紹介しているのですが、ここに載せいている計算式はたぶん合っていると思います。
(確認していたら、一部表記の間違いに気がついちゃいました...[これは修正済みです])
なので、楕円上の点を求めて、今、お使いの近似式が合っているかどうか確認されてみては
如何でしょうか?
【2008/12/11 03:46】 URL | Akira #- [ 編集]

確かに大変です。
私もプログラムを一から書くのであれば、
Xi^2の係数で 全体を割った式を使ったと思います。
実は、先輩(私は学生です)から引き継いだプログラムを使っておりまして、
その中では、Ax^2+Bxy+Cy^2+Dx+Ey-1=0
という式から長軸・短軸、中心、傾きを求めています。
その際の係数の使い方が合っているのかが気になったので、求め方をお尋ねしました。
ありがとうございました。
【2008/12/11 01:00】 URL | sus4 #- [ 編集]

記事を少し追記しました。
sus4さん、こんにちは。
少し話をしやすくするために、楕円の一般式を展開した式を追記しました。
だた、この記事にも書いてあるように、最後にX0,Y0,a,b,θを求める部分は
本当に自信がないので、教えられる程でもないのですが、sus4さんが書かれているように
Ax^2+Bxy+Cy^2+Dx+Ey-1=0
という置き方でも基本的には解けるんだと思います。たぶん...
(最初に書いてあるAx^2+Bxy+Cy^2+Dx+Ey+F=0を定数項のFで割っているのと
同じですから...)
しかし、今回追記した楕円の一般式を展開した式を見ても分かるように、定数項の部分は
ちょっと複雑な式になっているので、この定数項で式全体を割るよりも式全体に
sinθ^2、cosθ^2,1/a^2,1/b^2の項があるので、なんとなくXi^2の係数で
式全体を割って、
x^2+Axy+By^2+Cx+Dy+E=0
という形の式に置いて近似処理をしてみました。
この式で近似処理をしてA~Eの値を求めてX0,Y0,a,b,θを求める部分は、
一般式を展開した式を見ても、なんとなく大変そうな事が分かりますでしょ?!
【2008/12/06 16:01】 URL | Akira #- [ 編集]

遅くなってしまい申し訳ございません。
コメントへの返信が遅くなってしまい、申し訳ございません。

Ax^2+Bxy+Cy^2+Dx+Ey-1=0
という変数の置き方で、楕円の式を書いていたので、
その場合の楕円の各値(特に長軸,短軸)がどう表されるのか、計算しようとしていました。
「未知数5個に対して、式が5本あるので、X0,Y0,a,b,θを求めることができる。」
これはかなり複雑ですね。数式処理ソフトを使えば簡単に解けるかな・・
詳しく書いて頂き、ありがとうございました。
とりあえず計算してみます。
【2008/12/05 22:15】 URL | sus4 #- [ 編集]

まずは円の場合を参照願います。
このへんの解き方は円の最小二乗法の場合に関して、少し詳しく書いています。
http://imagingsolution.blog107.fc2.com/blog-entry-16.html
解き方の基本は円の場合とまったく同じです。
ただ、未知数が多いので、少し複雑になってしまいます。

そのページを参考にして今回、どうするか?というと、
A~Eまでの未知数をX0,Y0,a,b,sinθ,cosθであらわすと言っている部分は
円の場合の式の④~⑥式の部分に相当しています。
次にA~Eであらわされている式を2乗してΣをとります。(円の③式相当)
この式に関して、A~Eについて偏微分を行います。(円の⑦~⑧式相当)
これを行列であらわしたのが、このページにも書いてある行列の式になります。
すると、近似するデータ(Xi,Yi)を代入するとA~Eに関して解くことができます。
ここで、A~Eまでの値は求まったので、残す未知数は
X0,Y0,a,b,θ(sinθとcosθはθが求まると両方求まるので、1つとして数えます。)
の5個になります。
すると、未知数5個に対して、式が5本(A~Eまでの未知数をX0,Y0,a,b,sinθ,cosθであらわした式)成り立ちます。
(円の場合では未知数a、b、rの3個に対し、④~⑥の式が3本)
あとは、未知数5個に対して、式が5本あるので、連立方程式的にX0,Y0,a,b,θを求めることができるのですが、ここが少し解くのが難しかったです。
で、それを解いたつもりなのが、このページに書いてある、X0,Y0,a,b,θをA~Eであらわした式となっています。

って、書いていて聞かれている事が違う気がしてきた...
もし、単純に楕円の近似をプログラムで書きたいだけなら、

行列の最後の部分【(A B C D E) = になっている部分】のXi,Yiの部分に
近似するデータを代入します。(Σはただの合計です。)
あとは、逆行列の部分を解いて、行列の積を計算するとA~Eの値が求まります。
最後にX0,Y0,a,b,θをA~Eの値を使って計算れば、楕円の中心や傾きなどが求めることができます。

この処理をエクセルで解いたのが、このページでも紹介していますが、
http://imagingsolution.web.fc2.com/blog/LsmEllips.xls
のページになります。
これを見ると、プログラム的には参考になるのではないでしょうか?
【2008/11/22 10:13】 URL | Akira #- [ 編集]

ありがとうございます。
教えていただきありがとうございます。
楕円の一般式を展開し、A~Eまでの未知数をX0,Y0,a,b,sinθ,cosθで表す所まではわかるのですが、
その先の、X0,Y0,a,b,θに関して解くというのがわかりません。
実際にXi,Yiに値を入れて逆行列を計算すればよいのでしょうか。
【2008/11/22 01:20】 URL | sus4 #GKqFcwTE [ 編集]

もとい...
一般逆行列を使っても、X0,Y0,a,b,θを求める計算のところは同じでした...
(偏微分するところが、少し簡単にはなるんですけどね。)
【2008/11/20 02:25】 URL | Akira #- [ 編集]

ご覧頂きありがとうございます。
計算方法についてですが、この式を解いたときに、私の大学ノート3~4ページぐらいを使って解いた覚えがあるのですが、本当に簡単に説明すると、
一番上に書いてある楕円の一般式の2乗の部分などを展開して、Xi,Yiなどの項に関して =0となるように式を整理します。
次に、ここでの計算ではXi^2の項の係数で式全体を割るようにしています。
すると、この式に関して
Xi^2+AXiYi+Byi^2+CXi+DYi+E=0
と置いているので、A~Eまでの未知数をX0,Y0,a,b,Sinθ,Cosθであらわす事ができます。
さらに、A~Eまでの未知数は行列の式を解くと求めることができるので、X0,Y0,a,b,θに関して解くと、ここの記事に書いてある自信のない部分の式になるはず?です。
(もし違っていたら教えて下さい。)
この計算、やってみると分かって頂けると思うのですが、式が長くなって、意外と複雑です。
それよりも、いづれ記事にまとめようと思っているのですが、一般逆行列を使って解いた方がおそらく簡単なのですが、一般逆行列についてはご存知でしょうか?
実はその伏線を
行列で連立方程式を解く
http://imagingsolution.blog107.fc2.com/blog-entry-97.html
のページではっていたりもするのですが...
【2008/11/20 00:46】 URL | Akira #- [ 編集]

参考にさせて頂いております。
はじめまして。いつも利用させて頂いています。
最小二乗法によって求めた係数から、
楕円の中心x0、y0や傾きθ、長軸・短軸の長さを求めるプログラムを書いています。

「ここから先の計算が複雑で自信の無い部分なので、自己責任で参照願います。」以下の、
計算方法について、簡単にで構いませんので教えていただけないでしょうか。
【2008/11/19 21:13】 URL | sus4 #GKqFcwTE [ 編集]

そうですか...
やっぱり楕円体となると難しいですよね。
とりあえずはお役に立てているようなので、良かったです。
これからも、いろいろと書いていこうと思いますので、今後ともよろしくお願いいたします。
【2008/10/03 21:47】 URL | AKira #- [ 編集]

ありがとうございます
こんにちは。
今のところはHP掲載の2次元楕円近似を使用させていただいております。
楕円体は現状諦めています(笑)
2次元でさえ現状苦労しておりますので・・・
【2008/10/03 15:41】 URL | oni #- [ 編集]

なかなか難しそうですね。
oniさん。はじめまして。
楕円近似に関しては、昔かなり苦労したので、この記事を書いたのですが、楕円体だなんて本当に難しそうですね。楕円体の近似はやった事はありませんが、半径1の球の式
  X^2+Y^2+Z^2=1
を拡大縮小、回転、平行移動して一般式を作ればいいんでしょうかね?
ただ、三次元座標になると、求める未知数も多くなってしまうので、回転角など値が固定できるものがあるのであれば、固定した方がいいかもしれません。
もしくは、いづれ記事にしようと思っているのですが、擬似逆行列(一般化逆行列)を使って解けるかもしてません...
いづれにしても頑張って下さい。
【2008/10/02 23:03】 URL | AKira #- [ 編集]

使用させていただいております
はじめまして。
楕円体の近似を探していたところこちらのHPを見つけ利用させていただいております。
2次元の楕円近似でもなかなか見つからなかったため大変ありがたいです。
使用しているのは画像処理関係ではなく、楕円体形状基板の形状測定値から楕円近似を行っております。もちろん中心の直行する2ラインのみです。
現在難儀しているのが、計算結果の中心座標が測定座標より大きくなってしまうことです。(楕円の外に中心があるという結果)
解決するのにしばらく時間がかかりそうです。
【2008/10/02 19:21】 URL | oni #- [ 編集]

ありがとうございます。
T.S様、コメント頂きありがとうございます。
ただ、コメント欄に書いてある円上の点ですが、(sin ω、 cos ω)と書いてある点は
普通はωの角度の原点を右方向にとって、(cos ω、 sin ω)にすると思います。(何も気にせず書いていたので、間違えてしまったのですが、楕円を書くだけなら特に問題はないかと思います。)
この部分は最後のリンク(傾いた楕円の描画(楕円上の点の座標を求める方法))に、修正して書きましたので、そちらのページもよろしくお願い致します。
この楕円の書き方は、基本的に関数の拡大縮小や平行移動の考え方を使っているのですが、これも覚えておくと便利な場合がたまにあるので、地味にオススメです。
http://imagingsolution.blog107.fc2.com/blog-entry-46.html
【2008/09/08 00:24】 URL | Akira #- [ 編集]

同じく
私も以前,点の集合から楕円近似をする必要があったため,参考にさせていただきました.
ありがとうございました.
そして楕円近似した楕円を描く必要もあったため,A~Eが行列計算によって既知の値なので,本ページの上から3つ目の式に代入し,解の公式からy=の式を求め,楕円を描いていました.
しかし,管理人様のコメントのような美しい求め方もあるんだなとつくづく感心してしまいました.
【2008/09/07 22:08】 URL | T.S #- [ 編集]

しばし、お待ち下さい。
楕円の描画については傾きが無ければ普通に描画関数があるので簡単ですが、傾いている楕円の描画は私も少しはまりました。

別途コメントでもお世話になっていますので、いづれ詳細をまとめたいと思いますが、私は以下の方法で傾いた楕円を描画しています。

まず、中心が原点と一致した半径1の円を考えます。
すると円上の座標は
  (sinω、cosω)
となります。(ωは原点と円上の座標を結んだ線の傾き)
この座標をX軸方向にa倍、Y軸方向にb倍すると傾いていない楕円上の点が求まり
  (a×sinω、b×cosω)
となります。
この座標を原点周りにθ度回転し、
  X = a × sinω × cosθ - b × cosω × sinθ
  Y = a × sinω × sinθ + b × cosω × cosθ
さらにX0,Y0だけ平行移動すると
  X = a × sinω × cosθ - b × cosω × sinθ + X0
  Y = a × sinω × sinθ + b × cosω × cosθ + Y0
となり、求める傾いた楕円上の座標は求まります。
あとは例えばωを1度おきに計算して楕円を折れ線で書くと、それっぽい楕円になると思います。
お試し下さい。
【2008/08/07 00:38】 URL | Akira #- [ 編集]

とても参考になりました
画像処理の際に、楕円近似をしようと思っていたので、a,b,θ,X0,Y0の求め方がとても分かりやすく参考になりました。

できれば、求めたデータで理想?の楕円を描画したく、Xの値を代入して、Yの値を求めようと思ったのですが、楕円の公式がY=・・・にできないため、解けませんでした。
単純な数学の問題でしょうが、できれば、このあたりも載せていただけるとためになります。

ご検討ください。
【2008/08/06 22:26】 URL | かもめ #- [ 編集]


この記事に対するコメントの投稿














管理者にだけ表示を許可する


この記事に対するトラックバック
トラックバックURL
→http://imagingsolution.blog107.fc2.com/tb.php/20-b6081271
この記事にトラックバックする(FC2ブログユーザー)

現在の閲覧者数: / 合計