synfig-studio/trunk/src/synfigapp/blineconvert.cpp

   1 /* === S Y N F I G ========================================================= */
   2 /*!     \file blineconvert.cpp
   3 **      \brief Template File
   4 **
   5 **      $Id$
   6 **
   7 **      \legal
   8 **      Copyright (c) 2002-2005 Robert B. Quattlebaum Jr., Adrian Bentley
   9 **
  10 **      This package is free software; you can redistribute it and/or
  11 **      modify it under the terms of the GNU General Public License as
  12 **      published by the Free Software Foundation; either version 2 of
  13 **      the License, or (at your option) any later version.
  14 **
  15 **      This package is distributed in the hope that it will be useful,
  16 **      but WITHOUT ANY WARRANTY; without even the implied warranty of
  17 **      MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
  18 **      General Public License for more details.
  19 **      \endlegal
  20 */
  21 /* ========================================================================= */
  22
  23 /* === H E A D E R S ======================================================= */
  24
  25 #ifdef USING_PCH
  26 #       include "pch.h"
  27 #else
  28 #ifdef HAVE_CONFIG_H
  29 #       include <config.h>
  30 #endif
  31
  32 #include "blineconvert.h"
  33 #include <vector>
  34 #include <ETL/gaussian>
  35 #include <ETL/hermite>
  36 #include <ETL/clock>
  37 #include <float.h>
  38 #include <algorithm>
  39 #include <synfig/general.h>
  40 #include <cassert>
  41
  42 #include "general.h"
  43
  44 #endif
  45
  46 /* === U S I N G =========================================================== */
  47
  48 using namespace std;
  49 using namespace etl;
  50 using namespace synfig;
  51
  52 /* === M A C R O S ========================================================= */
  53
  54 #define EPSILON         (1e-10)
  55
  56 /* === G L O B A L S ======================================================= */
  57
  58 /* === P R O C E D U R E S ================================================= */
  59
  60 /* === M E T H O D S ======================================================= */
  61
  62
  63 //Derivative Functions for numerical approximation
  64
  65 //bias == 0 will get F' at f3, bias < 0 will get F' at f1, and bias > 0 will get F' at f5
  66 template < class T >
  67 inline void FivePointdt(T &df, const T &f1, const T &f2, const T &f3, const T &f4, const T &f5, int bias)
  68 {
  69         if (bias == 0)                          // middle
  70                 df = (f1 - f2*8 + f4*8 - f5)*(1/12.0f);
  71         else if (bias < 0)                      // left
  72                 df = (-f1*25 + f2*48 - f3*36 + f4*16 - f5*3)*(1/12.0f);
  73         else                                            // right
  74                 df = (f1*3 - f2*16 + f3*36 - f4*48 + f5*25)*(1/12.0f);
  75 }
  76
  77 template < class T >
  78 inline void ThreePointdt(T &df, const T &f1, const T &f2, const T &f3, int bias)
  79 {
  80         if (bias == 0)                          // middle
  81                 df = (-f1 + f3)*(1/2.0f);
  82         else if (bias < 0)                      // left
  83                 df = (-f1*3 + f2*4 - f3)*(1/2.0f);
  84         else                                            // right
  85                 df = (f1 - f2*4 + f3*3)*(1/2.0f);
  86 }
  87
  88 // template < class T >
  89 // inline void ThreePointddt(T &df, const T &f1, const T &f2, const T &f3, int bias)
  90 // {
  91 //      // a 3 point approximation pretends to have constant acceleration,
  92 //      // so only one algorithm needed for left, middle, or right
  93 //      df = (f1 -f2*2 + f3)*(1/2.0f);
  94 // }
  95 //
  96 // // WARNING -- totally broken
  97 // template < class T >
  98 // inline void FivePointddt(T &df, const T &f1, const T &f2, const T &f3, int bias)
  99 // {
 100 //      if(bias == 0)
 101 //      {
 102 //              assert(0); // !?
 103 //              //middle
 104 //              //df = (- f1 + f2*16 - f3*30 +  f4*16 - f5)*(1/12.0f);
 105 //      }/*else if(bias < 0)
 106 //      {
 107 //              //left
 108 //              df = (f1*7 - f2*26*4 + f3*19*6 - f4*14*4 + f5*11)*(1/12.0f);
 109 //      }else
 110 //      {
 111 //              //right
 112 //              df = (f1*3 - f2*16 + f3*36 - f4*48 + f5*25)*(1/12.0f);
 113 //      }*/
 114 //      //side ones don't work, use 3 point
 115 // }
 116 //
 117 // //implement an arbitrary derivative
 118 // //dumb algorithm
 119 // template < class T >
 120 // void DerivativeApprox(T &df, const T f[], const Real t[], int npoints, int indexval)
 121 // {
 122 //      /*
 123 //      Lj(x) = PI_i!=j (x - xi) / PI_i!=j (xj - xi)
 124 //
 125 //      so Lj'(x) = SUM_k PI_i!=j|k (x - xi) / PI_i!=j (xj - xi)
 126 //      */
 127 //
 128 //      unsigned int i,j,k,i0,i1;
 129 //
 130 //      Real Lpj,mult,div,tj;
 131 //      Real tval = t[indexval];
 132 //
 133 //      //sum k
 134 //      for(j=0;j<npoints;++j)
 135 //      {
 136 //              Lpj = 0;
 137 //              div = 1;
 138 //              tj = t[j];
 139 //
 140 //              for(k=0;k<npoints;++k)
 141 //              {
 142 //                      if(k != j) //because there is no summand for k == j, since that term is missing from the original equation
 143 //                      {
 144 //                              //summation for k
 145 //                              for(i=0;i<npoints;++i)
 146 //                              {
 147 //                                      if(i != k)
 148 //                                      {
 149 //                                              mult *= tval - t[i];
 150 //                                      }
 151 //                              }
 152 //
 153 //                              Lpj += mult; //add into the summation
 154 //
 155 //                              //since the ks follow the exact pattern we need for the divisor (use that too)
 156 //                              div *= tj - t[k];
 157 //                      }
 158 //              }
 159 //
 160 //              //get the actual coefficient
 161 //              Lpj /= div;
 162 //
 163 //              //add it in to the equation
 164 //              df += f[j]*Lpj;
 165 //      }
 166 // }
 167
 168 //END numerical derivatives
 169
 170 // template < class T >
 171 // inline int sign(T f, T tol)
 172 // {
 173 //      if(f < -tol) return -1;
 174 //      if(f > tol) return 1;
 175 //      return 0;
 176 // }
 177
 178 void GetFirstDerivatives(const std::vector<synfig::Point> &f, unsigned int left, unsigned int right, char *out, unsigned int dfstride)
 179 {
 180         unsigned int current = left;
 181
 182         if(right - left < 2)
 183                 return;
 184         else if(right - left == 2)
 185         {
 186                 synfig::Vector v = f[left+1] - f[left];
 187
 188                 //set both to the one we want
 189                 *(synfig::Vector*)out = v;
 190                 out += dfstride;
 191                 *(synfig::Vector*)out = v;
 192                 out += dfstride;
 193         }
 194         else if(right - left < 6/*5*/) //should use 3 point
 195         {
 196                 //left then middle then right
 197                 ThreePointdt(*(synfig::Vector*)out,f[left+0], f[left+1], f[left+2], -1);
 198                 current++;
 199                 out += dfstride;
 200
 201                 for(;current < right-1; current++, out += dfstride)
 202                         ThreePointdt(*(synfig::Vector*)out,f[current-1], f[current], f[current+1], 0);
 203
 204                 ThreePointdt(*(synfig::Vector*)out,f[right-3], f[right-2], f[right-1], 1);
 205                 current++;
 206                 out += dfstride;
 207
 208         }else //can use 5 point
 209         {
 210                 //left 2 then middle bunch then right two
 211                 //may want to use 3 point for inner edge ones
 212
 213                 FivePointdt(*(synfig::Vector*)out,f[left+0], f[left+1], f[left+2], f[left+3], f[left+4], -2);
 214                 out += dfstride;
 215                 FivePointdt(*(synfig::Vector*)out,f[left+1], f[left+2], f[left+3], f[left+4], f[left+5], -1);
 216                 out += dfstride;
 217                 current += 2;
 218
 219                 for(;current < right-2; current++, out += dfstride)
 220                         FivePointdt(*(synfig::Vector*)out,f[current-2], f[current-1], f[current], f[current+1], f[current+2], 0);
 221
 222                 FivePointdt(*(synfig::Vector*)out,f[right-6], f[right-5], f[right-4], f[right-3], f[right-2], 2);
 223                 out += dfstride;
 224                 FivePointdt(*(synfig::Vector*)out,f[right-5], f[right-4], f[right-3], f[right-2], f[right-1], 1);
 225                 out += dfstride;
 226                 current += 2;
 227         }
 228 }
 229
 230 void GetSimpleDerivatives(const std::vector<synfig::Point> &f, int left, int right,
 231                                                         std::vector<synfig::Point> &df, int outleft,
 232                                                         const std::vector<synfig::Real> &/*di*/)
 233 {
 234         int i1,i2,i;
 235         int offset = 2; //df = 1/2 (f[i+o]-f[i-o])
 236
 237         assert((int)df.size() >= right-left+outleft); //must be big enough
 238
 239         for(i = left; i < right; ++i)
 240         {
 241                 //right now indices (figure out distance later)
 242                 i1 = std::max(left,i-offset);
 243                 i2 = std::max(left,i+offset);
 244
 245                 df[outleft++] = (f[i2] - f[i1])*0.5f;
 246         }
 247 }
 248
 249 //get the curve error from the double sample list of work points (hopefully that's enough)
 250 Real CurveError(const synfig::Point *pts, unsigned int n, std::vector<synfig::Point> &work, int left, int right)
 251 {
 252         if(right-left < 2) return -1;
 253
 254         int i,j;
 255
 256         //get distances to each point
 257         Real d,dtemp,dsum;
 258         //synfig::Vector v,vt;
 259         //synfig::Point p1,p2;
 260         synfig::Point pi;
 261         std::vector<synfig::Point>::const_iterator it;//,end = work.begin()+right;
 262
 263         //unsigned int size = work.size();
 264
 265         //for each line, get distance
 266         d = 0; //starts at 0
 267         for(i = 0; i < (int)n; ++i)
 268         {
 269                 pi = pts[i];
 270
 271                 dsum = FLT_MAX;
 272
 273                 it = work.begin()+left;
 274                 //p2 = *it++; //put it at left+1
 275                 for(j = left/*+1*/; j < right; ++j,++it)
 276                 {
 277                         /*p1 = p2;
 278                         p2 = *it;
 279
 280                         v = p2 - p1;
 281                         vt = pi - p1;
 282
 283                         dtemp = v.mag_squared() > 1e-12 ? (vt*v)/v.mag_squared() : 0; //get the projected time value for the current line
 284
 285                         //get distance to line segment with the time value clamped 0-1
 286                         if(dtemp >= 1)  //use p+v
 287                         {
 288                                 vt += v; //makes it pp - (p+v)
 289                         }else if(dtemp > 0)     //use vt-proj
 290                         {
 291                                 vt -= v*dtemp; // vt - proj_v(vt)       //must normalize the projection vector to work
 292                         }
 293
 294                         //else use p
 295                         dtemp = vt.mag_squared();*/
 296
 297                         dtemp = (pi - *it).mag_squared();
 298                         if(dtemp < dsum)
 299                                 dsum = dtemp;
 300                 }
 301
 302                 //accumulate the points' min distance from the curve
 303                 d += sqrt(dsum);
 304         }
 305
 306         return d;
 307 }
 308
 309 typedef synfigapp::BLineConverter::cpindex cpindex;
 310
 311 //has the index data and the tangent scale data (relevant as it may be)
 312 int tessellate_curves(const std::vector<cpindex> &inds, const std::vector<Point> &f, const std::vector<synfig::Vector> &df, std::vector<Point> &work)
 313 {
 314         if(inds.size() < 2)
 315                 return 0;
 316
 317         etl::hermite<Point>     curve;
 318         int ntess = 0;
 319
 320         std::vector<cpindex>::const_iterator j = inds.begin(),j2, end = inds.end();
 321
 322         unsigned int ibase = inds[0].curind;
 323
 324         j2 = j++;
 325         for(; j != end; j2 = j++)
 326         {
 327                 //if this curve has invalid error (in j) then retessellate its work points (requires reparametrization, etc.)
 328                 if(j->error < 0)
 329                 {
 330                         //get the stepsize etc. for the number of points in here
 331                         unsigned int n = j->curind - j2->curind + 1; //thats the number of points in the span
 332                         unsigned int k, kend, i0, i3;
 333                         //so reset the right chunk
 334
 335                         Real t, dt = 1/(Real)(n*2-2); //assuming that they own only n points
 336
 337                         //start at first intermediate
 338                         t = 0;
 339
 340                         i0 = j2->curind; i3 = j->curind;
 341                         k = (i0-ibase)*2; //start on first intermediary point (2x+1)
 342                         kend = (i3-ibase)*2; //last point to set (not intermediary)
 343
 344                         //build hermite curve, it's easier
 345                         curve.p1() = f[i0];
 346                         curve.p2() = f[i3];
 347                         curve.t1() = df[i0-ibase] * (df[i0-ibase].mag_squared() > 1e-4
 348                                                                                  ? j2->tangentscale/df[i0-ibase].mag()
 349                                                                                  : j2->tangentscale);
 350                         curve.t2() = df[i3-ibase] * (df[i3-ibase].mag_squared() > 1e-4
 351                                                                                  ? j->tangentscale/df[i3-ibase].mag()
 352                                                                                  : j->tangentscale);
 353                         curve.sync();
 354
 355                         //MUST include the end point (since we are ignoring left one)
 356                         for(; k < kend; ++k, t += dt)
 357                         {
 358                                 work[k] = curve(t);
 359                         }
 360
 361                         work[k] = curve(1); //k == kend, t == 1 -> c(t) == p2
 362                         ++ntess;
 363                 }
 364         }
 365
 366         return ntess;
 367 }
 368
 369 synfigapp::BLineConverter::BLineConverter()
 370 {
 371         pixelwidth = 1;
 372         smoothness = 0.70f;
 373         width = 0;
 374 };
 375
 376 void
 377 synfigapp::BLineConverter::clear()
 378 {
 379         f.clear();
 380         f_w.clear();
 381         ftemp.clear();
 382         df.clear();
 383         cvt.clear();
 384         brk.clear();
 385         di.clear();
 386         d_i.clear();
 387         work.clear();
 388         curind.clear();
 389 }
 390
 391 void
 392 synfigapp::BLineConverter::operator () (std::list<synfig::BLinePoint> &out, const std::list<synfig::Point> &in,const std::list<synfig::Real> &in_w)
 393 {
 394         //Profiling information
 395         /*etl::clock::value_type initialprocess=0, curveval=0, breakeval=0, disteval=0;
 396         etl::clock::value_type preproceval=0, tesseval=0, erroreval=0, spliteval=0;
 397         unsigned int                    numpre=0, numtess=0, numerror=0, numsplit=0;
 398         etl::clock_realtime timer,total;*/
 399
 400         //total.reset();
 401         if(in.size()<=1)
 402                 return;
 403
 404         clear();
 405
 406         //removing digitization error harder than expected
 407
 408         //intended to fix little pen errors caused by the way people draw (tiny juts in opposite direction)
 409         //Different solutions
 410         //      Average at both end points (will probably eliminate many points at each end of the samples)
 411         //      Average after the break points are found (weird points would still affect the curve)
 412         //      Just always get rid of breaks at the beginning and end if they are a certain distance apart
 413         //              This is will be current approach so all we do now is try to remove duplicate points
 414
 415         //remove duplicate points - very bad for fitting
 416
 417         //timer.reset();
 418
 419         {
 420                 std::list<synfig::Point>::const_iterator i = in.begin(), end = in.end();
 421                 std::list<synfig::Real>::const_iterator iw = in_w.begin();
 422                 synfig::Point   c;
 423
 424                 if(in.size() == in_w.size())
 425                 {
 426                         for(;i != end; ++i,++iw)
 427                                 if(*i != c)             // eliminate duplicate points
 428                                 {
 429                                         f.push_back(c = *i);
 430                                         f_w.push_back(*iw);
 431                                 }
 432                 }
 433                 else
 434                 {
 435                         for(;i != end; ++i)
 436                                 if(*i != c)             // eliminate duplicate points
 437                                         f.push_back(c = *i);
 438                 }
 439         }
 440         //initialprocess = timer();
 441
 442         if(f.size()<=6)
 443                 return;
 444
 445         //get curvature information
 446         //timer.reset();
 447
 448         {
 449                 int i,i0,i1;
 450                 synfig::Vector v1,v2;
 451
 452                 cvt.resize(f.size());
 453
 454                 cvt.front() = 1;
 455                 cvt.back() = 1;
 456
 457                 for(i = 1; i < (int)f.size()-1; ++i)
 458                 {
 459                         i0 = std::max(0,i - 2);
 460                         i1 = std::min((int)(f.size()-1),i + 2);
 461
 462                         v1 = f[i] - f[i0];
 463                         v2 = f[i1] - f[i];
 464
 465                         cvt[i] = (v1*v2)/(v1.mag()*v2.mag());
 466                 }
 467         }
 468
 469         //curveval = timer();
 470         //synfig::info("calculated curvature");
 471
 472         //find corner points and interpolate inside those
 473         //timer.reset();
 474         {
 475                 //break at sharp derivative points
 476                 //TODO tolerance should be set based upon digitization resolution (length dependent index selection)
 477                 Real    tol = 0;                //break tolerance, for the cosine of the change in angle (really high curvature or something)
 478                 Real    fixdistsq = 4*width*width; //the distance to ignore breaks at the end points (for fixing stuff)
 479                 unsigned int i = 0;
 480
 481                 int             maxi = -1, last=0;
 482                 Real    minc = 1;
 483
 484                 brk.push_back(0);
 485
 486                 for(i = 1; i < cvt.size()-1; ++i)
 487                 {
 488                         //insert if too sharp (we need to break the tangents to insert onto the break list)
 489
 490                         if(cvt[i] < tol)
 491                         {
 492                                 if(cvt[i] < minc)
 493                                 {
 494                                         minc = cvt[i];
 495                                         maxi = i;
 496                                 }
 497                         }
 498                         else if(maxi >= 0)
 499                         {
 500                                 if(maxi >= last + 8)
 501                                 {
 502                                         //synfig::info("break: %d-%d",maxi+1,cvt.size());
 503                                         brk.push_back(maxi);
 504                                         last = maxi;
 505                                 }
 506                                 maxi = -1;
 507                                 minc = 1;
 508                         }
 509                 }
 510
 511                 brk.push_back(i);
 512
 513                 //postprocess for breaks too close to each other
 514                 Real d = 0;
 515                 Point p = f[brk.front()];
 516
 517                 //first set
 518                 for(i = 1; i < brk.size()-1; ++i) //do not want to include end point...
 519                 {
 520                         d = (f[brk[i]] - p).mag_squared();
 521                         if(d > fixdistsq) break; //don't want to group breaks if we ever get over the dist...
 522                 }
 523                 //want to erase all points before...
 524                 if(i != 1)
 525                         brk.erase(brk.begin(),brk.begin()+i-1);
 526
 527                 //end set
 528                 p = f[brk.back()];
 529                 for(i = brk.size()-2; i > 0; --i) //start at one in from the end
 530                 {
 531                         d = (f[brk[i]] - p).mag_squared();
 532                         if(d > fixdistsq) break; //don't want to group breaks if we ever get over the dist
 533                 }
 534                 if(i != brk.size()-2)
 535                         brk.erase(brk.begin()+i+2,brk.end()); //erase all points that we found... found none if i has not advanced
 536                 //must not include the one we ended up on
 537         }
 538         //breakeval = timer();
 539         //synfig::info("found break points: %d",brk.size());
 540
 541         //get the distance calculation of the entire curve (for tangent scaling)
 542
 543         //timer.reset();
 544         {
 545                 synfig::Point p1,p2;
 546
 547                 p1=p2=f[0];
 548
 549                 di.resize(f.size()); d_i.resize(f.size());
 550                 Real d = 0;
 551                 for(unsigned int i = 0; i < f.size();)
 552                 {
 553                         d += (d_i[i] = (p2-p1).mag());
 554                         di[i] = d;
 555
 556                         p1=p2;
 557                         p2=f[++i];
 558                 }
 559         }
 560         //disteval = timer();
 561         //synfig::info("calculated distance");
 562
 563         //now break at every point - calculate new derivatives each time
 564
 565         //TODO
 566         //must be sure that the break points are 3 or more apart
 567         //then must also store the breaks which are not smooth, etc.
 568         //and figure out tangents between there
 569
 570         //for each pair of break points (as long as they are far enough apart) recursively subdivide stuff
 571         //ignore the detected intermediate points
 572         {
 573                 unsigned int i0=0,i3=0,is=0;
 574                 int i=0,j=0;
 575
 576                 bool done = false;
 577
 578                 Real errortol = smoothness*pixelwidth; //???? what the hell should this value be
 579
 580                 BLinePoint a;
 581                 synfig::Vector v;
 582
 583                 //intemp = f; //don't want to smooth out the corners
 584
 585                 bool breaktan = false, setwidth;
 586                 a.set_split_tangent_flag(false);
 587                 //a.set_width(width);
 588                 a.set_width(1.0f);
 589
 590                 setwidth = (f.size() == f_w.size());
 591
 592                 for(j = 0; j < (int)brk.size() - 1; ++j)
 593                 {
 594                         //for b[j] to b[j+1] subdivide and stuff
 595                         i0 = brk[j];
 596                         i3 = brk[j+1];
 597
 598                         unsigned int size = i3-i0+1; //must include the end points
 599
 600                         //new derivatives
 601                         //timer.reset();
 602                         ftemp.assign(f.begin()+i0, f.begin()+i3+1);
 603                         for(i=0;i<20;++i)
 604                                 gaussian_blur_3(ftemp.begin(),ftemp.end(),false);
 605
 606                         df.resize(size);
 607
 608                         // Wondering whether the modification of the df vector
 609                         // using a char* pointer and pointer arithmetric was safe,
 610                         // I looked it up...
 611                         //
 612                         // http://www.open-std.org/jtc1/sc22/wg21/docs/papers/2007/n2369.pdf tells me:
 613                         //
 614                         //      23.2.5  Class template vector [vector]
 615                         //
 616                         //      [...] The elements of a vector are stored contiguously,
 617                         //      meaning that if v is a vector<T,Allocator> where T is
 618                         //      some type other than bool, then it obeys the identity
 619                         //      &v[n] == &v[0] + n for all 0 <= n < v.size().
 620                         //
 621                         GetFirstDerivatives(ftemp,0,size,(char*)&df[0],sizeof(df[0]));
 622
 623                         //GetSimpleDerivatives(ftemp,0,size,df,0,di);
 624                         //< don't have to worry about indexing stuff as it is all being taken care of right now
 625                         //preproceval += timer();
 626                         //numpre++;
 627
 628                         work.resize(size*2-1); //guarantee that all points will be tessellated correctly (one point in between every 2 adjacent points)
 629
 630                         //if size of work is size*2-1, the step size should be 1/(size*2 - 2)
 631                         //Real step = 1/(Real)(size*2 - 1);
 632
 633                         //start off with break points as indices
 634                         curind.clear();
 635                         curind.push_back(cpindex(i0,di[i3]-di[i0],0)); //0 error because no curve on the left
 636                         curind.push_back(cpindex(i3,di[i3]-di[i0],-1)); //error needs to be reevaluated
 637                         done = false; //we want to loop
 638
 639                         unsigned int dcount = 0;
 640
 641                         //while there are still enough points between us, and the error is too high subdivide (and invalidate neighbors that share tangents)
 642                         while(!done)
 643                         {
 644                                 //tessellate all curves with invalid error values
 645                                 work[0] = f[i0];
 646
 647                                 //timer.reset();
 648                                 /*numtess += */tessellate_curves(curind,f,df,work);
 649                                 //tesseval += timer();
 650
 651                                 //now get all error values
 652                                 //timer.reset();
 653                                 for(i = 1; i < (int)curind.size(); ++i)
 654                                 {
 655                                         if(curind[i].error < 0) //must have been retessellated, so now recalculate error value
 656                                         {
 657                                                 //evaluate error from points (starting at current index)
 658                                                 int size = curind[i].curind - curind[i-1].curind + 1;
 659                                                 curind[i].error = CurveError(&f[curind[i-1].curind], size,
 660                                                                                                          work,(curind[i-1].curind - i0)*2,(curind[i].curind - i0)*2+1);
 661
 662                                                 /*if(curind[i].error > 1.0e5)
 663                                                 {
 664                                                         synfig::info("Holy crap %d-%d error %f",curind[i-1].curind,curind[i].curind,curind[i].error);
 665                                                         curind[i].error = -1;
 666                                                         numtess += tessellate_curves(curind,f,df,work);
 667                                                         curind[i].error = CurveError(&f[curind[i-1].curind], size,
 668                                                                                                          work,0,work.size());//(curind[i-1].curind - i0)*2,(curind[i].curind - i0)*2+1);
 669                                                 }*/
 670                                                 //numerror++;
 671                                         }
 672                                 }
 673                                 //erroreval += timer();
 674
 675                                 //assume we're done
 676                                 done = true;
 677
 678                                 //check each error to see if it's too big, if so, then subdivide etc.
 679                                 int indsize = (int)curind.size();
 680                                 Real maxrelerror = 0;
 681                                 int maxi = -1;//, numpoints;
 682
 683                                 //timer.reset();
 684                                 //get the maximum error and split there
 685                                 for(i = 1; i < indsize; ++i)
 686                                 {
 687                                         //numpoints = curind[i].curind - curind[i-1].curind + 1;
 688
 689                                         if(curind[i].error > maxrelerror && curind[i].curind - curind[i-1].curind > 2) //only accept if it's valid
 690                                         {
 691                                                 maxrelerror = curind[i].error;
 692                                                 maxi = i;
 693                                         }
 694                                 }
 695
 696                                 //split if error is too great
 697                                 if(maxrelerror > errortol)
 698                                 {
 699                                         //add one to the left etc
 700                                         unsigned int    ibase = curind[maxi-1].curind, itop = curind[maxi].curind,
 701                                                                         ibreak = (ibase + itop)/2;
 702                                         Real scale, scale2;
 703
 704                                         assert(ibreak < f.size());
 705
 706                                         //synfig::info("Split %d -%d- %d, error: %f", ibase,ibreak,itop,maxrelerror);
 707
 708                                         if(ibase != itop)
 709                                         {
 710                                                 //invalidate current error of the changed tangents and add an extra segment
 711                                                 //enforce minimum tangents property
 712                                                 curind[maxi].error = -1;
 713                                                 curind[maxi-1].error = -1;
 714                                                 if(maxi+1 < indsize) curind[maxi+1].error = -1; //if there is a curve segment beyond this it will be effected as well
 715
 716                                                 scale = di[itop] - di[ibreak];
 717                                                 scale2 = maxi+1 < indsize ? di[curind[maxi+1].curind] - di[itop] : scale; //to the right valid?
 718                                                 curind[maxi].tangentscale = std::min(scale, scale2);
 719
 720                                                 scale = di[ibreak] - di[ibase];
 721                                                 scale2 = maxi >= 2 ? di[ibase] - di[curind[maxi-2].curind] : scale; // to the left valid -2 ?
 722                                                 curind[maxi-1].tangentscale = std::min(scale, scale2);
 723
 724                                                 scale = std::min(di[ibreak] - di[ibase], di[itop] - di[ibreak]);
 725
 726                                                 curind.insert(curind.begin()+maxi,cpindex(ibreak, scale, -1));
 727                                                 //curind.push_back(cpindex(ibreak, scale, -1));
 728                                                 //std::sort(curind.begin(), curind.end());
 729
 730                                                 done = false;
 731                                                 //numsplit++;
 732                                         }
 733                                 }
 734                                 //spliteval += timer();
 735
 736                                 dcount++;
 737                         }
 738
 739                         //insert the last point too (just set tangent for now
 740                         is = curind[0].curind;
 741
 742                         //first point inherits current tangent status
 743                         v = df[is - i0];
 744                         if(v.mag_squared() > EPSILON)
 745                                 v *= (curind[0].tangentscale/v.mag());
 746
 747                         if(!breaktan)
 748                                 a.set_tangent(v);
 749                         else a.set_tangent2(v);
 750
 751                         a.set_vertex(f[is]);
 752                         if(setwidth)a.set_width(f_w[is]);
 753
 754                         out.push_back(a);
 755                         a.set_split_tangent_flag(false); //won't need to break anymore
 756                         breaktan = false;
 757
 758                         for(i = 1; i < (int)curind.size()-1; ++i)
 759                         {
 760                                 is = curind[i].curind;
 761
 762                                 //first point inherits current tangent status
 763                                 v = df[is-i0];
 764                                 if(v.mag_squared() > EPSILON)
 765                                         v *= (curind[i].tangentscale/v.mag());
 766
 767                                 a.set_tangent(v); // always inside, so guaranteed to be smooth
 768                                 a.set_vertex(f[is]);
 769                                 if(setwidth)a.set_width(f_w[is]);
 770
 771                                 out.push_back(a);
 772                         }
 773
 774                         //set the last point's data
 775                         is = curind.back().curind; //should already be this
 776
 777                         v = df[is-i0];
 778                         if(v.mag_squared() > EPSILON)
 779                                 v *= (curind.back().tangentscale/v.mag());
 780
 781                         a.set_tangent1(v);
 782                         a.set_split_tangent_flag(true);
 783                         breaktan = true;
 784
 785                         //will get the vertex and tangent 2 from next round
 786                 }
 787
 788                 a.set_vertex(f[i3]);
 789                 a.set_split_tangent_flag(false);
 790                 if(setwidth)
 791                         a.set_width(f_w[i3]);
 792                 out.push_back(a);
 793
 794                 /*etl::clock::value_type totaltime = total(),
 795                                                            misctime = totaltime - initialprocess - curveval - breakeval - disteval
 796                                                                           - preproceval - tesseval - erroreval - spliteval;
 797
 798                 synfig::info(
 799                         "Curve Convert Profile:\n"
 800                         "\tInitial Preprocess:    %f\n"
 801                         "\tCurvature Calculation: %f\n"
 802                         "\tBreak Calculation:     %f\n"
 803                         "\tDistance Calculation:  %f\n"
 804                         "  Algorithm: (numtimes,totaltime)\n"
 805                         "\tPreprocess step:      (%d,%f)\n"
 806                         "\tTessellation step:    (%d,%f)\n"
 807                         "\tError step:           (%d,%f)\n"
 808                         "\tSplit step:           (%d,%f)\n"
 809                         "  Num Input: %d, Num Output: %d\n"
 810                         "  Total time: %f, Misc time: %f\n",
 811                         initialprocess, curveval,breakeval,disteval,
 812                         numpre,preproceval,numtess,tesseval,numerror,erroreval,numsplit,spliteval,
 813                         in.size(),out.size(),
 814                         totaltime,misctime);*/
 815
 816                 return;
 817         }
 818 }
 819
 820 void synfigapp::BLineConverter::EnforceMinWidth(std::list<synfig::BLinePoint> &bline, synfig::Real min_pressure)
 821 {
 822         std::list<synfig::BLinePoint>::iterator i = bline.begin(),
 823                                                                                         end = bline.end();
 824
 825         for(i = bline.begin(); i != end; ++i)
 826                 if(i->get_width() < min_pressure)
 827                         i->set_width(min_pressure);
 828 }