KiCad PCB EDA Suite
trace.cpp File Reference
#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include "auxiliary.h"
#include "curve.h"
#include "lists.h"
#include "potracelib.h"
#include "progress.h"
#include "trace.h"

Go to the source code of this file.

Classes

struct  opti_s
 

Macros

#define INFTY
 
#define COS179   -0.999847695156 /* the cosine of 179 degrees */
 
#define SAFE_CALLOC(var, n, typ)
 
#define TRY(x)
 

Typedefs

typedef double quadform_t[3][3]
 
typedef struct opti_s opti_t
 

Functions

static point_t dorth_infty (dpoint_t p0, dpoint_t p2)
 
static double dpara (dpoint_t p0, dpoint_t p1, dpoint_t p2)
 
static double ddenom (dpoint_t p0, dpoint_t p2)
 
static int cyclic (int a, int b, int c)
 
static void pointslope (privpath_t *pp, int i, int j, dpoint_t *ctr, dpoint_t *dir)
 
static double quadform (quadform_t Q, dpoint_t w)
 
static int xprod (point_t p1, point_t p2)
 
static double cprod (dpoint_t p0, dpoint_t p1, dpoint_t p2, dpoint_t p3)
 
static double iprod (dpoint_t p0, dpoint_t p1, dpoint_t p2)
 
static double iprod1 (dpoint_t p0, dpoint_t p1, dpoint_t p2, dpoint_t p3)
 
static double ddist (dpoint_t p, dpoint_t q)
 
static dpoint_t bezier (double t, dpoint_t p0, dpoint_t p1, dpoint_t p2, dpoint_t p3)
 
static double tangent (dpoint_t p0, dpoint_t p1, dpoint_t p2, dpoint_t p3, dpoint_t q0, dpoint_t q1)
 
static int calc_sums (privpath_t *pp)
 
static int calc_lon (privpath_t *pp)
 
static double penalty3 (privpath_t *pp, int i, int j)
 
static int bestpolygon (privpath_t *pp)
 
static int adjust_vertices (privpath_t *pp)
 
static void reverse (privcurve_t *curve)
 
static void smooth (privcurve_t *curve, double alphamax)
 
static int opti_penalty (privpath_t *pp, int i, int j, opti_t *res, double opttolerance, int *convc, double *areac)
 
static int opticurve (privpath_t *pp, double opttolerance)
 
int process_path (path_t *plist, const potrace_param_t *param, progress_t *progress)
 

Macro Definition Documentation

#define COS179   -0.999847695156 /* the cosine of 179 degrees */

Definition at line 26 of file trace.cpp.

Referenced by opti_penalty().

#define INFTY
Value:
10000000 /* it suffices that this is longer than any
* path; it need not be really infinite */

Definition at line 23 of file trace.cpp.

Referenced by calc_lon().

#define SAFE_CALLOC (   var,
  n,
  typ 
)
Value:
if( ( var = (typ*) calloc( n, sizeof( typ ) ) ) == NULL ) \
goto calloc_error

Definition at line 29 of file trace.cpp.

Referenced by adjust_vertices(), bestpolygon(), calc_lon(), calc_sums(), and opticurve().

#define TRY (   x)
Value:
if( x ) \
goto try_error

Definition at line 1494 of file trace.cpp.

Referenced by process_path().

Typedef Documentation

typedef struct opti_s opti_t

Definition at line 1117 of file trace.cpp.

typedef double quadform_t[3][3]

Definition at line 177 of file trace.cpp.

Function Documentation

static int adjust_vertices ( privpath_t pp)
static

Definition at line 772 of file trace.cpp.

References potrace_privpath_s::curve, potrace_privpath_s::len, potrace_privpath_s::m, min, mod, potrace_privpath_s::po, pointslope(), privcurve_init(), potrace_privpath_s::pt, quadform(), SAFE_CALLOC, sq, privcurve_s::vertex, point_s::x, potrace_dpoint_s::x, potrace_privpath_s::x0, point_s::y, potrace_dpoint_s::y, and potrace_privpath_s::y0.

Referenced by process_path().

773 {
774  int m = pp->m;
775  int* po = pp->po;
776  int n = pp->len;
777  point_t* pt = pp->pt;
778  int x0 = pp->x0;
779  int y0 = pp->y0;
780 
781  dpoint_t* ctr = NULL; /* ctr[m] */
782  dpoint_t* dir = NULL; /* dir[m] */
783  quadform_t* q = NULL; /* q[m] */
784  double v[3];
785  double d;
786  int i, j, k, l;
787  dpoint_t s;
788  int r;
789 
790  SAFE_CALLOC( ctr, m, dpoint_t );
791  SAFE_CALLOC( dir, m, dpoint_t );
792  SAFE_CALLOC( q, m, quadform_t );
793 
794  r = privcurve_init( &pp->curve, m );
795 
796  if( r )
797  {
798  goto calloc_error;
799  }
800 
801  /* calculate "optimal" point-slope representation for each line
802  * segment */
803  for( i = 0; i < m; i++ )
804  {
805  j = po[mod( i + 1, m )];
806  j = mod( j - po[i], n ) + po[i];
807  pointslope( pp, po[i], j, &ctr[i], &dir[i] );
808  }
809 
810  /* represent each line segment as a singular quadratic form; the
811  * distance of a point (x,y) from the line segment will be
812  * (x,y,1)Q(x,y,1)^t, where Q=q[i]. */
813  for( i = 0; i < m; i++ )
814  {
815  d = sq( dir[i].x ) + sq( dir[i].y );
816 
817  if( d == 0.0 )
818  {
819  for( j = 0; j < 3; j++ )
820  {
821  for( k = 0; k < 3; k++ )
822  {
823  q[i][j][k] = 0;
824  }
825  }
826  }
827  else
828  {
829  v[0] = dir[i].y;
830  v[1] = -dir[i].x;
831  v[2] = -v[1] * ctr[i].y - v[0] * ctr[i].x;
832 
833  for( l = 0; l < 3; l++ )
834  {
835  for( k = 0; k < 3; k++ )
836  {
837  q[i][l][k] = v[l] * v[k] / d;
838  }
839  }
840  }
841  }
842 
843  /* now calculate the "intersections" of consecutive segments.
844  * Instead of using the actual intersection, we find the point
845  * within a given unit square which minimizes the square distance to
846  * the two lines. */
847  for( i = 0; i < m; i++ )
848  {
849  quadform_t Q;
850  dpoint_t w;
851  double dx, dy;
852  double det;
853  double min, cand; /* minimum and candidate for minimum of quad. form */
854  double xmin, ymin; /* coordinates of minimum */
855  int z;
856 
857  /* let s be the vertex, in coordinates relative to x0/y0 */
858  s.x = pt[po[i]].x - x0;
859  s.y = pt[po[i]].y - y0;
860 
861  /* intersect segments i-1 and i */
862 
863  j = mod( i - 1, m );
864 
865  /* add quadratic forms */
866  for( l = 0; l < 3; l++ )
867  {
868  for( k = 0; k < 3; k++ )
869  {
870  Q[l][k] = q[j][l][k] + q[i][l][k];
871  }
872  }
873 
874  while( 1 )
875  {
876 /* minimize the quadratic form Q on the unit square */
877 /* find intersection */
878 
879 #ifdef HAVE_GCC_LOOP_BUG
880  /* work around gcc bug #12243 */
881  free( NULL );
882 #endif
883 
884  det = Q[0][0] * Q[1][1] - Q[0][1] * Q[1][0];
885 
886  if( det != 0.0 )
887  {
888  w.x = ( -Q[0][2] * Q[1][1] + Q[1][2] * Q[0][1] ) / det;
889  w.y = ( Q[0][2] * Q[1][0] - Q[1][2] * Q[0][0] ) / det;
890  break;
891  }
892 
893  /* matrix is singular - lines are parallel. Add another,
894  * orthogonal axis, through the center of the unit square */
895  if( Q[0][0] > Q[1][1] )
896  {
897  v[0] = -Q[0][1];
898  v[1] = Q[0][0];
899  }
900  else if( Q[1][1] )
901  {
902  v[0] = -Q[1][1];
903  v[1] = Q[1][0];
904  }
905  else
906  {
907  v[0] = 1;
908  v[1] = 0;
909  }
910 
911  d = sq( v[0] ) + sq( v[1] );
912  v[2] = -v[1] * s.y - v[0] * s.x;
913 
914  for( l = 0; l < 3; l++ )
915  {
916  for( k = 0; k < 3; k++ )
917  {
918  Q[l][k] += v[l] * v[k] / d;
919  }
920  }
921  }
922 
923  dx = fabs( w.x - s.x );
924  dy = fabs( w.y - s.y );
925 
926  if( dx <= .5 && dy <= .5 )
927  {
928  pp->curve.vertex[i].x = w.x + x0;
929  pp->curve.vertex[i].y = w.y + y0;
930  continue;
931  }
932 
933  /* the minimum was not in the unit square; now minimize quadratic
934  * on boundary of square */
935  min = quadform( Q, s );
936  xmin = s.x;
937  ymin = s.y;
938 
939  if( Q[0][0] == 0.0 )
940  {
941  goto fixx;
942  }
943 
944  for( z = 0; z < 2; z++ )
945  {
946  /* value of the y-coordinate */
947  w.y = s.y - 0.5 + z;
948  w.x = -( Q[0][1] * w.y + Q[0][2] ) / Q[0][0];
949  dx = fabs( w.x - s.x );
950  cand = quadform( Q, w );
951 
952  if( dx <= .5 && cand < min )
953  {
954  min = cand;
955  xmin = w.x;
956  ymin = w.y;
957  }
958  }
959 
960 fixx:
961 
962  if( Q[1][1] == 0.0 )
963  {
964  goto corners;
965  }
966 
967  for( z = 0; z < 2; z++ )
968  {
969  /* value of the x-coordinate */
970  w.x = s.x - 0.5 + z;
971  w.y = -( Q[1][0] * w.x + Q[1][2] ) / Q[1][1];
972  dy = fabs( w.y - s.y );
973  cand = quadform( Q, w );
974 
975  if( dy <= .5 && cand < min )
976  {
977  min = cand;
978  xmin = w.x;
979  ymin = w.y;
980  }
981  }
982 
983 corners:
984 
985  /* check four corners */
986  for( l = 0; l < 2; l++ )
987  {
988  for( k = 0; k < 2; k++ )
989  {
990  w.x = s.x - 0.5 + l;
991  w.y = s.y - 0.5 + k;
992  cand = quadform( Q, w );
993 
994  if( cand < min )
995  {
996  min = cand;
997  xmin = w.x;
998  ymin = w.y;
999  }
1000  }
1001  }
1002 
1003  pp->curve.vertex[i].x = xmin + x0;
1004  pp->curve.vertex[i].y = ymin + y0;
1005  continue;
1006  }
1007 
1008  free( ctr );
1009  free( dir );
1010  free( q );
1011  return 0;
1012 
1013 calloc_error:
1014  free( ctr );
1015  free( dir );
1016  free( q );
1017  return 1;
1018 }
long y
Definition: auxiliary.h:25
int privcurve_init(privcurve_t *curve, int n)
Definition: curve.cpp:98
static void pointslope(privpath_t *pp, int i, int j, dpoint_t *ctr, dpoint_t *dir)
Definition: trace.cpp:89
#define sq(a)
Definition: auxiliary.h:87
privcurve_t curve
Definition: curve.h:62
dpoint_t * vertex
Definition: curve.h:28
long x
Definition: auxiliary.h:24
point_t * pt
Definition: curve.h:53
#define SAFE_CALLOC(var, n, typ)
Definition: trace.cpp:29
static double quadform(quadform_t Q, dpoint_t w)
Definition: trace.cpp:180
#define mod(a, n)
Definition: greymap.cpp:24
double quadform_t[3][3]
Definition: trace.cpp:177
#define min(a, b)
Definition: auxiliary.h:85
static int bestpolygon ( privpath_t pp)
static

Definition at line 631 of file trace.cpp.

References potrace_privpath_s::len, potrace_privpath_s::lon, potrace_privpath_s::m, mod, penalty3(), potrace_privpath_s::po, and SAFE_CALLOC.

Referenced by process_path().

632 {
633  int i, j, m, k;
634  int n = pp->len;
635  double* pen = NULL; /* pen[n+1]: penalty vector */
636  int* prev = NULL; /* prev[n+1]: best path pointer vector */
637  int* clip0 = NULL; /* clip0[n]: longest segment pointer, non-cyclic */
638  int* clip1 = NULL; /* clip1[n+1]: backwards segment pointer, non-cyclic */
639  int* seg0 = NULL; /* seg0[m+1]: forward segment bounds, m<=n */
640  int* seg1 = NULL; /* seg1[m+1]: backward segment bounds, m<=n */
641  double thispen;
642  double best;
643  int c;
644 
645  SAFE_CALLOC( pen, n + 1, double );
646  SAFE_CALLOC( prev, n + 1, int );
647  SAFE_CALLOC( clip0, n, int );
648  SAFE_CALLOC( clip1, n + 1, int );
649  SAFE_CALLOC( seg0, n + 1, int );
650  SAFE_CALLOC( seg1, n + 1, int );
651 
652  /* calculate clipped paths */
653  for( i = 0; i < n; i++ )
654  {
655  c = mod( pp->lon[mod( i - 1, n )] - 1, n );
656 
657  if( c == i )
658  {
659  c = mod( i + 1, n );
660  }
661 
662  if( c < i )
663  {
664  clip0[i] = n;
665  }
666  else
667  {
668  clip0[i] = c;
669  }
670  }
671 
672  /* calculate backwards path clipping, non-cyclic. j <= clip0[i] iff
673  * clip1[j] <= i, for i,j=0..n. */
674  j = 1;
675 
676  for( i = 0; i < n; i++ )
677  {
678  while( j <= clip0[i] )
679  {
680  clip1[j] = i;
681  j++;
682  }
683  }
684 
685  /* calculate seg0[j] = longest path from 0 with j segments */
686  i = 0;
687 
688  for( j = 0; i < n; j++ )
689  {
690  seg0[j] = i;
691  i = clip0[i];
692  }
693 
694  seg0[j] = n;
695  m = j;
696 
697  /* calculate seg1[j] = longest path to n with m-j segments */
698  i = n;
699 
700  for( j = m; j > 0; j-- )
701  {
702  seg1[j] = i;
703  i = clip1[i];
704  }
705 
706  seg1[0] = 0;
707 
708  /* now find the shortest path with m segments, based on penalty3 */
709  /* note: the outer 2 loops jointly have at most n iterations, thus
710  * the worst-case behavior here is quadratic. In practice, it is
711  * close to linear since the inner loop tends to be short. */
712  pen[0] = 0;
713 
714  for( j = 1; j <= m; j++ )
715  {
716  for( i = seg1[j]; i <= seg0[j]; i++ )
717  {
718  best = -1;
719 
720  for( k = seg0[j - 1]; k >= clip1[i]; k-- )
721  {
722  thispen = penalty3( pp, k, i ) + pen[k];
723 
724  if( best < 0 || thispen < best )
725  {
726  prev[i] = k;
727  best = thispen;
728  }
729  }
730 
731  pen[i] = best;
732  }
733  }
734 
735  pp->m = m;
736  SAFE_CALLOC( pp->po, m, int );
737 
738  /* read off shortest path */
739  for( i = n, j = m - 1; i > 0; j-- )
740  {
741  i = prev[i];
742  pp->po[j] = i;
743  }
744 
745  free( pen );
746  free( prev );
747  free( clip0 );
748  free( clip1 );
749  free( seg0 );
750  free( seg1 );
751  return 0;
752 
753 calloc_error:
754  free( pen );
755  free( prev );
756  free( clip0 );
757  free( clip1 );
758  free( seg0 );
759  free( seg1 );
760  return 1;
761 }
static double penalty3(privpath_t *pp, int i, int j)
Definition: trace.cpp:573
#define SAFE_CALLOC(var, n, typ)
Definition: trace.cpp:29
#define mod(a, n)
Definition: greymap.cpp:24
static dpoint_t bezier ( double  t,
dpoint_t  p0,
dpoint_t  p1,
dpoint_t  p2,
dpoint_t  p3 
)
inlinestatic

Definition at line 260 of file trace.cpp.

References potrace_dpoint_s::x, and potrace_dpoint_s::y.

Referenced by SCH_LEGACY_PLUGIN_CACHE::loadBezier(), and opti_penalty().

261 {
262  double s = 1 - t;
263  dpoint_t res;
264 
265  /* Note: a good optimizing compiler (such as gcc-3) reduces the
266  * following to 16 multiplications, using common subexpression
267  * elimination. */
268 
269  res.x = s * s * s * p0.x + 3 * ( s * s * t ) * p1.x + 3 * ( t * t * s ) * p2.x
270  + t * t * t * p3.x;
271  res.y = s * s * s * p0.y + 3 * ( s * s * t ) * p1.y + 3 * ( t * t * s ) * p2.y
272  + t * t * t * p3.y;
273 
274  return res;
275 }
static int calc_lon ( privpath_t pp)
static

Definition at line 389 of file trace.cpp.

References abs, cyclic(), floordiv(), INFTY, potrace_privpath_s::len, potrace_privpath_s::lon, min, mod, potrace_privpath_s::pt, SAFE_CALLOC, sign(), point_s::x, xprod(), and point_s::y.

Referenced by process_path().

390 {
391  point_t* pt = pp->pt;
392  int n = pp->len;
393  int i, j, k, k1;
394  int ct[4], dir;
395  point_t constraint[2];
396  point_t cur;
397  point_t off;
398  int* pivk = NULL; /* pivk[n] */
399  int* nc = NULL; /* nc[n]: next corner */
400  point_t dk; /* direction of k-k1 */
401  int a, b, c, d;
402 
403  SAFE_CALLOC( pivk, n, int );
404  SAFE_CALLOC( nc, n, int );
405 
406  /* initialize the nc data structure. Point from each point to the
407  * furthest future point to which it is connected by a vertical or
408  * horizontal segment. We take advantage of the fact that there is
409  * always a direction change at 0 (due to the path decomposition
410  * algorithm). But even if this were not so, there is no harm, as
411  * in practice, correctness does not depend on the word "furthest"
412  * above. */
413  k = 0;
414 
415  for( i = n - 1; i >= 0; i-- )
416  {
417  if( pt[i].x != pt[k].x && pt[i].y != pt[k].y )
418  {
419  k = i + 1; /* necessarily i<n-1 in this case */
420  }
421 
422  nc[i] = k;
423  }
424 
425  SAFE_CALLOC( pp->lon, n, int );
426 
427  /* determine pivot points: for each i, let pivk[i] be the furthest k
428  * such that all j with i<j<k lie on a line connecting i,k. */
429 
430  for( i = n - 1; i >= 0; i-- )
431  {
432  ct[0] = ct[1] = ct[2] = ct[3] = 0;
433 
434  /* keep track of "directions" that have occurred */
435  dir = ( 3 + 3 * ( pt[mod( i + 1, n )].x - pt[i].x ) + ( pt[mod( i + 1, n )].y - pt[i].y ) )
436  / 2;
437  ct[dir]++;
438 
439  constraint[0].x = 0;
440  constraint[0].y = 0;
441  constraint[1].x = 0;
442  constraint[1].y = 0;
443 
444  /* find the next k such that no straight line from i to k */
445  k = nc[i];
446  k1 = i;
447 
448  while( 1 )
449  {
450  dir = ( 3 + 3 * sign( pt[k].x - pt[k1].x ) + sign( pt[k].y - pt[k1].y ) ) / 2;
451  ct[dir]++;
452 
453  /* if all four "directions" have occurred, cut this path */
454  if( ct[0] && ct[1] && ct[2] && ct[3] )
455  {
456  pivk[i] = k1;
457  goto foundk;
458  }
459 
460  cur.x = pt[k].x - pt[i].x;
461  cur.y = pt[k].y - pt[i].y;
462 
463  /* see if current constraint is violated */
464  if( xprod( constraint[0], cur ) < 0 || xprod( constraint[1], cur ) > 0 )
465  {
466  goto constraint_viol;
467  }
468 
469  /* else, update constraint */
470  if( abs( cur.x ) <= 1 && abs( cur.y ) <= 1 )
471  {
472  /* no constraint */
473  }
474  else
475  {
476  off.x = cur.x + ( ( cur.y >= 0 && ( cur.y > 0 || cur.x < 0 ) ) ? 1 : -1 );
477  off.y = cur.y + ( ( cur.x <= 0 && ( cur.x < 0 || cur.y < 0 ) ) ? 1 : -1 );
478 
479  if( xprod( constraint[0], off ) >= 0 )
480  {
481  constraint[0] = off;
482  }
483 
484  off.x = cur.x + ( ( cur.y <= 0 && ( cur.y < 0 || cur.x < 0 ) ) ? 1 : -1 );
485  off.y = cur.y + ( ( cur.x >= 0 && ( cur.x > 0 || cur.y < 0 ) ) ? 1 : -1 );
486 
487  if( xprod( constraint[1], off ) <= 0 )
488  {
489  constraint[1] = off;
490  }
491  }
492 
493  k1 = k;
494  k = nc[k1];
495 
496  if( !cyclic( k, i, k1 ) )
497  {
498  break;
499  }
500  }
501 
502 constraint_viol:
503  /* k1 was the last "corner" satisfying the current constraint, and
504  * k is the first one violating it. We now need to find the last
505  * point along k1..k which satisfied the constraint. */
506  dk.x = sign( pt[k].x - pt[k1].x );
507  dk.y = sign( pt[k].y - pt[k1].y );
508  cur.x = pt[k1].x - pt[i].x;
509  cur.y = pt[k1].y - pt[i].y;
510  /* find largest integer j such that xprod(constraint[0], cur+j*dk)
511  * >= 0 and xprod(constraint[1], cur+j*dk) <= 0. Use bilinearity
512  * of xprod. */
513  a = xprod( constraint[0], cur );
514  b = xprod( constraint[0], dk );
515  c = xprod( constraint[1], cur );
516  d = xprod( constraint[1], dk );
517  /* find largest integer j such that a+j*b>=0 and c+j*d<=0. This
518  * can be solved with integer arithmetic. */
519  j = INFTY;
520 
521  if( b < 0 )
522  {
523  j = floordiv( a, -b );
524  }
525 
526  if( d > 0 )
527  {
528  j = min( j, floordiv( -c, d ) );
529  }
530 
531  pivk[i] = mod( k1 + j, n );
532 foundk:;
533  } /* for i */
534 
535  /* clean up: for each i, let lon[i] be the largest k such that for
536  * all i' with i<=i'<k, i'<k<=pivk[i']. */
537 
538  j = pivk[n - 1];
539  pp->lon[n - 1] = j;
540 
541  for( i = n - 2; i >= 0; i-- )
542  {
543  if( cyclic( i + 1, pivk[i], j ) )
544  {
545  j = pivk[i];
546  }
547 
548  pp->lon[i] = j;
549  }
550 
551  for( i = n - 1; cyclic( mod( i + 1, n ), j, pp->lon[i] ); i-- )
552  {
553  pp->lon[i] = j;
554  }
555 
556  free( pivk );
557  free( nc );
558  return 0;
559 
560 calloc_error:
561  free( pivk );
562  free( nc );
563  return 1;
564 }
static int cyclic(int a, int b, int c)
Definition: trace.cpp:74
long y
Definition: auxiliary.h:25
#define INFTY
Definition: trace.cpp:23
static int xprod(point_t p1, point_t p2)
Definition: trace.cpp:204
#define abs(a)
Definition: auxiliary.h:84
long x
Definition: auxiliary.h:24
point_t * pt
Definition: curve.h:53
#define SAFE_CALLOC(var, n, typ)
Definition: trace.cpp:29
#define mod(a, n)
Definition: greymap.cpp:24
#define min(a, b)
Definition: auxiliary.h:85
int sign(T val)
Definition: math_util.h:44
static int floordiv(int a, int n)
Definition: auxiliary.h:70
static int calc_sums ( privpath_t pp)
static

Definition at line 327 of file trace.cpp.

References potrace_privpath_s::len, potrace_privpath_s::pt, SAFE_CALLOC, potrace_privpath_s::sums, point_s::x, sums_s::x, potrace_privpath_s::x0, sums_s::x2, sums_s::xy, point_s::y, sums_s::y, potrace_privpath_s::y0, and sums_s::y2.

Referenced by process_path().

328 {
329  int i, x, y;
330  int n = pp->len;
331 
332  SAFE_CALLOC( pp->sums, pp->len + 1, sums_t );
333 
334  /* origin */
335  pp->x0 = pp->pt[0].x;
336  pp->y0 = pp->pt[0].y;
337 
338  /* preparatory computation for later fast summing */
339  pp->sums[0].x2 = pp->sums[0].xy = pp->sums[0].y2 = pp->sums[0].x = pp->sums[0].y = 0;
340 
341  for( i = 0; i < n; i++ )
342  {
343  x = pp->pt[i].x - pp->x0;
344  y = pp->pt[i].y - pp->y0;
345  pp->sums[i + 1].x = pp->sums[i].x + x;
346  pp->sums[i + 1].y = pp->sums[i].y + y;
347  pp->sums[i + 1].x2 = pp->sums[i].x2 + (double) x * x;
348  pp->sums[i + 1].xy = pp->sums[i].xy + (double) x * y;
349  pp->sums[i + 1].y2 = pp->sums[i].y2 + (double) y * y;
350  }
351 
352  return 0;
353 
354 calloc_error:
355  return 1;
356 }
long y
Definition: auxiliary.h:25
long x
Definition: auxiliary.h:24
double y2
Definition: curve.h:41
double x
Definition: curve.h:37
double xy
Definition: curve.h:40
double y
Definition: curve.h:38
point_t * pt
Definition: curve.h:53
#define SAFE_CALLOC(var, n, typ)
Definition: trace.cpp:29
sums_t * sums
Definition: curve.h:57
Definition: curve.h:35
double x2
Definition: curve.h:39
static double cprod ( dpoint_t  p0,
dpoint_t  p1,
dpoint_t  p2,
dpoint_t  p3 
)
inlinestatic

Definition at line 211 of file trace.cpp.

References potrace_dpoint_s::x, and potrace_dpoint_s::y.

Referenced by opti_penalty(), and tangent().

212 {
213  double x1, y1, x2, y2;
214 
215  x1 = p1.x - p0.x;
216  y1 = p1.y - p0.y;
217  x2 = p3.x - p2.x;
218  y2 = p3.y - p2.y;
219 
220  return x1 * y2 - x2 * y1;
221 }
static int cyclic ( int  a,
int  b,
int  c 
)
inlinestatic

Definition at line 74 of file trace.cpp.

Referenced by calc_lon().

75 {
76  if( a <= c )
77  {
78  return a <= b && b < c;
79  }
80  else
81  {
82  return a <= b || b < c;
83  }
84 }
static double ddenom ( dpoint_t  p0,
dpoint_t  p2 
)
inlinestatic

Definition at line 65 of file trace.cpp.

References dorth_infty(), point_s::x, potrace_dpoint_s::x, point_s::y, and potrace_dpoint_s::y.

Referenced by smooth().

66 {
67  point_t r = dorth_infty( p0, p2 );
68 
69  return r.y * ( p2.x - p0.x ) - r.x * ( p2.y - p0.y );
70 }
long y
Definition: auxiliary.h:25
long x
Definition: auxiliary.h:24
static point_t dorth_infty(dpoint_t p0, dpoint_t p2)
Definition: trace.cpp:38
static double ddist ( dpoint_t  p,
dpoint_t  q 
)
inlinestatic

Definition at line 253 of file trace.cpp.

References sq, potrace_dpoint_s::x, and potrace_dpoint_s::y.

Referenced by opti_penalty().

254 {
255  return sqrt( sq( p.x - q.x ) + sq( p.y - q.y ) );
256 }
#define sq(a)
Definition: auxiliary.h:87
static point_t dorth_infty ( dpoint_t  p0,
dpoint_t  p2 
)
inlinestatic

Definition at line 38 of file trace.cpp.

References sign(), point_s::x, potrace_dpoint_s::x, point_s::y, and potrace_dpoint_s::y.

Referenced by ddenom().

39 {
40  point_t r;
41 
42  r.y = sign( p2.x - p0.x );
43  r.x = -sign( p2.y - p0.y );
44 
45  return r;
46 }
long y
Definition: auxiliary.h:25
long x
Definition: auxiliary.h:24
int sign(T val)
Definition: math_util.h:44
static double dpara ( dpoint_t  p0,
dpoint_t  p1,
dpoint_t  p2 
)
inlinestatic

Definition at line 50 of file trace.cpp.

References potrace_dpoint_s::x, and potrace_dpoint_s::y.

Referenced by opti_penalty(), opticurve(), and smooth().

51 {
52  double x1, y1, x2, y2;
53 
54  x1 = p1.x - p0.x;
55  y1 = p1.y - p0.y;
56  x2 = p2.x - p0.x;
57  y2 = p2.y - p0.y;
58 
59  return x1 * y2 - x2 * y1;
60 }
static double iprod ( dpoint_t  p0,
dpoint_t  p1,
dpoint_t  p2 
)
inlinestatic

Definition at line 225 of file trace.cpp.

References potrace_dpoint_s::x, and potrace_dpoint_s::y.

Referenced by opti_penalty().

226 {
227  double x1, y1, x2, y2;
228 
229  x1 = p1.x - p0.x;
230  y1 = p1.y - p0.y;
231  x2 = p2.x - p0.x;
232  y2 = p2.y - p0.y;
233 
234  return x1 * x2 + y1 * y2;
235 }
static double iprod1 ( dpoint_t  p0,
dpoint_t  p1,
dpoint_t  p2,
dpoint_t  p3 
)
inlinestatic

Definition at line 239 of file trace.cpp.

References potrace_dpoint_s::x, and potrace_dpoint_s::y.

Referenced by opti_penalty().

240 {
241  double x1, y1, x2, y2;
242 
243  x1 = p1.x - p0.x;
244  y1 = p1.y - p0.y;
245  x2 = p3.x - p2.x;
246  y2 = p3.y - p2.y;
247 
248  return x1 * x2 + y1 * y2;
249 }
static int opti_penalty ( privpath_t pp,
int  i,
int  j,
opti_t res,
double  opttolerance,
int *  convc,
double *  areac 
)
static

Definition at line 1122 of file trace.cpp.

References privcurve_s::alpha, opti_s::alpha, bezier(), privcurve_s::c, opti_s::c, COS179, cprod(), potrace_privpath_s::curve, ddist(), dpara(), interval(), iprod(), iprod1(), mod, privcurve_s::n, opti_s::pen, opti_s::s, sign(), sq, opti_s::t, tangent(), and privcurve_s::vertex.

Referenced by opticurve().

1129 {
1130  int m = pp->curve.n;
1131  int k, k1, k2, conv, i1;
1132  double area, alpha, d, d1, d2;
1133  dpoint_t p0, p1, p2, p3, pt;
1134  double A, R, A1, A2, A3, A4;
1135  double s, t;
1136 
1137  /* check convexity, corner-freeness, and maximum bend < 179 degrees */
1138 
1139  if( i == j )
1140  {
1141  /* sanity - a full loop can never be an opticurve */
1142  return 1;
1143  }
1144 
1145  k = i;
1146  i1 = mod( i + 1, m );
1147  k1 = mod( k + 1, m );
1148  conv = convc[k1];
1149 
1150  if( conv == 0 )
1151  {
1152  return 1;
1153  }
1154 
1155  d = ddist( pp->curve.vertex[i], pp->curve.vertex[i1] );
1156 
1157  for( k = k1; k != j; k = k1 )
1158  {
1159  k1 = mod( k + 1, m );
1160  k2 = mod( k + 2, m );
1161 
1162  if( convc[k1] != conv )
1163  {
1164  return 1;
1165  }
1166 
1167  if( sign( cprod( pp->curve.vertex[i], pp->curve.vertex[i1], pp->curve.vertex[k1],
1168  pp->curve.vertex[k2] ) )
1169  != conv )
1170  {
1171  return 1;
1172  }
1173 
1174  if( iprod1( pp->curve.vertex[i], pp->curve.vertex[i1], pp->curve.vertex[k1],
1175  pp->curve.vertex[k2] )
1176  < d * ddist( pp->curve.vertex[k1], pp->curve.vertex[k2] ) * COS179 )
1177  {
1178  return 1;
1179  }
1180  }
1181 
1182  /* the curve we're working in: */
1183  p0 = pp->curve.c[mod( i, m )][2];
1184  p1 = pp->curve.vertex[mod( i + 1, m )];
1185  p2 = pp->curve.vertex[mod( j, m )];
1186  p3 = pp->curve.c[mod( j, m )][2];
1187 
1188  /* determine its area */
1189  area = areac[j] - areac[i];
1190  area -= dpara( pp->curve.vertex[0], pp->curve.c[i][2], pp->curve.c[j][2] ) / 2;
1191 
1192  if( i >= j )
1193  {
1194  area += areac[m];
1195  }
1196 
1197  /* find intersection o of p0p1 and p2p3. Let t,s such that o =
1198  * interval(t,p0,p1) = interval(s,p3,p2). Let A be the area of the
1199  * triangle (p0,o,p3). */
1200 
1201  A1 = dpara( p0, p1, p2 );
1202  A2 = dpara( p0, p1, p3 );
1203  A3 = dpara( p0, p2, p3 );
1204  /* A4 = dpara(p1, p2, p3); */
1205  A4 = A1 + A3 - A2;
1206 
1207  if( A2 == A1 )
1208  {
1209  /* this should never happen */
1210  return 1;
1211  }
1212 
1213  t = A3 / ( A3 - A4 );
1214  s = A2 / ( A2 - A1 );
1215  A = A2 * t / 2.0;
1216 
1217  if( A == 0.0 )
1218  {
1219  /* this should never happen */
1220  return 1;
1221  }
1222 
1223  R = area / A; /* relative area */
1224  alpha = 2 - sqrt( 4 - R / 0.3 ); /* overall alpha for p0-o-p3 curve */
1225 
1226  res->c[0] = interval( t * alpha, p0, p1 );
1227  res->c[1] = interval( s * alpha, p3, p2 );
1228  res->alpha = alpha;
1229  res->t = t;
1230  res->s = s;
1231 
1232  p1 = res->c[0];
1233  p2 = res->c[1]; /* the proposed curve is now (p0,p1,p2,p3) */
1234 
1235  res->pen = 0;
1236 
1237  /* calculate penalty */
1238  /* check tangency with edges */
1239  for( k = mod( i + 1, m ); k != j; k = k1 )
1240  {
1241  k1 = mod( k + 1, m );
1242  t = tangent( p0, p1, p2, p3, pp->curve.vertex[k], pp->curve.vertex[k1] );
1243 
1244  if( t < -.5 )
1245  {
1246  return 1;
1247  }
1248 
1249  pt = bezier( t, p0, p1, p2, p3 );
1250  d = ddist( pp->curve.vertex[k], pp->curve.vertex[k1] );
1251 
1252  if( d == 0.0 )
1253  {
1254  /* this should never happen */
1255  return 1;
1256  }
1257 
1258  d1 = dpara( pp->curve.vertex[k], pp->curve.vertex[k1], pt ) / d;
1259 
1260  if( fabs( d1 ) > opttolerance )
1261  {
1262  return 1;
1263  }
1264 
1265  if( iprod( pp->curve.vertex[k], pp->curve.vertex[k1], pt ) < 0
1266  || iprod( pp->curve.vertex[k1], pp->curve.vertex[k], pt ) < 0 )
1267  {
1268  return 1;
1269  }
1270 
1271  res->pen += sq( d1 );
1272  }
1273 
1274  /* check corners */
1275  for( k = i; k != j; k = k1 )
1276  {
1277  k1 = mod( k + 1, m );
1278  t = tangent( p0, p1, p2, p3, pp->curve.c[k][2], pp->curve.c[k1][2] );
1279 
1280  if( t < -.5 )
1281  {
1282  return 1;
1283  }
1284 
1285  pt = bezier( t, p0, p1, p2, p3 );
1286  d = ddist( pp->curve.c[k][2], pp->curve.c[k1][2] );
1287 
1288  if( d == 0.0 )
1289  {
1290  /* this should never happen */
1291  return 1;
1292  }
1293 
1294  d1 = dpara( pp->curve.c[k][2], pp->curve.c[k1][2], pt ) / d;
1295  d2 = dpara( pp->curve.c[k][2], pp->curve.c[k1][2], pp->curve.vertex[k1] ) / d;
1296  d2 *= 0.75 * pp->curve.alpha[k1];
1297 
1298  if( d2 < 0 )
1299  {
1300  d1 = -d1;
1301  d2 = -d2;
1302  }
1303 
1304  if( d1 < d2 - opttolerance )
1305  {
1306  return 1;
1307  }
1308 
1309  if( d1 < d2 )
1310  {
1311  res->pen += sq( d1 - d2 );
1312  }
1313  }
1314 
1315  return 0;
1316 }
static double dpara(dpoint_t p0, dpoint_t p1, dpoint_t p2)
Definition: trace.cpp:50
static double iprod1(dpoint_t p0, dpoint_t p1, dpoint_t p2, dpoint_t p3)
Definition: trace.cpp:239
double * alpha
Definition: curve.h:29
#define sq(a)
Definition: auxiliary.h:87
static double ddist(dpoint_t p, dpoint_t q)
Definition: trace.cpp:253
static dpoint_t bezier(double t, dpoint_t p0, dpoint_t p1, dpoint_t p2, dpoint_t p3)
Definition: trace.cpp:260
privcurve_t curve
Definition: curve.h:62
dpoint_t * vertex
Definition: curve.h:28
static double iprod(dpoint_t p0, dpoint_t p1, dpoint_t p2)
Definition: trace.cpp:225
static dpoint_t interval(double lambda, dpoint_t a, dpoint_t b)
Definition: auxiliary.h:43
double alpha
Definition: trace.cpp:1115
static double tangent(dpoint_t p0, dpoint_t p1, dpoint_t p2, dpoint_t p3, dpoint_t q0, dpoint_t q1)
Definition: trace.cpp:281
dpoint_t c[2]
Definition: trace.cpp:1113
double s
Definition: trace.cpp:1114
static double cprod(dpoint_t p0, dpoint_t p1, dpoint_t p2, dpoint_t p3)
Definition: trace.cpp:211
double pen
Definition: trace.cpp:1112
#define COS179
Definition: trace.cpp:26
double t
Definition: trace.cpp:1114
int n
Definition: curve.h:20
#define mod(a, n)
Definition: greymap.cpp:24
int sign(T val)
Definition: math_util.h:44
dpoint_t(* c)[3]
Definition: curve.h:22
static int opticurve ( privpath_t pp,
double  opttolerance 
)
static

Definition at line 1322 of file trace.cpp.

References privcurve_s::alpha, opti_s::alpha, privcurve_s::alpha0, privcurve_s::alphacurve, privcurve_s::beta, privcurve_s::c, opti_s::c, potrace_privpath_s::curve, dpara(), interval(), mod, privcurve_s::n, potrace_privpath_s::ocurve, opti_penalty(), opti_s::pen, POTRACE_CURVETO, privcurve_init(), opti_s::s, SAFE_CALLOC, sign(), opti_s::t, privcurve_s::tag, and privcurve_s::vertex.

Referenced by process_path().

1323 {
1324  int m = pp->curve.n;
1325  int* pt = NULL; /* pt[m+1] */
1326  double* pen = NULL; /* pen[m+1] */
1327  int* len = NULL; /* len[m+1] */
1328  opti_t* opt = NULL; /* opt[m+1] */
1329  int om;
1330  int i, j, r;
1331  opti_t o;
1332  dpoint_t p0;
1333  int i1;
1334  double area;
1335  double alpha;
1336  double* s = NULL;
1337  double* t = NULL;
1338 
1339  int* convc = NULL; /* conv[m]: pre-computed convexities */
1340  double* areac = NULL; /* cumarea[m+1]: cache for fast area computation */
1341 
1342  SAFE_CALLOC( pt, m + 1, int );
1343  SAFE_CALLOC( pen, m + 1, double );
1344  SAFE_CALLOC( len, m + 1, int );
1345  SAFE_CALLOC( opt, m + 1, opti_t );
1346  SAFE_CALLOC( convc, m, int );
1347  SAFE_CALLOC( areac, m + 1, double );
1348 
1349  /* pre-calculate convexity: +1 = right turn, -1 = left turn, 0 = corner */
1350  for( i = 0; i < m; i++ )
1351  {
1352  if( pp->curve.tag[i] == POTRACE_CURVETO )
1353  {
1354  convc[i] = sign( dpara( pp->curve.vertex[mod( i - 1, m )], pp->curve.vertex[i],
1355  pp->curve.vertex[mod( i + 1, m )] ) );
1356  }
1357  else
1358  {
1359  convc[i] = 0;
1360  }
1361  }
1362 
1363  /* pre-calculate areas */
1364  area = 0.0;
1365  areac[0] = 0.0;
1366  p0 = pp->curve.vertex[0];
1367 
1368  for( i = 0; i < m; i++ )
1369  {
1370  i1 = mod( i + 1, m );
1371 
1372  if( pp->curve.tag[i1] == POTRACE_CURVETO )
1373  {
1374  alpha = pp->curve.alpha[i1];
1375  area += 0.3 * alpha * ( 4 - alpha )
1376  * dpara( pp->curve.c[i][2], pp->curve.vertex[i1], pp->curve.c[i1][2] ) / 2;
1377  area += dpara( p0, pp->curve.c[i][2], pp->curve.c[i1][2] ) / 2;
1378  }
1379 
1380  areac[i + 1] = area;
1381  }
1382 
1383  pt[0] = -1;
1384  pen[0] = 0;
1385  len[0] = 0;
1386 
1387  /* Fixme: we always start from a fixed point -- should find the best
1388  * curve cyclically */
1389 
1390  for( j = 1; j <= m; j++ )
1391  {
1392  /* calculate best path from 0 to j */
1393  pt[j] = j - 1;
1394  pen[j] = pen[j - 1];
1395  len[j] = len[j - 1] + 1;
1396 
1397  for( i = j - 2; i >= 0; i-- )
1398  {
1399  r = opti_penalty( pp, i, mod( j, m ), &o, opttolerance, convc, areac );
1400 
1401  if( r )
1402  {
1403  break;
1404  }
1405 
1406  if( len[j] > len[i] + 1 || ( len[j] == len[i] + 1 && pen[j] > pen[i] + o.pen ) )
1407  {
1408  pt[j] = i;
1409  pen[j] = pen[i] + o.pen;
1410  len[j] = len[i] + 1;
1411  opt[j] = o;
1412  }
1413  }
1414  }
1415 
1416  om = len[m];
1417  r = privcurve_init( &pp->ocurve, om );
1418 
1419  if( r )
1420  {
1421  goto calloc_error;
1422  }
1423 
1424  SAFE_CALLOC( s, om, double );
1425  SAFE_CALLOC( t, om, double );
1426 
1427  j = m;
1428 
1429  for( i = om - 1; i >= 0; i-- )
1430  {
1431  if( pt[j] == j - 1 )
1432  {
1433  pp->ocurve.tag[i] = pp->curve.tag[mod( j, m )];
1434  pp->ocurve.c[i][0] = pp->curve.c[mod( j, m )][0];
1435  pp->ocurve.c[i][1] = pp->curve.c[mod( j, m )][1];
1436  pp->ocurve.c[i][2] = pp->curve.c[mod( j, m )][2];
1437  pp->ocurve.vertex[i] = pp->curve.vertex[mod( j, m )];
1438  pp->ocurve.alpha[i] = pp->curve.alpha[mod( j, m )];
1439  pp->ocurve.alpha0[i] = pp->curve.alpha0[mod( j, m )];
1440  pp->ocurve.beta[i] = pp->curve.beta[mod( j, m )];
1441  s[i] = t[i] = 1.0;
1442  }
1443  else
1444  {
1445  pp->ocurve.tag[i] = POTRACE_CURVETO;
1446  pp->ocurve.c[i][0] = opt[j].c[0];
1447  pp->ocurve.c[i][1] = opt[j].c[1];
1448  pp->ocurve.c[i][2] = pp->curve.c[mod( j, m )][2];
1449  pp->ocurve.vertex[i] = interval(
1450  opt[j].s, pp->curve.c[mod( j, m )][2], pp->curve.vertex[mod( j, m )] );
1451  pp->ocurve.alpha[i] = opt[j].alpha;
1452  pp->ocurve.alpha0[i] = opt[j].alpha;
1453  s[i] = opt[j].s;
1454  t[i] = opt[j].t;
1455  }
1456 
1457  j = pt[j];
1458  }
1459 
1460  /* calculate beta parameters */
1461  for( i = 0; i < om; i++ )
1462  {
1463  i1 = mod( i + 1, om );
1464  pp->ocurve.beta[i] = s[i] / ( s[i] + t[i1] );
1465  }
1466 
1467  pp->ocurve.alphacurve = 1;
1468 
1469  free( pt );
1470  free( pen );
1471  free( len );
1472  free( opt );
1473  free( s );
1474  free( t );
1475  free( convc );
1476  free( areac );
1477  return 0;
1478 
1479 calloc_error:
1480  free( pt );
1481  free( pen );
1482  free( len );
1483  free( opt );
1484  free( s );
1485  free( t );
1486  free( convc );
1487  free( areac );
1488  return 1;
1489 }
static double dpara(dpoint_t p0, dpoint_t p1, dpoint_t p2)
Definition: trace.cpp:50
int privcurve_init(privcurve_t *curve, int n)
Definition: curve.cpp:98
double * alpha
Definition: curve.h:29
privcurve_t curve
Definition: curve.h:62
int alphacurve
Definition: curve.h:27
dpoint_t * vertex
Definition: curve.h:28
static dpoint_t interval(double lambda, dpoint_t a, dpoint_t b)
Definition: auxiliary.h:43
privcurve_t ocurve
Definition: curve.h:63
double alpha
Definition: trace.cpp:1115
dpoint_t c[2]
Definition: trace.cpp:1113
double s
Definition: trace.cpp:1114
double * beta
Definition: curve.h:31
#define POTRACE_CURVETO
Definition: potracelib.h:78
static int opti_penalty(privpath_t *pp, int i, int j, opti_t *res, double opttolerance, int *convc, double *areac)
Definition: trace.cpp:1122
#define SAFE_CALLOC(var, n, typ)
Definition: trace.cpp:29
double pen
Definition: trace.cpp:1112
double * alpha0
Definition: curve.h:30
int * tag
Definition: curve.h:21
double t
Definition: trace.cpp:1114
int n
Definition: curve.h:20
#define mod(a, n)
Definition: greymap.cpp:24
int sign(T val)
Definition: math_util.h:44
dpoint_t(* c)[3]
Definition: curve.h:22
static double penalty3 ( privpath_t pp,
int  i,
int  j 
)
static

Definition at line 573 of file trace.cpp.

References potrace_privpath_s::len, potrace_privpath_s::pt, potrace_privpath_s::sums, point_s::x, sums_s::x, sums_s::x2, sums_s::xy, point_s::y, sums_s::y, and sums_s::y2.

Referenced by bestpolygon().

574 {
575  int n = pp->len;
576  point_t* pt = pp->pt;
577  sums_t* sums = pp->sums;
578 
579  /* assume 0<=i<j<=n */
580  double x, y, x2, xy, y2;
581  double k;
582  double a, b, c, s;
583  double px, py, ex, ey;
584 
585  int r = 0; /* rotations from i to j */
586 
587  if( j >= n )
588  {
589  j -= n;
590  r = 1;
591  }
592 
593  /* critical inner loop: the "if" gives a 4.6 percent speedup */
594  if( r == 0 )
595  {
596  x = sums[j + 1].x - sums[i].x;
597  y = sums[j + 1].y - sums[i].y;
598  x2 = sums[j + 1].x2 - sums[i].x2;
599  xy = sums[j + 1].xy - sums[i].xy;
600  y2 = sums[j + 1].y2 - sums[i].y2;
601  k = j + 1 - i;
602  }
603  else
604  {
605  x = sums[j + 1].x - sums[i].x + sums[n].x;
606  y = sums[j + 1].y - sums[i].y + sums[n].y;
607  x2 = sums[j + 1].x2 - sums[i].x2 + sums[n].x2;
608  xy = sums[j + 1].xy - sums[i].xy + sums[n].xy;
609  y2 = sums[j + 1].y2 - sums[i].y2 + sums[n].y2;
610  k = j + 1 - i + n;
611  }
612 
613  px = ( pt[i].x + pt[j].x ) / 2.0 - pt[0].x;
614  py = ( pt[i].y + pt[j].y ) / 2.0 - pt[0].y;
615  ey = ( pt[j].x - pt[i].x );
616  ex = -( pt[j].y - pt[i].y );
617 
618  a = ( ( x2 - 2 * x * px ) / k + px * px );
619  b = ( ( xy - x * py - y * px ) / k + px * py );
620  c = ( ( y2 - 2 * y * py ) / k + py * py );
621 
622  s = ex * ex * a + 2 * ex * ey * b + ey * ey * c;
623 
624  return sqrt( s );
625 }
long y
Definition: auxiliary.h:25
long x
Definition: auxiliary.h:24
double y2
Definition: curve.h:41
double x
Definition: curve.h:37
double xy
Definition: curve.h:40
double y
Definition: curve.h:38
point_t * pt
Definition: curve.h:53
sums_t * sums
Definition: curve.h:57
Definition: curve.h:35
double x2
Definition: curve.h:39
static void pointslope ( privpath_t pp,
int  i,
int  j,
dpoint_t ctr,
dpoint_t dir 
)
static

Definition at line 89 of file trace.cpp.

References potrace_privpath_s::len, potrace_privpath_s::sums, sums_s::x, potrace_dpoint_s::x, sums_s::x2, sums_s::xy, sums_s::y, potrace_dpoint_s::y, and sums_s::y2.

Referenced by adjust_vertices().

90 {
91  /* assume i<j */
92 
93  int n = pp->len;
94  sums_t* sums = pp->sums;
95 
96  double x, y, x2, xy, y2;
97  double k;
98  double a, b, c, lambda2, l;
99  int r = 0; /* rotations from i to j */
100 
101  while( j >= n )
102  {
103  j -= n;
104  r += 1;
105  }
106 
107  while( i >= n )
108  {
109  i -= n;
110  r -= 1;
111  }
112 
113  while( j < 0 )
114  {
115  j += n;
116  r -= 1;
117  }
118 
119  while( i < 0 )
120  {
121  i += n;
122  r += 1;
123  }
124 
125  x = sums[j + 1].x - sums[i].x + r * sums[n].x;
126  y = sums[j + 1].y - sums[i].y + r * sums[n].y;
127  x2 = sums[j + 1].x2 - sums[i].x2 + r * sums[n].x2;
128  xy = sums[j + 1].xy - sums[i].xy + r * sums[n].xy;
129  y2 = sums[j + 1].y2 - sums[i].y2 + r * sums[n].y2;
130  k = j + 1 - i + r * n;
131 
132  ctr->x = x / k;
133  ctr->y = y / k;
134 
135  a = ( x2 - (double) x * x / k ) / k;
136  b = ( xy - (double) x * y / k ) / k;
137  c = ( y2 - (double) y * y / k ) / k;
138 
139  lambda2 = ( a + c + sqrt( ( a - c ) * ( a - c ) + 4 * b * b ) ) / 2; /* larger e.value */
140 
141  /* now find e.vector for lambda2 */
142  a -= lambda2;
143  c -= lambda2;
144 
145  if( fabs( a ) >= fabs( c ) )
146  {
147  l = sqrt( a * a + b * b );
148 
149  if( l != 0 )
150  {
151  dir->x = -b / l;
152  dir->y = a / l;
153  }
154  }
155  else
156  {
157  l = sqrt( c * c + b * b );
158 
159  if( l != 0 )
160  {
161  dir->x = -c / l;
162  dir->y = b / l;
163  }
164  }
165 
166  if( l == 0 )
167  {
168  dir->x = dir->y = 0; /* sometimes this can happen when k=4:
169  * the two eigenvalues coincide */
170  }
171 }
double y2
Definition: curve.h:41
double x
Definition: curve.h:37
double xy
Definition: curve.h:40
double y
Definition: curve.h:38
sums_t * sums
Definition: curve.h:57
Definition: curve.h:35
double x2
Definition: curve.h:39
int process_path ( path_t plist,
const potrace_param_t param,
progress_t progress 
)

Definition at line 1499 of file trace.cpp.

References adjust_vertices(), potrace_param_s::alphamax, bestpolygon(), calc_lon(), calc_sums(), progress_s::callback, potrace_privpath_s::curve, potrace_path_s::curve, potrace_privpath_s::fcurve, potrace_privpath_s::len, list_forall, potrace_privpath_s::ocurve, potrace_param_s::opticurve, opticurve(), potrace_param_s::opttolerance, potrace_path_s::priv, privcurve_to_curve(), progress_update(), reverse(), potrace_path_s::sign, smooth(), and TRY.

Referenced by potrace_trace().

1500 {
1501  path_t* p;
1502  double nn = 0, cn = 0;
1503 
1504  if( progress->callback )
1505  {
1506  /* precompute task size for progress estimates */
1507  nn = 0;
1508  list_forall( p, plist ) {
1509  nn += p->priv->len;
1510  }
1511  cn = 0;
1512  }
1513 
1514  /* call downstream function with each path */
1515  list_forall( p, plist ) {
1516  TRY( calc_sums( p->priv ) );
1517  TRY( calc_lon( p->priv ) );
1518  TRY( bestpolygon( p->priv ) );
1519  TRY( adjust_vertices( p->priv ) );
1520 
1521  if( p->sign == '-' )
1522  {
1523  /* reverse orientation of negative paths */
1524  reverse( &p->priv->curve );
1525  }
1526 
1527  smooth( &p->priv->curve, param->alphamax );
1528 
1529  if( param->opticurve )
1530  {
1531  TRY( opticurve( p->priv, param->opttolerance ) );
1532  p->priv->fcurve = &p->priv->ocurve;
1533  }
1534  else
1535  {
1536  p->priv->fcurve = &p->priv->curve;
1537  }
1538 
1539  privcurve_to_curve( p->priv->fcurve, &p->curve );
1540 
1541  if( progress->callback )
1542  {
1543  cn += p->priv->len;
1544  progress_update( cn / nn, progress );
1545  }
1546  }
1547 
1548  progress_update( 1.0, progress );
1549 
1550  return 0;
1551 
1552 try_error:
1553  return 1;
1554 }
static int calc_lon(privpath_t *pp)
Definition: trace.cpp:389
#define TRY(x)
Definition: trace.cpp:1494
static int bestpolygon(privpath_t *pp)
Definition: trace.cpp:631
double alphamax
Definition: potracelib.h:42
privcurve_t curve
Definition: curve.h:62
static int calc_sums(privpath_t *pp)
Definition: trace.cpp:327
privcurve_t ocurve
Definition: curve.h:63
static void smooth(privcurve_t *curve, double alphamax)
Definition: trace.cpp:1041
struct potrace_privpath_s * priv
Definition: potracelib.h:103
privcurve_t * fcurve
Definition: curve.h:64
static int opticurve(privpath_t *pp, double opttolerance)
Definition: trace.cpp:1322
void privcurve_to_curve(privcurve_t *pc, potrace_curve_t *c)
Definition: curve.cpp:122
#define list_forall(elt, list)
Definition: lists.h:44
double opttolerance
Definition: potracelib.h:44
static void progress_update(double d, progress_t *prog)
Definition: progress.h:29
static void reverse(privcurve_t *curve)
Definition: trace.cpp:1025
static int adjust_vertices(privpath_t *pp)
Definition: trace.cpp:772
potrace_curve_t curve
Definition: potracelib.h:96
void(* callback)(double progress, void *privdata)
Definition: progress.h:17
static double quadform ( quadform_t  Q,
dpoint_t  w 
)
inlinestatic

Definition at line 180 of file trace.cpp.

References potrace_dpoint_s::x, and potrace_dpoint_s::y.

Referenced by adjust_vertices().

181 {
182  double v[3];
183  int i, j;
184  double sum;
185 
186  v[0] = w.x;
187  v[1] = w.y;
188  v[2] = 1;
189  sum = 0.0;
190 
191  for( i = 0; i < 3; i++ )
192  {
193  for( j = 0; j < 3; j++ )
194  {
195  sum += v[i] * Q[i][j] * v[j];
196  }
197  }
198 
199  return sum;
200 }
static void reverse ( privcurve_t curve)
static

Definition at line 1025 of file trace.cpp.

References privcurve_s::n, and privcurve_s::vertex.

Referenced by AUTOPLACER::choose_side_for_fields(), process_path(), PNS::LINE::Reverse(), SHAPE_LINE_CHAIN::Reverse(), and PNS::SHOVE::runOptimizer().

1026 {
1027  int m = curve->n;
1028  int i, j;
1029  dpoint_t tmp;
1030 
1031  for( i = 0, j = m - 1; i < j; i++, j-- )
1032  {
1033  tmp = curve->vertex[i];
1034  curve->vertex[i] = curve->vertex[j];
1035  curve->vertex[j] = tmp;
1036  }
1037 }
dpoint_t * vertex
Definition: curve.h:28
int n
Definition: curve.h:20
static void smooth ( privcurve_t curve,
double  alphamax 
)
static

Definition at line 1041 of file trace.cpp.

References privcurve_s::alpha, privcurve_s::alpha0, privcurve_s::alphacurve, privcurve_s::beta, privcurve_s::c, ddenom(), dpara(), interval(), mod, privcurve_s::n, POTRACE_CORNER, POTRACE_CURVETO, privcurve_s::tag, and privcurve_s::vertex.

Referenced by process_path().

1042 {
1043  int m = curve->n;
1044 
1045  int i, j, k;
1046  double dd, denom, alpha;
1047  dpoint_t p2, p3, p4;
1048 
1049  /* examine each vertex and find its best fit */
1050  for( i = 0; i < m; i++ )
1051  {
1052  j = mod( i + 1, m );
1053  k = mod( i + 2, m );
1054  p4 = interval( 1 / 2.0, curve->vertex[k], curve->vertex[j] );
1055 
1056  denom = ddenom( curve->vertex[i], curve->vertex[k] );
1057 
1058  if( denom != 0.0 )
1059  {
1060  dd = dpara( curve->vertex[i], curve->vertex[j], curve->vertex[k] ) / denom;
1061  dd = fabs( dd );
1062  alpha = dd > 1 ? ( 1 - 1.0 / dd ) : 0;
1063  alpha = alpha / 0.75;
1064  }
1065  else
1066  {
1067  alpha = 4 / 3.0;
1068  }
1069 
1070  curve->alpha0[j] = alpha; /* remember "original" value of alpha */
1071 
1072  if( alpha >= alphamax )
1073  {
1074  /* pointed corner */
1075  curve->tag[j] = POTRACE_CORNER;
1076  curve->c[j][1] = curve->vertex[j];
1077  curve->c[j][2] = p4;
1078  }
1079  else
1080  {
1081  if( alpha < 0.55 )
1082  {
1083  alpha = 0.55;
1084  }
1085  else if( alpha > 1 )
1086  {
1087  alpha = 1;
1088  }
1089 
1090  p2 = interval( .5 + .5 * alpha, curve->vertex[i], curve->vertex[j] );
1091  p3 = interval( .5 + .5 * alpha, curve->vertex[k], curve->vertex[j] );
1092  curve->tag[j] = POTRACE_CURVETO;
1093  curve->c[j][0] = p2;
1094  curve->c[j][1] = p3;
1095  curve->c[j][2] = p4;
1096  }
1097 
1098  curve->alpha[j] = alpha; /* store the "cropped" value of alpha */
1099  curve->beta[j] = 0.5;
1100  }
1101 
1102  curve->alphacurve = 1;
1103 }
static double dpara(dpoint_t p0, dpoint_t p1, dpoint_t p2)
Definition: trace.cpp:50
double * alpha
Definition: curve.h:29
int alphacurve
Definition: curve.h:27
#define POTRACE_CORNER
Definition: potracelib.h:79
dpoint_t * vertex
Definition: curve.h:28
static dpoint_t interval(double lambda, dpoint_t a, dpoint_t b)
Definition: auxiliary.h:43
double * beta
Definition: curve.h:31
#define POTRACE_CURVETO
Definition: potracelib.h:78
static double ddenom(dpoint_t p0, dpoint_t p2)
Definition: trace.cpp:65
double * alpha0
Definition: curve.h:30
int * tag
Definition: curve.h:21
int n
Definition: curve.h:20
#define mod(a, n)
Definition: greymap.cpp:24
dpoint_t(* c)[3]
Definition: curve.h:22
static double tangent ( dpoint_t  p0,
dpoint_t  p1,
dpoint_t  p2,
dpoint_t  p3,
dpoint_t  q0,
dpoint_t  q1 
)
static

Definition at line 281 of file trace.cpp.

References cprod().

Referenced by opti_penalty().

283 {
284  double A, B, C; /* (1-t)^2 A + 2(1-t)t B + t^2 C = 0 */
285  double a, b, c; /* a t^2 + b t + c = 0 */
286  double d, s, r1, r2;
287 
288  A = cprod( p0, p1, q0, q1 );
289  B = cprod( p1, p2, q0, q1 );
290  C = cprod( p2, p3, q0, q1 );
291 
292  a = A - 2 * B + C;
293  b = -2 * A + 2 * B;
294  c = A;
295 
296  d = b * b - 4 * a * c;
297 
298  if( a == 0 || d < 0 )
299  {
300  return -1.0;
301  }
302 
303  s = sqrt( d );
304 
305  r1 = ( -b + s ) / ( 2 * a );
306  r2 = ( -b - s ) / ( 2 * a );
307 
308  if( r1 >= 0 && r1 <= 1 )
309  {
310  return r1;
311  }
312  else if( r2 >= 0 && r2 <= 1 )
313  {
314  return r2;
315  }
316  else
317  {
318  return -1.0;
319  }
320 }
static double cprod(dpoint_t p0, dpoint_t p1, dpoint_t p2, dpoint_t p3)
Definition: trace.cpp:211
static int xprod ( point_t  p1,
point_t  p2 
)
inlinestatic

Definition at line 204 of file trace.cpp.

References point_s::x, and point_s::y.

Referenced by calc_lon().

205 {
206  return p1.x * p2.y - p1.y * p2.x;
207 }
long y
Definition: auxiliary.h:25
long x
Definition: auxiliary.h:24