00001
00002 #include <iostream>
00003 #include <math.h>
00004 #include <vector>
00005 #include <list>
00006
00007 #include "SketchSpace.h"
00008 #include "Sketch.h"
00009 #include "Region.h"
00010 #include "visops.h"
00011
00012 #include "ShapeSpace.h"
00013 #include "ShapeRoot.h"
00014
00015 #include "PointData.h"
00016 #include "LineData.h"
00017 #include "ShapeLine.h"
00018
00019 #include "Crew/MapBuilder.h"
00020 #include "VRmixin.h"
00021
00022 #ifdef PLATFORM_APERIOS
00023
00024 template<typename T> static inline int signbit(T x) { return x<0; }
00025 #endif
00026
00027 using namespace std;
00028
00029 namespace DualCoding {
00030
00031 DATASTUFF_CC(LineData);
00032
00033 const float LineData::DEFAULT_MIN_LENGTH = 0.025f;
00034
00035
00036 const Point LineData::origin_pt = Point(0,0);
00037
00038 LineData::LineData(ShapeSpace& _space, const Point &p1, orientation_t orient)
00039 : BaseData(_space,getStaticType()), end1_pt(p1), end2_pt(),
00040 rho_norm(0), theta_norm(0), orientation(0), length(0) {
00041 int const width = space->getDualSpace().getWidth();
00042 int const height = space->getDualSpace().getHeight();
00043
00044
00045
00046
00047 float p2x=0, p2y=0;
00048 if ( fabs(orient-M_PI/2) < 0.001 ) {
00049 p2x = p1.coordX();
00050 p2y = p1.coordY() > height/2 ? 0 : height-1;
00051 } else {
00052 float slope = tan(orient);
00053 float intcpt = p1.coordY() - p1.coordX()*slope;
00054 p2x = p1.coordX() >= width/2 ? 0 : width-1;
00055 p2y = p2x * slope + intcpt;
00056 if ( p2y > height-1 ) {
00057 p2y = height-1;
00058 p2x = (p2y-intcpt) / slope;
00059 } else if ( p2y < 0 ) {
00060 p2y = 0;
00061 p2x = -intcpt/slope;
00062 }
00063 }
00064 end2_pt = Point(p2x,p2y);
00065 end1_pt.setValid(false);
00066 end1_pt.setActive(false);
00067 end2_pt.setValid(false);
00068 end2_pt.setActive(false);
00069 update_derived_properties();
00070 }
00071
00072 LineData::LineData(const LineData& other)
00073 : BaseData(other),
00074 end1_pt(other.end1_pt), end2_pt(other.end2_pt),
00075 rho_norm(other.rho_norm), theta_norm(other.theta_norm),
00076 orientation(other.orientation), length(other.length) {}
00077
00078 Point LineData::getCentroid() const { return (end1Pt()+end2Pt())*0.5f; }
00079
00080 void LineData::setInfinite(bool value) {
00081 end1_pt.setActive(!value);
00082 end2_pt.setActive(!value);
00083 }
00084
00085 #define NORMPOINT_MATCH_DISTSQ 400
00086 #define LINE_MATCH_OVERLAP -20
00087 #define LINE_ORIENTATION_MATCH M_PI/12
00088 #define LINE_PERPDIST 30
00089
00090
00091 bool LineData::isMatchFor(const ShapeRoot& other) const {
00092 if (!(isSameTypeAs(other) && isSameColorAs(other)))
00093 return false;
00094 else {
00095 const Shape<LineData>& other_ln = ShapeRootTypeConst(other,LineData);
00096 return isMatchFor(other_ln.getData());
00097 }
00098 }
00099
00100 bool LineData::isMatchFor(const LineData& other_line) const {
00101
00102
00103
00104
00105 AngPi theta_diff = float(theta_norm) - float(other_line.theta_norm);
00106 if ( (orientation_t)theta_diff > (orientation_t)M_PI/2 )
00107 theta_diff = (orientation_t)M_PI - theta_diff;
00108
00109
00110
00111
00112
00113 float const linePerpDist =
00114 space->getRefFrameType() == camcentric ? 20 : 100;
00115 AngPi const lineOrientationMatch =
00116 space->getRefFrameType() == camcentric ? M_PI/6 : M_PI/4;
00117
00118
00119
00120
00121
00122
00123 return
00124 theta_diff < lineOrientationMatch
00125 && (perpendicularDistanceFrom(other_line.getCentroid()) < linePerpDist || other_line.perpendicularDistanceFrom(getCentroid()) < linePerpDist)
00126 && isOverlappedWith(other_line,LINE_MATCH_OVERLAP);
00127 }
00128
00129 bool LineData::isOverlappedWith(const LineData& otherline, int amount) const {
00130 fmat::Column<3> thisVec = end2_pt.coords - end1_pt.coords;
00131 fmat::Column<2> v = fmat::SubMatrix<2,1>(thisVec);
00132 float theta = atan2(v[1],v[0]);
00133 fmat::Matrix<2,2> rot = fmat::rotation2D(-theta);
00134 v = rot * v;
00135 fmat::Column<3> otherPt1 = otherline.end1_pt.coords - end1_pt.coords;
00136 fmat::Column<3> otherPt2 = otherline.end2_pt.coords - end1_pt.coords;
00137 fmat::Column<2> op1 = rot * fmat::SubMatrix<2,1>(otherPt1);
00138 fmat::Column<2> op2 = rot * fmat::SubMatrix<2,1>(otherPt2);
00139 float thisMin = min(0.f, v[0]);
00140 float thisMax = max(0.f, v[0]);
00141 float otherMin = min(op1[0], op2[0]);
00142 float otherMax = max(op1[0], op2[0]);
00143 return (thisMin < otherMin) ? (otherMin+amount <= thisMax) : (thisMin+amount <= otherMax);
00144 }
00145
00146 void LineData::mergeWith(const ShapeRoot& other) {
00147 const Shape<LineData>& other_line = ShapeRootTypeConst(other,LineData);
00148 if (other_line->confidence <= 0)
00149 return;
00150
00151
00152
00153 float const totalLengthSq = length * length + other_line->length * other_line->length;
00154 EndPoint new_end1_pt = (end1_pt*length*length/totalLengthSq + other_line->end1Pt()*other_line->length*other_line->length/totalLengthSq);
00155 EndPoint new_end2_pt = (end2_pt*length*length/totalLengthSq + other_line->end2Pt()*other_line->length*other_line->length/totalLengthSq);
00156 new_end1_pt.valid = end1_pt.isValid() && other_line->end1Pt().isValid();
00157 new_end2_pt.valid = end2_pt.isValid() && other_line->end2Pt().isValid();
00158 end1_pt = new_end1_pt;
00159 end2_pt = new_end2_pt;
00160 update_derived_properties();
00161 }
00162
00163 bool LineData::isValidUpdate(coordinate_t c1_cur, coordinate_t c2_cur, coordinate_t c1_new, coordinate_t c2_new) {
00164 const float c1_noise = 10 + std::abs(c1_cur+c1_new) / 20.f;
00165 const float c2_noise = 10 + std::abs(c2_cur+c2_new) / 20.f;
00166 return (c1_new-c1_noise < c1_cur && c2_cur < c2_new+c2_noise);
00167 }
00168
00169
00170
00171 bool LineData::updateParams(const ShapeRoot& ground_root, bool force) {
00172 const Shape<LineData>& ground_line = ShapeRootTypeConst(ground_root,LineData);
00173 return updateParams(ground_line.getData(), force);
00174 }
00175
00176
00177 bool LineData::updateParams(const LineData &ground_line, bool force) {
00178
00179
00180
00181
00182
00183 const coordinate_t c1_cur = firstPtCoord();
00184 const coordinate_t c2_cur = secondPtCoord();
00185
00186 Point _end1_pt = firstPt();
00187 Point _end2_pt = secondPt();
00188
00189 updateLinePt(firstPt(), firstPtCoord(), firstPt(ground_line), firstPtCoord(ground_line), -1);
00190 updateLinePt(secondPt(), secondPtCoord(), secondPt(ground_line), secondPtCoord(ground_line), +1);
00191
00192
00193
00194 const coordinate_t c1_new = firstPtCoord();
00195 const coordinate_t c2_new = secondPtCoord();
00196
00197 if (isValidUpdate(c1_cur, c2_cur, c1_new, c2_new) || force){
00198
00199
00200 update_derived_properties();
00201 return true;
00202 }
00203
00204
00205 setEndPts(_end1_pt, _end2_pt);
00206 return false;
00207 }
00208
00209 void LineData::updateLinePt(EndPoint& localPt, coordinate_t local_coord,
00210 const EndPoint& groundPt, coordinate_t ground_coord,
00211 int sign) {
00212 if ( groundPt.isValid() ) {
00213 if ( localPt.isValid() )
00214 localPt.updateParams(groundPt);
00215 else
00216 localPt = groundPt;
00217 }
00218 else if ( (ground_coord - local_coord)*sign > 0 )
00219 localPt = groundPt;
00220 }
00221
00222 bool LineData::isAdmissible() const {
00223 if (end1Pt().isValid() && end2Pt().isValid())
00224 return length >= VRmixin::mapBuilder->curReq->minLineLength;
00225 else
00226 return length >= VRmixin::mapBuilder->curReq->minRayLength;
00227 }
00228
00229
00230 void LineData::printParams() const {
00231 cout << "Type = " << getTypeName() << " ID=" << getId() << " ParentID=" << getParentId() << endl;
00232 cout << " end1{" << end1Pt().coordX() << ", " << end1Pt().coordY() << "}"
00233 << " active=" << end1Pt().isActive()
00234 << " valid=" << end1Pt().isValid() << endl;
00235
00236 cout << " end2{" << end2Pt().coordX() << ", " << end2Pt().coordY() << "}"
00237 << " active=" << end2Pt().isActive()
00238 << " valid=" << end2Pt().isValid() << std::endl;
00239
00240 cout << " rho_norm=" << rho_norm
00241 << ", theta_norm=" << theta_norm
00242 << ", orientation=" << getOrientation()
00243 << ", length=" << getLength() << endl;
00244
00245 printf(" color = %d %d %d\n",getColor().red,getColor().green,getColor().blue);
00246
00247 cout << " mobile=" << getMobile()
00248 << ", viewable=" << isViewable() << endl;
00249
00250 vector<float> abc = lineEquation_abc();
00251 printf(" equ = %f %f %f\n",abc[0],abc[1],abc[2]);
00252 }
00253
00254 void LineData::printEnds() const {
00255 cout << " end1{" << end1Pt().coordX() << ", " << end1Pt().coordY() << "}";
00256 cout << " active=" << end1Pt().isActive() << ", valid=" << end1Pt().isValid() << endl;
00257 cout << " end2{" << end2Pt().coordX() << ", " << end2Pt().coordY() << "}";
00258 cout << " active=" << end2Pt().isActive() << ", valid=" << end2Pt().isValid() << endl;
00259
00260 }
00261
00262
00263
00264
00265
00266
00267 void LineData::applyTransform(const fmat::Transform& Tmat, const ReferenceFrameType_t newref) {
00268 end1Pt().applyTransform(Tmat,newref);
00269 end2Pt().applyTransform(Tmat,newref);
00270 update_derived_properties();
00271 }
00272
00273 void LineData::projectToGround(const fmat::Transform& camToBase, const PlaneEquation& groundplane) {
00274 end1Pt().projectToGround(camToBase,groundplane);
00275 end2Pt().projectToGround(camToBase,groundplane);
00276 update_derived_properties();
00277 }
00278
00279
00280
00281
00282 EndPoint& LineData::leftPt() { return end1Pt().isLeftOf(end2Pt()) ? end1_pt : end2_pt; }
00283 const EndPoint& LineData::leftPt() const { return end1Pt().isLeftOf(end2Pt()) ? end1_pt : end2_pt; }
00284 EndPoint& LineData::rightPt() { return end1Pt().isLeftOf(end2Pt()) ? end2_pt : end1_pt; }
00285 const EndPoint& LineData::rightPt() const { return end1Pt().isLeftOf(end2Pt()) ? end2_pt : end1_pt; }
00286 EndPoint& LineData::topPt() { return end1Pt().isAbove(end2Pt()) ? end1_pt : end2_pt; }
00287 const EndPoint& LineData::topPt() const { return end1Pt().isAbove(end2Pt()) ? end1_pt : end2_pt; }
00288 EndPoint& LineData::bottomPt() { return end1Pt().isAbove(end2Pt()) ? end2_pt : end1_pt; }
00289 const EndPoint& LineData::bottomPt() const { return end1Pt().isAbove(end2Pt()) ? end2_pt : end1_pt; }
00290
00291 Shape<PointData> LineData::leftPtShape() {
00292 Shape<PointData> result(new PointData(*space, leftPt()));
00293 result->setName("leftPt");
00294 result->inheritFrom(*this);
00295 result->setViewable(false);
00296 return result;
00297 }
00298
00299 Shape<PointData> LineData::rightPtShape() {
00300 Shape<PointData> result(new PointData(*space, rightPt()));
00301 result->setName("rightPt");
00302 result->inheritFrom(*this);
00303 result->setViewable(false);
00304 return result;
00305 }
00306
00307 Shape<PointData> LineData::topPtShape() {
00308 Shape<PointData> result(new PointData(*space, topPt()));
00309 result->setName("topPt");
00310 result->inheritFrom(*this);
00311 result->setViewable(false);
00312 return result;
00313 }
00314
00315 Shape<PointData> LineData::bottomPtShape() {
00316 Shape<PointData> result(new PointData(*space, bottomPt()));
00317 result->setName("bottomPt");
00318 result->inheritFrom(*this);
00319 result->setViewable(false);
00320 return result;
00321 }
00322
00323 EndPoint& LineData::firstPt() {
00324 if ( isNotVertical() )
00325 if ( end1Pt().coordX() < end2Pt().coordX() )
00326 return end1Pt();
00327 else return end2Pt();
00328 else
00329 if ( end1Pt().coordY() < end2Pt().coordY() )
00330 return end1Pt();
00331 else return end2Pt();
00332 }
00333
00334 const EndPoint& LineData::firstPt() const {
00335 if ( isNotVertical() )
00336 if ( end1Pt().coordX() < end2Pt().coordX() )
00337 return end1Pt();
00338 else return end2Pt();
00339 else
00340 if ( end1Pt().coordY() < end2Pt().coordY() )
00341 return end1Pt();
00342 else return end2Pt();
00343 }
00344
00345 const EndPoint& LineData::firstPt(const LineData &otherline) const {
00346 if ( isNotVertical() )
00347 if ( otherline.end1Pt().coordX() < otherline.end2Pt().coordX() )
00348 return otherline.end1Pt();
00349 else return otherline.end2Pt();
00350 else
00351 if ( otherline.end1Pt().coordY() < otherline.end2Pt().coordY() )
00352 return otherline.end1Pt();
00353 else return otherline.end2Pt();
00354 }
00355
00356 EndPoint& LineData::secondPt() {
00357 if ( isNotVertical() )
00358 if ( end1Pt().coordX() > end2Pt().coordX() )
00359 return end1Pt();
00360 else return end2Pt();
00361 else
00362 if ( end1Pt().coordY() > end2Pt().coordY() )
00363 return end1Pt();
00364 else return end2Pt();
00365 }
00366
00367 const EndPoint& LineData::secondPt() const {
00368 if ( isNotVertical() )
00369 if ( end1Pt().coordX() > end2Pt().coordX() )
00370 return end1Pt();
00371 else return end2Pt();
00372 else
00373 if ( end1Pt().coordY() > end2Pt().coordY() )
00374 return end1Pt();
00375 else return end2Pt();
00376 }
00377
00378 const EndPoint& LineData::secondPt(const LineData &otherline) const {
00379 if ( isNotVertical() )
00380 if ( otherline.end1Pt().coordX() > otherline.end2Pt().coordX() )
00381 return otherline.end1Pt();
00382 else return otherline.end2Pt();
00383 else
00384 if ( otherline.end1Pt().coordY() > otherline.end2Pt().coordY() )
00385 return otherline.end1Pt();
00386 else return otherline.end2Pt();
00387 }
00388
00389 Shape<PointData> LineData::firstPtShape() {
00390 Shape<PointData> result(new PointData(*space, firstPt()));
00391 result->setName("firstPt");
00392 result->inheritFrom(*this);
00393 result->setViewable(false);
00394 return result;
00395 }
00396
00397 Shape<PointData> LineData::secondPtShape() {
00398 Shape<PointData> result(new PointData(*space, secondPt()));
00399 result->setName("secondPt");
00400 result->inheritFrom(*this);
00401 result->setViewable(false);
00402 return result;
00403 }
00404
00405 coordinate_t LineData::firstPtCoord() const {
00406 return isNotVertical() ? firstPt().coordX() : firstPt().coordY();
00407 }
00408
00409 coordinate_t LineData::firstPtCoord(const LineData &otherline) const {
00410 return isNotVertical() ?
00411 firstPt(otherline).coordX() :
00412 firstPt(otherline).coordY();
00413 }
00414
00415 coordinate_t LineData::secondPtCoord() const {
00416 return isNotVertical() ? secondPt().coordX() : secondPt().coordY();
00417 }
00418
00419 coordinate_t LineData::secondPtCoord(const LineData &otherline) const {
00420 return isNotVertical() ?
00421 secondPt(otherline).coordX() :
00422 secondPt(otherline).coordY();
00423 }
00424
00425
00426
00427
00428 void LineData::setEndPts(const EndPoint& _end1_pt, const EndPoint& _end2_pt) {
00429 end1_pt.setCoords(_end1_pt);
00430 end1_pt.setActive(_end1_pt.isActive());
00431 end1_pt.setValid(_end1_pt.isValid());
00432 end1_pt.setNumUpdates(_end1_pt.numUpdates());
00433
00434 end2_pt.setCoords(_end2_pt);
00435 end2_pt.setActive(_end2_pt.isActive());
00436 end2_pt.setValid(_end2_pt.isValid());
00437 end2_pt.setNumUpdates(_end2_pt.numUpdates());
00438
00439 update_derived_properties();
00440 }
00441
00442
00443
00444
00445 std::pair<float,float> LineData::lineEquation_mb() const {
00446 float m;
00447 if ((fabs(end2Pt().coordX() - end1Pt().coordX()) * BIG_SLOPE)
00448 <= fabs(end2Pt().coordY() - end2Pt().coordY()))
00449 m = BIG_SLOPE;
00450 else
00451 m = (end2Pt().coordY() - end1Pt().coordY())/(end2Pt().coordX() - end1Pt().coordX());
00452 float b = end1Pt().coordY() - m * end1Pt().coordX();
00453 return pair<float,float>(m,b);
00454 }
00455
00456
00457
00458 std::vector<float> LineData::lineEquation_abc_xz() const {
00459 float dx = end2Pt().coordX() - end1Pt().coordX();
00460 float dz = end2Pt().coordZ() - end1Pt().coordZ();
00461
00462 std::vector<float> abc(3,1);
00463 float& a = abc[0];
00464 float& b = abc[1];
00465 float& c = abc[2];
00466
00467
00468 if((dx == 0)
00469 || (dz/dx > BIG_SLOPE)) {
00470 a = 1;
00471 b = 0;
00472 c = end1Pt().coordX();
00473 }
00474
00475
00476 else if((dz == 0)
00477 || (dx/dz > BIG_SLOPE)) {
00478 a = 0;
00479 b = 1;
00480 c = end1Pt().coordZ();
00481 }
00482
00483
00484 else {
00485 a = 1;
00486 b = (end1Pt().coordX() - end2Pt().coordX())
00487 / (end2Pt().coordZ() - end1Pt().coordZ());
00488 c = end1Pt().coordX() + b*end1Pt().coordZ();
00489 }
00490
00491 return(abc);
00492
00493 }
00494
00495
00496 std::vector<float> LineData::lineEquation_abc() const {
00497 float dx = end2Pt().coordX() - end1Pt().coordX();
00498 float dy = end2Pt().coordY() - end1Pt().coordY();
00499
00500 std::vector<float> abc(3,1);
00501 float& a = abc[0];
00502 float& b = abc[1];
00503 float& c = abc[2];
00504
00505
00506 if( std::abs(dx) < 1.0e-6f || Slope(dy/dx) > BIG_SLOPE) {
00507 a = 1;
00508 b = 0;
00509 c = end1Pt().coordX();
00510 }
00511
00512
00513 else if ( std::abs(dy) < 1.0e-6f || Slope(dx/dy) > BIG_SLOPE ) {
00514 a = 0;
00515 b = 1;
00516 c = end1Pt().coordY();
00517 }
00518
00519
00520 else {
00521 a = 1;
00522
00523
00524 b = -dx / dy;
00525 c = end1Pt().coordX() + b*end1Pt().coordY();
00526 }
00527
00528 return(abc);
00529 }
00530
00531
00532
00533
00534 void LineData::update_derived_properties() {
00535 rho_norm = perpendicularDistanceFrom(origin_pt);
00536 const vector<float> abc = lineEquation_abc();
00537 const float& a1 = abc[0];
00538 const float& b1 = abc[1];
00539 const float& c1 = abc[2];
00540 const float c1sign = (c1 >= 0) ? 1 : -1;
00541 theta_norm = atan2(b1*c1sign, a1*c1sign);
00542 orientation = theta_norm + AngPi((orientation_t)M_PI/2);
00543 length = end1Pt().distanceFrom(end2Pt());
00544 const ReferenceFrameType_t ref = getRefFrameType();
00545 end1_pt.setRefFrameType(ref);
00546 end2_pt.setRefFrameType(ref);
00547 deleteRendering();
00548 }
00549
00550 bool LineData::isNotVertical() const {
00551 const AngPi threshold = (orientation_t)M_PI / 3;
00552 const AngPi orient = getOrientation();
00553 return (orient <= threshold) || (orient >= (orientation_t)M_PI - threshold);
00554 }
00555
00556
00557
00558
00559
00560
00561
00562
00563
00564
00565
00566
00567
00568
00569
00570
00571
00572
00573
00574
00575
00576
00577
00578
00579
00580
00581
00582
00583
00584
00585 bool LineData::pointIsLeftOf(const Point& pt) const {
00586 const EndPoint &p1 =
00587 ( getRefFrameType() == camcentric )
00588 ? ((end1_pt.coordY() < end2_pt.coordY()) ? end1_pt : end2_pt)
00589 : ((end1_pt.coordX() < end2_pt.coordX()) ? end1_pt : end2_pt);
00590 const EndPoint &p2 = (&p1 == &end1_pt) ? end2_pt : end1_pt;
00591 if ( (getRefFrameType()==camcentric) ? (p1.coordY() == p2.coordY()) : (p1.coordX() == p2.coordX()) )
00592 return false;
00593 else
00594 return ( (p2.coordX() - p1.coordX()) * (pt.coordY() - p1.coordY()) -
00595 (p2.coordY() - p1.coordY()) * (pt.coordX() - p1.coordX()) ) > 0;
00596 }
00597
00598 bool LineData::pointIsRightOf(const Point& pt) const {
00599 const EndPoint &p1 =
00600 ( getRefFrameType() == camcentric )
00601 ? ((end1_pt.coordY() < end2_pt.coordY()) ? end1_pt : end2_pt)
00602 : ((end1_pt.coordX() < end2_pt.coordX()) ? end1_pt : end2_pt);
00603 const EndPoint &p2 = (&p1 == &end1_pt) ? end2_pt : end1_pt;
00604 if ( (getRefFrameType()==camcentric) ? (p1.coordY() == p2.coordY()) : (p1.coordX() == p2.coordX()) )
00605 return false;
00606 else
00607 return ( (p2.coordX() - p1.coordX()) * (pt.coordY() - p1.coordY()) -
00608 (p2.coordY() - p1.coordY()) * (pt.coordX() - p1.coordX()) ) < 0;
00609 }
00610
00611 bool LineData::pointIsAbove(const Point& pt) const {
00612 const EndPoint &p1 =
00613 ( getRefFrameType() == camcentric )
00614 ? ((end1_pt.coordX() < end2_pt.coordX()) ? end1_pt : end2_pt)
00615 : ((end1_pt.coordY() < end2_pt.coordY()) ? end1_pt : end2_pt);
00616 const EndPoint &p2 = (&p1 == &end1_pt) ? end2_pt : end1_pt;
00617 if ( (getRefFrameType()==camcentric) ? (p1.coordX() == p2.coordX()) : (p1.coordY() == p2.coordY()) )
00618 return false;
00619 else
00620 return ( (p2.coordX() - p1.coordX()) * (pt.coordY() - p1.coordY()) -
00621 (p2.coordY() - p1.coordY()) * (pt.coordX() - p1.coordX()) ) > 0;
00622 }
00623
00624 bool LineData::pointIsBelow(const Point& pt) const {
00625 const EndPoint &p1 =
00626 ( getRefFrameType() == camcentric )
00627 ? ((end1_pt.coordX() < end2_pt.coordX()) ? end1_pt : end2_pt)
00628 : ((end1_pt.coordY() < end2_pt.coordY()) ? end1_pt : end2_pt);
00629 const EndPoint &p2 = (&p1 == &end1_pt) ? end2_pt : end1_pt;
00630 if ( (getRefFrameType()==camcentric) ? (p1.coordX() == p2.coordX()) : (p1.coordY() == p2.coordY()) )
00631 return false;
00632 else
00633 return ( (p2.coordX() - p1.coordX()) * (pt.coordY() - p1.coordY()) -
00634 (p2.coordY() - p1.coordY()) * (pt.coordX() - p1.coordX()) ) < 0;
00635 }
00636
00637
00638
00639 bool LineData::isLongerThan(const Shape<LineData>& other) const {
00640 return length > other->length; }
00641
00642 bool LineData::isLongerThan(float ref_length) const {
00643 return length > ref_length; }
00644
00645 bool LineData::isShorterThan(const Shape<LineData>& other) const {
00646 return length < other->length; }
00647
00648 bool LineData::isShorterThan(float ref_length) const {
00649 return length < ref_length; }
00650
00651 bool LineData::isBetween(const Point &p, const LineData &other) const {
00652 if (getOrientation() == other.getOrientation()) {
00653 float dl = perpendicularDistanceFrom(other.end1Pt());
00654 return (perpendicularDistanceFrom(p) <= dl && other.perpendicularDistanceFrom(p) <= dl);
00655 }
00656 else {
00657 bool b;
00658 const LineData p_line (*space, p,
00659 intersectionWithLine(other, b, b));
00660 const AngPi theta_pline = p_line.getOrientation();
00661 const AngPi theta_greater =
00662 (getOrientation() > other.getOrientation()) ? getOrientation() : other.getOrientation();
00663 const AngPi theta_smaller =
00664 (getOrientation() < other.getOrientation()) ? getOrientation() : other.getOrientation();
00665 if (theta_greater - theta_smaller > M_PI/2)
00666 return (theta_pline >= theta_greater || theta_pline <= theta_smaller);
00667 else
00668 return (theta_pline <= theta_greater && theta_pline >= theta_smaller);
00669 }
00670 }
00671
00672
00673
00674
00675 bool
00676 LineData::intersectsLine(const Shape<LineData>& other) const {
00677 return intersectsLine(other.getData());
00678 }
00679
00680 bool
00681 LineData::intersectsLine(const LineData& other) const {
00682
00683 pair<float,float> F = lineEquation_mb();
00684
00685
00686 pair<float,float> G = other.lineEquation_mb();
00687
00688
00689
00690
00691
00692
00693
00694
00695
00696
00697
00698
00699
00700
00701
00702
00703
00704
00705
00706 float r3 = F.first * other.end1Pt().coordX() + F.second - other.end1Pt().coordY();
00707
00708
00709 float r4 = F.first * other.end2Pt().coordX() + F.second - other.end2Pt().coordY();
00710
00711
00712
00713
00714
00715
00716 if((r3 != 0)
00717 && (r4 != 0)
00718 && (signbit(r3) == signbit(r4))
00719 )
00720 return false;
00721
00722
00723
00724
00725
00726 float r1 = G.first * end1Pt().coordX() + G.second - end1Pt().coordY();
00727
00728
00729 float r2 = G.first * end2Pt().coordX() + G.second - end2Pt().coordY();
00730
00731
00732
00733
00734
00735
00736 if((r1 != 0)
00737 && (r2 != 0)
00738 && (signbit(r1) == signbit(r2))
00739 )
00740 return false;
00741
00742
00743 return true;
00744 }
00745
00746
00747 Point LineData::intersectionWithLine(const Shape<LineData>& other,
00748 bool& intersection_on_this,
00749 bool& intersection_on_other) const {
00750 return intersectionWithLine(other.getData(), intersection_on_this,intersection_on_other);
00751 }
00752
00753 Point
00754 LineData::intersectionWithLine(const LineData& other,
00755 bool& intersection_on_this, bool& intersection_on_other) const {
00756
00757
00758
00759
00760
00761
00762 float x1, x2, x3, x4, y1, y2, y3, y4;
00763 x1 = end1Pt().coordX();
00764 x2 = end2Pt().coordX();
00765 x3 = other.end1Pt().coordX();
00766 x4 = other.end2Pt().coordX();
00767 y1 = end1Pt().coordY();
00768 y2 = end2Pt().coordY();
00769 y3 = other.end1Pt().coordY();
00770 y4 = other.end2Pt().coordY();
00771
00772
00773
00774
00775
00776
00777
00778
00779
00780
00781
00782
00783 float denom = ((y4-y3)*(x2-x1) - (x4-x3)*(y2-y1));
00784 float u_a_numerator = ((x4-x3)*(y1-y3) - (y4-y3)*(x1-x3));
00785 float u_b_numerator = ((x2-x1)*(y1-y3) - (y2-y1)*(x1-x3));
00786
00787
00788 if(denom == 0.0) {
00789 if (u_a_numerator == 0.0 && u_b_numerator == 0.0) {
00790
00791
00792
00793 return(end1Pt());
00794 }
00795 else {
00796 cout << x1 << " " << x2 << " " << x3 << " " << x4 << " "
00797 << y1 << " " << y2 << " " << y3 << " " << y4 << endl;
00798 cout << "this theta; " << getOrientation() << ", other theta: " << other.getOrientation() << endl;
00799 cerr << "ERROR in intersectionWithLine: lines are parallel!\n";
00800 return(Point(-9999.0f,-99999.0f));
00801 }
00802 }
00803
00804 else {
00805 float u_a = u_a_numerator / denom;
00806 float u_b = u_b_numerator / denom;
00807
00808
00809 if(0<=u_a && u_a<=1)
00810 intersection_on_this = true;
00811 else
00812 intersection_on_this = false;
00813
00814
00815 if(0<=u_b && u_b<=1)
00816 intersection_on_other = true;
00817 else
00818 intersection_on_other = false;
00819
00820 return(Point((x1+u_a*(x2-x1)),
00821 (y1+u_a*(y2-y1))));
00822 }
00823 }
00824
00825
00826
00827
00828 float LineData::perpendicularDistanceFrom(const Point& otherPt) const {
00829
00830
00831
00832
00833 vector<float> abc = lineEquation_abc();
00834 float& a = abc[0];
00835 float& b = abc[1];
00836 float& c = abc[2];
00837
00838
00839
00840
00841
00842 float d = fabs((a * otherPt.coordX() + b * otherPt.coordY() - c)/sqrt(a*a + b*b));
00843
00844 return(d);
00845 }
00846
00847
00848
00849
00850
00851
00852
00853
00854 #define EXTRACT_LINE_MIN_AREA 60
00855
00856
00857
00858 #define EXTRACT_LINE_MIN_SCORE 15 // was 25
00859 #define BIG_SLOPE_CS 5000.0
00860
00861 bool LineData::pointsOnPerimeterOfWindow(Sketch<bool> const& sketch, double t, double r, Point& pt1, Point& pt2) {
00862 const int width = sketch->getWidth(), height = sketch->getHeight();
00863 r -= RSIZE / 2;
00864 t *= M_PI/180;
00865
00866 pt1.setCoords(cos(t), sin(t));
00867 pt1 *= r;
00868 t += M_PI/2.0;
00869 pt2.setCoords(cos(t), sin(t));
00870 pt2 *= 10;
00871 pt2 += pt1;
00872
00873
00874 if (fabs(pt2.coordX() - pt1.coordX()) <= 0.001) {
00875 pt2.setCoords(pt1.coordX(), height);
00876 return true;
00877 }
00878
00879 double m = (pt2.coordY() - pt1.coordY()) / (pt2.coordX() - pt1.coordX());
00880 double b = pt1.coordY() - m*pt1.coordX();
00881
00882 static vector<Point> accOfPoints;
00883 accOfPoints.clear();
00884
00885
00886 if (b >= 0 && b <= height)
00887 accOfPoints.push_back(Point(0, b));
00888
00889
00890 if ((m*width + b) >= 0 && (m*width + b) <= height)
00891 accOfPoints.push_back(Point( width, m*width + b));
00892
00893
00894 if ((height-b)/m >= 0 && (height-b)/m <= width)
00895 accOfPoints.push_back(Point((height-b)/m, height));
00896
00897
00898 if (-b/m >= 0 && -b/m <= width)
00899 accOfPoints.push_back(Point( -b/m, 0));
00900
00901
00902 if (accOfPoints.size() != 2)
00903 return false;
00904
00905 pt1 = accOfPoints[0];
00906 pt2 = accOfPoints[1];
00907 return true;
00908 }
00909
00910 double LineData::getScore(int t, int r, int height, short hough[TSIZE][RSIZE], float ptsArray[TSIZE][RSIZE]) {
00911 int tally = hough[t][r];
00912
00913
00914
00915
00916
00917 if (ptsArray[t][r] == 0)
00918 return 0;
00919 double scale = min((height / ptsArray[t][r]), 1.0f);
00920 if (tally * scale > 15)
00921 return tally * scale;
00922 return 0;
00923 }
00924
00925 vector<Shape<LineData> > LineData::houghExtractLines(Sketch<bool> const& sketch,
00926 Sketch<bool> const& occluders,
00927 const int num_lines) {
00928 const int width = sketch->getWidth(), height = sketch->getHeight();
00929 SketchSpace &SkS = sketch->getSpace();
00930 ShapeSpace &ShS = SkS.getDualSpace();
00931 NEW_SKETCH_N(skelDist,uint,visops::mdist(sketch));
00932 NEW_SKETCH_N(occ,bool,occluders)
00933 short hough[TSIZE][RSIZE];
00934 float ptsArray[TSIZE][RSIZE];
00935
00936
00937 memset(hough, 0, TSIZE*RSIZE*sizeof(short));
00938 memset(ptsArray, 0, TSIZE*RSIZE*sizeof(float));
00939 Point pt1, pt2;
00940 for (int t = 0; t < TSIZE; t++)
00941 for (int r = 0; r < RSIZE; r++)
00942 if (pointsOnPerimeterOfWindow(sketch, t, r, pt1, pt2))
00943 ptsArray[t][r] = (pt1-pt2).xyNorm();
00944
00945
00946 NEW_SKETCH_N(edges, bool, visops::non_bounds(sketch, 2));
00947
00948
00949 for (int x = 0; x < width; x++) {
00950 for (int y = 0; y < height; y++) {
00951 if ( edges(x,y) ) {
00952 for ( int t = 0; t < TSIZE; t++ ) {
00953 double theta = M_PI/180.0 * t;
00954 int r = int(x*cos(theta) + y*sin(theta) + RSIZE/2);
00955 hough[t][r]++;
00956 }
00957 }
00958 }
00959 }
00960
00961
00962 for(int t = 0; t < TSIZE; t++) {
00963 for(int r = 0; r < RSIZE; r++) {
00964 hough[t][r] = (short)getScore(t, r, height, hough, ptsArray);
00965 }
00966 }
00967
00968 vector<Shape<LineData> > lines_vec;
00969 while ( (int)lines_vec.size() < num_lines ) {
00970 int maxScore = 0;
00971 int maxT = 0;
00972 int maxR = 0;
00973
00974 for (int i = 0; i < TSIZE; i++)
00975 for(int j = 0; j < RSIZE; j++)
00976 if (hough[i][j] > maxScore) {
00977 maxScore = hough[i][j];
00978 maxT = i;
00979 maxR = j;
00980 }
00981
00982 if (maxScore < EXTRACT_LINE_MIN_SCORE) break;
00983 if (! pointsOnPerimeterOfWindow(sketch, maxT, maxR, pt1, pt2)) {
00984 std::cout << "Bad line " << endl;
00985 continue;
00986 }
00987
00988
00989 float x_normpoint = (maxR - RSIZE/2)*cos(maxT * M_PI/180.0);
00990 float y_normpoint = (maxR - RSIZE/2)*sin(maxT * M_PI/180.0);
00991
00992 float m = std::max(std::min(std::tan((maxT * M_PI/180.0)+AngPi((orientation_t)M_PI/2)), BIG_SLOPE_CS), -BIG_SLOPE_CS);
00993 float b = y_normpoint - m*x_normpoint;
00994 Point pt(0, b);
00995 AngPi orientation = atan(m);
00996 int count = 0;
00997 for ( int startval = 0;; ) {
00998 if ( ++count > 6 ) break;
00999 Shape<LineData> extracted_line(ShS, pt, orientation);
01000 if ( std::abs(m) <= 1 )
01001 startval = extracted_line->scanHorizForEndPts(skelDist,occluders,m,b,startval);
01002 else
01003 startval = extracted_line->scanVertForEndPts(skelDist,occluders,m,b,startval);
01004 if ( startval < 0 ) {
01005 extracted_line.deleteShape();
01006 break;
01007 }
01008 extracted_line->setColor(sketch->getColor());
01009 extracted_line->setParentId(sketch->getViewableId());
01010 extracted_line->firstPt().checkValidity(width,height,beg_dist_thresh);
01011 extracted_line->secondPt().checkValidity(width,height,beg_dist_thresh);
01012
01013 SkS.applyTmatinv(extracted_line->firstPt().getCoords());
01014 SkS.applyTmatinv(extracted_line->secondPt().getCoords());
01015 extracted_line->update_derived_properties();
01016
01017 bool matched = false;
01018 for ( vector<Shape<LineData> >::iterator ln_it = lines_vec.begin(); ln_it != lines_vec.end(); ln_it++ ) {
01019 Shape<LineData> &ln = *ln_it;
01020 if ( extracted_line->isMatchFor(ln) ) {
01021 matched = true;
01022
01023 if ( extracted_line->getLength() > 3 + ln->getLength() ) {
01024
01025
01026
01027
01028
01029 ln_it->deleteShape();
01030 *ln_it = extracted_line;
01031 }
01032 else {
01033 extracted_line.deleteShape();
01034
01035
01036
01037
01038
01039 }
01040 break;
01041 }
01042 }
01043 if (matched == false)
01044 lines_vec.push_back(extracted_line);
01045 }
01046
01047
01048 const int tRange = 3;
01049 const int rRange = 3;
01050 for (int i = maxT-tRange; i <= maxT+tRange; i++ )
01051 for (int j = maxR-rRange; j <= maxR+rRange; j++ )
01052 if (i >=0 && i < TSIZE && j >=0 && j < RSIZE) {
01053
01054 hough[i][j] = 0;
01055 }
01056 }
01057 return lines_vec;
01058 }
01059
01060 Shape<LineData> LineData::extractLine(Sketch<bool>& sketch) {
01061
01062 NEW_SKETCH_N(occluders, bool, visops::copy(sketch));
01063 return houghExtractLines(sketch, occluders, 1)[0];
01064
01065
01066 }
01067
01068 Shape<LineData> LineData::extractLine(Sketch<bool>& sketch, const Sketch<bool>& occlusions) {
01069
01070
01071
01072
01073 return houghExtractLines(sketch, occlusions, 1)[0];
01074 }
01075
01076 vector<Shape<LineData> > LineData::extractLines(Sketch<bool> const& sketch, const int num_lines) {
01077
01078 NEW_SKETCH_N(occluders, bool, visops::copy(sketch));
01079 return houghExtractLines(sketch, occluders, num_lines);
01080
01081
01082 }
01083
01084 vector<Shape<LineData> > LineData::extractLines(Sketch<bool> const& sketch, Sketch<bool> const& occluders, const int num_lines) {
01085
01086
01087
01088
01089 return houghExtractLines(sketch, occluders, num_lines);
01090 }
01091
01092
01093
01094
01095
01096
01097
01098
01099
01100
01101
01102
01103
01104
01105
01106
01107
01108
01109
01110
01111
01112
01113
01114
01115
01116
01117
01118
01119
01120
01121
01122
01123
01124
01125
01126
01127
01128
01129
01130
01131
01132
01133
01134
01135
01136
01137
01138
01139
01140
01141
01142
01143
01144
01145
01146
01147
01148
01149
01150
01151
01152
01153
01154
01155
01156
01157
01158
01159
01160
01161
01162
01163
01164
01165
01166
01167
01168
01169
01170
01171
01172
01173
01174
01175
01176
01177
01178
01179
01180
01181
01182
01183
01184
01185
01186
01187
01188
01189
01190
01191
01192
01193
01194
01195
01196
01197
01198
01199
01200
01201
01202
01203
01204
01205
01206
01207
01208
01209
01210
01211
01212
01213
01214
01215
01216
01217
01218
01219
01220
01221
01222
01223
01224
01225
01226
01227
01228
01229
01230
01231
01232
01233
01234 Shape<LineData> LineData::splitLine(ShapeSpace &ShS, Region &skelchunk,
01235 Sketch<bool> &skeleton, const Sketch<bool> &occlusions) {
01236
01237 Point bounds[4] = {skelchunk.findTopBoundPoint(), skelchunk.findLeftBoundPoint(),
01238 skelchunk.findRightBoundPoint(), skelchunk.findBotBoundPoint()};
01239
01240
01241
01242
01243 for (int i = 0; i < 4; i++) {
01244 for (int j = i+1; j < 4; j++) {
01245 if (bounds[i].distanceFrom(bounds[j]) > 20 && ! skelchunk.isContained((bounds[i]+bounds[j])/2.f, 3)) {
01246
01247 LineData ln(ShS,bounds[i],bounds[j]);
01248 Point most_distant(skelchunk.mostDistantPtFrom(ln));
01249 ShS.getDualSpace().applyTmat(most_distant.getCoords());
01250 int const clear_dist = 15;
01251 int x1 = max(0, int(most_distant.coordX() - clear_dist));
01252 int x2 = min(skeleton.width-1, int(most_distant.coordX() + clear_dist));
01253 int y1 = max(0, int(most_distant.coordY() - clear_dist));
01254 int y2 = min(skeleton.height-1, int(most_distant.coordY() + clear_dist));
01255 for (int x = x1; x <= x2; x++)
01256 for (int y = y1; y <= y2; y++)
01257 skeleton(x,y) = false;
01258 return extractLine(skeleton, occlusions);
01259 }
01260 }
01261 }
01262 ShapeRoot invalid;
01263 return ShapeRootType(invalid,LineData);
01264 }
01265
01266
01267
01268
01269
01270
01271
01272
01273
01274
01275
01276
01277
01278
01279
01280
01281
01282
01283
01284
01285
01286
01287
01288
01289
01290
01291
01292
01293
01294
01295
01296 int LineData::scanHorizForEndPts(const Sketch<uint>& skelDist,
01297 const Sketch<bool>& occlusions,
01298 float m, float b, int xstart) {
01299
01300
01301
01302 bool on_line = false;
01303 int curxstart = -1;
01304 int possxstop = -1;
01305 int bestlength = -1;
01306 int goodpixels;
01307 int const xstop = skelDist.height - line_min_length;
01308 while ( xstart < xstop ) {
01309 possxstop = -1;
01310 for (int x = xstart, y, dist=0; x < skelDist.width; x++) {
01311 y = int(m*x+b);
01312 if (y < 0 || y >= skelDist.height) continue;
01313 dist = skelDist(x,y);
01314 if (on_line == false) {
01315 if (dist <= beg_dist_thresh) {
01316 on_line = true;
01317 goodpixels = 0;
01318 curxstart = x - dist;
01319
01320 int curystart;
01321 bool backupflag = false;
01322 while ( curxstart >= 0 && (curystart=int(m*curxstart+b)) >= 0 && curystart < skelDist.height ) {
01323 if ( occlusions(curxstart,curystart) || skelDist(curxstart,curystart) == 0 ) {
01324 --curxstart;
01325 backupflag = true;
01326 } else {
01327 curxstart += backupflag;
01328 break;
01329 }
01330 }
01331 if ( curxstart < 0)
01332 curxstart = 0;
01333 }
01334 }
01335 else {
01336 if ( dist <= end_dist_thresh || occlusions(x,y) ) {
01337 ++goodpixels;
01338 possxstop = x;
01339 }
01340 else if ( dist > extractorGapTolerance ) {
01341
01342 on_line = false;
01343 bestlength = possxstop - curxstart;
01344 if ( bestlength >= line_min_length )
01345 break;
01346 }
01347 }
01348 }
01349 bestlength = possxstop - curxstart;
01350
01351 if ( bestlength <= 0 )
01352 return -1;
01353 if ( bestlength >= line_min_length && float(goodpixels)/bestlength >= 0.8 ) {
01354 float y1 = m*curxstart + b;
01355 float y2 = m*possxstop + b;
01356 setEndPts(Point(curxstart,y1), Point(possxstop,y2));
01357
01358 return possxstop + end_dist_thresh;
01359 }
01360 xstart = max(xstart+1,possxstop) + end_dist_thresh;
01361 }
01362
01363 return -1;
01364 }
01365
01366 int LineData::scanVertForEndPts(const Sketch<uint>& skelDist,
01367 const Sketch<bool>& occlusions,
01368 float m, float b, int ystart) {
01369
01370
01371
01372 bool on_line = false;
01373 int curystart = -1;
01374 int possystop = -1;
01375 int bestlength = -1;
01376 int goodpixels = 0;
01377 int const ystop = skelDist.height - line_min_length;
01378 while ( ystart < ystop ) {
01379
01380 possystop = -1;
01381 for (int x, y = ystart, dist=0; y < skelDist.height; y++) {
01382 x = int((y-b)/m);
01383 if (x < 0 || x >= skelDist.width) continue;
01384 dist = skelDist(x,y);
01385 if (on_line == false) {
01386 if (dist <= beg_dist_thresh) {
01387 on_line = true;
01388 goodpixels = 0;
01389 curystart = y - dist;
01390
01391 int curxstart;
01392 bool backupflag = false;
01393 while ( curystart >= 0 && (curxstart=int((curystart-b)/m)) >= 0 && curxstart < skelDist.width )
01394 if ( occlusions(curxstart,curystart) || skelDist(curxstart,curystart) == 0 ) {
01395 --curystart;
01396 backupflag = true;
01397 } else {
01398 curystart += backupflag;
01399 break;
01400 }
01401 if ( curystart < 0)
01402 curystart = 0;
01403 }
01404 }
01405 else {
01406 if ( dist <= end_dist_thresh || occlusions(x,y) ) {
01407 ++goodpixels;
01408 possystop = y;
01409 }
01410 else if ( dist > extractorGapTolerance ) {
01411
01412 on_line = false;
01413 bestlength = possystop - curystart;
01414 if ( bestlength >= line_min_length )
01415 break;
01416 }
01417 }
01418 }
01419 bestlength = possystop - curystart;
01420
01421 if ( bestlength <= 0 )
01422 return -1;
01423 if ( bestlength >= line_min_length && float(goodpixels)/bestlength >= 0.8 ) {
01424 float x1 = (curystart-b)/m;
01425 float x2 = (possystop-b)/m;
01426 setEndPts(Point(x1,curystart), Point(x2,possystop));
01427 return possystop + end_dist_thresh;
01428 }
01429
01430
01431 ystart = max(ystart+1,possystop) + end_dist_thresh;
01432
01433 }
01434
01435 return -1;
01436 }
01437
01438
01439 void LineData::clearLine(Sketch<bool>& sketch) {
01440 const Sketch<bool>& line_rendering = getRendering();
01441 uint const clear_dist = 5;
01442 Sketch<bool> not_too_close = (visops::mdist(line_rendering) >= clear_dist);
01443 sketch &= not_too_close;
01444 }
01445
01446
01447
01448
01449
01450
01451
01452 Sketch<bool>& LineData::getRendering() {
01453 if ( ! end1Pt().rendering_valid || ! end2Pt().rendering_valid )
01454 deleteRendering();
01455 if ( rendering_sketch == NULL )
01456 rendering_sketch = render();
01457 return *rendering_sketch;
01458 }
01459
01460
01461
01462 Sketch<bool>* LineData::render() const {
01463 SketchSpace &renderspace = space->getDualSpace();
01464 Sketch<bool>* draw_result =
01465 new Sketch<bool>(renderspace, "render("+getName()+")");
01466
01467
01468 (*draw_result)->inheritFrom(*this);
01469 *draw_result = 0;
01470 renderOnTo(draw_result);
01471 return draw_result;
01472 }
01473
01474 void LineData::renderOnTo(Sketch<bool>* draw_result) const {
01475 SketchSpace &renderspace = space->getDualSpace();
01476 int const width = renderspace.getWidth();
01477 int const height = renderspace.getHeight();
01478 float x1,y1,x2,y2;
01479 setDrawCoords(x1, y1, x2, y2, width, height);
01480
01481 drawline2d(*draw_result, (int)x1, (int)y1, (int)x2, (int)y2);
01482 end1Pt().rendering_valid = true;
01483 end2Pt().rendering_valid = true;
01484 }
01485
01486
01487 void LineData::setDrawCoords(float& x1,float& y1, float& x2, float& y2,
01488 const int width, const int height) const {
01489 EndPoint e1(end1Pt());
01490 EndPoint e2(end2Pt());
01491
01492 space->getDualSpace().applyTmat(e1.getCoords());
01493 space->getDualSpace().applyTmat(e2.getCoords());
01494 const EndPoint &left_point = e1.coordX() <= e2.coordX() ? e1 : e2;
01495 const EndPoint &right_point = e1.coordX() <= e2.coordX() ? e2 : e1;
01496
01497
01498
01499 if( left_point.coordY() == right_point.coordY()) {
01500 if (!left_point.isActive())
01501 x1 = 0;
01502 else x1 = max(0.0f,min(width-1.0f,left_point.coordX()));
01503
01504 if (!right_point.isActive())
01505 x2 = width-1;
01506 else x2 = max(0.0f,min(width-1.0f,right_point.coordX()));
01507
01508 y1 = left_point.coordY();
01509 y2 = y1;
01510 }
01511 else if ( left_point.coordX() == right_point.coordX()) {
01512
01513 const EndPoint &top_point = e1.coordY() <= e2.coordY() ? e1 : e2;
01514 const EndPoint &bottom_point = e1.coordY() <= e2.coordY() ? e2 : e1;
01515
01516
01517 if(!top_point.isActive())
01518 y1 = 0;
01519 else y1 = max(0.0f,min(height-1.0f,top_point.coordY()));
01520
01521 if (!bottom_point.isActive())
01522 y2 = height-1;
01523 else y2 = max(0.f,min(height-1.0f,bottom_point.coordY()));
01524
01525 x1 = left_point.coordX();
01526 x2 = x1;
01527 }
01528
01529 else {
01530 float const m = (right_point.coordY()-left_point.coordY())/(right_point.coordX()-left_point.coordX());
01531 float const b = left_point.coordY() - m*left_point.coordX();
01532
01533
01534 int const i0x = (int)((0-b)/m);
01535 int const ihx = (int)(((height-1)-b)/m);
01536 int const i0y = (int)(m*0+b);
01537 int const iwy = (int)(m*(width-1)+b);
01538
01539
01540 if ( left_point.isActive() ) {
01541 x1 = left_point.coordX();
01542 y1 = left_point.coordY();
01543 }
01544
01545
01546 else {
01547
01548
01549 if ( i0y >= 0 && i0y < height ) {
01550 x1 = 0;
01551 y1 = i0y;
01552 }
01553
01554
01555 else {
01556
01557
01558 if ( i0x < ihx ) {
01559 x1 = i0x;
01560 y1 = 0;
01561 }
01562
01563
01564 else {
01565 x1 = ihx;
01566 y1 = height-1;
01567 }
01568 }
01569 }
01570
01571
01572 if ( right_point.isActive() ) {
01573 x2 = right_point.coordX();
01574 y2 = right_point.coordY();
01575 }
01576 else {
01577 if ( iwy >= 0 && iwy < height ) {
01578 x2 = width-1;
01579 y2 = iwy;
01580 }
01581 else {
01582 if ( i0x > ihx ) {
01583 x2 = i0x;
01584 y2 = 0;
01585 }
01586 else {
01587 x2 = ihx;
01588 y2 = height-1;
01589 }
01590 }
01591 }
01592 }
01593 }
01594
01595 void LineData::drawline2d(Sketch<bool>& canvas, int x0, int y0, int x1, int y1) {
01596 int width = canvas->getWidth();
01597 int height = canvas->getHeight();
01598
01599
01600
01601
01602 int dx = x1-x0, ax = abs(dx)<<1, sx = signbit(dx) ? -1 : 1;
01603 int dy = y1-y0, ay = abs(dy)<<1, sy = signbit(dy) ? -1 : 1;
01604 int d, x=x0, y=y0;
01605
01606 if ( ax > ay ) {
01607
01608 d = ay-(ax>>1);
01609 for (;;) {
01610 if (x >= 0 && y >= 0 && x < width && y < height)
01611 canvas(x,y) = true;
01612
01613 if (x==x1)
01614 return;
01615
01616 if ( d >= 0 ) {
01617 y += sy;
01618 d -= ax;
01619 }
01620
01621 x += sx;
01622 d += ay;
01623 }
01624 }
01625 else {
01626
01627 d = ax-(ay>>1);
01628 for (;;) {
01629 if (x >= 0 && y >= 0 && x < width && y < height)
01630 canvas(x,y) = true;
01631
01632 if ( y == y1 )
01633 return;
01634
01635 if ( d >= 0 ) {
01636 x += sx;
01637 d -= ay;
01638 }
01639
01640 y += sy;
01641 d += ax;
01642 }
01643 }
01644 }
01645
01646
01647
01648
01649 std::vector<Shape<LineData> > LineData::houghTransform(const Sketch<bool>& fat,
01650 const Sketch<bool>& skinny,
01651 const Sketch<bool>& occlusions,
01652 const size_t num_lines, int minLength)
01653 {
01654 if ( minLength == -1 )
01655 minLength = static_cast<int>(fat->getWidth() * DEFAULT_MIN_LENGTH);
01656 std::vector<Shape<LineData> > result;
01657 ShapeSpace& ShS = fat->getSpace().getDualSpace();
01658
01659 const int width = fat->getWidth(), height = fat->getHeight();
01660 const int numR = 120, numTheta = 120;
01661 const int numRf = 40, numThetaf = 40;
01662 int hsize = numR*numTheta, hsizef = numRf * numThetaf;
01663 int htab[hsize], hfine[hsizef];
01664 float minTheta = 0, maxTheta = 2*(float)M_PI;
01665 float minR = 0, maxR = std::sqrt((float)width*(float)width+(float)height*(float)height);
01666 float thetaSpread = maxTheta - minTheta;
01667 float rSpread = maxR - minR;
01668
01669
01670 for (int i = 0; i < hsize; i++)
01671 htab[i] = 0;
01672
01673
01674 float theta, r;
01675 int ridx;
01676 for (int x = 0; x < width; x++) {
01677 for (int y = 0; y < height; y++) {
01678 if (fat(x,y) == true) {
01679 for (int tidx = 0; tidx < numTheta ; tidx++) {
01680 theta = minTheta + tidx*thetaSpread/numTheta;
01681 r = (float)x*cos(theta)+(float)y*sin(theta);
01682 ridx = (int)((r-minR)*(float)numR/rSpread);
01683 if (ridx >= 0 && ridx < numR)
01684 htab[ridx+tidx*numR]++;
01685 }
01686 }
01687 }
01688 }
01689
01690
01691
01692
01693
01694
01695
01696
01697
01698
01699
01700
01701
01702
01703
01704
01705
01706
01707
01708
01709
01710 int max_val = -1, max_i = 0;
01711 while(result.size() < num_lines) {
01712
01713
01714 max_val = -1;
01715 max_i = 0;
01716 for (int i = 0; i < numR*numTheta; i++) {
01717 if (max_val < htab[i]) {
01718 max_val = htab[i];
01719 max_i = i;
01720 }
01721 }
01722
01723
01724 if (max_val < minLength)
01725 break;
01726
01727 float bestR = (max_i%numR)*rSpread/numR + minR;
01728 float bestTheta = (max_i/numR)*thetaSpread/numTheta + minTheta;
01729
01730
01731 const float fthetaSpread = (float)M_PI/40.f;
01732 const float frSpread = maxR/40.f;
01733 float fminTheta = bestTheta-fthetaSpread/2.f;
01734 float fminR = bestR-frSpread/2.f;
01735
01736 for (int i = 0; i < hsizef; i++)
01737 hfine[i] = 0;
01738
01739
01740 float thetaf, rf;
01741 int ridxf;
01742 for (int x = 0; x < width; x++) {
01743 for (int y = 0; y < height; y++) {
01744 if (skinny(x,y) == true) {
01745 for (int tidx = 0; tidx < numThetaf ; tidx++) {
01746 thetaf = fminTheta + tidx*fthetaSpread/numThetaf;
01747 rf = (float)x*cos(thetaf)+(float)y*sin(thetaf);
01748 ridxf = (int)((rf-fminR)*(float)numRf/frSpread);
01749 if (ridxf >= 0 && ridxf < numRf)
01750 hfine[ridxf+tidx*numRf]++;
01751 }
01752 }
01753 }
01754 }
01755
01756
01757 int max_val_f = -1, max_i_f = 0;
01758 for (int i = 0; i < numRf*numThetaf; i++) {
01759 if (max_val_f < hfine[i]) {
01760 max_val_f = hfine[i];
01761 max_i_f = i;
01762 }
01763 }
01764 float bestRf = (max_i_f%numRf)*frSpread/numRf + fminR;
01765 float bestThetaf = (max_i_f/numRf)*fthetaSpread/numThetaf + fminTheta;
01766
01767
01768
01769 float lineLen = 500;
01770
01771 const float x1 = bestRf*cos(bestThetaf), y1 = bestRf*sin(bestThetaf);
01772 const float x2 = x1+lineLen*std::cos(bestThetaf+(float)M_PI_2), y2 = y1+lineLen*std::sin(bestThetaf+(float)M_PI_2);
01773 const float x3 = x1-lineLen*std::cos(bestThetaf+(float)M_PI_2), y3 = y1-lineLen*std::sin(bestThetaf+(float)M_PI_2);
01774
01775 Point e1(x3,y3), e2(x2,y2);
01776 Shape<LineData> line(ShS,e1,e2);
01777 line->setColor(skinny->getColor());
01778 line->setParentId(skinny->getViewableId());
01779
01780
01781 float m = (y1 != 0) ? -(x1/y1) : BIG_SLOPE;
01782 float b = y1 - m*x1;
01783
01784 NEW_SKETCH_N(skelDist,uint,visops::mdist(skinny));
01785 if ( abs(m) <= 1 )
01786 line->scanHorizForEndPts(skelDist,occlusions,m,b,0);
01787 else
01788 line->scanVertForEndPts(skelDist,occlusions,m,b,0);
01789
01790
01791
01792 line->end1_pt.checkValidity(width,height,beg_dist_thresh);
01793 line->end2_pt.checkValidity(width,height,beg_dist_thresh);
01794
01795
01796
01797
01798 float len = line->getLength();
01799 EndPoint start = line->end1_pt;
01800 EndPoint dx = (line->end2_pt - start)/len;
01801 float onCount = 0;
01802 for (float i=0; i<len; i++) {
01803 EndPoint cur = start+dx*i;
01804 if (skinny((int)(cur.coordX()), (int)(cur.coordY()))) {
01805 onCount++;
01806 }
01807 }
01808
01809 std::cout << "Line with " << onCount << " / " << len << " pixels";
01810 if (line->isLongerThan(minLength) && (onCount / len) > .5 ) {
01811 std::cout<<" accepted";
01812 result.push_back(line);
01813 }
01814 std::cout<<std::endl;
01815
01816
01817
01818 for(int row = -5; row <= 5; row++) {
01819 for (int col = -(5-abs(row)); col <=5-abs(row); col++) {
01820 int index = (max_i - (max_i%numR) + ((max_i+col)%numR)+row*numR + numR*numTheta)%(numR*numTheta);
01821 std::cout << "Empty hough " << (index%numR) << "," << (index/numR) << " " << std::endl;
01822 htab[(max_i - (max_i%numR) + ((max_i+col)%numR)+row*numR + numR*numTheta)%(numR*numTheta)] = 0;
01823 }
01824
01825 }
01826
01827 }
01828
01829 return result;
01830 }
01831
01832 bool LineData::linesParallel(Shape<LineData> l1, Shape<LineData>l2)
01833 {
01834 const float maxDTheta = .15f, minDThetaR = 10.0f;
01835 float dTheta = l1->getOrientation() - l2->getOrientation();
01836 if (dTheta > (float)M_PI_2)
01837 dTheta = dTheta - (float)M_PI;
01838 if (std::abs(dTheta) < maxDTheta)
01839 {
01840 if (std::abs (100*dTheta + (l1->rho_norm - l2->rho_norm)) > minDThetaR)
01841 return true;
01842 }
01843 return false;
01844 }
01845
01846 LineData& LineData::operator=(const LineData& other) {
01847 if (&other == this)
01848 return *this;
01849 BaseData::operator=(other);
01850 end1_pt = other.end1_pt;
01851 end2_pt = other.end2_pt;
01852 rho_norm = other.rho_norm;
01853 theta_norm = other.theta_norm;
01854 orientation = other.orientation;
01855 length = other.length;
01856 return *this;
01857 }
01858
01859
01860
01861 bool LineData::pointsOnSameSide(const Point& p1, const Point& p2)
01862 {
01863 float dx = end2_pt.coordX() - end1_pt.coordX();
01864 float dy = end2_pt.coordY() - end1_pt.coordY();
01865
01866 float p1val = (p1.coordY() - end1_pt.coordY())*dx - (p1.coordX() - end1_pt.coordX())*dy;
01867 float p2val = (p2.coordY() - end1_pt.coordY())*dx - (p2.coordX() - end1_pt.coordX())*dy;
01868
01869 return (p1val>0) == (p2val>0);
01870 }
01871
01872
01873
01874 bool LineData::pointOnLine(const Point& p)
01875 {
01876 const float BOUNDS_EXTEND = 1;
01877 float dx = end2_pt.coordX() - end1_pt.coordX();
01878 float dy = end2_pt.coordY() - end1_pt.coordY();
01879 float val = (p.coordY() - end1_pt.coordY())*dx - (p.coordX() - end1_pt.coordX())*dy;
01880 val /= length;
01881 bool inBounds = (p.coordX() >= (leftPt().coordX() - BOUNDS_EXTEND)) &&
01882 (p.coordX() <= (rightPt().coordX() + BOUNDS_EXTEND)) &&
01883 (p.coordY() >= (topPt().coordY() - BOUNDS_EXTEND)) &&
01884 (p.coordY() <= (bottomPt().coordY() + BOUNDS_EXTEND));
01885 return (val > -1) && (val < 1) && inBounds;
01886 }
01887
01888
01889
01890
01891 bool LineData::LineDataLengthLessThan::operator() (const LineData &line1, const LineData &line2) const {
01892 return line1.getLength() < line2.getLength();
01893 }
01894
01895 bool LineData::LengthLessThan::operator() (const Shape<LineData> &line1, const Shape<LineData> &line2) const {
01896 return line1->getLength() < line2->getLength();
01897 }
01898
01899 bool LineData::ParallelTest::operator() (const Shape<LineData> &line1, const Shape<LineData> &line2) const {
01900 return angdist(line1->getOrientation(),line2->getOrientation()) <= tolerance;
01901 }
01902
01903 bool LineData::PerpendicularTest::operator() (const Shape<LineData> &line1, const Shape<LineData> &line2) const {
01904 return angdist(angdist(line1->getOrientation(),line2->getOrientation()), (orientation_t)M_PI/2) <= tolerance;
01905 }
01906
01907 bool LineData::ColinearTest::operator() (const Shape<LineData> &line1, const Shape<LineData> &line2) const {
01908 return ParallelTest(ang_tol)(line1,line2) &&
01909 abs( line1->getRhoNorm() - line2->getRhoNorm() ) <= dist_tol;
01910 }
01911
01912 bool LineData::IsHorizontal::operator() (const Shape<LineData> &line) const {
01913 const AngPi orient = line->getOrientation();
01914 return (orient <= threshold) || (orient >= M_PI - threshold);
01915 }
01916
01917 bool LineData::IsVertical::operator() (const Shape<LineData> &line) const {
01918 const AngPi orient = line->getOrientation();
01919 return (orient >= threshold) && (orient <= M_PI - threshold);
01920 }
01921
01922 }