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

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