27 #include <Base/Geometry/HomgLine2D.hh>
29 #include <MathAlgo/PolynomialSolve.hh>
30 #include <MathAlgo/SVD.hh>
31 #include <Base/Math/Matrix2x2.hh>
32 #include <Base/Geometry/HomgPoint2D.hh>
33 #include <Base/Image/Image.hh>
34 #include <Geometry/Quadric3D.hh>
35 #include <Geometry/PMatrix.hh>
36 #include <Geometry/FMatrixEstimation.hh>
38 #define SVD_ZERO_THRESH DBL_EPSILON * 16.0
45 return ((*
this)[0][0]*(*
this)[1][1]-(*
this)[0][1]*(*
this)[0][1]);
50 if (fabs(GetDeterminant())>DBL_EPSILON) {
55 if (GetSignature() != 0) {
65 if (fabs(GetDeterminant())>DBL_EPSILON) {
70 if (GetSignature() != 1) {
80 if (fabs(GetDeterminant())>DBL_EPSILON) {
85 if (GetSignature() != 2) {
95 if (fabs(GetDeterminant())<=DBL_EPSILON) {
101 if (GetSignature() != 1) {
108 if (GetDiscriminant()<-DBL_EPSILON) {
112 if ((*
this)[0][0]+(*
this)[1][1]>0.0) {
123 if (fabs(GetDeterminant())<=DBL_EPSILON) {
129 if (GetSignature() != 1) {
135 if (fabs(GetDiscriminant())<DBL_EPSILON) {
144 if (fabs(GetDeterminant())<=DBL_EPSILON) {
150 if (GetSignature() != 1) {
157 if (GetDiscriminant()>=-DBL_EPSILON) {
166 if (fabs(GetDeterminant())<=DBL_EPSILON) {
171 if (GetSignature() != 3) {
180 double norm = this->NormFrobenius();
181 if (fabs(norm)<1e-9)
return false;
183 return (fabs(GetDeterminant() / (norm*norm*norm))>1e-9);
187 const unsigned int Signature = GetSignature();
201 const double disc = GetDiscriminant();
212 BIASASSERT(IsPoint());
216 if (fabs(thepoint[2])>1e-10) thepoint.
Homogenize();
224 BIASERR(
"Better use ImageDraw::Ellipse for Conic drawing.");
228 const unsigned int uiHeight = img.
GetHeight();
229 const unsigned int uiWidth = img.
GetWidth();
234 register int LastSign = 0;
239 for (
register unsigned int y=0;y<uiHeight; y++) {
241 dist = LocatePoint(p);
242 if (dist>0.0) LastSign = 1;
243 else if (dist<0.0) LastSign = -1;
245 for (
register unsigned int x=0; x<uiWidth; x++) {
247 dist = LocatePoint(p);
250 ppI[y][x] = (ppI[y][x]>127)?0:255;
254 }
else if (dist>0.0) {
256 ppI[y][x] = (ppI[y][x]>127)?0:255;
262 ppI[y][x] = (ppI[y][x]>127)?0:255;
268 for (
register unsigned int x=0;x<uiWidth; x++) {
270 dist = LocatePoint(p);
271 if (dist>0.0) LastSign = 1;
272 else if (dist<0.0) LastSign = -1;
274 for (
register unsigned int y=0; y<uiHeight; y++) {
276 dist = LocatePoint(p);
280 ppI[y][x] = (ppI[y][x]>127)?0:255;
284 }
else if (dist>0.0) {
286 ppI[y][x] = (ppI[y][x]>127)?0:255;
292 ppI[y][x] = (ppI[y][x]>127)?0:255;
300 const double& radius_a,
const double& radius_b) {
303 (*this)[0][0] = 1.0 / (radius_a*radius_a);
304 (*this)[1][1] = 1.0 / (radius_b*radius_b);
311 Transform[0][2] = -Center[0];
312 Transform[1][2] = -Center[1];
313 Transform[2][2] = Center[2];
317 Rotation[1][1] = Rotation[0][0] = cos(dAngle);
318 Rotation[0][1] = -(Rotation[1][0] = sin(dAngle));
321 Transform = Rotation * Transform;
333 const double& dScale)
338 CONIC2D_TYPE cx=center[0]/center[2], cy=center[1]/center[2];
344 BIASASSERT(cov.
Invert(icov)==0);
346 const double iscale=1.0/(dScale*dScale);
347 (*this)[0][0] = icov[0][0] * iscale;
348 (*this)[0][1] = icov[0][1] * iscale;
349 (*this)[1][0] = icov[1][0] * iscale;
350 (*this)[1][1] = icov[1][1] * iscale;
360 (*this)[0][2]=(*this)[2][0]=-(cx*icov[0][0]+cy*icov[0][1]);
361 (*this)[1][2]=(*this)[2][1]=-(cx*icov[1][0]+cy*icov[1][1]);
363 (*this)[2][2]=-cx*(*this)[0][2]-cy*(*this)[1][2]-1.0;
371 const double dDisc = GetDiscriminant();
372 if (fabs(dDisc) <= DBL_EPSILON) {
377 const double a = (*this)[0][0];
378 const double b = 0.5*((*this)[1][0]+(*this)[0][1]);
379 const double c = (*this)[1][1];
380 const double d = 0.5*((*this)[2][0]+(*this)[0][2]);
381 const double e = 0.5*((*this)[2][1]+(*this)[1][2]);
384 center[0] = (b * e - c * d) / dDisc;
385 center[1] = (b * d - a * e) / dDisc;
389 cov[0][0]=(*this)[1][1]/dDisc;
390 cov[1][0]=-(*this)[1][0]/dDisc;
391 cov[0][1]=-(*this)[0][1]/dDisc;
392 cov[1][1]=(*this)[0][0]/dDisc;
404 if (!GetPointAndCovariance(c, cov)){
415 cout <<
"va "<<va[0]<<
", "<<va[1]<<
"\tla :"<<la
416 <<
"\nvb "<<vb[0]<<
", "<<vb[1]<<
"\tlb :"<<lb<<endl;
427 cout <<
"center "<<center[0]<<
", "<<center[1]<<
"\ta: "<<a[0]<<
", "<<a[1]
428 <<
"\tb: "<<b[0]<<
", "<<b[1]<<endl;
435 double& radius_b)
const {
449 const double& a = (*this)[0][0];
450 const double b = 0.5*((*this)[1][0]+(*this)[0][1]);
451 const double& c = (*this)[1][1];
452 const double d = 0.5*((*this)[2][0]+(*this)[0][2]);
453 const double e = 0.5*((*this)[2][1]+(*this)[1][2]);
455 const double dDisc = GetDiscriminant();
456 if (fabs(dDisc) <= DBL_EPSILON) {
462 Center[0] = (b * e - c * d) / dDisc;
463 Center[1] = (b * d - a * e) / dDisc;
466 const double q0 = LocatePoint(Center);
468 double theta, lambda1, lambda2;
469 const double p = (a + c) / 2.0;
470 const double q = (a - c) / 2.0;
471 const double r = sqrt(q * q + b * b);
472 if ((fabs(b) <= DBL_EPSILON) && (fabs(q) <=DBL_EPSILON))
475 theta = 0.5 * atan2(b, q);
483 if ((fabs(lambda1) <= DBL_EPSILON) ||
484 (fabs(lambda2) <= DBL_EPSILON) ||
485 (fabs(q0) <= DBL_EPSILON)) {
493 radius_a = 1.0 / sqrt(fabs(lambda1));
494 radius_b = 1.0 / sqrt(fabs(lambda2));
502 double& NormalizationFactor){
504 const double b = 0.5*((*this)[1][0]+(*this)[0][1]);
505 const double d = 0.5*((*this)[2][0]+(*this)[0][2]);
506 const double e = 0.5*((*this)[2][1]+(*this)[1][2]);
508 const double dDisc = GetDiscriminant();
509 if (fabs(dDisc) <= DBL_EPSILON) {
515 Center[0] = (b * e - (*this)[1][1] * d) / dDisc;
516 Center[1] = (b * d - (*this)[0][0] * e) / dDisc;
519 const double q0 = -LocatePoint(Center);
526 const double N = sqrt(GetDiscriminant()) / M_PI;
528 const double scale = log(ProbBorder) / (-1.0*q0*(log(N)+1.0));
530 NormalizationFactor = sqrt(GetDiscriminant()) / M_PI;
541 DualConic = P * DualQuadric * P.
Transpose();
544 (*this) = DualConic.GetDualConic(UseSVD);
561 for (
unsigned int i=0; i<3; i++) {
564 eigenimage = (*this) * eigenvector;
567 S[i] = eigenimage * eigenvector;
568 if (S[i]<-SVD_ZERO_THRESH) Signature--;
569 else if (S[i]>SVD_ZERO_THRESH) Signature++;
574 if (Signature<0)
return (
unsigned int)(-Signature);
575 return (
unsigned int)Signature;
583 vector<double> Coefficients;
585 const double r = Radius;
586 const double u = C[0]/C[2];
587 const double v = C[1]/C[2];
589 const double a = (*this)[0][0];
590 const double b = (*this)[1][1];
591 const double c = (*this)[0][1] + (*this)[1][0];
592 const double d = (*this)[0][2] + (*this)[2][0];
593 const double e = (*this)[1][2] + (*this)[2][1];
594 const double f = (*this)[2][2];
597 Coefficients.resize(5);
600 Coefficients[4] = -2*a*b + b*b + a*a + c*c;
602 Coefficients[3] = 2*u*c*a + 2*u*c*b - 4*a*a*v + 2*d*c + 2*e*b + 4*b*a*v
605 Coefficients[2] = 2*f*b + 2*u*c*e - 2*a*a*r*r + 2*a*a*u*u + v*v*c*c
606 - r*r*c*c + u*u*c*c + e*e + 2*u*d*a + d*d + 2*b*a*u*u
607 + 4*e*a*v - 2*a*f - 4*v*c*d - 4*u*c*a*v
608 + 6*a*a*v*v + 2*u*d*b - 2*b*a*v*v + 2*b*a*r*r;
610 Coefficients[1] = - 2*r*r*c*d - 4*a*a*u*u*v + 2*u*c*f - 2*e*a*v*v
611 + 2*e*a*r*r + 2*v*v*c*d + 2*f*e - 4*u*d*a*v + 4*f*a*v
612 - 2*v*d*d + 2*u*u*u*c*a + 2*u*c*a*v*v - 2*u*c*a*r*r
613 + 2*e*a*u*u + 2*u*d*e + 4*a*a*v*r*r - 4*a*a*v*v*v
616 Coefficients[0] = 2*u*d*f + 2*a*a*u*u*v*v - 2*a*a*u*u*r*r - 2*a*a*v*v*r*r
617 + 2*u*u*u*d*a + 2*f*a*u*u - 2*f*a*v*v + 2*f*a*r*r + f*f
618 + a*a*u*u*u*u + a*a*v*v*v*v + a*a*r*r*r*r + u*u*d*d
619 + v*v*d*d - r*r*d*d + 2*u*d*a*v*v - 2*u*d*a*r*r;
674 if (!IsProper()) UseSVD =
true;
682 SVD csvd(m, SVD_ZERO_THRESH,
false);
685 if (NullspaceDim==0)
return (csvd.
Invert());
689 for (
unsigned int i=NullspaceDim; i>0; i--) {
699 if (thenorm< SVD_ZERO_THRESH)
return Conic2D(SM);
706 const CONIC2D_TYPE& a = (*this)[0][0];
707 const CONIC2D_TYPE& b = (*this)[0][1];
708 const CONIC2D_TYPE& c = (*this)[0][2];
709 const CONIC2D_TYPE& e = (*this)[1][1];
710 const CONIC2D_TYPE& f = (*this)[1][2];
711 const CONIC2D_TYPE& g = (*this)[2][2];
714 Dual[0][0] = e*g - f*f;
715 Dual[0][1] = Dual[1][0] = c*f - b*g;
716 Dual[0][2] = Dual[2][0] = b*f - c*e;
717 Dual[1][1] = a*g - c*c;
718 Dual[1][2] = Dual[2][1] = b*c - a*f;
719 Dual[2][2] = a*e - b*b;
727 BIASCDOUT(D_CONIC2D_TAN,
"conic: "<<*
this<<endl);
728 BIASCDOUT(D_CONIC2D_TAN,
"x: "<<x<<endl);
739 BIASCDOUT(D_CONIC2D_TAN,
"polar: "<<polar<<endl);
745 Mult(polar_cross, dc);
746 dc = (polar_cross * dc);
748 BIASCDOUT(D_CONIC2D_TAN,
"dc: "<<dc<<endl);
768 BIASCDOUT(D_CONIC2D_TAN,
"normed dc: "<<dc<<endl);
775 double a=dc[0][0], b=dc[0][1], c=dc[1][1], d=dc[0][2], e=dc[1][2];
783 double y1a, y1b, y2a, y2b, y1, y2, min;
785 double tmp=e*e*0.25 - c*0.5, tmp2;
792 BIASERR(
"no solution for y2: tmp<0 "<<tmp);
795 BIASCDOUT(D_CONIC2D_TAN,
"y2a: "<<2.0*(e-y2a)*y2a-c
796 <<
"\ny2b: "<<2.0*(e-y2b)*y2b-c<<endl);
799 tmp=d*d*0.25 - a*0.5;
806 BIASERR(
"no solution for y1: tmp<0 "<<tmp);
809 BIASCDOUT(D_CONIC2D_TAN,
"y1a: "<<2*d*y1a-2*y1a*y1a-a
810 <<
"\ny1b: "<<2*d*y1b-2*y1b*y1b-a<<endl);
816 BIASCDOUT(D_CONIC2D_TAN,
"y1a, y2a ("<<y1a<<
", "<<y2a<<
")"<<endl);
818 min = fabs((d-y1)*y2+(e-y2)*y1 - b);
819 BIASCDOUT(D_CONIC2D_TAN,
"y1a, y2a = "<<min<<endl);
821 BIASCDOUT(D_CONIC2D_TAN,
"y1a, y2b ("<<y1a<<
", "<<y2b<<
")"<<endl);
822 tmp = fabs((d-y1a)*y2b+(e-y2b)*y1a - b);
827 BIASCDOUT(D_CONIC2D_TAN,
"y1a, y2b = "<<min<<endl);
829 BIASCDOUT(D_CONIC2D_TAN,
"not taking y1a, y2b = "<<tmp<<endl);
832 BIASCDOUT(D_CONIC2D_TAN,
"y1b, y2a ("<<y1b<<
", "<<y2a<<
")"<<endl);
833 tmp = fabs((d-y1b)*y2a+(e-y2a)*y1b - b);
838 BIASCDOUT(D_CONIC2D_TAN,
"y1b, y2a = "<<min<<endl);
840 BIASCDOUT(D_CONIC2D_TAN,
"not taking y1b, y2a = "<<tmp<<endl);
844 BIASCDOUT(D_CONIC2D_TAN,
"y1b, y2b ("<<y1b<<
", "<<y2b<<
")"<<endl);
845 tmp = fabs((d-y1b)*y2b+(e-y2b)*y1b - b);
850 BIASCDOUT(D_CONIC2D_TAN,
"y1b, y2b = "<<min<<endl);
852 BIASCDOUT(D_CONIC2D_TAN,
"not taking y1b, y2b = "<<tmp<<endl);
856 p2.
Set(d-y1, e-y2, 1.0);
858 BIASCDOUT(D_CONIC2D_TAN,
"tangent points: "<<p1<<
"\t"<<p2<<endl);
870 if (thetype==DoubleLine) {
877 for (
unsigned int i=0; i<3; i++) {
878 for (
unsigned int j=0; j<3; j++) {
880 (*
this)[(i+1)%3][(j+2)%3],
881 (*
this)[(i+2)%3][(j+1)%3],
882 (*
this)[(i+2)%3][(j+2)%3]);
883 cofactors[i][j] = submatrix.
det();
893 for (
unsigned int i=0; i<3; i++) {
894 if (fabs(Dual[i][i])>themax) {
895 themax = fabs(Dual[i][i]);
900 BIASERR(
"didnt find positive diagonal element");
903 if (Dual[index][index]<0.0) {
907 double beta = sqrt(themax);
911 Dual[index][2]/beta);
925 for (
unsigned int i=0; i<3; i++) {
926 for (
unsigned int j=0; j<3; j++) {
927 if (fabs(C[i][j])>1e-10) {
943 }
else if (thetype==SingleLine) {
945 if ((A[0][0]<=0.0) && (A[1][1]<=0.0) && (A[1][1]<=0.0)) {
948 for (
unsigned int i=0; i<3; i++) {
949 if (A[i][i]>=0.0) g[i] = sqrt(A[i][i]);
950 else if (A[i][i]>-1e-5) g[i] = 0;
953 BIASERR(
"line has complex coefficients, need sqrt from "<<A[i][i]);
970 bool signabequal = (Aabsolute[0][1]>0.0);
971 bool signbcequal = (Aabsolute[1][2]>0.0);
972 bool signacequal = (Aabsolute[0][2]>0.0);
979 BIASERR(
"inconsistent data: "<<Aabsolute[0][1]<<
" "
980 <<Aabsolute[1][2]<<
" "<<Aabsolute[0][2]);
984 BIASERR(
"inconsistent data: "<<Aabsolute[0][1]<<
" "
985 <<Aabsolute[1][2]<<
" "<<Aabsolute[0][2]);
993 BIASERR(
"inconsistent data: "<<Aabsolute[0][1]<<
" "
994 <<Aabsolute[1][2]<<
" "<<Aabsolute[0][2]);
1002 BIASERR(
"inconsistent data: "<<Aabsolute[0][1]<<
" "
1003 <<Aabsolute[1][2]<<
" "<<Aabsolute[0][2]);
1011 BIASERR(
"GetLines called for Conic which does not consist of one or two lines! Type is "<<thetype);
1020 std::vector<HomgPoint2D>& intersectionPoints)
const {
1021 intersectionPoints.clear();
1024 if (otherconic.
IsEmpty())
return 0;
1025 switch (GetConicType()) {
1029 if (otherconic.
GetPoint(solution)!=0)
return 1;
1030 intersectionPoints.push_back(solution);
1037 vector<HomgLine2D> thelines;
1040 vector<HomgPoint2D> thepoints;
1041 for (
unsigned int i=0; i<thelines.size(); i++) {
1043 for (
unsigned int j=0; j<thepoints.size(); j++) {
1044 thepoints[j].Normalize();
1045 bool isNewPoint =
true;
1046 for (
unsigned int k=0; k<intersectionPoints.size(); k++)
1047 if (
Equal(thepoints[j][0], intersectionPoints[k][0]) &&
1048 Equal(thepoints[j][1], intersectionPoints[k][1]) &&
1049 Equal(thepoints[j][2], intersectionPoints[k][2]))
1051 if (isNewPoint) intersectionPoints.push_back(thepoints[j]);
1060 vector<HomgLine2D> otherlines;
1063 for (
unsigned int i=0; i<thelines.size(); i++) {
1064 for (
unsigned int j=0; j<otherlines.size(); j++) {
1066 thelines[i].CrossProduct(otherlines[j]);
1067 const double norm = intersection.
NormL2();
1068 if (norm<1e-10)
continue;
1069 intersection /= norm;
1070 bool isNewPoint =
true;
1071 for (
unsigned int k=0; k<intersectionPoints.size(); k++)
1072 if (
Equal(intersection[0], intersectionPoints[k][0]) &&
1073 Equal(intersection[1], intersectionPoints[k][1]) &&
1074 Equal(intersection[2], intersectionPoints[k][2]))
1076 if (isNewPoint) intersectionPoints.push_back(intersection);
1083 if (otherconic.
GetPoint(point)!=0)
return -1;
1084 if (fabs(LocatePoint(point))>1e-9)
return 0;
1085 intersectionPoints.push_back(point);
1094 if (GetPoint(point)!=0)
return -1;
1095 if (fabs(otherconic.
LocatePoint(point))>1e-9)
return 0;
1096 intersectionPoints.push_back(point);
1099 case Empty:
return 0;
1102 }
else if (!otherconic.
IsProper()) {
1107 return IntersectConicProper(otherconic, intersectionPoints);
1113 std::vector<HomgPoint2D>& intersectionPoints)
const {
1115 BIASASSERT(IsProper());
1121 intersectionPoints.clear();
1138 for (
unsigned int i=0; i<9; i++) {
1139 f1[i] = (*this)[i/3][i%3];
1140 f2[i] = otherconic[i/3][i%3];
1142 vector<double> polynomial;
1145 std::vector<double > solutions;
1150 if (PS.
Solve(polynomial, solutions)<1) {
1151 BIASERR(
"PS failed.");
1154 if (solutions.empty()) {
1162 Conic2D degenerated(solutions[0]*(*
this) + (1.0-solutions[0])*otherconic);
1164 if (degeneratednorm<1e-10)
return 1;
1165 degenerated /= degeneratednorm;
1180 vector<HomgLine2D> lines;
1182 for (
unsigned int i=0; i<lines.size(); i++) {
1183 vector<HomgPoint2D> points;
1184 this->IntersectLine(lines[i],points);
1190 for (
unsigned int j=0 ;j<points.size(); j++) {
1192 points[j].Normalize();
1193 bool isNewPoint =
true;
1194 for (
unsigned int k=0; k<intersectionPoints.size(); k++)
1195 if (
Equal(points[j][0], intersectionPoints[k][0]) &&
1196 Equal(points[j][1], intersectionPoints[k][1]) &&
1197 Equal(points[j][2], intersectionPoints[k][2])) isNewPoint =
false;
1198 if (isNewPoint) intersectionPoints.push_back(points[j]);
1207 std::vector<HomgPoint2D>& intersectionPoints)
const {
1208 intersectionPoints.clear();
1209 if (theline.
NormL2()<1e-10)
return -1;
1212 myline /= myline.
NormL2();
1218 for (
unsigned int i=0; i<3; i++) {
1219 if (fabs(myline[i])>themax) {
1220 themax = fabs(myline[i]);
1227 BIASERR(
"Supplied invalid line 000");
1235 B[(theindex+1)%3][(theindex+2)%3],
1236 B[(theindex+2)%3][(theindex+1)%3],
1237 B[(theindex+2)%3][(theindex+2)%3]);
1239 double det = submatrix.det();
1246 double alpha = sqrt(det)/myline[theindex];
1252 bool foundit =
false;
1253 for (
unsigned int i=0; i<3; i++) {
1254 for (
unsigned int j=0; j<3; j++) {
1255 if (fabs(C[i][j])>1e-10) {
1268 if (!foundit)
return -2;
1273 intersectionPoints.push_back(p);
1275 intersectionPoints.push_back(q);
1282 theline /= theline.
NormL2();
1290 theline1 /= theline1.
NormL2();
1291 theline2 /= theline2.
NormL2();
1300 BIASERR(
"untested");
1304 const double norm = n.
NormL2();
1305 if (norm<1e-9)
return -1;
1308 pointcov.
Mult(n, distance);
1319 BIASERR(
"untested");
1321 if (GetLinearizedOffset(point, distance, pointcov)!=0){
1322 BIASERR(
"could not estimate distance");
1325 double l = LocatePoint(point);
1326 return (l*l / (4.0 *
1327 point.
ScalarProduct((*
this) * pointcov * (*
this) * point)));
void SetTwoLines(HomgLine2D theline1, HomgLine2D theline2)
create a degenerate point conic consisting of two lines
Vector< double > GetNullvector(const int last_offset=0)
return one of the nullvectors.
void GetRow(const unsigned int row, Vector3< T > &r) const
extract one row ('Zeile') from ths matrix (for convenience)
unsigned int GetSignature() const
computes signature of matrix, (no.
int Invert(Matrix2x2< T > &result) const
analyticaly inverts matrix
Conic2D GetDualConic(bool UseSVD=false) const
returns the dual of this conic, which is the inverse
class HomgPoint2D describes a point with 2 degrees of freedom in projective coordinates.
class representing the covariance matrix of a homogenous point 2D
T det() const
calculate the determinante
computes and holds the singular value decomposition of a rectangular (not necessarily quadratic) Matr...
void CheckCoefficients(std::vector< POLYNOMIALSOLVE_TYPE > &coeff, double eps=DBL_EPSILON)
removes leading zeros in coefficients
int IntersectConicProper(const Conic2D &otherconic, std::vector< HomgPoint2D > &intersectionPoints) const
compute intersection of 2 non-degenerate (full-rank) conics
bool GetEllipseParameters(HomgPoint2D &Center, double &dAngle, double &radius_a, double &radius_b) const
computes explicit ellipse representation
bool IsEllipse() const
check whether this conic is an ellipse
is a 'fixed size' quadratic matrix of dim.
Matrix< T > Transpose() const
transpose function, storing data destination matrix
bool IsDoubleLine() const
check whether this conic is a double line
void ScalarProduct(const Vector3< T > &argvec, T &result) const
scalar product (=inner product) of two vectors, storing the result in result
HomgPoint2D & Homogenize()
homogenize class data member elements to W==1 by divison by W
void OuterProduct(const Vector3< T > &v, Matrix3x3< T > &res) const
outer product, constructs a matrix.
bool IsEmpty() const
check whether this conic contains no points
double GetDiscriminant() const
void SetEllipse(const HomgPoint2D &Center, const double &dAngle, const double &radius_a, const double &radius_b)
construct an ellipse with explicit parameters
bool IsParabola() const
check whether this conic is a parabola
unsigned int GetWidth() const
base class for solving polynomial equations
int GetLines(std::vector< HomgLine2D > &lines) const
retrieve one or two lines from a degenerate point conic
Matrix< T > GetSkewSymmetricMatrix() const
constructs a skew symmetric 3x3 matrix from (*this), which can be used instead of the cross product ...
A quadric as a 4x4 matrix describing a surface in 3d projective space defined by a quadratic equation...
double GetPointDistance(const HomgPoint2D &point, HomgPoint2D &distance, const HomgPoint2DCov &pointcov=HomgPoint2DCov(MatrixIdentity)) const
return the squared mahalanobis distance of the point to the conic
const int RelNullspaceDim() const
compare singular values against greatest, not absolute
double NormFrobenius() const
void Mult(const Vector3< T > &argvec, Vector3< T > &destvec) const
matrix - vector multiplicate this matrix with Vector3, storing the result in destvec calculates: dest...
double LocatePoint(const HomgPoint2D &point2D) const
determines if the point is inside/on/outside the conic
int IntersectConic(const Conic2D &otherconic, std::vector< HomgPoint2D > &intersectionPoints) const
compute intersection of 2 conics
bool IsAtInfinity() const
bool IsHyperbola() const
check whether this conic is a hyperbola
bool HasRealSolution(const std::vector< POLYNOMIALSOLVE_TYPE > &coeff)
determine whether the polynomial has a solution in IR
int IntersectLine(const HomgLine2D &theline, std::vector< HomgPoint2D > &intersectionPoints) const
compute intersection of a conic with a line
Matrix< double > GetV() const
return V
const Matrix< double > & GetVT() const
return VT (=transposed(V))
functions for estimating a fundamental matrix (FMatrix) given a set of 2d-2d correspondences (no outl...
unsigned int GetHeight() const
const Vector< double > & GetS() const
return S which is a vector of the singular values of A in descending order.
a line l = (a b c)^T is a form of the implicit straight line equation 0 = a*x + b*y + c if homogenize...
void SetQuadricProjection(const Quadric3D &Q, const PMatrix &P, bool UseSVD=false)
constructs a conic that is a projection of the quadric Q using camera matrix P
bool IsProper() const
is this a proper conic (parabola/hyperbola/ellipse) with rank3 or some of the degenerate cases ...
bool GetPointAndCovariance(HomgPoint2D ¢er, Matrix2x2< CONIC2D_TYPE > &cov) const
reconstructs center and covariance matrix from conic
int GetTangentPoints(const HomgPoint2D &x, HomgPoint2D &p1, HomgPoint2D &p2) const
Returns the two points p1 and p2 on the conic subject to the constraint that the two tangents at the ...
void SetPointAndCovariance(const HomgPoint2D ¢er, const Matrix2x2< CONIC2D_TYPE > &cov, const double &dScale=GAUSS2D_CONFIDENCE_39_PERCENT)
constructs a conic from a point and a covariance matrix
void GetDetPolynomial(const BIAS::Vector< double > &f1, const BIAS::Vector< double > &f2, std::vector< POLYNOMIALSOLVE_TYPE > &v)
set up a polynomial in alpha given Det(alpha*F1+(1-alpha)*F2)=0, where the f vectors are interpreted ...
bool Equal(const T left, const T right, const T eps)
comparison function for floating point values See http://www.boost.org/libs/test/doc/components/test_...
int EigenvalueDecomposition(T &value1, Vector2< T > &vector1, T &value2, Vector2< T > &vector2) const
Eigenvalue decomposition.
bool IntersectsCircle(const HomgPoint2D &C, double Radius) const
returns true if the given Conic intersects with the circle of radius Radius around point C...
matrix class with arbitrary size, indexing is row major.
ConicType GetConicType() const
return the type of point conic of (*this)
Matrix3x3< T > Transpose() const
returns transposed matrix tested 12.06.2002
bool IsSingleLine() const
check whether this conic is a single line
int Solve(const std::vector< POLYNOMIALSOLVE_TYPE > &coeff, std::vector< POLYNOMIALSOLVE_TYPE > &sol)
solve polynomial of arbitrary order n coeff[n]x^n+...+coeff[2]x^2+coeff[1]x+coeff[0]=0 coeff[n]!=0...
const Matrix< double > & GetU() const
return U U is a m x m orthogonal matrix
Matrix< double > Invert()
returns pseudoinverse of A = U * S * V^T A^+ = V * S^+ * U^T
void Draw(Image< unsigned char > &img) const
draws conic into an image using a brute force method
bool EllipsoidSilhouetteToGauss(const double &ProbBorder, HomgPoint2D &Center, double &NormalizationFactor)
assign a probability of a gaussian distribution to a conic to measure probabilities in the image plan...
void SetSingleLine(HomgLine2D theline)
create a degenerate point conic consisting only of one line
describes a projective 3D -> 2D mapping in homogenous coordinates
void SetAsCrossProductMatrix(const Vector3< T > &vec)
Sets matrix from vector as cross product matrix of this vector.
Quadric3D GetDualQuadric(bool UseSVD=false) const
computes the dual quadric (points <-> planes)
int GetLinearizedOffset(const HomgPoint2D &point, HomgPoint2D &distance, const HomgPoint2DCov &pointcov=HomgPoint2DCov(MatrixIdentity)) const
compute approximate offset of a measured point to a known conic
void Set(const HOMGPOINT2D_TYPE &x, const HOMGPOINT2D_TYPE &y)
set elementwise with given 2 euclidean scalar values.
A 3x3 matrix representing a conic (cone/plane intersection)
Vector3< T > & Normalize()
normalize this vector to length 1
double NormL2() const
the L2 norm sqrt(a^2 + b^2 + c^2)
bool IsPoint() const
check whether this conic is a point
void SetIdentity()
set the elements of this matrix to the identity matrix (possibly overriding the inherited method) ...
int GetPoint(HomgPoint2D &thepoint) const
if the conic consists of a single point (degenerate case, rank=2) retrieve this point ...
const StorageType ** GetImageDataArray() const
overloaded GetImageDataArray() from ImageBase