OSDN Git Service

0ce8d52bcf40da93a3c2ddebd43bb0784ac2b721
[nyartoolkit-and/nyartoolkit-and.git] / src / jp / nyatla / nyartoolkit / core / transmat / NyARTransMat.java
1 /* \r
2  * PROJECT: NyARToolkit (Extension)\r
3  * --------------------------------------------------------------------------------\r
4  * This work is based on the original ARToolKit developed by\r
5  *   Hirokazu Kato\r
6  *   Mark Billinghurst\r
7  *   HITLab, University of Washington, Seattle\r
8  * http://www.hitl.washington.edu/artoolkit/\r
9  *\r
10  * The NyARToolkit is Java edition ARToolKit class library.\r
11  * Copyright (C)2008-2009 Ryo Iizuka\r
12  *\r
13  * This program is free software: you can redistribute it and/or modify\r
14  * it under the terms of the GNU General Public License as published by\r
15  * the Free Software Foundation, either version 3 of the License, or\r
16  * (at your option) any later version.\r
17  * \r
18  * This program is distributed in the hope that it will be useful,\r
19  * but WITHOUT ANY WARRANTY; without even the implied warranty of\r
20  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the\r
21  * GNU General Public License for more details.\r
22  *\r
23  * You should have received a copy of the GNU General Public License\r
24  * along with this program.  If not, see <http://www.gnu.org/licenses/>.\r
25  * \r
26  * For further information please contact.\r
27  *      http://nyatla.jp/nyatoolkit/\r
28  *      <airmail(at)ebony.plala.or.jp> or <nyatla(at)nyatla.jp>\r
29  * \r
30  */\r
31 package jp.nyatla.nyartoolkit.core.transmat;\r
32 \r
33 import jp.nyatla.nyartoolkit.NyARException;\r
34 import jp.nyatla.nyartoolkit.core.param.*;\r
35 import jp.nyatla.nyartoolkit.core.squaredetect.NyARSquare;\r
36 import jp.nyatla.nyartoolkit.core.transmat.solver.*;\r
37 import jp.nyatla.nyartoolkit.core.transmat.optimize.*;\r
38 import jp.nyatla.nyartoolkit.core.transmat.rotmatrix.*;\r
39 import jp.nyatla.nyartoolkit.core.types.*;\r
40 import jp.nyatla.nyartoolkit.core.types.matrix.*;\r
41 \r
42 /**\r
43  * This class calculates ARMatrix from square information and holds it. --\r
44  * 変換行列を計算して、結果を保持するクラス。\r
45  * \r
46  */\r
47 public class NyARTransMat implements INyARTransMat\r
48 {\r
49         private final static double FIT_DIFF_THRESHOLD_CONT_OPT = 1.0;\r
50         private final static double FIT_DIFF_THRESHOLD_CONT = 0.1;\r
51         private final static double FIT_DIFF_THRESHOLD_OPT = 1.0;\r
52         \r
53 \r
54         private final NyARDoublePoint2d _center=new NyARDoublePoint2d(0,0);\r
55         private final NyARTransOffset _offset=new NyARTransOffset();\r
56         private NyARPerspectiveProjectionMatrix _projection_mat_ref;\r
57         protected NyARRotMatrix _rotmatrix;\r
58         protected INyARTransportVectorSolver _transsolver;\r
59         protected NyARPartialDifferentiationOptimize _mat_optimize;\r
60 \r
61 \r
62         private NyARCameraDistortionFactor _ref_dist_factor;\r
63 \r
64         /**\r
65          * 派生クラスで自分でメンバオブジェクトを指定したい場合はこちらを使う。\r
66          *\r
67          */\r
68         protected NyARTransMat()\r
69         {\r
70                 //_calculator,_rotmatrix,_mat_optimizeをコンストラクタの終了後に\r
71                 //作成して割り当ててください。\r
72                 return;\r
73         }\r
74         public NyARTransMat(NyARParam i_param) throws NyARException\r
75         {\r
76                 final NyARCameraDistortionFactor dist=i_param.getDistortionFactor();\r
77                 final NyARPerspectiveProjectionMatrix pmat=i_param.getPerspectiveProjectionMatrix();\r
78                 this._transsolver=new NyARTransportVectorSolver(pmat,4);\r
79                 //互換性が重要な時は、NyARRotMatrix_ARToolKitを使うこと。\r
80                 //理屈はNyARRotMatrix_NyARToolKitもNyARRotMatrix_ARToolKitも同じだけど、少しだけ値がずれる。\r
81                 this._rotmatrix = new NyARRotMatrix(pmat);\r
82                 this._mat_optimize=new NyARPartialDifferentiationOptimize(pmat);\r
83                 this._ref_dist_factor=dist;\r
84                 this._projection_mat_ref=pmat;\r
85         }\r
86 \r
87         public void setCenter(double i_x, double i_y)\r
88         {\r
89                 this._center.x= i_x;\r
90                 this._center.y= i_y;\r
91         }\r
92 \r
93 \r
94 \r
95 \r
96         /**\r
97          * 頂点順序をi_directionに対応して並べ替えます。\r
98          * @param i_square\r
99          * @param i_direction\r
100          * @param o_sqvertex_ref\r
101          * @param o_liner_ref\r
102          */\r
103         private final void initVertexOrder(NyARSquare i_square, int i_direction, NyARDoublePoint2d[] o_sqvertex_ref, NyARLinear[] o_liner_ref)\r
104         {\r
105                 //頂点順序を考慮した矩形の頂点情報\r
106                 o_sqvertex_ref[0]= i_square.sqvertex[(4 - i_direction) % 4];\r
107                 o_sqvertex_ref[1]= i_square.sqvertex[(5 - i_direction) % 4];\r
108                 o_sqvertex_ref[2]= i_square.sqvertex[(6 - i_direction) % 4];\r
109                 o_sqvertex_ref[3]= i_square.sqvertex[(7 - i_direction) % 4];    \r
110                 o_liner_ref[0]=i_square.line[(4 - i_direction) % 4];\r
111                 o_liner_ref[1]=i_square.line[(5 - i_direction) % 4];\r
112                 o_liner_ref[2]=i_square.line[(6 - i_direction) % 4];\r
113                 o_liner_ref[3]=i_square.line[(7 - i_direction) % 4];\r
114                 return;\r
115         }\r
116 \r
117 \r
118         private final NyARDoublePoint2d[] __transMat_sqvertex_ref = new NyARDoublePoint2d[4];\r
119         private final NyARDoublePoint2d[] __transMat_vertex_2d = NyARDoublePoint2d.createArray(4);\r
120         private final NyARDoublePoint3d[] __transMat_vertex_3d = NyARDoublePoint3d.createArray(4);\r
121         private final NyARLinear[] __transMat_linear_ref=new NyARLinear[4];\r
122         private final NyARDoublePoint3d __transMat_trans=new NyARDoublePoint3d();\r
123         /**\r
124          * double arGetTransMat( ARMarkerInfo *marker_info,double center[2], double width, double conv[3][4] )\r
125          * \r
126          * @param i_square\r
127          * 計算対象のNyARSquareオブジェクト\r
128          * @param i_direction\r
129          * @param i_width\r
130          * @return\r
131          * @throws NyARException\r
132          */\r
133         public void transMat(final NyARSquare i_square, int i_direction, double i_width, NyARTransMatResult o_result_conv) throws NyARException\r
134         {\r
135                 final NyARDoublePoint2d[] sqvertex_ref = __transMat_sqvertex_ref;\r
136                 final NyARLinear[] linear_ref=__transMat_linear_ref;\r
137                 final NyARDoublePoint3d trans=this.__transMat_trans;\r
138                 \r
139                 //計算用に頂点情報を初期化(順番調整)\r
140                 initVertexOrder(i_square, i_direction, sqvertex_ref,linear_ref);\r
141                 \r
142                 //平行移動量計算機に、2D座標系をセット\r
143                 NyARDoublePoint2d[] vertex_2d=this.__transMat_vertex_2d;\r
144                 NyARDoublePoint3d[] vertex_3d=this.__transMat_vertex_3d;\r
145                 this._ref_dist_factor.ideal2ObservBatch(sqvertex_ref, vertex_2d,4);             \r
146                 this._transsolver.set2dVertex(vertex_2d,4);\r
147                 \r
148                 //基準矩形の3D座標系を作成\r
149                 this._offset.setSquare(i_width,this._center);\r
150 \r
151                 //回転行列を計算\r
152                 this._rotmatrix.initRotBySquare(linear_ref,sqvertex_ref);\r
153                 \r
154                 //回転後の3D座標系から、平行移動量を計算\r
155                 this._rotmatrix.getPoint3dBatch(this._offset.vertex,vertex_3d,4);\r
156                 this._transsolver.solveTransportVector(vertex_3d,trans);\r
157                 \r
158                 //計算結果の最適化(平行移動量と回転行列の最適化)\r
159                 o_result_conv.error=this.optimize(this._rotmatrix, trans, this._transsolver,this._offset.vertex, vertex_2d);\r
160                 \r
161                 // マトリクスの保存\r
162                 this.updateMatrixValue(this._rotmatrix, this._offset.point, trans,o_result_conv);\r
163                 return;\r
164         }\r
165 \r
166         /*\r
167          * (non-Javadoc)\r
168          * @see jp.nyatla.nyartoolkit.core.transmat.INyARTransMat#transMatContinue(jp.nyatla.nyartoolkit.core.NyARSquare, int, double, jp.nyatla.nyartoolkit.core.transmat.NyARTransMatResult)\r
169          */\r
170         public void transMatContinue(NyARSquare i_square, int i_direction, double i_width, NyARTransMatResult o_result_conv) throws NyARException\r
171         {\r
172                 final NyARDoublePoint2d[] sqvertex_ref = __transMat_sqvertex_ref;\r
173                 final NyARLinear[] linear_ref=__transMat_linear_ref;\r
174                 final NyARDoublePoint3d trans=this.__transMat_trans;\r
175 \r
176                 // io_result_convが初期値なら、transMatで計算する。\r
177                 if (!o_result_conv.has_value) {\r
178                         this.transMat(i_square, i_direction, i_width, o_result_conv);\r
179                         return;\r
180                 }\r
181 \r
182                 //計算用に頂点情報を初期化(順番調整)\r
183                 initVertexOrder(i_square, i_direction, sqvertex_ref,linear_ref);\r
184 \r
185                 \r
186                 //平行移動量計算機に、2D座標系をセット\r
187                 NyARDoublePoint2d[] vertex_2d=this.__transMat_vertex_2d;\r
188                 NyARDoublePoint3d[] vertex_3d=this.__transMat_vertex_3d;\r
189                 this._ref_dist_factor.ideal2ObservBatch(sqvertex_ref, vertex_2d,4);             \r
190                 this._transsolver.set2dVertex(vertex_2d,4);\r
191                 \r
192                 //基準矩形の3D座標系を作成\r
193                 this._offset.setSquare(i_width,this._center);\r
194 \r
195                 //回転行列を計算\r
196                 this._rotmatrix.initRotByPrevResult(o_result_conv);\r
197                 \r
198                 //回転後の3D座標系から、平行移動量を計算\r
199                 this._rotmatrix.getPoint3dBatch(this._offset.vertex,vertex_3d,4);\r
200                 this._transsolver.solveTransportVector(vertex_3d,trans);\r
201 \r
202                 //現在のエラーレートを計算しておく\r
203                 double min_err=errRate(this._rotmatrix,trans, this._offset.vertex, vertex_2d,4,vertex_3d);\r
204                 NyARDoubleMatrix33 rot=this.__rot;\r
205                 //エラーレートが前回のエラー値より閾値分大きかったらアゲイン\r
206                 if(min_err<o_result_conv.error+FIT_DIFF_THRESHOLD_CONT){\r
207                         rot.setValue(this._rotmatrix);\r
208                         //最適化してみる。\r
209                         for (int i = 0;i<5; i++) {\r
210                                 //変換行列の最適化\r
211                                 this._mat_optimize.modifyMatrix(rot, trans, this._offset.vertex, vertex_2d, 4);\r
212                                 double err=errRate(rot,trans,this._offset.vertex, vertex_2d,4,vertex_3d);\r
213                                 //System.out.println("E:"+err);\r
214                                 if(min_err-err<FIT_DIFF_THRESHOLD_CONT_OPT){\r
215                                         //System.out.println("BREAK");\r
216                                         break;\r
217                                 }\r
218                                 this._transsolver.solveTransportVector(vertex_3d, trans);\r
219                                 this._rotmatrix.setValue(rot);\r
220                                 min_err=err;\r
221                         }\r
222                         this.updateMatrixValue(this._rotmatrix, this._offset.point, trans,o_result_conv);\r
223                 }else{\r
224                         //回転行列を計算\r
225                         this._rotmatrix.initRotBySquare(linear_ref,sqvertex_ref);\r
226                         \r
227                         //回転後の3D座標系から、平行移動量を計算\r
228                         this._rotmatrix.getPoint3dBatch(this._offset.vertex,vertex_3d,4);\r
229                         this._transsolver.solveTransportVector(vertex_3d,trans);\r
230                         \r
231                         //計算結果の最適化(平行移動量と回転行列の最適化)\r
232                         min_err=this.optimize(this._rotmatrix, trans, this._transsolver,this._offset.vertex, vertex_2d);\r
233                         this.updateMatrixValue(this._rotmatrix, this._offset.point, trans,o_result_conv);\r
234                 }\r
235                 o_result_conv.error=min_err;\r
236                 return;\r
237         }\r
238         private NyARDoubleMatrix33 __rot=new NyARDoubleMatrix33();\r
239         private double optimize(NyARRotMatrix io_rotmat,NyARDoublePoint3d io_transvec,INyARTransportVectorSolver i_solver,NyARDoublePoint3d[] i_offset_3d,NyARDoublePoint2d[] i_2d_vertex) throws NyARException\r
240         {\r
241                 //System.out.println("START");\r
242                 NyARDoublePoint3d[] vertex_3d=this.__transMat_vertex_3d;\r
243                 //初期のエラー値を計算\r
244                 double min_err=errRate(io_rotmat, io_transvec, i_offset_3d, i_2d_vertex,4,vertex_3d);\r
245                 NyARDoubleMatrix33 rot=this.__rot;\r
246                 rot.setValue(io_rotmat);\r
247                 for (int i = 0;i<5; i++) {\r
248                         //変換行列の最適化\r
249                         this._mat_optimize.modifyMatrix(rot, io_transvec, i_offset_3d, i_2d_vertex, 4);\r
250                         double err=errRate(rot,io_transvec, i_offset_3d, i_2d_vertex,4,vertex_3d);\r
251                         //System.out.println("E:"+err);\r
252                         if(min_err-err<FIT_DIFF_THRESHOLD_OPT){\r
253                                 //System.out.println("BREAK");\r
254                                 break;\r
255                         }\r
256                         i_solver.solveTransportVector(vertex_3d, io_transvec);\r
257                         io_rotmat.setValue(rot);\r
258                         min_err=err;\r
259                 }\r
260                 //System.out.println("END");\r
261                 return min_err;\r
262         }\r
263         \r
264         //エラーレート計算機\r
265         public double errRate(NyARDoubleMatrix33 io_rot,NyARDoublePoint3d i_trans, NyARDoublePoint3d[] i_vertex3d, NyARDoublePoint2d[] i_vertex2d,int i_number_of_vertex,NyARDoublePoint3d[] o_rot_vertex) throws NyARException\r
266         {\r
267                 NyARPerspectiveProjectionMatrix cp = this._projection_mat_ref;\r
268                 final double cp00=cp.m00;\r
269                 final double cp01=cp.m01;\r
270                 final double cp02=cp.m02;\r
271                 final double cp11=cp.m11;\r
272                 final double cp12=cp.m12;\r
273 \r
274                 double err=0;\r
275                 for(int i=0;i<i_number_of_vertex;i++){\r
276                         double x3d,y3d,z3d;\r
277                         o_rot_vertex[i].x=x3d=io_rot.m00*i_vertex3d[i].x+io_rot.m01*i_vertex3d[i].y+io_rot.m02*i_vertex3d[i].z;\r
278                         o_rot_vertex[i].y=y3d=io_rot.m10*i_vertex3d[i].x+io_rot.m11*i_vertex3d[i].y+io_rot.m12*i_vertex3d[i].z;\r
279                         o_rot_vertex[i].z=z3d=io_rot.m20*i_vertex3d[i].x+io_rot.m21*i_vertex3d[i].y+io_rot.m22*i_vertex3d[i].z;\r
280                         x3d+=i_trans.x;\r
281                         y3d+=i_trans.y;\r
282                         z3d+=i_trans.z;\r
283                         \r
284                         //射影変換\r
285                         double x2d=x3d*cp00+y3d*cp01+z3d*cp02;\r
286                         double y2d=y3d*cp11+z3d*cp12;\r
287                         double h2d=z3d;\r
288                         \r
289                         //エラーレート計算\r
290                         double t1=i_vertex2d[i].x-x2d/h2d;\r
291                         double t2=i_vertex2d[i].y-y2d/h2d;\r
292                         err+=t1*t1+t2*t2;\r
293                         \r
294                 }\r
295                 return err/i_number_of_vertex;\r
296         }               \r
297         \r
298         \r
299         \r
300         /**\r
301          * パラメータで変換行列を更新します。\r
302          * \r
303          * @param i_rot\r
304          * @param i_off\r
305          * @param i_trans\r
306          */\r
307         public void updateMatrixValue(NyARRotMatrix i_rot, NyARDoublePoint3d i_off, NyARDoublePoint3d i_trans,NyARTransMatResult o_result)\r
308         {\r
309                 o_result.m00=i_rot.m00;\r
310                 o_result.m01=i_rot.m01;\r
311                 o_result.m02=i_rot.m02;\r
312                 o_result.m03=i_rot.m00 * i_off.x + i_rot.m01 * i_off.y + i_rot.m02 * i_off.z + i_trans.x;\r
313 \r
314                 o_result.m10 = i_rot.m10;\r
315                 o_result.m11 = i_rot.m11;\r
316                 o_result.m12 = i_rot.m12;\r
317                 o_result.m13 = i_rot.m10 * i_off.x + i_rot.m11 * i_off.y + i_rot.m12 * i_off.z + i_trans.y;\r
318 \r
319                 o_result.m20 = i_rot.m20;\r
320                 o_result.m21 = i_rot.m21;\r
321                 o_result.m22 = i_rot.m22;\r
322                 o_result.m23 = i_rot.m20 * i_off.x + i_rot.m21 * i_off.y + i_rot.m22 * i_off.z + i_trans.z;\r
323 \r
324                 o_result.has_value = true;\r
325                 return;\r
326         }       \r
327 }\r