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