Shamusworld >> Repos - architektonas/blob - src/base/rs_spline.cpp

   1 // rs_spline.cpp
   2 //
   3 // Part of the Architektonas Project
   4 // Originally part of QCad Community Edition by Andrew Mustun
   5 // Extensively rewritten and refactored by James L. Hammons
   6 // (C) 2010 Underground Software
   7 //
   8 // JLH = James L. Hammons <jlhamm@acm.org>
   9 //
  10 // Who  When        What
  11 // ---  ----------  -----------------------------------------------------------
  12 // JLH  06/02/2010  Added this text. :-)
  13 //
  14
  15 #include "rs_spline.h"
  16
  17 #include "rs_debug.h"
  18 #include "rs_graphicview.h"
  19 #include "drawing.h"
  20 #include "paintintf.h"
  21
  22 /**
  23  * Constructor.
  24  */
  25 RS_Spline::RS_Spline(RS_EntityContainer * parent, const RS_SplineData & d):
  26         RS_EntityContainer(parent), data(d)
  27 {
  28         calculateBorders();
  29 }
  30
  31 /**
  32  * Destructor.
  33  */
  34 RS_Spline::~RS_Spline()
  35 {
  36 }
  37
  38 RS_Entity * RS_Spline::clone()
  39 {
  40         RS_Spline * l = new RS_Spline(*this);
  41 #warning "!!! Need to deal with setAutoDelete() Qt3->Qt4 !!!"
  42 //      l->entities.setAutoDelete(entities.autoDelete());
  43         l->initId();
  44         l->detach();
  45         return l;
  46 }
  47
  48 /**     @return RS2::EntitySpline */
  49 /*virtual*/ RS2::EntityType RS_Spline::rtti() const
  50 {
  51         return RS2::EntitySpline;
  52 }
  53
  54 /** @return false */
  55 /*virtual*/ bool RS_Spline::isEdge() const
  56 {
  57         return false;
  58 }
  59
  60 /** @return Copy of data that defines the spline. */
  61 RS_SplineData RS_Spline::getData() const
  62 {
  63         return data;
  64 }
  65
  66 /** Sets the splines degree (1-3). */
  67 void RS_Spline::setDegree(int deg)
  68 {
  69         if (deg >= 1 && deg <= 3)
  70                 data.degree = deg;
  71 }
  72
  73 /** @return Degree of this spline curve (1-3).*/
  74 int RS_Spline::getDegree()
  75 {
  76         return data.degree;
  77 }
  78
  79 /** @return 0. */
  80 int RS_Spline::getNumberOfKnots()
  81 {
  82         return 0;
  83 }
  84
  85 /** @return Number of control points. */
  86 int RS_Spline::getNumberOfControlPoints()
  87 {
  88         return data.controlPoints.count();
  89 }
  90
  91 /**
  92  * @retval true if the spline is closed.
  93  * @retval false otherwise.
  94  */
  95 bool RS_Spline::isClosed()
  96 {
  97         return data.closed;
  98 }
  99
 100 /**
 101  * Sets the closed falg of this spline.
 102  */
 103 void RS_Spline::setClosed(bool c)
 104 {
 105         data.closed = c;
 106         update();
 107 }
 108
 109 void RS_Spline::calculateBorders()
 110 {
 111     /*minV = Vector::minimum(data.startpoint, data.endpoint);
 112     maxV = Vector::maximum(data.startpoint, data.endpoint);
 113
 114     Q3ValueList<Vector>::iterator it;
 115     for (it = data.controlPoints.begin();
 116     it!=data.controlPoints.end(); ++it) {
 117
 118     minV = Vector::minimum(*it, minV);
 119     maxV = Vector::maximum(*it, maxV);
 120 }
 121     */
 122 }
 123
 124 VectorSolutions RS_Spline::getRefPoints()
 125 {
 126         VectorSolutions ret(data.controlPoints.count());
 127
 128         int i = 0;
 129         QList<Vector>::iterator it;
 130
 131         for(it=data.controlPoints.begin(); it!=data.controlPoints.end(); ++it, ++i)
 132         {
 133                 ret.set(i, (*it));
 134         }
 135
 136         return ret;
 137 }
 138
 139 Vector RS_Spline::getNearestRef(const Vector & coord, double * dist)
 140 {
 141         //return getRefPoints().getClosest(coord, dist);
 142         return RS_Entity::getNearestRef(coord, dist);
 143 }
 144
 145 Vector RS_Spline::getNearestSelectedRef(const Vector & coord, double * dist)
 146 {
 147         //return getRefPoints().getClosest(coord, dist);
 148         return RS_Entity::getNearestSelectedRef(coord, dist);
 149 }
 150
 151 /**
 152  * Updates the internal polygon of this spline. Called when the
 153  * spline or it's data, position, .. changes.
 154  */
 155 void RS_Spline::update()
 156 {
 157         RS_DEBUG->print("RS_Spline::update");
 158
 159         clear();
 160
 161         if (isUndone())
 162                 return;
 163
 164         if (data.degree < 1 || data.degree > 3)
 165         {
 166                 RS_DEBUG->print("RS_Spline::update: invalid degree: %d", data.degree);
 167                 return;
 168         }
 169
 170         if (data.controlPoints.count() < (uint)data.degree + 1)
 171         {
 172                 RS_DEBUG->print("RS_Spline::update: not enough control points");
 173                 return;
 174         }
 175
 176         resetBorders();
 177
 178 //      Q3ValueList<Vector> tControlPoints = data.controlPoints;
 179         QList<Vector> tControlPoints = data.controlPoints;
 180
 181         if (data.closed)
 182         {
 183                 for(int i=0; i<data.degree; ++i)
 184                         tControlPoints.append(data.controlPoints[i]);
 185         }
 186
 187         int i;
 188         int npts = tControlPoints.count();
 189         // order:
 190         int k = data.degree + 1;
 191         // resolution:
 192         int p1 = getGraphicVariableInt("$SPLINESEGS", 8) * npts;
 193
 194         double * b = new double[npts * 3 + 1];
 195         double * h = new double[npts + 1];
 196         double * p = new double[p1 * 3 + 1];
 197
 198 //      Q3ValueList<Vector>::iterator it;
 199         QList<Vector>::iterator it;
 200         i = 1;
 201
 202         for(it=tControlPoints.begin(); it!=tControlPoints.end(); ++it)
 203         {
 204                 b[i + 0] = (*it).x;
 205                 b[i + 1] = (*it).y;
 206                 b[i + 2] = 0.0;
 207
 208                 RS_DEBUG->print("RS_Spline::update: b[%d]: %f/%f", i, b[i], b[i + 1]);
 209                 i += 3;
 210         }
 211
 212         // set all homogeneous weighting factors to 1.0
 213         for(i=1; i<=npts; i++)
 214                 h[i] = 1.0;
 215
 216         for(i=1; i<=3*p1; i++)
 217                 p[i] = 0.0;
 218
 219         if (data.closed)
 220                 rbsplinu(npts, k, p1, b, h, p);
 221         else
 222                 rbspline(npts, k, p1, b, h, p);
 223
 224         Vector prev(false);
 225
 226         for(i=1; i<=3*p1; i=i+3)
 227         {
 228                 if (prev.valid)
 229                 {
 230                         RS_Line * line = new RS_Line(this, RS_LineData(prev, Vector(p[i], p[i + 1])));
 231                         line->setLayer(NULL);
 232                         line->setPen(RS_Pen(RS2::FlagInvalid));
 233                         addEntity(line);
 234                 }
 235
 236                 prev = Vector(p[i], p[i+1]);
 237                 minV = Vector::minimum(prev, minV);
 238                 maxV = Vector::maximum(prev, maxV);
 239         }
 240
 241         delete[] b;
 242         delete[] h;
 243         delete[] p;
 244 }
 245
 246 Vector RS_Spline::getNearestEndpoint(const Vector & coord, double * dist)
 247 {
 248     double minDist = RS_MAXDOUBLE;
 249     double d;
 250     Vector ret(false);
 251
 252     for (uint i=0; i<data.controlPoints.count(); i++) {
 253         d = data.controlPoints[i].distanceTo(coord);
 254
 255         if (d<minDist) {
 256             minDist = d;
 257             ret = data.controlPoints[i];
 258         }
 259     }
 260     if (dist!=NULL) {
 261         *dist = minDist;
 262     }
 263         return ret;
 264 }
 265
 266 /*
 267 // The default implementation of RS_EntityContainer is inaccurate but
 268 //   has to do for now..
 269 Vector RS_Spline::getNearestPointOnEntity(const Vector& coord,
 270         bool onEntity, double* dist, RS_Entity** entity) {
 271 }
 272 */
 273
 274 Vector RS_Spline::getNearestCenter(const Vector & /*coord*/, double * dist)
 275 {
 276         if (dist != NULL)
 277                 *dist = RS_MAXDOUBLE;
 278
 279         return Vector(false);
 280 }
 281
 282 Vector RS_Spline::getNearestMiddle(const Vector & /*coord*/, double * dist)
 283 {
 284         if (dist!=NULL)
 285                 *dist = RS_MAXDOUBLE;
 286
 287         return Vector(false);
 288 }
 289
 290 Vector RS_Spline::getNearestDist(double /*distance*/, const Vector & /*coord*/, double * dist)
 291 {
 292         if (dist != NULL)
 293                 *dist = RS_MAXDOUBLE;
 294
 295         return Vector(false);
 296 }
 297
 298 void RS_Spline::move(Vector offset)
 299 {
 300 //      Q3ValueList<Vector>::iterator it;
 301         QList<Vector>::iterator it;
 302
 303         for(it=data.controlPoints.begin(); it!=data.controlPoints.end(); ++it)
 304                 (*it).move(offset);
 305
 306         update();
 307 }
 308
 309 void RS_Spline::rotate(Vector center, double angle)
 310 {
 311 //      Q3ValueList<Vector>::iterator it;
 312         QList<Vector>::iterator it;
 313
 314         for(it=data.controlPoints.begin(); it!=data.controlPoints.end(); ++it)
 315                 (*it).rotate(center, angle);
 316
 317         update();
 318 }
 319
 320 void RS_Spline::scale(Vector center, Vector factor)
 321 {
 322 //      Q3ValueList<Vector>::iterator it;
 323         QList<Vector>::iterator it;
 324
 325         for(it=data.controlPoints.begin(); it!=data.controlPoints.end(); ++it)
 326                 (*it).scale(center, factor);
 327
 328         update();
 329 }
 330
 331 void RS_Spline::mirror(Vector axisPoint1, Vector axisPoint2)
 332 {
 333 //      Q3ValueList<Vector>::iterator it;
 334         QList<Vector>::iterator it;
 335
 336         for(it=data.controlPoints.begin(); it!=data.controlPoints.end(); ++it)
 337                 (*it).mirror(axisPoint1, axisPoint2);
 338
 339         update();
 340 }
 341
 342 void RS_Spline::moveRef(const Vector & ref, const Vector & offset)
 343 {
 344 //      Q3ValueList<Vector>::iterator it;
 345         QList<Vector>::iterator it;
 346
 347         for(it=data.controlPoints.begin(); it!=data.controlPoints.end(); ++it)
 348         {
 349                 if (ref.distanceTo(*it) < 1.0e-4)
 350                         (*it).move(offset);
 351         }
 352
 353         update();
 354 }
 355
 356 //void RS_Spline::draw(RS_Painter* painter, RS_GraphicView* view, double /*patternOffset*/)
 357 void RS_Spline::draw(PaintInterface * painter, RS_GraphicView * view, double /*patternOffset*/)
 358 {
 359         if (painter == NULL || view == NULL)
 360                 return;
 361
 362         RS_Entity * e = firstEntity(RS2::ResolveNone);
 363         double offset = 0.0;
 364
 365         if (e != NULL)
 366         {
 367                 view->drawEntity(e);
 368                 offset += e->getLength();
 369                 //RS_DEBUG->print("offset: %f\nlength was: %f", offset, e->getLength());
 370         }
 371
 372         for (RS_Entity * e=nextEntity(RS2::ResolveNone); e!=NULL; e = nextEntity(RS2::ResolveNone))
 373         {
 374                 view->drawEntityPlain(e, -offset);
 375                 offset += e->getLength();
 376                 //RS_DEBUG->print("offset: %f\nlength was: %f", offset, e->getLength());
 377         }
 378 }
 379
 380 /**
 381  * Todo: draw the spline, user patterns.
 382  */
 383 /*
 384 void RS_Spline::draw(RS_Painter* painter, RS_GraphicView* view) {
 385    if (painter==NULL || view==NULL) {
 386        return;
 387    }
 388
 389    / *
 390       if (data.controlPoints.count()>0) {
 391           Vector prev(false);
 392           Q3ValueList<Vector>::iterator it;
 393           for (it = data.controlPoints.begin(); it!=data.controlPoints.end(); ++it) {
 394               if (prev.valid) {
 395                   painter->drawLine(view->toGui(prev),
 396                                     view->toGui(*it));
 397               }
 398               prev = (*it);
 399           }
 400       }
 401    * /
 402
 403    int i;
 404    int npts = data.controlPoints.count();
 405    // order:
 406    int k = 4;
 407    // resolution:
 408    int p1 = 100;
 409
 410    double* b = new double[npts*3+1];
 411    double* h = new double[npts+1];
 412    double* p = new double[p1*3+1];
 413
 414    Q3ValueList<Vector>::iterator it;
 415    i = 1;
 416    for (it = data.controlPoints.begin(); it!=data.controlPoints.end(); ++it) {
 417        b[i] = (*it).x;
 418        b[i+1] = (*it).y;
 419        b[i+2] = 0.0;
 420
 421         RS_DEBUG->print("RS_Spline::draw: b[%d]: %f/%f", i, b[i], b[i+1]);
 422         i+=3;
 423    }
 424
 425    // set all homogeneous weighting factors to 1.0
 426    for (i=1; i <= npts; i++) {
 427        h[i] = 1.0;
 428    }
 429
 430    //
 431    for (i = 1; i <= 3*p1; i++) {
 432        p[i] = 0.0;
 433    }
 434
 435    rbspline(npts,k,p1,b,h,p);
 436
 437    Vector prev(false);
 438    for (i = 1; i <= 3*p1; i=i+3) {
 439        if (prev.valid) {
 440            painter->drawLine(view->toGui(prev),
 441                              view->toGui(Vector(p[i], p[i+1])));
 442        }
 443        prev = Vector(p[i], p[i+1]);
 444    }
 445 }
 446 */
 447
 448 /**
 449  * @return The reference points of the spline.
 450  */
 451 //Q3ValueList<Vector> RS_Spline::getControlPoints()
 452 QList<Vector> RS_Spline::getControlPoints()
 453 {
 454         return data.controlPoints;
 455 }
 456
 457 /**
 458  * Appends the given point to the control points.
 459  */
 460 void RS_Spline::addControlPoint(const Vector & v)
 461 {
 462         data.controlPoints.append(v);
 463 }
 464
 465 /**
 466  * Removes the control point that was last added.
 467  */
 468 void RS_Spline::removeLastControlPoint()
 469 {
 470         data.controlPoints.pop_back();
 471 }
 472
 473 /**
 474  * Generates B-Spline open knot vector with multiplicity
 475  * equal to the order at the ends.
 476  */
 477 void RS_Spline::knot(int num, int order, int knotVector[])
 478 {
 479         knotVector[1] = 0;
 480
 481         for(int i=2; i<=num+order; i++)
 482         {
 483                 if ((i > order) && (i < num + 2))
 484                 {
 485                         knotVector[i] = knotVector[i - 1] + 1;
 486                 }
 487                 else
 488                 {
 489                         knotVector[i] = knotVector[i - 1];
 490                 }
 491         }
 492 }
 493
 494 /**
 495  * Generates rational B-spline basis functions for an open knot vector.
 496  */
 497 void RS_Spline::rbasis(int c, double t, int npts, int x[], double h[], double r[])
 498 {
 499         int nplusc;
 500         int i, k;
 501         double d, e;
 502         double sum;
 503         //double temp[36];
 504
 505         nplusc = npts + c;
 506
 507         double * temp = new double[nplusc + 1];
 508
 509         // calculate the first order nonrational basis functions n[i]
 510         for(i=1; i<=nplusc-1; i++)
 511         {
 512                 if ((t >= x[i]) && (t < x[i + 1]))
 513                         temp[i] = 1;
 514                 else
 515                         temp[i] = 0;
 516         }
 517
 518         /* calculate the higher order nonrational basis functions */
 519
 520         for(k=2; k<=c; k++)
 521         {
 522                 for(i=1; i<=nplusc-k; i++)
 523                 {
 524                         // if the lower order basis function is zero skip the calculation
 525                         if (temp[i] != 0)
 526                                 d = ((t - x[i]) * temp[i]) / (x[i + k - 1] - x[i]);
 527                         else
 528                                 d = 0;
 529
 530                         // if the lower order basis function is zero skip the calculation
 531                         if (temp[i + 1] != 0)
 532                                 e = ((x[i + k] - t) * temp[i + 1]) / (x[i + k] - x[i + 1]);
 533                         else
 534                                 e = 0;
 535
 536                         temp[i] = d + e;
 537                 }
 538         }
 539
 540         // pick up last point
 541         if (t == (double)x[nplusc])
 542                 temp[npts] = 1;
 543
 544         // calculate sum for denominator of rational basis functions
 545         sum = 0.;
 546
 547         for(i=1; i<=npts; i++)
 548                 sum = sum + temp[i] * h[i];
 549
 550         // form rational basis functions and put in r vector
 551         for(i=1; i<=npts; i++)
 552         {
 553                 if (sum != 0)
 554                         r[i] = (temp[i] * h[i]) / (sum);
 555                 else
 556                         r[i] = 0;
 557         }
 558
 559         delete[] temp;
 560 }
 561
 562 /**
 563  * Generates a rational B-spline curve using a uniform open knot vector.
 564  */
 565 void RS_Spline::rbspline(int npts, int k, int p1, double b[], double h[], double p[])
 566 {
 567         int i, j, icount, jcount;
 568         int i1;
 569         //int x[30]; /* allows for 20 data points with basis function of order 5 */
 570         int nplusc;
 571
 572         double step;
 573         double t;
 574         //double nbasis[20];
 575         double temp;
 576
 577         nplusc = npts + k;
 578
 579         int * x = new int[nplusc + 1];
 580         double * nbasis = new double[npts + 1];
 581
 582         // zero and redimension the knot vector and the basis array
 583
 584         for(i = 0; i <= npts; i++)
 585                 nbasis[i] = 0.0;
 586
 587         for(i = 0; i <= nplusc; i++)
 588                 x[i] = 0;
 589
 590         // generate the uniform open knot vector
 591         knot(npts, k, x);
 592
 593         icount = 0;
 594
 595         // calculate the points on the rational B-spline curve
 596         t = 0;
 597         step = ((double)x[nplusc]) / ((double)(p1 - 1));
 598
 599         for(i1=1; i1<= p1; i1++)
 600         {
 601                 if ((double)x[nplusc] - t < 5e-6)
 602                         t = (double)x[nplusc];
 603
 604                 // generate the basis function for this value of t
 605                 rbasis(k, t, npts, x, h, nbasis);
 606
 607                 // generate a point on the curve
 608                 for(j=1; j<=3; j++)
 609                 {
 610                         jcount = j;
 611                         p[icount + j] = 0.0;
 612
 613                         // Do local matrix multiplication
 614                         for(i=1; i<=npts; i++)
 615                         {
 616                                 temp = nbasis[i] * b[jcount];
 617                                 p[icount + j] = p[icount + j] + temp;
 618                                 jcount = jcount + 3;
 619                         }
 620                 }
 621
 622                 icount = icount + 3;
 623                 t = t + step;
 624         }
 625
 626         delete[] x;
 627         delete[] nbasis;
 628 }
 629
 630 void RS_Spline::knotu(int num, int order, int knotVector[])
 631 {
 632         int nplusc = num + order;
 633         int nplus2 = num + 2;
 634         knotVector[1] = 0;
 635
 636         for(int i=2; i<=nplusc; i++)
 637                 knotVector[i] = i - 1;
 638 }
 639
 640 void RS_Spline::rbsplinu(int npts, int k, int p1, double b[], double h[], double p[])
 641 {
 642         int i, j, icount, jcount;
 643         int i1;
 644         //int x[30];            /* allows for 20 data points with basis function of order 5 */
 645         int nplusc;
 646
 647         double step;
 648         double t;
 649         //double nbasis[20];
 650         double temp;
 651
 652         nplusc = npts + k;
 653
 654         int * x = new int[nplusc + 1];
 655         double * nbasis = new double[npts + 1];
 656
 657         /*  zero and redimension the knot vector and the basis array */
 658
 659         for(i=0; i<=npts; i++)
 660                 nbasis[i] = 0.0;
 661
 662         for(i=0; i<=nplusc; i++)
 663                 x[i] = 0;
 664
 665         /* generate the uniform periodic knot vector */
 666
 667         knotu(npts, k, x);
 668
 669         /*
 670                 printf("The knot vector is ");
 671                 for (i = 1; i <= nplusc; i++){
 672                         printf(" %d ", x[i]);
 673                 }
 674                 printf("\n");
 675
 676                 printf("The usable parameter range is ");
 677                 for (i = k; i <= npts+1; i++){
 678                         printf(" %d ", x[i]);
 679                 }
 680                 printf("\n");
 681         */
 682
 683         icount = 0;
 684
 685         /*    calculate the points on the rational B-spline curve */
 686
 687         t = k - 1;
 688         step = ((double)((npts) - (k - 1))) / ((double)(p1 - 1));
 689
 690         for(i1=1; i1<=p1; i1++)
 691         {
 692                 if ((double)x[nplusc] - t < 5e-6)
 693                         t = (double)x[nplusc];
 694
 695                 rbasis(k, t, npts, x, h, nbasis);      /* generate the basis function for this value of t */
 696                 /*
 697                                 printf("t = %f \n",t);
 698                                 printf("nbasis = ");
 699                                 for (i = 1; i <= npts; i++){
 700                                         printf("%f  ",nbasis[i]);
 701                                 }
 702                                 printf("\n");
 703                 */
 704                 for(j=1; j<=3; j++)      /* generate a point on the curve */
 705                 {
 706                         jcount = j;
 707                         p[icount + j] = 0.0;
 708
 709                         for(i=1; i<=npts; i++)  /* Do local matrix multiplication */
 710                         {
 711                                 temp = nbasis[i] * b[jcount];
 712                                 p[icount + j] = p[icount + j] + temp;
 713                                 /*
 714                                                                 printf("jcount,nbasis,b,nbasis*b,p = %d %f %f %f %f\n",jcount,nbasis[i],b[jcount],temp,p[icount+j]);
 715                                 */
 716                                 jcount = jcount + 3;
 717                         }
 718                 }
 719                 /*
 720                                 printf("icount, p %d %f %f %f \n",icount,p[icount+1],p[icount+2],p[icount+3]);
 721                 */
 722                 icount = icount + 3;
 723                 t = t + step;
 724         }
 725
 726         delete[] x;
 727         delete[] nbasis;
 728 }
 729
 730 /**
 731  * Dumps the spline's data to stdout.
 732  */
 733 std::ostream & operator<<(std::ostream & os, const RS_Spline & l)
 734 {
 735         os << " Spline: " << l.getData() << "\n";
 736         return os;
 737 }