Commit 6950425
Create a library of the Graphics Gems code
Changed files (17)
modifieddodo.pymodifiedinstall.pydeletedlibbirdfont/GraphicsGems/GGVecLib.cdeletedlibbirdfont/GraphicsGems/GraphicsGems.hdeletedlibbirdfont/GraphicsGems/fit_cubic.cmodifiedlibbirdfont/Path.valamodifiedlibbirdfont/StrokeTool.valaaddedlibbirdgems/FitCubic.vala~addedlibbirdgems/GGVecLib.caddedlibbirdgems/GraphicsGems.haddedlibbirdgems/fit_cubic.caddedlibbirdgems/fit_cubic.c~modifiedscripts/bavala.pymodifiedscripts/build.pymodifiedscripts/linux_build.pymodifiedscripts/version.pymodifiedscripts/windows_build.py
Diff
@@ -1,5 +1,5 @@ """-Copyright (C) 2012, 2013, 2014 Eduardo Naufel Schettino and Johan Mattsson+Copyright (C) 2012, 2013, 2014 2015 Eduardo Naufel Schettino and Johan Mattsson This program is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by@@ -30,6 +30,7 @@ from scripts import config DOIT_CONFIG = { 'default_tasks': [ 'build',+ 'libbirdgems', 'libbirdxml', 'libbirdfont', 'birdfont',@@ -63,6 +64,10 @@ else: 'posixtypes' ] +LIBBIRD_LIBS = [+ 'glib-2.0'+]+ def task_build (): if not os.path.exists ("build/configured"): print ("Project is not configured")@@ -108,41 +113,46 @@ def task_libbirdxml(): yield libbirdxml.gen_ln() -libbird = Vala(src='libbirdfont', build='build', library='birdfont', so_version=version.SO_VERSION, pkg_libs=LIBS, vala_deps=[libbirdxml])+libbirdgems = Vala(src='libbirdgems', build='build', library='birdgems', so_version=version.LIBBIRDGEMS_SO_VERSION, pkg_libs=LIBBIRD_LIBS, vala_deps=[])+def task_libbirdgems():+ yield libbirdgems.gen_c(valac_options)+ yield libbirdgems.gen_o([])+ yield libbirdgems.gen_so('-L ./build -l m')+ yield libbirdgems.gen_ln()+++libbird = Vala(src='libbirdfont', build='build', library='birdfont', so_version=version.SO_VERSION, pkg_libs=LIBS, vala_deps=[libbirdgems, libbirdxml]) def task_libbirdfont(): yield libbird.gen_c(valac_options) yield libbird.gen_o(['-fPIC -I./build/', """-D 'GETTEXT_PACKAGE="birdfont"'"""])- yield libbird.gen_so('-L ./build -l birdxml')+ yield libbird.gen_so('-L ./build -l birdxml -L ./build -l birdgems') yield libbird.gen_ln() def task_birdfont():- bird = Vala(src='birdfont', build='build', pkg_libs=LIBS, vala_deps=[libbird, libbirdxml])+ bird = Vala(src='birdfont', build='build', pkg_libs=LIBS, vala_deps=[libbird, libbirdxml, libbirdgems]) yield bird.gen_c(valac_options) yield bird.gen_bin(["""-D 'GETTEXT_PACKAGE="birdfont"' """]) def task_birdfont_autotrace():- exp = Vala(src='birdfont-autotrace', build='build', pkg_libs=LIBS, vala_deps=[libbird, libbirdxml])+ exp = Vala(src='birdfont-autotrace', build='build', pkg_libs=LIBS, vala_deps=[libbird, libbirdxml, libbirdgems]) yield exp.gen_c(valac_options) yield exp.gen_bin(["""-D 'GETTEXT_PACKAGE="birdfont"' """]) def task_birdfont_export():- exp = Vala(src='birdfont-export', build='build', pkg_libs=LIBS, vala_deps=[libbird, libbirdxml])+ exp = Vala(src='birdfont-export', build='build', pkg_libs=LIBS, vala_deps=[libbird, libbirdxml, libbirdgems]) yield exp.gen_c(valac_options) yield exp.gen_bin(["""-D 'GETTEXT_PACKAGE="birdfont"' """]) def task_birdfont_import():- exp = Vala(src='birdfont-import', build='build', pkg_libs=LIBS, vala_deps=[libbird, libbirdxml])+ exp = Vala(src='birdfont-import', build='build', pkg_libs=LIBS, vala_deps=[libbird, libbirdxml, libbirdgems]) yield exp.gen_c(valac_options) yield exp.gen_bin(["""-D 'GETTEXT_PACKAGE="birdfont"' """]) --- def task_compile_translations (): """translate po files""" return {@@ -1,6 +1,6 @@ #!/usr/bin/python """-Copyright (C) 2013 2014 Johan Mattsson+Copyright (C) 2013 2014 2015 Johan Mattsson This program is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by@@ -168,6 +168,24 @@ elif os.path.isfile ('build/bin/libbirdxml.' + version.LIBBIRDXML_SO_VERSION + ' link (libdir, 'libbirdxml.' + version.LIBBIRDXML_SO_VERSION + '.dylib', ' libbirdxml.dylib') else: print ("Can't find libbirdxml.")+ exit (1)+++if os.path.isfile ('build/bin/libbirdgems.so.' + version.LIBBIRDGEMS_SO_VERSION):+ install ('build/bin/libbirdgems.so.' + version.LIBBIRDGEMS_SO_VERSION, libdir, 644)+ link (libdir, 'libbirdgems.so.' + version.LIBBIRDGEMS_SO_VERSION, ' libbirdgems.so.' + version.LIBBIRDXML_SO_VERSION_MAJOR)+ link (libdir, 'libbirdgems.so.' + version.LIBBIRDGEMS_SO_VERSION, ' libbirdgems.so')+elif os.path.isfile ('build/libbirdgems.so.' + version.LIBBIRDGEMS_SO_VERSION):+ install ('build/libbirdgems.so.' + version.LIBBIRDGEMS_SO_VERSION, libdir, 644)+ link (libdir, 'libbirdgems.so.' + version.LIBBIRDGEMS_SO_VERSION, ' libbirdgems.so.' + version.LIBBIRDGEMS_SO_VERSION_MAJOR)+ link (libdir, 'libbirdgems.so.' + version.LIBBIRDGEMS_SO_VERSION, ' libbirdgems.so')+elif os.path.isfile ('build/bin/libbirdgems.' + version.LIBBIRDGEMS_SO_VERSION + '.dylib'):+ install ('build/bin/libbirdgems.' + version.LIBBIRDGEMS_SO_VERSION + '.dylib', libdir, 644)+ link (libdir, 'libbirdgems.' + version.LIBBIRDGEMS_SO_VERSION + '.dylib', ' libbirdgems.dylib.' + version.LIBBIRDGEMS_SO_VERSION_MAJOR)+ link (libdir, 'libbirdgems.' + version.LIBBIRDGEMS_SO_VERSION + '.dylib', ' libbirdgems.dylib')+else:+ print ("Can't find libbirdgems.")+ exit (1) #manpagesdeleted file mode 100644@@ -1,458 +0,0 @@-/* -2d and 3d Vector C Library -by Andrew Glassner-from "Graphics Gems", Academic Press, 1990-*/--/*- * EULA: The Graphics Gems code is copyright-protected. In other words, - * you cannot claim the text of the code as your own and resell it. - * Using the code is permitted in any program, product, or library, - * non-commercial or commercial. Giving credit is not required, though - * is a nice gesture. The code comes as-is, and if there are any flaws - * or problems with any Gems code, nobody involved with Gems - authors,- * editors, publishers, or webmasters - are to be held responsible. - * Basically, don't be a jerk, and remember that anything free comes - * with no guarantee. - */--#include <math.h>-#include "GraphicsGems.h"--/******************/-/* 2d Library */-/******************/--/* returns squared length of input vector */ -double V2SquaredLength(a) -Vector2 *a;-{ return((a->x * a->x)+(a->y * a->y));- }- -/* returns length of input vector */-double V2Length(a) -Vector2 *a;-{- return(sqrt(V2SquaredLength(a)));- }- -/* negates the input vector and returns it */-Vector2 *V2Negate(v) -Vector2 *v;-{- v->x = -v->x; v->y = -v->y;- return(v);- }--/* normalizes the input vector and returns it */-Vector2 *V2Normalize(v) -Vector2 *v;-{-double len = V2Length(v);- if (len != 0.0) { v->x /= len; v->y /= len; }- return(v);- }---/* scales the input vector to the new length and returns it */-Vector2 *V2Scale(v, newlen) -Vector2 *v;-double newlen;-{-double len = V2Length(v);- if (len != 0.0) { v->x *= newlen/len; v->y *= newlen/len; }- return(v);- }--/* return vector sum c = a+b */-Vector2 *V2Add(a, b, c)-Vector2 *a, *b, *c;-{- c->x = a->x+b->x; c->y = a->y+b->y;- return(c);- }- -/* return vector difference c = a-b */-Vector2 *V2Sub(a, b, c)-Vector2 *a, *b, *c;-{- c->x = a->x-b->x; c->y = a->y-b->y;- return(c);- }--/* return the dot product of vectors a and b */-double V2Dot(a, b) -Vector2 *a, *b; -{- return((a->x*b->x)+(a->y*b->y));- }--/* linearly interpolate between vectors by an amount alpha */-/* and return the resulting vector. */-/* When alpha=0, result=lo. When alpha=1, result=hi. */-Vector2 *V2Lerp(lo, hi, alpha, result) -Vector2 *lo, *hi, *result; -double alpha;-{- result->x = LERP(alpha, lo->x, hi->x);- result->y = LERP(alpha, lo->y, hi->y);- return(result);- }---/* make a linear combination of two vectors and return the result. */-/* result = (a * ascl) + (b * bscl) */-Vector2 *V2Combine (a, b, result, ascl, bscl) -Vector2 *a, *b, *result;-double ascl, bscl;-{- result->x = (ascl * a->x) + (bscl * b->x);- result->y = (ascl * a->y) + (bscl * b->y);- return(result);- }--/* multiply two vectors together component-wise */-Vector2 *V2Mul (a, b, result) -Vector2 *a, *b, *result;-{- result->x = a->x * b->x;- result->y = a->y * b->y;- return(result);- }--/* return the distance between two points */-double V2DistanceBetween2Points(a, b)-Point2 *a, *b;-{-double dx = a->x - b->x;-double dy = a->y - b->y;- return(sqrt((dx*dx)+(dy*dy)));- }--/* return the vector perpendicular to the input vector a */-Vector2 *V2MakePerpendicular(a, ap)-Vector2 *a, *ap;-{- ap->x = -a->y;- ap->y = a->x;- return(ap);- }--/* create, initialize, and return a new vector */-Vector2 *V2New(x, y)-double x, y;-{-Vector2 *v = NEWTYPE(Vector2);- v->x = x; v->y = y; - return(v);- }- --/* create, initialize, and return a duplicate vector */-Vector2 *V2Duplicate(a)-Vector2 *a;-{-Vector2 *v = NEWTYPE(Vector2);- v->x = a->x; v->y = a->y; - return(v);- }- -/* multiply a point by a matrix and return the transformed point */-Point2 *V2MulPointByMatrix(p, m)-Point2 *p;-Matrix3 *m;-{-double w;-Point2 ptmp;- ptmp.x = (p->x * m->element[0][0]) + - (p->y * m->element[1][0]) + m->element[2][0];- ptmp.y = (p->x * m->element[0][1]) + - (p->y * m->element[1][1]) + m->element[2][1];- w = (p->x * m->element[0][2]) + - (p->y * m->element[1][2]) + m->element[2][2];- if (w != 0.0) { ptmp.x /= w; ptmp.y /= w; }- *p = ptmp;- return(p);- }--/* multiply together matrices c = ab */-/* note that c must not point to either of the input matrices */-Matrix3 *V2MatMul(a, b, c)-Matrix3 *a, *b, *c;-{-int i, j, k;- for (i=0; i<3; i++) {- for (j=0; j<3; j++) {- c->element[i][j] = 0;- for (k=0; k<3; k++) c->element[i][j] += - a->element[i][k] * b->element[k][j];- }- }- return(c);- }-----/******************/-/* 3d Library */-/******************/- -/* returns squared length of input vector */ -double V3SquaredLength(a) -Vector3 *a;-{- return((a->x * a->x)+(a->y * a->y)+(a->z * a->z));- }--/* returns length of input vector */-double V3Length(a) -Vector3 *a;-{- return(sqrt(V3SquaredLength(a)));- }--/* negates the input vector and returns it */-Vector3 *V3Negate(v) -Vector3 *v;-{- v->x = -v->x; v->y = -v->y; v->z = -v->z;- return(v);- }--/* normalizes the input vector and returns it */-Vector3 *V3Normalize(v) -Vector3 *v;-{-double len = V3Length(v);- if (len != 0.0) { v->x /= len; v->y /= len; v->z /= len; }- return(v);- }--/* scales the input vector to the new length and returns it */-Vector3 *V3Scale(v, newlen) -Vector3 *v;-double newlen;-{-double len = V3Length(v);- if (len != 0.0) {- v->x *= newlen/len; v->y *= newlen/len; v->z *= newlen/len;- }- return(v);- }---/* return vector sum c = a+b */-Vector3 *V3Add(a, b, c)-Vector3 *a, *b, *c;-{- c->x = a->x+b->x; c->y = a->y+b->y; c->z = a->z+b->z;- return(c);- }- -/* return vector difference c = a-b */-Vector3 *V3Sub(a, b, c)-Vector3 *a, *b, *c;-{- c->x = a->x-b->x; c->y = a->y-b->y; c->z = a->z-b->z;- return(c);- }--/* return the dot product of vectors a and b */-double V3Dot(a, b) -Vector3 *a, *b; -{- return((a->x*b->x)+(a->y*b->y)+(a->z*b->z));- }--/* linearly interpolate between vectors by an amount alpha */-/* and return the resulting vector. */-/* When alpha=0, result=lo. When alpha=1, result=hi. */-Vector3 *V3Lerp(lo, hi, alpha, result) -Vector3 *lo, *hi, *result; -double alpha;-{- result->x = LERP(alpha, lo->x, hi->x);- result->y = LERP(alpha, lo->y, hi->y);- result->z = LERP(alpha, lo->z, hi->z);- return(result);- }--/* make a linear combination of two vectors and return the result. */-/* result = (a * ascl) + (b * bscl) */-Vector3 *V3Combine (a, b, result, ascl, bscl) -Vector3 *a, *b, *result;-double ascl, bscl;-{- result->x = (ascl * a->x) + (bscl * b->x);- result->y = (ascl * a->y) + (bscl * b->y);- result->z = (ascl * a->z) + (bscl * b->z);- return(result);- }---/* multiply two vectors together component-wise and return the result */-Vector3 *V3Mul (a, b, result) -Vector3 *a, *b, *result;-{- result->x = a->x * b->x;- result->y = a->y * b->y;- result->z = a->z * b->z;- return(result);- }--/* return the distance between two points */-double V3DistanceBetween2Points(a, b)-Point3 *a, *b;-{-double dx = a->x - b->x;-double dy = a->y - b->y;-double dz = a->z - b->z;- return(sqrt((dx*dx)+(dy*dy)+(dz*dz)));- }--/* return the cross product c = a cross b */-Vector3 *V3Cross(a, b, c)-Vector3 *a, *b, *c;-{- c->x = (a->y*b->z) - (a->z*b->y);- c->y = (a->z*b->x) - (a->x*b->z);- c->z = (a->x*b->y) - (a->y*b->x);- return(c);- }--/* create, initialize, and return a new vector */-Vector3 *V3New(x, y, z)-double x, y, z;-{-Vector3 *v = NEWTYPE(Vector3);- v->x = x; v->y = y; v->z = z;- return(v);- }--/* create, initialize, and return a duplicate vector */-Vector3 *V3Duplicate(a)-Vector3 *a;-{-Vector3 *v = NEWTYPE(Vector3);- v->x = a->x; v->y = a->y; v->z = a->z;- return(v);- }-- -/* multiply a point by a matrix and return the transformed point */-Point3 *V3MulPointByMatrix(p, m)-Point3 *p;-Matrix4 *m;-{-double w;-Point3 ptmp;- ptmp.x = (p->x * m->element[0][0]) + (p->y * m->element[1][0]) + - (p->z * m->element[2][0]) + m->element[3][0];- ptmp.y = (p->x * m->element[0][1]) + (p->y * m->element[1][1]) + - (p->z * m->element[2][1]) + m->element[3][1];- ptmp.z = (p->x * m->element[0][2]) + (p->y * m->element[1][2]) + - (p->z * m->element[2][2]) + m->element[3][2];- w = (p->x * m->element[0][3]) + (p->y * m->element[1][3]) + - (p->z * m->element[2][3]) + m->element[3][3];- if (w != 0.0) { ptmp.x /= w; ptmp.y /= w; ptmp.z /= w; }- *p = ptmp;- return(p);- }--/* multiply together matrices c = ab */-/* note that c must not point to either of the input matrices */-Matrix4 *V3MatMul(a, b, c)-Matrix4 *a, *b, *c;-{-int i, j, k;- for (i=0; i<4; i++) {- for (j=0; j<4; j++) {- c->element[i][j] = 0;- for (k=0; k<4; k++) c->element[i][j] += - a->element[i][k] * b->element[k][j];- }- }- return(c);- }--/* binary greatest common divisor by Silver and Terzian. See Knuth */-/* both inputs must be >= 0 */-gcd(u, v)-int u, v;-{-int t, f;- if ((u<0) || (v<0)) return(1); /* error if u<0 or v<0 */- f = 1;- while ((0 == (u%2)) && (0 == (v%2))) {- u>>=1; v>>=1, f*=2;- }- if (u&01) { t = -v; goto B4; } else { t = u; }- B3: if (t > 0) { t >>= 1; } else { t = -((-t) >> 1); }- B4: if (0 == (t%2)) goto B3;-- if (t > 0) u = t; else v = -t;- if (0 != (t = u - v)) goto B3;- return(u*f);- } --/***********************/-/* Useful Routines */-/***********************/--/* return roots of ax^2+bx+c */-/* stable algebra derived from Numerical Recipes by Press et al.*/-int quadraticRoots(a, b, c, roots)-double a, b, c, *roots;-{-double d, q;-int count = 0;- d = (b*b)-(4*a*c);- if (d < 0.0) { *roots = *(roots+1) = 0.0; return(0); }- q = -0.5 * (b + (SGN(b)*sqrt(d)));- if (a != 0.0) { *roots++ = q/a; count++; }- if (q != 0.0) { *roots++ = c/q; count++; }- return(count);- }---/* generic 1d regula-falsi step. f is function to evaluate */-/* interval known to contain root is given in left, right */-/* returns new estimate */-double RegulaFalsi(f, left, right)-double (*f)(), left, right;-{-double d = (*f)(right) - (*f)(left);- if (d != 0.0) return (right - (*f)(right)*(right-left)/d);- return((left+right)/2.0);- }--/* generic 1d Newton-Raphson step. f is function, df is derivative */-/* x is current best guess for root location. Returns new estimate */-double NewtonRaphson(f, df, x)-double (*f)(), (*df)(), x;-{-double d = (*df)(x);- if (d != 0.0) return (x-((*f)(x)/d));- return(x-1.0);- }---/* hybrid 1d Newton-Raphson/Regula Falsi root finder. */-/* input function f and its derivative df, an interval */-/* left, right known to contain the root, and an error tolerance */-/* Based on Blinn */-double findroot(left, right, tolerance, f, df)-double left, right, tolerance;-double (*f)(), (*df)();-{-double newx = left;- while (ABS((*f)(newx)) > tolerance) {- newx = NewtonRaphson(f, df, newx);- if (newx < left || newx > right) - newx = RegulaFalsi(f, left, right);- if ((*f)(newx) * (*f)(left) <= 0.0) right = newx; - else left = newx;- }- return(newx);- } deleted file mode 100644@@ -1,170 +0,0 @@-/* - * GraphicsGems.h - * Version 1.0 - Andrew Glassner- * from "Graphics Gems", Academic Press, 1990- */--/*- * EULA: The Graphics Gems code is copyright-protected. In other words, - * you cannot claim the text of the code as your own and resell it. - * Using the code is permitted in any program, product, or library, - * non-commercial or commercial. Giving credit is not required, though - * is a nice gesture. The code comes as-is, and if there are any flaws - * or problems with any Gems code, nobody involved with Gems - authors,- * editors, publishers, or webmasters - are to be held responsible. - * Basically, don't be a jerk, and remember that anything free comes - * with no guarantee. - */- -#ifndef GG_H--#define GG_H 1--/*********************/-/* 2d geometry types */-/*********************/--typedef struct Point2Struct { /* 2d point */- double x, y;- } Point2;-typedef Point2 Vector2;--typedef struct IntPoint2Struct { /* 2d integer point */- int x, y;- } IntPoint2;--typedef struct Matrix3Struct { /* 3-by-3 matrix */- double element[3][3];- } Matrix3;--typedef struct Box2dStruct { /* 2d box */- Point2 min, max;- } Box2;- --/*********************/-/* 3d geometry types */-/*********************/--typedef struct Point3Struct { /* 3d point */- double x, y, z;- } Point3;-typedef Point3 Vector3;--typedef struct IntPoint3Struct { /* 3d integer point */- int x, y, z;- } IntPoint3;---typedef struct Matrix4Struct { /* 4-by-4 matrix */- double element[4][4];- } Matrix4;--typedef struct Box3dStruct { /* 3d box */- Point3 min, max;- } Box3;----/***********************/-/* one-argument macros */-/***********************/--/* absolute value of a */-#define ABS(a) (((a)<0) ? -(a) : (a))--/* round a to nearest int */-#define ROUND(a) ((a)>0 ? (int)((a)+0.5) : -(int)(0.5-(a)))--/* take sign of a, either -1, 0, or 1 */-#define ZSGN(a) (((a)<0) ? -1 : (a)>0 ? 1 : 0) --/* take binary sign of a, either -1, or 1 if >= 0 */-#define SGN(a) (((a)<0) ? -1 : 1)--/* shout if something that should be true isn't */-#define ASSERT(x) \-if (!(x)) fprintf(stderr," Assert failed: x\n");--/* square a */-#define SQR(a) ((a)*(a)) ---/***********************/-/* two-argument macros */-/***********************/--/* find minimum of a and b */-#define MIN(a,b) (((a)<(b))?(a):(b)) --/* find maximum of a and b */-#define MAX(a,b) (((a)>(b))?(a):(b)) --/* swap a and b (see Gem by Wyvill) */-#define SWAP(a,b) { a^=b; b^=a; a^=b; }--/* linear interpolation from l (when a=0) to h (when a=1)*/-/* (equal to (a*h)+((1-a)*l) */-#define LERP(a,l,h) ((l)+(((h)-(l))*(a)))--/* clamp the input to the specified range */-#define CLAMP(v,l,h) ((v)<(l) ? (l) : (v) > (h) ? (h) : v)---/****************************/-/* memory allocation macros */-/****************************/--/* create a new instance of a structure (see Gem by Hultquist) */-#define NEWSTRUCT(x) (struct x *)(malloc((unsigned)sizeof(struct x)))--/* create a new instance of a type */-#define NEWTYPE(x) (x *)(malloc((unsigned)sizeof(x)))---/********************/-/* useful constants */-/********************/--#define PI 3.141592 /* the venerable pi */-#define PITIMES2 6.283185 /* 2 * pi */-#define PIOVER2 1.570796 /* pi / 2 */-#define E 2.718282 /* the venerable e */-#define SQRT2 1.414214 /* sqrt(2) */-#define SQRT3 1.732051 /* sqrt(3) */-#define GOLDEN 1.618034 /* the golden ratio */-#define DTOR 0.017453 /* convert degrees to radians */-#define RTOD 57.29578 /* convert radians to degrees */---/************/-/* booleans */-/************/--#define TRUE 1-#define FALSE 0-#define ON 1-#define OFF 0-typedef int boolean; /* boolean data type */-typedef boolean flag; /* flag data type */--extern double V2SquaredLength(), V2Length();-extern double V2Dot(), V2DistanceBetween2Points(); -extern Vector2 *V2Negate(), *V2Normalize(), *V2Scale(), *V2Add(), *V2Sub();-extern Vector2 *V2Lerp(), *V2Combine(), *V2Mul(), *V2MakePerpendicular();-extern Vector2 *V2New(), *V2Duplicate();-extern Point2 *V2MulPointByMatrix();-extern Matrix3 *V2MatMul();--extern double V3SquaredLength(), V3Length();-extern double V3Dot(), V3DistanceBetween2Points();-extern Vector3 *V3Normalize(), *V3Scale(), *V3Add(), *V3Sub();-extern Vector3 *V3Lerp(), *V3Combine(), *V3Mul(), *V3Cross();-extern Vector3 *V3New(), *V3Duplicate();-extern Point3 *V3MulPointByMatrix();-extern Matrix4 *V3MatMul();--extern double RegulaFalsi(), NewtonRaphson(), findroot();--#endif-deleted file mode 100644@@ -1,570 +0,0 @@-/*-An Algorithm for Automatically Fitting Digitized Curves-by Philip J. Schneider-from "Graphics Gems", Academic Press, 1990-Adapted to BirdFont by Johan Mattsson 2015-*/--/*- * EULA: The Graphics Gems code is copyright-protected. In other words, - * you cannot claim the text of the code as your own and resell it. - * Using the code is permitted in any program, product, or library, - * non-commercial or commercial. Giving credit is not required, though - * is a nice gesture. The code comes as-is, and if there are any flaws - * or problems with any Gems code, nobody involved with Gems - authors,- * editors, publishers, or webmasters - are to be held responsible. - * Basically, don't be a jerk, and remember that anything free comes - * with no guarantee. - */--/* fit_cubic.c */ -/* Piecewise cubic fitting code */--#include "GraphicsGems.h" -#include "birdfont.h"-#include <stdio.h>-#include <malloc.h>-#include <math.h>-#include <stdio.h>--typedef Point2 *BezierCurve;--/* Forward declarations */-void FitCurve();-static void FitCubic();-static double *Reparameterize();-static double NewtonRaphsonRootFind();-static Point2 BezierII();-static double B0(), B1(), B2(), B3();-static Vector2 ComputeLeftTangent();-static Vector2 ComputeRightTangent();-static Vector2 ComputeCenterTangent();-static double ComputeMaxError();-static double *ChordLengthParameterize();-static BezierCurve GenerateBezier();-static Vector2 V2AddII();-static Vector2 V2ScaleIII();-static Vector2 V2SubII();--#define MAXPOINTS 1000 /* The most points you can have */--BirdFontPath* simplified_path = NULL;--BirdFontPath* fit_bezier_curve_to_line (BirdFontPath* path, int start, int stop, double error) {- int i, j, n;- Point2* points;- BirdFontPath* result; -- if (simplified_path != NULL) {- printf ("Simplifier is not thread safe.");- return bird_font_path_new ();- }-- n = stop - start + 1;- - if (n > bird_font_path_size (path) || n < 0) {- printf ("Index of of bounds in simplifier");- return bird_font_path_new ();- }- - points = malloc (n * sizeof (Point2)); - - j = 0;- for (i = start; i <= stop; i++) {- points[j].x = bird_font_path_get_x_at (path, i);- points[j].y = bird_font_path_get_y_at (path, i);- j++;- }- - simplified_path = bird_font_path_new ();- FitCurve(points, n, error);- free (points);-- result = simplified_path;- simplified_path = NULL;- return result;-}--DrawBezierCurve(int n, Point2* curve) {- if (n == 3) {- bird_font_path_add_cubic_bezier_points (simplified_path,- curve[0].x, curve[0].y, - curve[1].x, curve[1].y, - curve[2].x, curve[2].y,- curve[3].x, curve[3].y);- } else {- printf ("Expecting three points\n");- }-}--/*- * FitCurve :- * Fit a Bezier curve to a set of digitized points - */-void FitCurve(d, nPts, error)- Point2 *d; /* Array of digitized points */- int nPts; /* Number of digitized points */- double error; /* User-defined error squared */-{- Vector2 tHat1, tHat2; /* Unit tangent vectors at endpoints */-- tHat1 = ComputeLeftTangent(d, 0);- tHat2 = ComputeRightTangent(d, nPts - 1);- FitCubic(d, 0, nPts - 1, tHat1, tHat2, error);-}--/*- * FitCubic :- * Fit a Bezier curve to a (sub)set of digitized points- */-static void FitCubic(d, first, last, tHat1, tHat2, error)- Point2 *d; /* Array of digitized points */- int first, last; /* Indices of first and last pts in region */- Vector2 tHat1, tHat2; /* Unit tangent vectors at endpoints */- double error; /* User-defined error squared */-{- BezierCurve bezCurve; /*Control points of fitted Bezier curve*/- double *u; /* Parameter values for point */- double *uPrime; /* Improved parameter values */- double maxError; /* Maximum fitting error */- int splitPoint; /* Point to split point set at */- int nPts; /* Number of points in subset */- double iterationError; /*Error below which you try iterating */- int maxIterations = 4; /* Max times to try iterating */- Vector2 tHatCenter; /* Unit tangent vector at splitPoint */- int i; -- iterationError = error * error;- nPts = last - first + 1;-- /* Use heuristic if region only has two points in it */- if (nPts == 2) {- double dist = V2DistanceBetween2Points(&d[last], &d[first]) / 3.0;-- bezCurve = (Point2 *)malloc(4 * sizeof(Point2));- bezCurve[0] = d[first];- bezCurve[3] = d[last];- V2Add(&bezCurve[0], V2Scale(&tHat1, dist), &bezCurve[1]);- V2Add(&bezCurve[3], V2Scale(&tHat2, dist), &bezCurve[2]);- DrawBezierCurve(3, bezCurve);- free((void *)bezCurve);- return;- }-- /* Parameterize points, and attempt to fit curve */- u = ChordLengthParameterize(d, first, last);- bezCurve = GenerateBezier(d, first, last, u, tHat1, tHat2);-- /* Find max deviation of points to fitted curve */- maxError = ComputeMaxError(d, first, last, bezCurve, u, &splitPoint);- if (maxError < error) {- DrawBezierCurve(3, bezCurve);- free((void *)u);- free((void *)bezCurve);- return;- }--- /* If error not too large, try some reparameterization */- /* and iteration */- if (maxError < iterationError) {- for (i = 0; i < maxIterations; i++) {- uPrime = Reparameterize(d, first, last, u, bezCurve);- free((void *)bezCurve);- bezCurve = GenerateBezier(d, first, last, uPrime, tHat1, tHat2);- maxError = ComputeMaxError(d, first, last,- bezCurve, uPrime, &splitPoint);- if (maxError < error) {- DrawBezierCurve(3, bezCurve);- free((void *)u);- free((void *)bezCurve);- free((void *)uPrime);- return;- }- free((void *)u);- u = uPrime;- }- }-- /* Fitting failed -- split at max error point and fit recursively */- free((void *)u);- free((void *)bezCurve);- tHatCenter = ComputeCenterTangent(d, splitPoint);- FitCubic(d, first, splitPoint, tHat1, tHatCenter, error);- V2Negate(&tHatCenter);- FitCubic(d, splitPoint, last, tHatCenter, tHat2, error);-}---/*- * GenerateBezier :- * Use least-squares method to find Bezier control points for region.- *- */-static BezierCurve GenerateBezier(d, first, last, uPrime, tHat1, tHat2)- Point2 *d; /* Array of digitized points */- int first, last; /* Indices defining region */- double *uPrime; /* Parameter values for region */- Vector2 tHat1, tHat2; /* Unit tangents at endpoints */-{- int i;- Vector2 A[MAXPOINTS][2]; /* Precomputed rhs for eqn */- int nPts; /* Number of pts in sub-curve */- double C[2][2]; /* Matrix C */- double X[2]; /* Matrix X */- double det_C0_C1, /* Determinants of matrices */- det_C0_X,- det_X_C1;- double alpha_l, /* Alpha values, left and right */- alpha_r;- Vector2 tmp; /* Utility variable */- BezierCurve bezCurve; /* RETURN bezier curve ctl pts */-- bezCurve = (Point2 *)malloc(4 * sizeof(Point2));- nPts = last - first + 1;-- - /* Compute the A's */- for (i = 0; i < nPts; i++) {- Vector2 v1, v2;- v1 = tHat1;- v2 = tHat2;- V2Scale(&v1, B1(uPrime[i]));- V2Scale(&v2, B2(uPrime[i]));- A[i][0] = v1;- A[i][1] = v2;- }-- /* Create the C and X matrices */- C[0][0] = 0.0;- C[0][1] = 0.0;- C[1][0] = 0.0;- C[1][1] = 0.0;- X[0] = 0.0;- X[1] = 0.0;-- for (i = 0; i < nPts; i++) {- C[0][0] += V2Dot(&A[i][0], &A[i][0]);- C[0][1] += V2Dot(&A[i][0], &A[i][1]);-/* C[1][0] += V2Dot(&A[i][0], &A[i][1]);*/ - C[1][0] = C[0][1];- C[1][1] += V2Dot(&A[i][1], &A[i][1]);-- tmp = V2SubII(d[first + i],- V2AddII(- V2ScaleIII(d[first], B0(uPrime[i])),- V2AddII(- V2ScaleIII(d[first], B1(uPrime[i])),- V2AddII(- V2ScaleIII(d[last], B2(uPrime[i])),- V2ScaleIII(d[last], B3(uPrime[i]))))));- -- X[0] += V2Dot(&A[i][0], &tmp);- X[1] += V2Dot(&A[i][1], &tmp);- }-- /* Compute the determinants of C and X */- det_C0_C1 = C[0][0] * C[1][1] - C[1][0] * C[0][1];- det_C0_X = C[0][0] * X[1] - C[1][0] * X[0];- det_X_C1 = X[0] * C[1][1] - X[1] * C[0][1];-- /* Finally, derive alpha values */- alpha_l = (det_C0_C1 == 0) ? 0.0 : det_X_C1 / det_C0_C1;- alpha_r = (det_C0_C1 == 0) ? 0.0 : det_C0_X / det_C0_C1;-- /* If alpha negative, use the Wu/Barsky heuristic (see text) */- /* (if alpha is 0, you get coincident control points that lead to- * divide by zero in any subsequent NewtonRaphsonRootFind() call. */- double segLength = V2DistanceBetween2Points(&d[last], &d[first]);- double epsilon = 1.0e-6 * segLength;- if (alpha_l < epsilon || alpha_r < epsilon)- {- /* fall back on standard (probably inaccurate) formula, and subdivide further if needed. */- double dist = segLength / 3.0;- bezCurve[0] = d[first];- bezCurve[3] = d[last];- V2Add(&bezCurve[0], V2Scale(&tHat1, dist), &bezCurve[1]);- V2Add(&bezCurve[3], V2Scale(&tHat2, dist), &bezCurve[2]);- return (bezCurve);- }-- /* First and last control points of the Bezier curve are */- /* positioned exactly at the first and last data points */- /* Control points 1 and 2 are positioned an alpha distance out */- /* on the tangent vectors, left and right, respectively */- bezCurve[0] = d[first];- bezCurve[3] = d[last];- V2Add(&bezCurve[0], V2Scale(&tHat1, alpha_l), &bezCurve[1]);- V2Add(&bezCurve[3], V2Scale(&tHat2, alpha_r), &bezCurve[2]);- return (bezCurve);-}---/*- * Reparameterize:- * Given set of points and their parameterization, try to find- * a better parameterization.- *- */-static double *Reparameterize(d, first, last, u, bezCurve)- Point2 *d; /* Array of digitized points */- int first, last; /* Indices defining region */- double *u; /* Current parameter values */- BezierCurve bezCurve; /* Current fitted curve */-{- int nPts = last-first+1; - int i;- double *uPrime; /* New parameter values */-- uPrime = (double *)malloc(nPts * sizeof(double));- for (i = first; i <= last; i++) {- uPrime[i-first] = NewtonRaphsonRootFind(bezCurve, d[i], u[i-- first]);- }- return (uPrime);-}----/*- * NewtonRaphsonRootFind :- * Use Newton-Raphson iteration to find better root.- */-static double NewtonRaphsonRootFind(Q, P, u)- BezierCurve Q; /* Current fitted curve */- Point2 P; /* Digitized point */- double u; /* Parameter value for "P" */-{- double numerator, denominator;- Point2 Q1[3], Q2[2]; /* Q' and Q'' */- Point2 Q_u, Q1_u, Q2_u; /*u evaluated at Q, Q', & Q'' */- double uPrime; /* Improved u */- int i;- - /* Compute Q(u) */- Q_u = BezierII(3, Q, u);- - /* Generate control vertices for Q' */- for (i = 0; i <= 2; i++) {- Q1[i].x = (Q[i+1].x - Q[i].x) * 3.0;- Q1[i].y = (Q[i+1].y - Q[i].y) * 3.0;- }- - /* Generate control vertices for Q'' */- for (i = 0; i <= 1; i++) {- Q2[i].x = (Q1[i+1].x - Q1[i].x) * 2.0;- Q2[i].y = (Q1[i+1].y - Q1[i].y) * 2.0;- }- - /* Compute Q'(u) and Q''(u) */- Q1_u = BezierII(2, Q1, u);- Q2_u = BezierII(1, Q2, u);- - /* Compute f(u)/f'(u) */- numerator = (Q_u.x - P.x) * (Q1_u.x) + (Q_u.y - P.y) * (Q1_u.y);- denominator = (Q1_u.x) * (Q1_u.x) + (Q1_u.y) * (Q1_u.y) +- (Q_u.x - P.x) * (Q2_u.x) + (Q_u.y - P.y) * (Q2_u.y);- if (denominator == 0.0f) return u;-- /* u = u - f(u)/f'(u) */- uPrime = u - (numerator/denominator);- return (uPrime);-}-- - -/*- * Bezier :- * Evaluate a Bezier curve at a particular parameter value- * - */-static Point2 BezierII(degree, V, t)- int degree; /* The degree of the bezier curve */- Point2 *V; /* Array of control points */- double t; /* Parametric value to find point for */-{- int i, j; - Point2 Q; /* Point on curve at parameter t */- Point2 *Vtemp; /* Local copy of control points */-- /* Copy array */- Vtemp = (Point2 *)malloc((unsigned)((degree+1) - * sizeof (Point2)));- for (i = 0; i <= degree; i++) {- Vtemp[i] = V[i];- }-- /* Triangle computation */- for (i = 1; i <= degree; i++) { - for (j = 0; j <= degree-i; j++) {- Vtemp[j].x = (1.0 - t) * Vtemp[j].x + t * Vtemp[j+1].x;- Vtemp[j].y = (1.0 - t) * Vtemp[j].y + t * Vtemp[j+1].y;- }- }-- Q = Vtemp[0];- free((void *)Vtemp);- return Q;-}---/*- * B0, B1, B2, B3 :- * Bezier multipliers- */-static double B0(u)- double u;-{- double tmp = 1.0 - u;- return (tmp * tmp * tmp);-}---static double B1(u)- double u;-{- double tmp = 1.0 - u;- return (3 * u * (tmp * tmp));-}--static double B2(u)- double u;-{- double tmp = 1.0 - u;- return (3 * u * u * tmp);-}--static double B3(u)- double u;-{- return (u * u * u);-}----/*- * ComputeLeftTangent, ComputeRightTangent, ComputeCenterTangent :- *Approximate unit tangents at endpoints and "center" of digitized curve- */-static Vector2 ComputeLeftTangent(d, end)- Point2 *d; /* Digitized points*/- int end; /* Index to "left" end of region */-{- Vector2 tHat1;- tHat1 = V2SubII(d[end+1], d[end]);- tHat1 = *V2Normalize(&tHat1);- return tHat1;-}--static Vector2 ComputeRightTangent(d, end)- Point2 *d; /* Digitized points */- int end; /* Index to "right" end of region */-{- Vector2 tHat2;- tHat2 = V2SubII(d[end-1], d[end]);- tHat2 = *V2Normalize(&tHat2);- return tHat2;-}---static Vector2 ComputeCenterTangent(d, center)- Point2 *d; /* Digitized points */- int center; /* Index to point inside region */-{- Vector2 V1, V2, tHatCenter;-- V1 = V2SubII(d[center-1], d[center]);- V2 = V2SubII(d[center], d[center+1]);- tHatCenter.x = (V1.x + V2.x)/2.0;- tHatCenter.y = (V1.y + V2.y)/2.0;- tHatCenter = *V2Normalize(&tHatCenter);- return tHatCenter;-}---/*- * ChordLengthParameterize :- * Assign parameter values to digitized points - * using relative distances between points.- */-static double *ChordLengthParameterize(d, first, last)- Point2 *d; /* Array of digitized points */- int first, last; /* Indices defining region */-{- int i; - double *u; /* Parameterization */-- u = (double *)malloc((unsigned)(last-first+1) * sizeof(double));-- u[0] = 0.0;- for (i = first+1; i <= last; i++) {- u[i-first] = u[i-first-1] +- V2DistanceBetween2Points(&d[i], &d[i-1]);- }-- for (i = first + 1; i <= last; i++) {- u[i-first] = u[i-first] / u[last-first];- }-- return(u);-}-----/*- * ComputeMaxError :- * Find the maximum squared distance of digitized points- * to fitted curve.-*/-static double ComputeMaxError(d, first, last, bezCurve, u, splitPoint)- Point2 *d; /* Array of digitized points */- int first, last; /* Indices defining region */- BezierCurve bezCurve; /* Fitted Bezier curve */- double *u; /* Parameterization of points */- int *splitPoint; /* Point of maximum error */-{- int i;- double maxDist; /* Maximum error */- double dist; /* Current error */- Point2 P; /* Point on curve */- Vector2 v; /* Vector from point to curve */-- *splitPoint = (last - first + 1)/2;- maxDist = 0.0;- for (i = first + 1; i < last; i++) {- P = BezierII(3, bezCurve, u[i-first]);- v = V2SubII(P, d[i]);- dist = V2SquaredLength(&v);- if (dist >= maxDist) {- maxDist = dist;- *splitPoint = i;- }- }- return (maxDist);-}-static Vector2 V2AddII(a, b)- Vector2 a, b;-{- Vector2 c;- c.x = a.x + b.x; c.y = a.y + b.y;- return (c);-}-static Vector2 V2ScaleIII(v, s)- Vector2 v;- double s;-{- Vector2 result;- result.x = v.x * s; result.y = v.y * s;- return (result);-}--static Vector2 V2SubII(a, b)- Vector2 a, b;-{- Vector2 c;- c.x = a.x - b.x; c.y = a.y - b.y;- return (c);-}-@@ -17,9 +17,6 @@ using Math; namespace BirdFont { -[CCode (cname = "fit_bezier_curve_to_line")]-public extern static Path fit_bezier_curve_to_line (Path lines, int start, int stop, double error);- public enum Direction { CLOCKWISE, COUNTER_CLOCKWISE@@ -2370,20 +2367,6 @@ public class Path { start.get_right_handle ().move_to_coordinate (x1, y1); end.get_left_handle ().move_to_coordinate (x2, y2); }- - public int size () {- return points.size;- }- - public double get_x_at (int i) {- return_val_if_fail (0 <= i < points.size, 0);- return points.get (i).x;- }-- public double get_y_at (int i) {- return_val_if_fail (0 <= i < points.size, 0);- return points.get (i).y;- } } }@@ -157,10 +157,10 @@ public class StrokeTool : Tool { stop = p.points.size - 1; } - ep_start = p.points.get (start);+ ep_start = p.points.get (start); ep = p.points.get (stop); - segment = fit_bezier_curve_to_line (p, start, stop, 0.0007);+ segment = fit_bezier_path (p, start, stop, 0.0007); if (stroke_selected) { // FIXME: DELETE ((!) BirdFont.get_current_font ().get_glyph ("l")).add_path (segment.copy ());@@ -230,6 +230,49 @@ public class StrokeTool : Tool { return simplified; } + static Path fit_bezier_path (Path p, int start, int stop, double error) {+ int index, size;+ Path simplified;+ double[] lines;+ double[] result;+ EditPoint ep;+ + simplified = new Path ();++ return_val_if_fail (0 <= start < p.points.size, simplified);+ return_val_if_fail (0 <= stop < p.points.size, simplified);++ size = stop - start + 1;+ lines = new double[2 * size];+ + index = 0;+ + for (int i = start; i <= stop; i++) {+ ep = p.points.get (i);+ lines[index] = ep.x;+ index++;+ + lines[index] = ep.y;+ index++;+ }+ + return_val_if_fail (2 * size == index, new Path ());+ + Gems.fit_bezier_curve_to_line (lines, error, out result);++ return_val_if_fail (!is_null (result), simplified);+ + for (int i = 0; i + 7 < result.length; i += 8) {+ simplified.add_cubic_bezier_points (+ result[i], result[i + 1],+ result[i + 2], result[i + 3],+ result[i + 4], result[i + 5],+ result[i + 6], result[i + 7]);+ }+ + return simplified;+ }+ static bool segment_is_line (Path p, int start, int stop) { EditPoint p0, p1, p2; new file mode 100644@@ -0,0 +1,3 @@+class FitCubic {++}new file mode 100644@@ -0,0 +1,458 @@+/* +2d and 3d Vector C Library +by Andrew Glassner+from "Graphics Gems", Academic Press, 1990+*/++/*+ * EULA: The Graphics Gems code is copyright-protected. In other words, + * you cannot claim the text of the code as your own and resell it. + * Using the code is permitted in any program, product, or library, + * non-commercial or commercial. Giving credit is not required, though + * is a nice gesture. The code comes as-is, and if there are any flaws + * or problems with any Gems code, nobody involved with Gems - authors,+ * editors, publishers, or webmasters - are to be held responsible. + * Basically, don't be a jerk, and remember that anything free comes + * with no guarantee. + */++#include <math.h>+#include "GraphicsGems.h"++/******************/+/* 2d Library */+/******************/++/* returns squared length of input vector */ +double V2SquaredLength(a) +Vector2 *a;+{ return((a->x * a->x)+(a->y * a->y));+ }+ +/* returns length of input vector */+double V2Length(a) +Vector2 *a;+{+ return(sqrt(V2SquaredLength(a)));+ }+ +/* negates the input vector and returns it */+Vector2 *V2Negate(v) +Vector2 *v;+{+ v->x = -v->x; v->y = -v->y;+ return(v);+ }++/* normalizes the input vector and returns it */+Vector2 *V2Normalize(v) +Vector2 *v;+{+double len = V2Length(v);+ if (len != 0.0) { v->x /= len; v->y /= len; }+ return(v);+ }+++/* scales the input vector to the new length and returns it */+Vector2 *V2Scale(v, newlen) +Vector2 *v;+double newlen;+{+double len = V2Length(v);+ if (len != 0.0) { v->x *= newlen/len; v->y *= newlen/len; }+ return(v);+ }++/* return vector sum c = a+b */+Vector2 *V2Add(a, b, c)+Vector2 *a, *b, *c;+{+ c->x = a->x+b->x; c->y = a->y+b->y;+ return(c);+ }+ +/* return vector difference c = a-b */+Vector2 *V2Sub(a, b, c)+Vector2 *a, *b, *c;+{+ c->x = a->x-b->x; c->y = a->y-b->y;+ return(c);+ }++/* return the dot product of vectors a and b */+double V2Dot(a, b) +Vector2 *a, *b; +{+ return((a->x*b->x)+(a->y*b->y));+ }++/* linearly interpolate between vectors by an amount alpha */+/* and return the resulting vector. */+/* When alpha=0, result=lo. When alpha=1, result=hi. */+Vector2 *V2Lerp(lo, hi, alpha, result) +Vector2 *lo, *hi, *result; +double alpha;+{+ result->x = LERP(alpha, lo->x, hi->x);+ result->y = LERP(alpha, lo->y, hi->y);+ return(result);+ }+++/* make a linear combination of two vectors and return the result. */+/* result = (a * ascl) + (b * bscl) */+Vector2 *V2Combine (a, b, result, ascl, bscl) +Vector2 *a, *b, *result;+double ascl, bscl;+{+ result->x = (ascl * a->x) + (bscl * b->x);+ result->y = (ascl * a->y) + (bscl * b->y);+ return(result);+ }++/* multiply two vectors together component-wise */+Vector2 *V2Mul (a, b, result) +Vector2 *a, *b, *result;+{+ result->x = a->x * b->x;+ result->y = a->y * b->y;+ return(result);+ }++/* return the distance between two points */+double V2DistanceBetween2Points(a, b)+Point2 *a, *b;+{+double dx = a->x - b->x;+double dy = a->y - b->y;+ return(sqrt((dx*dx)+(dy*dy)));+ }++/* return the vector perpendicular to the input vector a */+Vector2 *V2MakePerpendicular(a, ap)+Vector2 *a, *ap;+{+ ap->x = -a->y;+ ap->y = a->x;+ return(ap);+ }++/* create, initialize, and return a new vector */+Vector2 *V2New(x, y)+double x, y;+{+Vector2 *v = NEWTYPE(Vector2);+ v->x = x; v->y = y; + return(v);+ }+ ++/* create, initialize, and return a duplicate vector */+Vector2 *V2Duplicate(a)+Vector2 *a;+{+Vector2 *v = NEWTYPE(Vector2);+ v->x = a->x; v->y = a->y; + return(v);+ }+ +/* multiply a point by a matrix and return the transformed point */+Point2 *V2MulPointByMatrix(p, m)+Point2 *p;+Matrix3 *m;+{+double w;+Point2 ptmp;+ ptmp.x = (p->x * m->element[0][0]) + + (p->y * m->element[1][0]) + m->element[2][0];+ ptmp.y = (p->x * m->element[0][1]) + + (p->y * m->element[1][1]) + m->element[2][1];+ w = (p->x * m->element[0][2]) + + (p->y * m->element[1][2]) + m->element[2][2];+ if (w != 0.0) { ptmp.x /= w; ptmp.y /= w; }+ *p = ptmp;+ return(p);+ }++/* multiply together matrices c = ab */+/* note that c must not point to either of the input matrices */+Matrix3 *V2MatMul(a, b, c)+Matrix3 *a, *b, *c;+{+int i, j, k;+ for (i=0; i<3; i++) {+ for (j=0; j<3; j++) {+ c->element[i][j] = 0;+ for (k=0; k<3; k++) c->element[i][j] += + a->element[i][k] * b->element[k][j];+ }+ }+ return(c);+ }+++++/******************/+/* 3d Library */+/******************/+ +/* returns squared length of input vector */ +double V3SquaredLength(a) +Vector3 *a;+{+ return((a->x * a->x)+(a->y * a->y)+(a->z * a->z));+ }++/* returns length of input vector */+double V3Length(a) +Vector3 *a;+{+ return(sqrt(V3SquaredLength(a)));+ }++/* negates the input vector and returns it */+Vector3 *V3Negate(v) +Vector3 *v;+{+ v->x = -v->x; v->y = -v->y; v->z = -v->z;+ return(v);+ }++/* normalizes the input vector and returns it */+Vector3 *V3Normalize(v) +Vector3 *v;+{+double len = V3Length(v);+ if (len != 0.0) { v->x /= len; v->y /= len; v->z /= len; }+ return(v);+ }++/* scales the input vector to the new length and returns it */+Vector3 *V3Scale(v, newlen) +Vector3 *v;+double newlen;+{+double len = V3Length(v);+ if (len != 0.0) {+ v->x *= newlen/len; v->y *= newlen/len; v->z *= newlen/len;+ }+ return(v);+ }+++/* return vector sum c = a+b */+Vector3 *V3Add(a, b, c)+Vector3 *a, *b, *c;+{+ c->x = a->x+b->x; c->y = a->y+b->y; c->z = a->z+b->z;+ return(c);+ }+ +/* return vector difference c = a-b */+Vector3 *V3Sub(a, b, c)+Vector3 *a, *b, *c;+{+ c->x = a->x-b->x; c->y = a->y-b->y; c->z = a->z-b->z;+ return(c);+ }++/* return the dot product of vectors a and b */+double V3Dot(a, b) +Vector3 *a, *b; +{+ return((a->x*b->x)+(a->y*b->y)+(a->z*b->z));+ }++/* linearly interpolate between vectors by an amount alpha */+/* and return the resulting vector. */+/* When alpha=0, result=lo. When alpha=1, result=hi. */+Vector3 *V3Lerp(lo, hi, alpha, result) +Vector3 *lo, *hi, *result; +double alpha;+{+ result->x = LERP(alpha, lo->x, hi->x);+ result->y = LERP(alpha, lo->y, hi->y);+ result->z = LERP(alpha, lo->z, hi->z);+ return(result);+ }++/* make a linear combination of two vectors and return the result. */+/* result = (a * ascl) + (b * bscl) */+Vector3 *V3Combine (a, b, result, ascl, bscl) +Vector3 *a, *b, *result;+double ascl, bscl;+{+ result->x = (ascl * a->x) + (bscl * b->x);+ result->y = (ascl * a->y) + (bscl * b->y);+ result->z = (ascl * a->z) + (bscl * b->z);+ return(result);+ }+++/* multiply two vectors together component-wise and return the result */+Vector3 *V3Mul (a, b, result) +Vector3 *a, *b, *result;+{+ result->x = a->x * b->x;+ result->y = a->y * b->y;+ result->z = a->z * b->z;+ return(result);+ }++/* return the distance between two points */+double V3DistanceBetween2Points(a, b)+Point3 *a, *b;+{+double dx = a->x - b->x;+double dy = a->y - b->y;+double dz = a->z - b->z;+ return(sqrt((dx*dx)+(dy*dy)+(dz*dz)));+ }++/* return the cross product c = a cross b */+Vector3 *V3Cross(a, b, c)+Vector3 *a, *b, *c;+{+ c->x = (a->y*b->z) - (a->z*b->y);+ c->y = (a->z*b->x) - (a->x*b->z);+ c->z = (a->x*b->y) - (a->y*b->x);+ return(c);+ }++/* create, initialize, and return a new vector */+Vector3 *V3New(x, y, z)+double x, y, z;+{+Vector3 *v = NEWTYPE(Vector3);+ v->x = x; v->y = y; v->z = z;+ return(v);+ }++/* create, initialize, and return a duplicate vector */+Vector3 *V3Duplicate(a)+Vector3 *a;+{+Vector3 *v = NEWTYPE(Vector3);+ v->x = a->x; v->y = a->y; v->z = a->z;+ return(v);+ }++ +/* multiply a point by a matrix and return the transformed point */+Point3 *V3MulPointByMatrix(p, m)+Point3 *p;+Matrix4 *m;+{+double w;+Point3 ptmp;+ ptmp.x = (p->x * m->element[0][0]) + (p->y * m->element[1][0]) + + (p->z * m->element[2][0]) + m->element[3][0];+ ptmp.y = (p->x * m->element[0][1]) + (p->y * m->element[1][1]) + + (p->z * m->element[2][1]) + m->element[3][1];+ ptmp.z = (p->x * m->element[0][2]) + (p->y * m->element[1][2]) + + (p->z * m->element[2][2]) + m->element[3][2];+ w = (p->x * m->element[0][3]) + (p->y * m->element[1][3]) + + (p->z * m->element[2][3]) + m->element[3][3];+ if (w != 0.0) { ptmp.x /= w; ptmp.y /= w; ptmp.z /= w; }+ *p = ptmp;+ return(p);+ }++/* multiply together matrices c = ab */+/* note that c must not point to either of the input matrices */+Matrix4 *V3MatMul(a, b, c)+Matrix4 *a, *b, *c;+{+int i, j, k;+ for (i=0; i<4; i++) {+ for (j=0; j<4; j++) {+ c->element[i][j] = 0;+ for (k=0; k<4; k++) c->element[i][j] += + a->element[i][k] * b->element[k][j];+ }+ }+ return(c);+ }++/* binary greatest common divisor by Silver and Terzian. See Knuth */+/* both inputs must be >= 0 */+gcd(u, v)+int u, v;+{+int t, f;+ if ((u<0) || (v<0)) return(1); /* error if u<0 or v<0 */+ f = 1;+ while ((0 == (u%2)) && (0 == (v%2))) {+ u>>=1; v>>=1, f*=2;+ }+ if (u&01) { t = -v; goto B4; } else { t = u; }+ B3: if (t > 0) { t >>= 1; } else { t = -((-t) >> 1); }+ B4: if (0 == (t%2)) goto B3;++ if (t > 0) u = t; else v = -t;+ if (0 != (t = u - v)) goto B3;+ return(u*f);+ } ++/***********************/+/* Useful Routines */+/***********************/++/* return roots of ax^2+bx+c */+/* stable algebra derived from Numerical Recipes by Press et al.*/+int quadraticRoots(a, b, c, roots)+double a, b, c, *roots;+{+double d, q;+int count = 0;+ d = (b*b)-(4*a*c);+ if (d < 0.0) { *roots = *(roots+1) = 0.0; return(0); }+ q = -0.5 * (b + (SGN(b)*sqrt(d)));+ if (a != 0.0) { *roots++ = q/a; count++; }+ if (q != 0.0) { *roots++ = c/q; count++; }+ return(count);+ }+++/* generic 1d regula-falsi step. f is function to evaluate */+/* interval known to contain root is given in left, right */+/* returns new estimate */+double RegulaFalsi(f, left, right)+double (*f)(), left, right;+{+double d = (*f)(right) - (*f)(left);+ if (d != 0.0) return (right - (*f)(right)*(right-left)/d);+ return((left+right)/2.0);+ }++/* generic 1d Newton-Raphson step. f is function, df is derivative */+/* x is current best guess for root location. Returns new estimate */+double NewtonRaphson(f, df, x)+double (*f)(), (*df)(), x;+{+double d = (*df)(x);+ if (d != 0.0) return (x-((*f)(x)/d));+ return(x-1.0);+ }+++/* hybrid 1d Newton-Raphson/Regula Falsi root finder. */+/* input function f and its derivative df, an interval */+/* left, right known to contain the root, and an error tolerance */+/* Based on Blinn */+double findroot(left, right, tolerance, f, df)+double left, right, tolerance;+double (*f)(), (*df)();+{+double newx = left;+ while (ABS((*f)(newx)) > tolerance) {+ newx = NewtonRaphson(f, df, newx);+ if (newx < left || newx > right) + newx = RegulaFalsi(f, left, right);+ if ((*f)(newx) * (*f)(left) <= 0.0) right = newx; + else left = newx;+ }+ return(newx);+ } new file mode 100644@@ -0,0 +1,170 @@+/* + * GraphicsGems.h + * Version 1.0 - Andrew Glassner+ * from "Graphics Gems", Academic Press, 1990+ */++/*+ * EULA: The Graphics Gems code is copyright-protected. In other words, + * you cannot claim the text of the code as your own and resell it. + * Using the code is permitted in any program, product, or library, + * non-commercial or commercial. Giving credit is not required, though + * is a nice gesture. The code comes as-is, and if there are any flaws + * or problems with any Gems code, nobody involved with Gems - authors,+ * editors, publishers, or webmasters - are to be held responsible. + * Basically, don't be a jerk, and remember that anything free comes + * with no guarantee. + */+ +#ifndef GG_H++#define GG_H 1++/*********************/+/* 2d geometry types */+/*********************/++typedef struct Point2Struct { /* 2d point */+ double x, y;+ } Point2;+typedef Point2 Vector2;++typedef struct IntPoint2Struct { /* 2d integer point */+ int x, y;+ } IntPoint2;++typedef struct Matrix3Struct { /* 3-by-3 matrix */+ double element[3][3];+ } Matrix3;++typedef struct Box2dStruct { /* 2d box */+ Point2 min, max;+ } Box2;+ ++/*********************/+/* 3d geometry types */+/*********************/++typedef struct Point3Struct { /* 3d point */+ double x, y, z;+ } Point3;+typedef Point3 Vector3;++typedef struct IntPoint3Struct { /* 3d integer point */+ int x, y, z;+ } IntPoint3;+++typedef struct Matrix4Struct { /* 4-by-4 matrix */+ double element[4][4];+ } Matrix4;++typedef struct Box3dStruct { /* 3d box */+ Point3 min, max;+ } Box3;++++/***********************/+/* one-argument macros */+/***********************/++/* absolute value of a */+#define ABS(a) (((a)<0) ? -(a) : (a))++/* round a to nearest int */+#define ROUND(a) ((a)>0 ? (int)((a)+0.5) : -(int)(0.5-(a)))++/* take sign of a, either -1, 0, or 1 */+#define ZSGN(a) (((a)<0) ? -1 : (a)>0 ? 1 : 0) ++/* take binary sign of a, either -1, or 1 if >= 0 */+#define SGN(a) (((a)<0) ? -1 : 1)++/* shout if something that should be true isn't */+#define ASSERT(x) \+if (!(x)) fprintf(stderr," Assert failed: x\n");++/* square a */+#define SQR(a) ((a)*(a)) +++/***********************/+/* two-argument macros */+/***********************/++/* find minimum of a and b */+#define MIN(a,b) (((a)<(b))?(a):(b)) ++/* find maximum of a and b */+#define MAX(a,b) (((a)>(b))?(a):(b)) ++/* swap a and b (see Gem by Wyvill) */+#define SWAP(a,b) { a^=b; b^=a; a^=b; }++/* linear interpolation from l (when a=0) to h (when a=1)*/+/* (equal to (a*h)+((1-a)*l) */+#define LERP(a,l,h) ((l)+(((h)-(l))*(a)))++/* clamp the input to the specified range */+#define CLAMP(v,l,h) ((v)<(l) ? (l) : (v) > (h) ? (h) : v)+++/****************************/+/* memory allocation macros */+/****************************/++/* create a new instance of a structure (see Gem by Hultquist) */+#define NEWSTRUCT(x) (struct x *)(malloc((unsigned)sizeof(struct x)))++/* create a new instance of a type */+#define NEWTYPE(x) (x *)(malloc((unsigned)sizeof(x)))+++/********************/+/* useful constants */+/********************/++#define PI 3.141592 /* the venerable pi */+#define PITIMES2 6.283185 /* 2 * pi */+#define PIOVER2 1.570796 /* pi / 2 */+#define E 2.718282 /* the venerable e */+#define SQRT2 1.414214 /* sqrt(2) */+#define SQRT3 1.732051 /* sqrt(3) */+#define GOLDEN 1.618034 /* the golden ratio */+#define DTOR 0.017453 /* convert degrees to radians */+#define RTOD 57.29578 /* convert radians to degrees */+++/************/+/* booleans */+/************/++#define TRUE 1+#define FALSE 0+#define ON 1+#define OFF 0+typedef int boolean; /* boolean data type */+typedef boolean flag; /* flag data type */++extern double V2SquaredLength(), V2Length();+extern double V2Dot(), V2DistanceBetween2Points(); +extern Vector2 *V2Negate(), *V2Normalize(), *V2Scale(), *V2Add(), *V2Sub();+extern Vector2 *V2Lerp(), *V2Combine(), *V2Mul(), *V2MakePerpendicular();+extern Vector2 *V2New(), *V2Duplicate();+extern Point2 *V2MulPointByMatrix();+extern Matrix3 *V2MatMul();++extern double V3SquaredLength(), V3Length();+extern double V3Dot(), V3DistanceBetween2Points();+extern Vector3 *V3Normalize(), *V3Scale(), *V3Add(), *V3Sub();+extern Vector3 *V3Lerp(), *V3Combine(), *V3Mul(), *V3Cross();+extern Vector3 *V3New(), *V3Duplicate();+extern Point3 *V3MulPointByMatrix();+extern Matrix4 *V3MatMul();++extern double RegulaFalsi(), NewtonRaphson(), findroot();++#endif+new file mode 100644@@ -0,0 +1,617 @@+/*+An Algorithm for Automatically Fitting Digitized Curves+by Philip J. Schneider+from "Graphics Gems", Academic Press, 1990+Adapted to BirdFont by Johan Mattsson 2015+*/++/*+ * EULA: The Graphics Gems code is copyright-protected. In other words, + * you cannot claim the text of the code as your own and resell it. + * Using the code is permitted in any program, product, or library, + * non-commercial or commercial. Giving credit is not required, though + * is a nice gesture. The code comes as-is, and if there are any flaws + * or problems with any Gems code, nobody involved with Gems - authors,+ * editors, publishers, or webmasters - are to be held responsible. + * Basically, don't be a jerk, and remember that anything free comes + * with no guarantee. + */++/* fit_cubic.c */ +/* Piecewise cubic fitting code */++#include "GraphicsGems.h"+#include <stdio.h>+#include <malloc.h>+#include <math.h>+#include <stdio.h>++typedef Point2 *BezierCurve;++/* Forward declarations */+void FitCurve();+static void FitCubic();+static double *Reparameterize();+static double NewtonRaphsonRootFind();+static Point2 BezierII();+static double B0(), B1(), B2(), B3();+static Vector2 ComputeLeftTangent();+static Vector2 ComputeRightTangent();+static Vector2 ComputeCenterTangent();+static double ComputeMaxError();+static double *ChordLengthParameterize();+static BezierCurve GenerateBezier();+static Vector2 V2AddII();+static Vector2 V2ScaleIII();+static Vector2 V2SubII();++#define MAXPOINTS 1000 /* The most points you can have */++int simplified_path_buffer_size = 0;+int simplified_path_size = 0;+double* simplified_path = NULL;++/** Generates an a bezier path and returns the length of the output array. */ +void fit_bezier_curve_to_line (+ double* lines,+ int lines_size,+ double error,+ double** bezier_path, /** generated Bezier curves */+ int* bezier_path_size)+{+ int i, j;+ Point2* points;+ int npoints;++ if (npoints % 2 != 0) {+ fprintf (stderr, "Odd number of coordinates in fit_bezier_curve_to_line."); + return;+ }++ if (lines == NULL || lines_size == 0) {+ fprintf (stderr, "No lines in fit_bezier_curve_to_line.");+ return;+ }++ if (simplified_path != NULL) {+ fprintf (stderr, "Path simplification is alreading running.");+ return; + }+ + if (bezier_path == NULL) {+ fprintf (stderr, "No destination for output buffer in fit_bezier_curve_to_line");+ return; + }++ if (bezier_path_size == NULL) {+ fprintf (stderr, "No destination for bezier_path_size in fit_bezier_curve_to_line");+ return; + }+ + npoints = lines_size / 2;+ + printf ("npoints: %d alloc: %d\n", npoints, npoints * sizeof (Point2));+ + points = malloc (npoints * sizeof (Point2));++ j = 0;+ for (i = 0; i < npoints; i++) {+ points[i].x = lines[j];+ points[i].y = lines[j + 1];+ j += 2;+ }+ + printf ("simplified_path\n");+ simplified_path = malloc (8 * npoints * sizeof (double));+ simplified_path_buffer_size = 8 * npoints;+ simplified_path_size = 0;+ + FitCurve(points, npoints, error);+ + *bezier_path = simplified_path;+ *bezier_path_size = simplified_path_size;+ + simplified_path = NULL;+ + free (points);+}++DrawBezierCurve(int n, Point2* curve)+{+ int i;++ if (simplified_path_size + 8 > simplified_path_buffer_size) {+ fprintf (stderr, "The bezier buffer is full (%d).\n", simplified_path_buffer_size);+ return;+ }++ if (n != 3) {+ fprintf (stderr, "Expecting three points\n");+ return;+ }+ + i = simplified_path_size;++ simplified_path[i + 0] = curve[0].x;+ simplified_path[i + 1] = curve[0].y;+ simplified_path[i + 2] = curve[1].x;+ simplified_path[i + 3] = curve[1].y;+ simplified_path[i + 4] = curve[2].x;+ simplified_path[i + 5] = curve[2].y;+ simplified_path[i + 6] = curve[3].x;+ simplified_path[i + 7] = curve[3].y;++ simplified_path_size = i + 8;+}++/*+ * FitCurve :+ * Fit a Bezier curve to a set of digitized points + */+void FitCurve(d, nPts, error)+ Point2 *d; /* Array of digitized points */+ int nPts; /* Number of digitized points */+ double error; /* User-defined error squared */+{+ Vector2 tHat1, tHat2; /* Unit tangent vectors at endpoints */++ tHat1 = ComputeLeftTangent(d, 0);+ tHat2 = ComputeRightTangent(d, nPts - 1);+ FitCubic(d, 0, nPts - 1, tHat1, tHat2, error);+}++/*+ * FitCubic :+ * Fit a Bezier curve to a (sub)set of digitized points+ */+static void FitCubic(d, first, last, tHat1, tHat2, error)+ Point2 *d; /* Array of digitized points */+ int first, last; /* Indices of first and last pts in region */+ Vector2 tHat1, tHat2; /* Unit tangent vectors at endpoints */+ double error; /* User-defined error squared */+{+ BezierCurve bezCurve; /*Control points of fitted Bezier curve*/+ double *u; /* Parameter values for point */+ double *uPrime; /* Improved parameter values */+ double maxError; /* Maximum fitting error */+ int splitPoint; /* Point to split point set at */+ int nPts; /* Number of points in subset */+ double iterationError; /*Error below which you try iterating */+ int maxIterations = 4; /* Max times to try iterating */+ Vector2 tHatCenter; /* Unit tangent vector at splitPoint */+ int i; ++ iterationError = error * error;+ nPts = last - first + 1;++ /* Use heuristic if region only has two points in it */+ if (nPts == 2) {+ double dist = V2DistanceBetween2Points(&d[last], &d[first]) / 3.0;++ bezCurve = (Point2 *)malloc(4 * sizeof(Point2));+ bezCurve[0] = d[first];+ bezCurve[3] = d[last];+ V2Add(&bezCurve[0], V2Scale(&tHat1, dist), &bezCurve[1]);+ V2Add(&bezCurve[3], V2Scale(&tHat2, dist), &bezCurve[2]);+ DrawBezierCurve(3, bezCurve);+ free((void *)bezCurve);+ return;+ }++ /* Parameterize points, and attempt to fit curve */+ u = ChordLengthParameterize(d, first, last);+ bezCurve = GenerateBezier(d, first, last, u, tHat1, tHat2);++ /* Find max deviation of points to fitted curve */+ maxError = ComputeMaxError(d, first, last, bezCurve, u, &splitPoint);+ if (maxError < error) {+ DrawBezierCurve(3, bezCurve);+ free((void *)u);+ free((void *)bezCurve);+ return;+ }+++ /* If error not too large, try some reparameterization */+ /* and iteration */+ if (maxError < iterationError) {+ for (i = 0; i < maxIterations; i++) {+ uPrime = Reparameterize(d, first, last, u, bezCurve);+ free((void *)bezCurve);+ bezCurve = GenerateBezier(d, first, last, uPrime, tHat1, tHat2);+ maxError = ComputeMaxError(d, first, last,+ bezCurve, uPrime, &splitPoint);+ if (maxError < error) {+ DrawBezierCurve(3, bezCurve);+ free((void *)u);+ free((void *)bezCurve);+ free((void *)uPrime);+ return;+ }+ free((void *)u);+ u = uPrime;+ }+ }++ /* Fitting failed -- split at max error point and fit recursively */+ free((void *)u);+ free((void *)bezCurve);+ tHatCenter = ComputeCenterTangent(d, splitPoint);+ FitCubic(d, first, splitPoint, tHat1, tHatCenter, error);+ V2Negate(&tHatCenter);+ FitCubic(d, splitPoint, last, tHatCenter, tHat2, error);+}+++/*+ * GenerateBezier :+ * Use least-squares method to find Bezier control points for region.+ *+ */+static BezierCurve GenerateBezier(d, first, last, uPrime, tHat1, tHat2)+ Point2 *d; /* Array of digitized points */+ int first, last; /* Indices defining region */+ double *uPrime; /* Parameter values for region */+ Vector2 tHat1, tHat2; /* Unit tangents at endpoints */+{+ int i;+ Vector2 A[MAXPOINTS][2]; /* Precomputed rhs for eqn */+ int nPts; /* Number of pts in sub-curve */+ double C[2][2]; /* Matrix C */+ double X[2]; /* Matrix X */+ double det_C0_C1, /* Determinants of matrices */+ det_C0_X,+ det_X_C1;+ double alpha_l, /* Alpha values, left and right */+ alpha_r;+ Vector2 tmp; /* Utility variable */+ BezierCurve bezCurve; /* RETURN bezier curve ctl pts */++ bezCurve = (Point2 *)malloc(4 * sizeof(Point2));+ nPts = last - first + 1;++ + /* Compute the A's */+ for (i = 0; i < nPts; i++) {+ Vector2 v1, v2;+ v1 = tHat1;+ v2 = tHat2;+ V2Scale(&v1, B1(uPrime[i]));+ V2Scale(&v2, B2(uPrime[i]));+ A[i][0] = v1;+ A[i][1] = v2;+ }++ /* Create the C and X matrices */+ C[0][0] = 0.0;+ C[0][1] = 0.0;+ C[1][0] = 0.0;+ C[1][1] = 0.0;+ X[0] = 0.0;+ X[1] = 0.0;++ for (i = 0; i < nPts; i++) {+ C[0][0] += V2Dot(&A[i][0], &A[i][0]);+ C[0][1] += V2Dot(&A[i][0], &A[i][1]);+/* C[1][0] += V2Dot(&A[i][0], &A[i][1]);*/ + C[1][0] = C[0][1];+ C[1][1] += V2Dot(&A[i][1], &A[i][1]);++ tmp = V2SubII(d[first + i],+ V2AddII(+ V2ScaleIII(d[first], B0(uPrime[i])),+ V2AddII(+ V2ScaleIII(d[first], B1(uPrime[i])),+ V2AddII(+ V2ScaleIII(d[last], B2(uPrime[i])),+ V2ScaleIII(d[last], B3(uPrime[i]))))));+ ++ X[0] += V2Dot(&A[i][0], &tmp);+ X[1] += V2Dot(&A[i][1], &tmp);+ }++ /* Compute the determinants of C and X */+ det_C0_C1 = C[0][0] * C[1][1] - C[1][0] * C[0][1];+ det_C0_X = C[0][0] * X[1] - C[1][0] * X[0];+ det_X_C1 = X[0] * C[1][1] - X[1] * C[0][1];++ /* Finally, derive alpha values */+ alpha_l = (det_C0_C1 == 0) ? 0.0 : det_X_C1 / det_C0_C1;+ alpha_r = (det_C0_C1 == 0) ? 0.0 : det_C0_X / det_C0_C1;++ /* If alpha negative, use the Wu/Barsky heuristic (see text) */+ /* (if alpha is 0, you get coincident control points that lead to+ * divide by zero in any subsequent NewtonRaphsonRootFind() call. */+ double segLength = V2DistanceBetween2Points(&d[last], &d[first]);+ double epsilon = 1.0e-6 * segLength;+ if (alpha_l < epsilon || alpha_r < epsilon)+ {+ /* fall back on standard (probably inaccurate) formula, and subdivide further if needed. */+ double dist = segLength / 3.0;+ bezCurve[0] = d[first];+ bezCurve[3] = d[last];+ V2Add(&bezCurve[0], V2Scale(&tHat1, dist), &bezCurve[1]);+ V2Add(&bezCurve[3], V2Scale(&tHat2, dist), &bezCurve[2]);+ return (bezCurve);+ }++ /* First and last control points of the Bezier curve are */+ /* positioned exactly at the first and last data points */+ /* Control points 1 and 2 are positioned an alpha distance out */+ /* on the tangent vectors, left and right, respectively */+ bezCurve[0] = d[first];+ bezCurve[3] = d[last];+ V2Add(&bezCurve[0], V2Scale(&tHat1, alpha_l), &bezCurve[1]);+ V2Add(&bezCurve[3], V2Scale(&tHat2, alpha_r), &bezCurve[2]);+ return (bezCurve);+}+++/*+ * Reparameterize:+ * Given set of points and their parameterization, try to find+ * a better parameterization.+ *+ */+static double *Reparameterize(d, first, last, u, bezCurve)+ Point2 *d; /* Array of digitized points */+ int first, last; /* Indices defining region */+ double *u; /* Current parameter values */+ BezierCurve bezCurve; /* Current fitted curve */+{+ int nPts = last-first+1; + int i;+ double *uPrime; /* New parameter values */++ uPrime = (double *)malloc(nPts * sizeof(double));+ for (i = first; i <= last; i++) {+ uPrime[i-first] = NewtonRaphsonRootFind(bezCurve, d[i], u[i-+ first]);+ }+ return (uPrime);+}++++/*+ * NewtonRaphsonRootFind :+ * Use Newton-Raphson iteration to find better root.+ */+static double NewtonRaphsonRootFind(Q, P, u)+ BezierCurve Q; /* Current fitted curve */+ Point2 P; /* Digitized point */+ double u; /* Parameter value for "P" */+{+ double numerator, denominator;+ Point2 Q1[3], Q2[2]; /* Q' and Q'' */+ Point2 Q_u, Q1_u, Q2_u; /*u evaluated at Q, Q', & Q'' */+ double uPrime; /* Improved u */+ int i;+ + /* Compute Q(u) */+ Q_u = BezierII(3, Q, u);+ + /* Generate control vertices for Q' */+ for (i = 0; i <= 2; i++) {+ Q1[i].x = (Q[i+1].x - Q[i].x) * 3.0;+ Q1[i].y = (Q[i+1].y - Q[i].y) * 3.0;+ }+ + /* Generate control vertices for Q'' */+ for (i = 0; i <= 1; i++) {+ Q2[i].x = (Q1[i+1].x - Q1[i].x) * 2.0;+ Q2[i].y = (Q1[i+1].y - Q1[i].y) * 2.0;+ }+ + /* Compute Q'(u) and Q''(u) */+ Q1_u = BezierII(2, Q1, u);+ Q2_u = BezierII(1, Q2, u);+ + /* Compute f(u)/f'(u) */+ numerator = (Q_u.x - P.x) * (Q1_u.x) + (Q_u.y - P.y) * (Q1_u.y);+ denominator = (Q1_u.x) * (Q1_u.x) + (Q1_u.y) * (Q1_u.y) ++ (Q_u.x - P.x) * (Q2_u.x) + (Q_u.y - P.y) * (Q2_u.y);+ if (denominator == 0.0f) return u;++ /* u = u - f(u)/f'(u) */+ uPrime = u - (numerator/denominator);+ return (uPrime);+}++ + +/*+ * Bezier :+ * Evaluate a Bezier curve at a particular parameter value+ * + */+static Point2 BezierII(degree, V, t)+ int degree; /* The degree of the bezier curve */+ Point2 *V; /* Array of control points */+ double t; /* Parametric value to find point for */+{+ int i, j; + Point2 Q; /* Point on curve at parameter t */+ Point2 *Vtemp; /* Local copy of control points */++ /* Copy array */+ Vtemp = (Point2 *)malloc((unsigned)((degree+1) + * sizeof (Point2)));+ for (i = 0; i <= degree; i++) {+ Vtemp[i] = V[i];+ }++ /* Triangle computation */+ for (i = 1; i <= degree; i++) { + for (j = 0; j <= degree-i; j++) {+ Vtemp[j].x = (1.0 - t) * Vtemp[j].x + t * Vtemp[j+1].x;+ Vtemp[j].y = (1.0 - t) * Vtemp[j].y + t * Vtemp[j+1].y;+ }+ }++ Q = Vtemp[0];+ free((void *)Vtemp);+ return Q;+}+++/*+ * B0, B1, B2, B3 :+ * Bezier multipliers+ */+static double B0(u)+ double u;+{+ double tmp = 1.0 - u;+ return (tmp * tmp * tmp);+}+++static double B1(u)+ double u;+{+ double tmp = 1.0 - u;+ return (3 * u * (tmp * tmp));+}++static double B2(u)+ double u;+{+ double tmp = 1.0 - u;+ return (3 * u * u * tmp);+}++static double B3(u)+ double u;+{+ return (u * u * u);+}++++/*+ * ComputeLeftTangent, ComputeRightTangent, ComputeCenterTangent :+ *Approximate unit tangents at endpoints and "center" of digitized curve+ */+static Vector2 ComputeLeftTangent(d, end)+ Point2 *d; /* Digitized points*/+ int end; /* Index to "left" end of region */+{+ Vector2 tHat1;+ tHat1 = V2SubII(d[end+1], d[end]);+ tHat1 = *V2Normalize(&tHat1);+ return tHat1;+}++static Vector2 ComputeRightTangent(d, end)+ Point2 *d; /* Digitized points */+ int end; /* Index to "right" end of region */+{+ Vector2 tHat2;+ tHat2 = V2SubII(d[end-1], d[end]);+ tHat2 = *V2Normalize(&tHat2);+ return tHat2;+}+++static Vector2 ComputeCenterTangent(d, center)+ Point2 *d; /* Digitized points */+ int center; /* Index to point inside region */+{+ Vector2 V1, V2, tHatCenter;++ V1 = V2SubII(d[center-1], d[center]);+ V2 = V2SubII(d[center], d[center+1]);+ tHatCenter.x = (V1.x + V2.x)/2.0;+ tHatCenter.y = (V1.y + V2.y)/2.0;+ tHatCenter = *V2Normalize(&tHatCenter);+ return tHatCenter;+}+++/*+ * ChordLengthParameterize :+ * Assign parameter values to digitized points + * using relative distances between points.+ */+static double *ChordLengthParameterize(d, first, last)+ Point2 *d; /* Array of digitized points */+ int first, last; /* Indices defining region */+{+ int i; + double *u; /* Parameterization */++ u = (double *)malloc((unsigned)(last-first+1) * sizeof(double));++ u[0] = 0.0;+ for (i = first+1; i <= last; i++) {+ u[i-first] = u[i-first-1] ++ V2DistanceBetween2Points(&d[i], &d[i-1]);+ }++ for (i = first + 1; i <= last; i++) {+ u[i-first] = u[i-first] / u[last-first];+ }++ return(u);+}+++++/*+ * ComputeMaxError :+ * Find the maximum squared distance of digitized points+ * to fitted curve.+*/+static double ComputeMaxError(d, first, last, bezCurve, u, splitPoint)+ Point2 *d; /* Array of digitized points */+ int first, last; /* Indices defining region */+ BezierCurve bezCurve; /* Fitted Bezier curve */+ double *u; /* Parameterization of points */+ int *splitPoint; /* Point of maximum error */+{+ int i;+ double maxDist; /* Maximum error */+ double dist; /* Current error */+ Point2 P; /* Point on curve */+ Vector2 v; /* Vector from point to curve */++ *splitPoint = (last - first + 1)/2;+ maxDist = 0.0;+ for (i = first + 1; i < last; i++) {+ P = BezierII(3, bezCurve, u[i-first]);+ v = V2SubII(P, d[i]);+ dist = V2SquaredLength(&v);+ if (dist >= maxDist) {+ maxDist = dist;+ *splitPoint = i;+ }+ }+ return (maxDist);+}+static Vector2 V2AddII(a, b)+ Vector2 a, b;+{+ Vector2 c;+ c.x = a.x + b.x; c.y = a.y + b.y;+ return (c);+}+static Vector2 V2ScaleIII(v, s)+ Vector2 v;+ double s;+{+ Vector2 result;+ result.x = v.x * s; result.y = v.y * s;+ return (result);+}++static Vector2 V2SubII(a, b)+ Vector2 a, b;+{+ Vector2 c;+ c.x = a.x - b.x; c.y = a.y - b.y;+ return (c);+}+new file mode 100644@@ -0,0 +1,601 @@+/*+An Algorithm for Automatically Fitting Digitized Curves+by Philip J. Schneider+from "Graphics Gems", Academic Press, 1990+Adapted to BirdFont by Johan Mattsson 2015+*/++/*+ * EULA: The Graphics Gems code is copyright-protected. In other words, + * you cannot claim the text of the code as your own and resell it. + * Using the code is permitted in any program, product, or library, + * non-commercial or commercial. Giving credit is not required, though + * is a nice gesture. The code comes as-is, and if there are any flaws + * or problems with any Gems code, nobody involved with Gems - authors,+ * editors, publishers, or webmasters - are to be held responsible. + * Basically, don't be a jerk, and remember that anything free comes + * with no guarantee. + */++/* fit_cubic.c */ +/* Piecewise cubic fitting code */++#include "GraphicsGems.h"+#include <stdio.h>+#include <malloc.h>+#include <math.h>+#include <stdio.h>++typedef Point2 *BezierCurve;++/* Forward declarations */+void FitCurve();+static void FitCubic();+static double *Reparameterize();+static double NewtonRaphsonRootFind();+static Point2 BezierII();+static double B0(), B1(), B2(), B3();+static Vector2 ComputeLeftTangent();+static Vector2 ComputeRightTangent();+static Vector2 ComputeCenterTangent();+static double ComputeMaxError();+static double *ChordLengthParameterize();+static BezierCurve GenerateBezier();+static Vector2 V2AddII();+static Vector2 V2ScaleIII();+static Vector2 V2SubII();++#define MAXPOINTS 1000 /* The most points you can have */++int simplified_path_buffer_size = 0;+int simplified_path_size = 0;+double* simplified_path = NULL;++/** Generates an a bezier path and returns the length of the output array. */ +int fit_bezier_curve_to_line (lines, lines_size, error, bezier_path)+ double* lines,+ int lines_size,+ double error,+ double** bezier_path /** generated Bezier curves */+{+ int i, j, n;+ Point2* points+ int npoints;++ npoints = lines_size / 2;+ bezier_path = NULL;++ if (lines_size % 2 != 0) {+ printf ("Odd number of coordinates."); + return 0;+ }++ if (lines == NULL || lines_size == 0) {+ printf ("No lines."); + return 0;+ }++ if (simplified_path != NULL) {+ printf ("Path simplification is alreading running."); + return 0; + }++ points = malloc (npoints * sizeof (Point2));++ j = 0;+ for (i = 0; i <= npoints; i++) {+ points[i].x = lines[j];+ points[i].y = lines[j + 1];+ j += 2;+ }+ + simplified_path = malloc (4 * npoints * sizeof (double));+ simplified_path_buffer_size = 4 * npoints;+ + FitCurve(points, n, error);+ free (points);++ bezier_path = simplified_path;+ simplified_path = NULL;++ return simplified_path_size;+}++DrawBezierCurve(int n, Point2* curve)+{+ int i;++ if (simplified_path_size + 4 > simplified_path_buffer_size) {+ printf ("Out of memmory\n");+ return;+ }++ if (n != 3) {+ printf ("Expecting three points\n");+ }+ + i = simplified_path_size;++ simplified_path[i++] = curve[0].x;+ simplified_path[i++] = curve[0].y;+ simplified_path[i++] = curve[1].x;+ simplified_path[i++] = curve[1].y;+ simplified_path[i++] = curve[2].x;+ simplified_path[i++] = curve[2].y;+ simplified_path[i++] = curve[3].x;+ simplified_path[i++] = curve[3].y;++ simplified_path_size = i;+}++/*+ * FitCurve :+ * Fit a Bezier curve to a set of digitized points + */+void FitCurve(d, nPts, error)+ Point2 *d; /* Array of digitized points */+ int nPts; /* Number of digitized points */+ double error; /* User-defined error squared */+{+ Vector2 tHat1, tHat2; /* Unit tangent vectors at endpoints */++ tHat1 = ComputeLeftTangent(d, 0);+ tHat2 = ComputeRightTangent(d, nPts - 1);+ FitCubic(d, 0, nPts - 1, tHat1, tHat2, error);+}++/*+ * FitCubic :+ * Fit a Bezier curve to a (sub)set of digitized points+ */+static void FitCubic(d, first, last, tHat1, tHat2, error)+ Point2 *d; /* Array of digitized points */+ int first, last; /* Indices of first and last pts in region */+ Vector2 tHat1, tHat2; /* Unit tangent vectors at endpoints */+ double error; /* User-defined error squared */+{+ BezierCurve bezCurve; /*Control points of fitted Bezier curve*/+ double *u; /* Parameter values for point */+ double *uPrime; /* Improved parameter values */+ double maxError; /* Maximum fitting error */+ int splitPoint; /* Point to split point set at */+ int nPts; /* Number of points in subset */+ double iterationError; /*Error below which you try iterating */+ int maxIterations = 4; /* Max times to try iterating */+ Vector2 tHatCenter; /* Unit tangent vector at splitPoint */+ int i; ++ iterationError = error * error;+ nPts = last - first + 1;++ /* Use heuristic if region only has two points in it */+ if (nPts == 2) {+ double dist = V2DistanceBetween2Points(&d[last], &d[first]) / 3.0;++ bezCurve = (Point2 *)malloc(4 * sizeof(Point2));+ bezCurve[0] = d[first];+ bezCurve[3] = d[last];+ V2Add(&bezCurve[0], V2Scale(&tHat1, dist), &bezCurve[1]);+ V2Add(&bezCurve[3], V2Scale(&tHat2, dist), &bezCurve[2]);+ DrawBezierCurve(3, bezCurve);+ free((void *)bezCurve);+ return;+ }++ /* Parameterize points, and attempt to fit curve */+ u = ChordLengthParameterize(d, first, last);+ bezCurve = GenerateBezier(d, first, last, u, tHat1, tHat2);++ /* Find max deviation of points to fitted curve */+ maxError = ComputeMaxError(d, first, last, bezCurve, u, &splitPoint);+ if (maxError < error) {+ DrawBezierCurve(3, bezCurve);+ free((void *)u);+ free((void *)bezCurve);+ return;+ }+++ /* If error not too large, try some reparameterization */+ /* and iteration */+ if (maxError < iterationError) {+ for (i = 0; i < maxIterations; i++) {+ uPrime = Reparameterize(d, first, last, u, bezCurve);+ free((void *)bezCurve);+ bezCurve = GenerateBezier(d, first, last, uPrime, tHat1, tHat2);+ maxError = ComputeMaxError(d, first, last,+ bezCurve, uPrime, &splitPoint);+ if (maxError < error) {+ DrawBezierCurve(3, bezCurve);+ free((void *)u);+ free((void *)bezCurve);+ free((void *)uPrime);+ return;+ }+ free((void *)u);+ u = uPrime;+ }+ }++ /* Fitting failed -- split at max error point and fit recursively */+ free((void *)u);+ free((void *)bezCurve);+ tHatCenter = ComputeCenterTangent(d, splitPoint);+ FitCubic(d, first, splitPoint, tHat1, tHatCenter, error);+ V2Negate(&tHatCenter);+ FitCubic(d, splitPoint, last, tHatCenter, tHat2, error);+}+++/*+ * GenerateBezier :+ * Use least-squares method to find Bezier control points for region.+ *+ */+static BezierCurve GenerateBezier(d, first, last, uPrime, tHat1, tHat2)+ Point2 *d; /* Array of digitized points */+ int first, last; /* Indices defining region */+ double *uPrime; /* Parameter values for region */+ Vector2 tHat1, tHat2; /* Unit tangents at endpoints */+{+ int i;+ Vector2 A[MAXPOINTS][2]; /* Precomputed rhs for eqn */+ int nPts; /* Number of pts in sub-curve */+ double C[2][2]; /* Matrix C */+ double X[2]; /* Matrix X */+ double det_C0_C1, /* Determinants of matrices */+ det_C0_X,+ det_X_C1;+ double alpha_l, /* Alpha values, left and right */+ alpha_r;+ Vector2 tmp; /* Utility variable */+ BezierCurve bezCurve; /* RETURN bezier curve ctl pts */++ bezCurve = (Point2 *)malloc(4 * sizeof(Point2));+ nPts = last - first + 1;++ + /* Compute the A's */+ for (i = 0; i < nPts; i++) {+ Vector2 v1, v2;+ v1 = tHat1;+ v2 = tHat2;+ V2Scale(&v1, B1(uPrime[i]));+ V2Scale(&v2, B2(uPrime[i]));+ A[i][0] = v1;+ A[i][1] = v2;+ }++ /* Create the C and X matrices */+ C[0][0] = 0.0;+ C[0][1] = 0.0;+ C[1][0] = 0.0;+ C[1][1] = 0.0;+ X[0] = 0.0;+ X[1] = 0.0;++ for (i = 0; i < nPts; i++) {+ C[0][0] += V2Dot(&A[i][0], &A[i][0]);+ C[0][1] += V2Dot(&A[i][0], &A[i][1]);+/* C[1][0] += V2Dot(&A[i][0], &A[i][1]);*/ + C[1][0] = C[0][1];+ C[1][1] += V2Dot(&A[i][1], &A[i][1]);++ tmp = V2SubII(d[first + i],+ V2AddII(+ V2ScaleIII(d[first], B0(uPrime[i])),+ V2AddII(+ V2ScaleIII(d[first], B1(uPrime[i])),+ V2AddII(+ V2ScaleIII(d[last], B2(uPrime[i])),+ V2ScaleIII(d[last], B3(uPrime[i]))))));+ ++ X[0] += V2Dot(&A[i][0], &tmp);+ X[1] += V2Dot(&A[i][1], &tmp);+ }++ /* Compute the determinants of C and X */+ det_C0_C1 = C[0][0] * C[1][1] - C[1][0] * C[0][1];+ det_C0_X = C[0][0] * X[1] - C[1][0] * X[0];+ det_X_C1 = X[0] * C[1][1] - X[1] * C[0][1];++ /* Finally, derive alpha values */+ alpha_l = (det_C0_C1 == 0) ? 0.0 : det_X_C1 / det_C0_C1;+ alpha_r = (det_C0_C1 == 0) ? 0.0 : det_C0_X / det_C0_C1;++ /* If alpha negative, use the Wu/Barsky heuristic (see text) */+ /* (if alpha is 0, you get coincident control points that lead to+ * divide by zero in any subsequent NewtonRaphsonRootFind() call. */+ double segLength = V2DistanceBetween2Points(&d[last], &d[first]);+ double epsilon = 1.0e-6 * segLength;+ if (alpha_l < epsilon || alpha_r < epsilon)+ {+ /* fall back on standard (probably inaccurate) formula, and subdivide further if needed. */+ double dist = segLength / 3.0;+ bezCurve[0] = d[first];+ bezCurve[3] = d[last];+ V2Add(&bezCurve[0], V2Scale(&tHat1, dist), &bezCurve[1]);+ V2Add(&bezCurve[3], V2Scale(&tHat2, dist), &bezCurve[2]);+ return (bezCurve);+ }++ /* First and last control points of the Bezier curve are */+ /* positioned exactly at the first and last data points */+ /* Control points 1 and 2 are positioned an alpha distance out */+ /* on the tangent vectors, left and right, respectively */+ bezCurve[0] = d[first];+ bezCurve[3] = d[last];+ V2Add(&bezCurve[0], V2Scale(&tHat1, alpha_l), &bezCurve[1]);+ V2Add(&bezCurve[3], V2Scale(&tHat2, alpha_r), &bezCurve[2]);+ return (bezCurve);+}+++/*+ * Reparameterize:+ * Given set of points and their parameterization, try to find+ * a better parameterization.+ *+ */+static double *Reparameterize(d, first, last, u, bezCurve)+ Point2 *d; /* Array of digitized points */+ int first, last; /* Indices defining region */+ double *u; /* Current parameter values */+ BezierCurve bezCurve; /* Current fitted curve */+{+ int nPts = last-first+1; + int i;+ double *uPrime; /* New parameter values */++ uPrime = (double *)malloc(nPts * sizeof(double));+ for (i = first; i <= last; i++) {+ uPrime[i-first] = NewtonRaphsonRootFind(bezCurve, d[i], u[i-+ first]);+ }+ return (uPrime);+}++++/*+ * NewtonRaphsonRootFind :+ * Use Newton-Raphson iteration to find better root.+ */+static double NewtonRaphsonRootFind(Q, P, u)+ BezierCurve Q; /* Current fitted curve */+ Point2 P; /* Digitized point */+ double u; /* Parameter value for "P" */+{+ double numerator, denominator;+ Point2 Q1[3], Q2[2]; /* Q' and Q'' */+ Point2 Q_u, Q1_u, Q2_u; /*u evaluated at Q, Q', & Q'' */+ double uPrime; /* Improved u */+ int i;+ + /* Compute Q(u) */+ Q_u = BezierII(3, Q, u);+ + /* Generate control vertices for Q' */+ for (i = 0; i <= 2; i++) {+ Q1[i].x = (Q[i+1].x - Q[i].x) * 3.0;+ Q1[i].y = (Q[i+1].y - Q[i].y) * 3.0;+ }+ + /* Generate control vertices for Q'' */+ for (i = 0; i <= 1; i++) {+ Q2[i].x = (Q1[i+1].x - Q1[i].x) * 2.0;+ Q2[i].y = (Q1[i+1].y - Q1[i].y) * 2.0;+ }+ + /* Compute Q'(u) and Q''(u) */+ Q1_u = BezierII(2, Q1, u);+ Q2_u = BezierII(1, Q2, u);+ + /* Compute f(u)/f'(u) */+ numerator = (Q_u.x - P.x) * (Q1_u.x) + (Q_u.y - P.y) * (Q1_u.y);+ denominator = (Q1_u.x) * (Q1_u.x) + (Q1_u.y) * (Q1_u.y) ++ (Q_u.x - P.x) * (Q2_u.x) + (Q_u.y - P.y) * (Q2_u.y);+ if (denominator == 0.0f) return u;++ /* u = u - f(u)/f'(u) */+ uPrime = u - (numerator/denominator);+ return (uPrime);+}++ + +/*+ * Bezier :+ * Evaluate a Bezier curve at a particular parameter value+ * + */+static Point2 BezierII(degree, V, t)+ int degree; /* The degree of the bezier curve */+ Point2 *V; /* Array of control points */+ double t; /* Parametric value to find point for */+{+ int i, j; + Point2 Q; /* Point on curve at parameter t */+ Point2 *Vtemp; /* Local copy of control points */++ /* Copy array */+ Vtemp = (Point2 *)malloc((unsigned)((degree+1) + * sizeof (Point2)));+ for (i = 0; i <= degree; i++) {+ Vtemp[i] = V[i];+ }++ /* Triangle computation */+ for (i = 1; i <= degree; i++) { + for (j = 0; j <= degree-i; j++) {+ Vtemp[j].x = (1.0 - t) * Vtemp[j].x + t * Vtemp[j+1].x;+ Vtemp[j].y = (1.0 - t) * Vtemp[j].y + t * Vtemp[j+1].y;+ }+ }++ Q = Vtemp[0];+ free((void *)Vtemp);+ return Q;+}+++/*+ * B0, B1, B2, B3 :+ * Bezier multipliers+ */+static double B0(u)+ double u;+{+ double tmp = 1.0 - u;+ return (tmp * tmp * tmp);+}+++static double B1(u)+ double u;+{+ double tmp = 1.0 - u;+ return (3 * u * (tmp * tmp));+}++static double B2(u)+ double u;+{+ double tmp = 1.0 - u;+ return (3 * u * u * tmp);+}++static double B3(u)+ double u;+{+ return (u * u * u);+}++++/*+ * ComputeLeftTangent, ComputeRightTangent, ComputeCenterTangent :+ *Approximate unit tangents at endpoints and "center" of digitized curve+ */+static Vector2 ComputeLeftTangent(d, end)+ Point2 *d; /* Digitized points*/+ int end; /* Index to "left" end of region */+{+ Vector2 tHat1;+ tHat1 = V2SubII(d[end+1], d[end]);+ tHat1 = *V2Normalize(&tHat1);+ return tHat1;+}++static Vector2 ComputeRightTangent(d, end)+ Point2 *d; /* Digitized points */+ int end; /* Index to "right" end of region */+{+ Vector2 tHat2;+ tHat2 = V2SubII(d[end-1], d[end]);+ tHat2 = *V2Normalize(&tHat2);+ return tHat2;+}+++static Vector2 ComputeCenterTangent(d, center)+ Point2 *d; /* Digitized points */+ int center; /* Index to point inside region */+{+ Vector2 V1, V2, tHatCenter;++ V1 = V2SubII(d[center-1], d[center]);+ V2 = V2SubII(d[center], d[center+1]);+ tHatCenter.x = (V1.x + V2.x)/2.0;+ tHatCenter.y = (V1.y + V2.y)/2.0;+ tHatCenter = *V2Normalize(&tHatCenter);+ return tHatCenter;+}+++/*+ * ChordLengthParameterize :+ * Assign parameter values to digitized points + * using relative distances between points.+ */+static double *ChordLengthParameterize(d, first, last)+ Point2 *d; /* Array of digitized points */+ int first, last; /* Indices defining region */+{+ int i; + double *u; /* Parameterization */++ u = (double *)malloc((unsigned)(last-first+1) * sizeof(double));++ u[0] = 0.0;+ for (i = first+1; i <= last; i++) {+ u[i-first] = u[i-first-1] ++ V2DistanceBetween2Points(&d[i], &d[i-1]);+ }++ for (i = first + 1; i <= last; i++) {+ u[i-first] = u[i-first] / u[last-first];+ }++ return(u);+}+++++/*+ * ComputeMaxError :+ * Find the maximum squared distance of digitized points+ * to fitted curve.+*/+static double ComputeMaxError(d, first, last, bezCurve, u, splitPoint)+ Point2 *d; /* Array of digitized points */+ int first, last; /* Indices defining region */+ BezierCurve bezCurve; /* Fitted Bezier curve */+ double *u; /* Parameterization of points */+ int *splitPoint; /* Point of maximum error */+{+ int i;+ double maxDist; /* Maximum error */+ double dist; /* Current error */+ Point2 P; /* Point on curve */+ Vector2 v; /* Vector from point to curve */++ *splitPoint = (last - first + 1)/2;+ maxDist = 0.0;+ for (i = first + 1; i < last; i++) {+ P = BezierII(3, bezCurve, u[i-first]);+ v = V2SubII(P, d[i]);+ dist = V2SquaredLength(&v);+ if (dist >= maxDist) {+ maxDist = dist;+ *splitPoint = i;+ }+ }+ return (maxDist);+}+static Vector2 V2AddII(a, b)+ Vector2 a, b;+{+ Vector2 c;+ c.x = a.x + b.x; c.y = a.y + b.y;+ return (c);+}+static Vector2 V2ScaleIII(v, s)+ Vector2 v;+ double s;+{+ Vector2 result;+ result.x = v.x * s; result.y = v.y * s;+ return (result);+}++static Vector2 V2SubII(a, b)+ Vector2 a, b;+{+ Vector2 c;+ c.x = a.x - b.x; c.y = a.y - b.y;+ return (c);+}+@@ -82,13 +82,15 @@ class Vala(object): if library: self.header = join(build, library) + '.h'- self.vapi = join(build, library) + '.vapi'+ self.vapi = join(build, library) + '.vapi' # generated vapi file+ self.other_vapi_files = get_sources_path (src, '*.vapi') # other vapi files self.so = join(build, src) + '.so.' + so_version self.so_link = join(build, src) + '.so' self.so_link_name = src + '.so' self.so_version = so_version self.so_name = 'lib' + library + '.so.' + so_version-+ else:+ self.other_vapi_files = [] def gen_c(self, opts): """translate code from vala to C and create .vapi"""@@ -107,9 +109,10 @@ class Vala(object): dep_vapi = [d.vapi for d in self.vala_deps] action = cmd('valac', options, params, dep_vapi, self.vala) targets = self.cc[:]+ if self.library: targets += [self.header, self.vapi]- + for f in self.c: yield { 'name': 'copy_c_' + f,@@ -118,14 +121,25 @@ class Vala(object): 'cp ' + f + ' ' + self.build + '/' + self.src + '/' ], }+ + for f in self.other_vapi_files:+ yield {+ 'name': 'vapi_files_' + f,+ 'actions': [ + 'mkdir -p '+ self.build + '/', + 'cp ' + f + ' ' + self.build + '/'+ ],+ } - print (action) - yield {- 'name': 'compile_c',- 'actions': [ action ],- 'file_dep': self.vala + dep_vapi,- 'targets': targets,- }+ print (action)+ + if not self.vala == []:+ yield {+ 'name': 'compile_c',+ 'actions': [ action ],+ 'file_dep': self.vala + dep_vapi,+ 'targets': targets,+ } def gen_o(self, opts):@@ -1,6 +1,6 @@ #!/usr/bin/python """-Copyright (C) 2013, 2014 Johan Mattsson+Copyright (C) 2013, 2014 2015 Johan Mattsson This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public License as @@ -164,6 +164,61 @@ def libbirdxml(prefix, cc, cflags, ldflags, valac, valaflags, library, nonNull = elif library.find ('.dylib') > -1: run ("""cd build/bin && ln -sf """ + library + " libbirdxml.dylib") ++def libbirdgems(prefix, cc, cflags, ldflags, valac, valaflags, library, nonNull = True):+ #libbirdfont+ run("mkdir -p build/libbirdgems")+ run("mkdir -p build/bin")++ experimentalNonNull = ""+ if nonNull:+ experimentalNonNull = "--enable-experimental-non-null"++ run(valac + """\+ -C \+ """ + valaflags + """ \+ --pkg posix \+ --vapidir=./ \+ --basedir build/libbirdgems/ \+ """ + experimentalNonNull + """ \+ --enable-experimental \+ --library libbirdgems \+ -H build/libbirdxml/birdxml.h \+ libbirdxml/*.vala \+ """)+ + if cc == "":+ print ("Skipping compilation");+ else:+ run(cc + " " + cflags + """ \+ -c build/libbirdgems/*.c \+ """)+ + run("mv ./*.o build/libbirdgems/ ")++ if library.endswith (".dylib"):+ sonameparam = "" # gcc on mac os does not have the soname parameter+ else:+ sonameparam = "-Wl,-soname," + library+ + run(cc + " " + ldflags + """ \+ -shared \+ """ + sonameparam + """ \+ build/libbirdgems/*.o \+ $(pkg-config --libs glib-2.0) \+ $(pkg-config --libs gobject-2.0) \+ -o """ + library)+ run("mv " + library + " build/bin/")+ + if os.path.exists("build/bin/libbirdgems.so"):+ run ("cd build/bin && unlink libbirdgems.so")++ # create link to the versioned library+ if library.find ('.so') > -1:+ run ("""cd build/bin && ln -sf """ + library + " libbirdgems.so")+ elif library.find ('.dylib') > -1:+ run ("""cd build/bin && ln -sf """ + library + " libbirdgems.dylib")+ def birdfont_export(prefix, cc, cflags, ldflags, valac, valaflags, nonNull = True): # birdfont-export@@ -1,6 +1,6 @@ #!/usr/bin/python """-Copyright (C) 2013 2014 Johan Mattsson+Copyright (C) 2013 2014 2015 Johan Mattsson This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public License as @@ -53,11 +53,15 @@ ldflags = options.ldflags + " " + "$(pkg-config --libs gdk-pixbuf-2.0)" library_cflags = options.cflags library_ldflags= options.ldflags + " -Wl,-soname," + "libbirdfont.so." + version.LIBBIRDXML_SO_VERSION +birdgemslibrary_cflags = options.cflags+birdgemslibrary_ldflags= options.ldflags + " -Wl,-soname," + "birdgems.so." + version.LIBBIRDGEMS_SO_VERSION+ xmllibrary_cflags = options.cflags xmllibrary_ldflags= options.ldflags + " -Wl,-soname," + "libbirdxml.so." + version.SO_VERSION configfile.write_config (prefix) compile_translations()+build.libbirdgems(prefix, cc, birdgemslibrary_cflags, birdgemslibrary_ldflags, valac, valaflags, "birdgems.so." + version.LIBBIRDGEMS_SO_VERSION, False) build.libbirdxml(prefix, cc, xmllibrary_cflags, xmllibrary_ldflags, valac, valaflags, "libbirdxml.so." + version.LIBBIRDXML_SO_VERSION, False) build.libbirdfont(prefix, cc, library_cflags, library_ldflags, valac, valaflags, "libbirdfont.so." + version.SO_VERSION, False) build.birdfont_autotrace(prefix, cc, cflags, ldflags, valac, valaflags, False)@@ -21,3 +21,7 @@ SO_VERSION = SO_VERSION_MAJOR + '.' + SO_VERSION_MINOR LIBBIRDXML_SO_VERSION_MAJOR = '0' LIBBIRDXML_SO_VERSION_MINOR = '0' LIBBIRDXML_SO_VERSION = LIBBIRDXML_SO_VERSION_MAJOR + '.' + LIBBIRDXML_SO_VERSION_MINOR;++LIBBIRDGEMS_SO_VERSION_MAJOR = '0'+LIBBIRDGEMS_SO_VERSION_MINOR = '0'+LIBBIRDGEMS_SO_VERSION = LIBBIRDGEMS_SO_VERSION_MAJOR + '.' + LIBBIRDGEMS_SO_VERSION_MINOR;@@ -30,6 +30,7 @@ from translations import compile_translations from run import run compile_translations()+build.libbirdgems(prefix, cc, cflags, library_ldflags, valac, valaflags, "libbirdgems.dll") build.libbirdxml(prefix, cc, cflags, library_ldflags, valac, valaflags, "libbirdxml.dll") build.libbirdfont(prefix, cc, cflags, library_ldflags, valac, valaflags, "libbirdfont.dll")