Orange branch commit.

This commit adds new underlying uv unwrapper code, intended to be
more extensible. At the moment this has a re-implementation of LSCM.
This has not been activated yet, since it doesn't add anything new.

What's new is the stretch minimize tool from tuhopuu. It works by
selecting some some uv's in the uv editor window, and then pressing
ctrl+V. The uv's on the boundary stay fixed.

More stuff will follow as I port it over & fix it.

1 files changed, 1785 insertions, 0 deletions

diff --git a/source/blender/src/parametrizer.c b/source/blender/src/parametrizer.c
new file mode 100644
index 00000000000..ff5d92e4604
--- /dev/null
+++ b/source/blender/src/parametrizer.c
@@ -0,0 +1,1785 @@
+
+#include "MEM_guardedalloc.h"
+
+#include "BLI_memarena.h"
+#include "BLI_arithb.h"
+#include "BLI_rand.h"
+
+#include "BKE_utildefines.h"
+
+#include "BIF_editsima.h"
+#include "BIF_toolbox.h"
+
+#include "ONL_opennl.h"
+
+#include "parametrizer.h"
+#include "parametrizer_intern.h"
+
+#include <math.h>
+#include <stdlib.h>
+#include <stdio.h>
+#include <string.h>
+
+#if defined(_WIN32)
+#define M_PI 3.14159265358979323846
+#endif
+
+/* Hash */
+
+static int PHashSizes[] = {
+	1, 3, 5, 11, 17, 37, 67, 131, 257, 521, 1031, 2053, 4099, 8209, 
+	16411, 32771, 65537, 131101, 262147, 524309, 1048583, 2097169, 
+	4194319, 8388617, 16777259, 33554467, 67108879, 134217757, 268435459
+};
+
+#define PHASH_hash(ph, item) (((unsigned long) (item))%((unsigned int) (ph)->cursize))
+
+PHash *phash_new(int sizehint)
+{
+	PHash *ph = (PHash*)MEM_callocN(sizeof(PHash), "PHash");
+	ph->size = 0;
+	ph->cursize_id = 0;
+	ph->first = NULL;
+
+	while (PHashSizes[ph->cursize_id] < sizehint)
+		ph->cursize_id++;
+
+	ph->cursize = PHashSizes[ph->cursize_id];
+	ph->buckets = (PHashLink**)MEM_callocN(ph->cursize*sizeof(*ph->buckets), "PHashBuckets");
+
+	return ph;
+}
+
+void phash_delete(PHash *ph)
+{
+	MEM_freeN(ph->buckets);
+	MEM_freeN(ph);
+}
+
+void phash_delete_with_links(PHash *ph)
+{
+	PHashLink *link, *next=NULL;
+
+	for (link = ph->first; link; link = next) {
+		next = link->next;
+		MEM_freeN(link);
+	}
+
+	phash_delete(ph);
+}
+
+int phash_size(PHash *ph)
+{
+	return ph->size;
+}
+
+void phash_insert(PHash *ph, PHashLink *link)
+{
+	int size = ph->cursize;
+	int hash = PHASH_hash(ph, link->key);
+	PHashLink *lookup = ph->buckets[hash];
+
+	if (lookup == NULL) {
+		/* insert in front of the list */
+		ph->buckets[hash] = link;
+		link->next = ph->first;
+		ph->first = link;
+	}
+	else {
+		/* insert after existing element */
+		link->next = lookup->next;
+		lookup->next = link;
+	}
+		
+	ph->size++;
+
+	if (ph->size > (size*3)) {
+		PHashLink *next = NULL, *first = ph->first;
+
+		ph->cursize = PHashSizes[++ph->cursize_id];
+		MEM_freeN(ph->buckets);
+		ph->buckets = (PHashLink**)MEM_callocN(ph->cursize*sizeof(*ph->buckets), "PHashBuckets");
+		ph->size = 0;
+		ph->first = NULL;
+
+		for (link = first; link; link = next) {
+			next = link->next;
+			phash_insert(ph, link);
+		}
+	}
+}
+
+PHashLink *phash_lookup(PHash *ph, PHashKey key)
+{
+	PHashLink *link;
+	int hash = PHASH_hash(ph, key);
+
+	for (link = ph->buckets[hash]; link; link = link->next)
+		if (link->key == key)
+			return link;
+		else if (PHASH_hash(ph, link->key) != hash)
+			return NULL;
+	
+	return link;
+}
+
+PHashLink *phash_next(PHash *ph, PHashKey key, PHashLink *link)
+{
+	int hash = PHASH_hash(ph, key);
+
+	for (link = link->next; link; link = link->next)
+		if (link->key == key)
+			return link;
+		else if (PHASH_hash(ph, link->key) != hash)
+			return NULL;
+	
+	return link;
+}
+
+/* Heap */
+
+#define PHEAP_PARENT(i) ((i-1)>>1)
+#define PHEAP_LEFT(i)   ((i<<1)+1)
+#define PHEAP_RIGHT(i)  ((i<<1)+2)
+#define PHEAP_COMPARE(a, b) (a->value < b->value)
+#define PHEAP_EQUALS(a, b) (a->value == b->value)
+#define PHEAP_SWAP(heap, i, j) \
+	{ SWAP(int, heap->tree[i]->index, heap->tree[j]->index); \
+	  SWAP(PHeapLink*, heap->tree[i], heap->tree[j]);  }
+
+static void pheap_down(PHeap *heap, int i)
+{
+	while (P_TRUE) {
+		int size = heap->size, smallest;
+		int l = PHEAP_LEFT(i);
+		int r = PHEAP_RIGHT(i);
+
+		smallest = ((l < size) && PHEAP_COMPARE(heap->tree[l], heap->tree[i]))? l: i;
+
+		if ((r < size) && PHEAP_COMPARE(heap->tree[r], heap->tree[smallest]))
+			smallest = r;
+		
+		if (smallest == i)
+			break;
+
+		PHEAP_SWAP(heap, i, smallest);
+		i = smallest;
+	}
+}
+
+static void pheap_up(PHeap *heap, int i)
+{
+	while (i > 0) {
+		int p = PHEAP_PARENT(i);
+
+		if (PHEAP_COMPARE(heap->tree[p], heap->tree[i]))
+			break;
+
+		PHEAP_SWAP(heap, p, i);
+		i = p;
+	}
+}
+
+PHeap *pheap_new()
+{
+	/* TODO: replace mallocN with something faster */
+
+	PHeap *heap = (PHeap*)MEM_callocN(sizeof(PHeap), "PHeap");
+	heap->bufsize = 1;
+	heap->tree = (PHeapLink**)MEM_mallocN(sizeof(PHeapLink*), "PHeapTree");
+
+	return heap;
+}
+
+void pheap_delete(PHeap *heap)
+{
+	MEM_freeN(heap->tree);
+	MEM_freeN(heap);
+}
+
+PHeapLink *pheap_insert(PHeap *heap, float value, void *ptr)
+{
+	PHeapLink *link;
+
+	if ((heap->size + 1) > heap->bufsize) {
+		int newsize = heap->bufsize*2;
+
+		PHeapLink **ntree = (PHeapLink**)MEM_mallocN(newsize*sizeof(PHeapLink*), "PHeapTree");
+		memcpy(ntree, heap->tree, sizeof(PHeapLink*)*heap->size);
+		MEM_freeN(heap->tree);
+
+		heap->tree = ntree;
+		heap->bufsize = newsize;
+	}
+
+	param_assert(heap->size < heap->bufsize);
+
+	link = MEM_mallocN(sizeof *link, "PHeapLink");
+	link->value = value;
+	link->ptr = ptr;
+	link->index = heap->size;
+
+	heap->tree[link->index] = link;
+
+	heap->size++;
+
+	pheap_up(heap, heap->size-1);
+
+	return link;
+}
+
+int pheap_empty(PHeap *heap)
+{
+	return (heap->size == 0);
+}
+
+int pheap_size(PHeap *heap)
+{
+	return heap->size;
+}
+
+void *pheap_min(PHeap *heap)
+{
+	return heap->tree[0]->ptr;
+}
+
+void *pheap_popmin(PHeap *heap)
+{
+	void *ptr = heap->tree[0]->ptr;
+
+	MEM_freeN(heap->tree[0]);
+
+	if (heap->size == 1)
+		heap->size--;
+	else {
+		PHEAP_SWAP(heap, 0, heap->size-1);
+		heap->size--;
+
+		pheap_down(heap, 0);
+	}
+
+	return ptr;
+}
+
+static void pheap_remove(PHeap *heap, PHeapLink *link)
+{
+	int i = link->index;
+
+	while (i > 0) {
+		int p = PHEAP_PARENT(i);
+
+		PHEAP_SWAP(heap, p, i);
+		i = p;
+	}
+
+	pheap_popmin(heap);
+}
+
+/* Construction */
+
+PEdge *p_wheel_edge_next(PEdge *e)
+{
+	return e->next->next->pair;
+}
+
+PEdge *p_wheel_edge_prev(PEdge *e)
+{   
+	return (e->pair)? e->pair->next: NULL;
+}
+
+static PVert *p_vert_add(PChart *chart, PHashKey key, float *co, PEdge *e)
+{
+	PVert *v = (PVert*)BLI_memarena_alloc(chart->handle->arena, sizeof *v);
+	v->co = co;
+	v->link.key = key;
+	v->edge = e;
+
+	phash_insert(chart->verts, (PHashLink*)v);
+
+	return v;
+}
+
+static PVert *p_vert_lookup(PChart *chart, PHashKey key, float *co, PEdge *e)
+{
+	PVert *v = (PVert*)phash_lookup(chart->verts, key);
+
+	if (v)
+		return v;
+	else
+		return p_vert_add(chart, key, co, e);
+}
+
+static PVert *p_vert_copy(PChart *chart, PVert *v)
+{
+	PVert *nv = (PVert*)BLI_memarena_alloc(chart->handle->arena, sizeof *nv);
+	nv->co = v->co;
+	nv->uv[0] = v->uv[0];
+	nv->uv[1] = v->uv[1];
+	nv->link.key = v->link.key;
+	nv->edge = v->edge;
+
+	phash_insert(chart->verts, (PHashLink*)nv);
+
+	return nv;
+}
+
+static PEdge *p_edge_lookup(PChart *chart, PHashKey *vkeys)
+{
+	PHashKey key = vkeys[0]^vkeys[1];
+	PEdge *e = (PEdge*)phash_lookup(chart->edges, key);
+
+	while (e) {
+		if ((e->vert->link.key == vkeys[0]) && (e->next->vert->link.key == vkeys[1]))
+			return e;
+		else if ((e->vert->link.key == vkeys[1]) && (e->next->vert->link.key == vkeys[0]))
+			return e;
+
+		e = (PEdge*)phash_next(chart->edges, key, (PHashLink*)e);
+	}
+
+	return NULL;
+}
+
+static void p_face_flip(PFace *f)
+{
+	PEdge *e1 = f->edge, *e2 = e1->next, *e3 = e2->next;
+	PVert *v1 = e1->vert, *v2 = e2->vert, *v3 = e3->vert;
+	int f1 = e1->flag, f2 = e2->flag, f3 = e3->flag;
+
+	e1->vert = v2;
+	e1->next = e3;
+	e1->flag = (f1 & ~PEDGE_VERTEX_FLAGS) | (f2 & PEDGE_VERTEX_FLAGS);
+
+	e2->vert = v3;
+	e2->next = e1;
+	e2->flag = (f2 & ~PEDGE_VERTEX_FLAGS) | (f3 & PEDGE_VERTEX_FLAGS);
+
+	e3->vert = v1;
+	e3->next = e2;
+	e3->flag = (f3 & ~PEDGE_VERTEX_FLAGS) | (f1 & PEDGE_VERTEX_FLAGS);
+}
+
+static void p_vert_load_pin_uvs(PVert *v)
+{
+	PEdge *e;
+	int nedges = 0;
+
+	v->uv[0] = v->uv[1] = 0.0f;
+	nedges = 0;
+	e = v->edge;
+	do {
+		if (e->orig_uv && (e->flag & PEDGE_PIN)) {
+			v->flag |= PVERT_PIN;
+			v->uv[0] += e->orig_uv[0];
+			v->uv[1] += e->orig_uv[1];
+			nedges++;
+		}
+
+		e = p_wheel_edge_next(e);
+	} while (e && e != (v->edge));
+
+	if (nedges > 0) {
+		v->uv[0] /= nedges;
+		v->uv[1] /= nedges;
+	}
+}
+
+static void p_vert_load_select_uvs(PVert *v)
+{
+	PEdge *e;
+	int nedges = 0;
+
+	v->uv[0] = v->uv[1] = 0.0f;
+	nedges = 0;
+	e = v->edge;
+	do {
+		if (e->orig_uv && (e->flag & PEDGE_SELECT))
+			v->flag |= PVERT_SELECT;
+
+		v->uv[0] += e->orig_uv[0];
+		v->uv[1] += e->orig_uv[1];
+		nedges++;
+
+		e = p_wheel_edge_next(e);
+	} while (e && e != (v->edge));
+
+	if (nedges > 0) {
+		v->uv[0] /= nedges;
+		v->uv[1] /= nedges;
+	}
+}
+
+static void p_extrema_verts(PChart *chart, PVert **v1, PVert **v2)
+{
+	float minv[3], maxv[3], dirlen;
+	PVert *v, *minvert[3], *maxvert[3];
+	int i, dir;
+
+	/* find minimum and maximum verts over x/y/z axes */
+	minv[0] = minv[1] = minv[2] = 1e20;
+	maxv[0] = maxv[1] = maxv[2] = -1e20;
+
+	minvert[0] = minvert[1] = minvert[2] = NULL;
+	maxvert[0] = maxvert[1] = maxvert[2] = NULL;
+
+	for (v = (PVert*)chart->verts->first; v; v=v->link.next) {
+		for (i = 0; i < 3; i++) {
+			if (v->co[i] < minv[i]) {
+				minv[i] = v->co[i];
+				minvert[i] = v;
+			}
+			if (v->co[i] > maxv[i]) {
+				maxv[i] = v->co[i];
+				maxvert[i] = v;
+			}
+		}
+	}
+
+	/* find axes with longest distance */
+	dir = 0;
+	dirlen = -1.0;
+
+	for (i = 0; i < 3; i++) {
+		if (maxv[i] - minv[i] > dirlen) {
+			dir = i;
+			dirlen = maxv[i] - minv[i];
+		}
+	}
+
+	if (minvert[dir] == maxvert[dir]) {
+		/* degenerate case */
+		PFace *f = (PFace*)chart->faces->first;
+		*v1 = f->edge->vert;
+		*v2 = f->edge->next->vert;
+	}
+	else {
+		*v1 = minvert[dir];
+		*v2 = maxvert[dir];
+	}
+}
+
+static float p_vec_normalise(float *v)
+{
+   float d;
+   
+    d = sqrt(v[0]*v[0] + v[1]*v[1] + v[2]*v[2]);
+
+    if(d != 0.0f) {
+		d = 1.0f/d;
+
+		v[0] *= d;
+		v[1] *= d;
+		v[2] *= d;
+	}
+
+	return d;
+}
+
+static float p_vec_angle_cos(float *v1, float *v2, float *v3)
+{
+	float d1[3], d2[3];
+
+	d1[0] = v1[0] - v2[0];
+	d1[1] = v1[1] - v2[1];
+	d1[2] = v1[2] - v2[2];
+
+	d2[0] = v3[0] - v2[0];
+	d2[1] = v3[1] - v2[1];
+	d2[2] = v3[2] - v2[2];
+
+	p_vec_normalise(d1);
+	p_vec_normalise(d2);
+
+	return d1[0]*d2[0] + d1[1]*d2[1] + d1[2]*d2[2];
+}
+
+static float p_vec_angle(float *v1, float *v2, float *v3)
+{
+	float dot = p_vec_angle_cos(v1, v2, v3);
+
+	if (dot <= -1.0f)
+		return (float)M_PI;
+	else if (dot >= 1.0f)
+		return 0.0f;
+	else
+		return (float)acos(dot);
+}
+
+static void p_face_angles(PFace *f, float *a1, float *a2, float *a3)
+{
+	PEdge *e1 = f->edge, *e2 = e1->next, *e3 = e2->next;
+	PVert *v1 = e1->vert, *v2 = e2->vert, *v3 = e3->vert;
+
+	*a1 = p_vec_angle(v3->co, v1->co, v2->co);
+	*a2 = p_vec_angle(v1->co, v2->co, v3->co);
+	*a3 = M_PI - *a2 - *a1;
+}
+
+static float p_face_area(PFace *f)
+{
+	PEdge *e1 = f->edge, *e2 = e1->next, *e3 = e2->next;
+	PVert *v1 = e1->vert, *v2 = e2->vert, *v3 = e3->vert;
+
+	return AreaT3Dfl(v1->co, v2->co, v3->co);
+}
+
+static float p_face_uv_area_signed(PFace *f)
+{
+	PEdge *e1 = f->edge, *e2 = e1->next, *e3 = e2->next;
+	PVert *v1 = e1->vert, *v2 = e2->vert, *v3 = e3->vert;
+
+	return 0.5f*(((v2->uv[0]-v1->uv[0]) * (v3->uv[1]-v1->uv[1])) - 
+	            ((v3->uv[0]-v1->uv[0]) * (v2->uv[1]-v1->uv[1])));
+}
+
+static float p_face_uv_area(PFace *f)
+{
+	return fabs(p_face_uv_area_signed(f));
+}
+
+static void p_chart_area(PChart *chart, float *uv_area, float *area)
+{
+	PFace *f;
+
+	*uv_area = *area = 0.0f;
+
+	for (f=(PFace*)chart->faces->first; f; f=f->link.next) {
+		*uv_area += p_face_uv_area(f);
+		*area += p_face_area(f);
+	}
+}
+
+static PChart *p_chart_new(PHandle *handle)
+{
+	PChart *chart = (PChart*)MEM_callocN(sizeof*chart, "PChart");
+	chart->verts = phash_new(1);
+	chart->edges = phash_new(1);
+	chart->faces = phash_new(1);
+	chart->handle = handle;
+
+	return chart;
+}
+
+static void p_chart_delete(PChart *chart)
+{
+	/* the actual links are free by memarena */
+	phash_delete(chart->verts);
+	phash_delete(chart->edges);
+	phash_delete(chart->faces);
+
+	MEM_freeN(chart);
+}
+
+static PBool p_edge_implicit_seam(PEdge *e, PEdge *ep)
+{
+	float *uv1, *uv2, *uvp1, *uvp2;
+	float limit[2];
+
+	uv1 = e->orig_uv;
+	uv2 = e->next->orig_uv;
+
+	if (e->vert->link.key == ep->vert->link.key) {
+		uvp1 = ep->orig_uv;
+		uvp2 = ep->next->orig_uv;
+	}
+	else {
+		uvp1 = ep->next->orig_uv;
+		uvp2 = ep->orig_uv;
+	}
+
+	get_connected_limit_tface_uv(limit);
+
+	if((fabs(uv1[0]-uvp1[0]) > limit[0]) && (fabs(uv1[1]-uvp1[1]) > limit[1])) {
+		e->flag |= PEDGE_SEAM;
+		ep->flag |= PEDGE_SEAM;
+		return P_TRUE;
+	}
+	if((fabs(uv2[0]-uvp2[0]) > limit[0]) && (fabs(uv2[1]-uvp2[1]) > limit[1])) {
+		e->flag |= PEDGE_SEAM;
+		ep->flag |= PEDGE_SEAM;
+		return P_TRUE;
+	}
+	
+	return P_FALSE;
+}
+
+static PBool p_edge_has_pair(PChart *chart, PEdge *e, PEdge **pair, PBool impl)
+{
+	PHashKey key;
+	PEdge *pe;
+	PVert *v1, *v2;
+	PHashKey key1 = e->vert->link.key;
+	PHashKey key2 = e->next->vert->link.key;
+
+	if (e->flag & PEDGE_SEAM)
+		return P_FALSE;
+	
+	key = key1 ^ key2;
+	pe = (PEdge*)phash_lookup(chart->edges, key);
+	*pair = NULL;
+
+	while (pe) {
+		if (pe != e) {
+			v1 = pe->vert;
+			v2 = pe->next->vert;
+
+			if (((v1->link.key == key1) && (v2->link.key == key2)) ||
+				((v1->link.key == key2) && (v2->link.key == key1))) {
+
+				/* don't connect seams and t-junctions */
+				if ((pe->flag & PEDGE_SEAM) || *pair ||
+				    (impl && p_edge_implicit_seam(e, pe))) {
+					*pair = NULL;
+					return P_FALSE;
+				}
+
+				*pair = pe;
+			}
+		}
+
+		pe = (PEdge*)phash_next(chart->edges, key, (PHashLink*)pe);
+	}
+
+	if (*pair && (e->vert == (*pair)->vert)) {
+		if ((*pair)->next->pair || (*pair)->next->next->pair) {
+			/* non unfoldable, maybe mobius ring or klein bottle */
+			*pair = NULL;
+			return P_FALSE;
+		}
+	}
+
+	return (*pair != NULL);
+}
+
+static PBool p_edge_connect_pair(PChart *chart, PEdge *e, PEdge ***stack, PBool impl)
+{
+	PEdge *pair;
+
+	if(!e->pair && p_edge_has_pair(chart, e, &pair, impl)) {
+		if (e->vert == pair->vert)
+			p_face_flip(pair->face);
+
+		e->pair = pair;
+		pair->pair = e;
+
+		if (!(pair->face->flag & PFACE_CONNECTED)) {
+			**stack = pair;
+			(*stack)++;
+		}
+	}
+
+	return (e->pair != NULL);
+}
+
+static int p_connect_pairs(PChart *chart, PBool impl)
+{
+	PEdge **stackbase = MEM_mallocN(sizeof*stackbase * phash_size(chart->faces), "Pstackbase");
+	PEdge **stack = stackbase;
+	PFace *f, *first;
+	PEdge *e, *e1, *e2;
+	int ncharts = 0;
+
+	/* connect pairs, count edges, set vertex-edge pointer to a pairless edge */
+	for (first=(PFace*)chart->faces->first; first; first=first->link.next) {
+		if (first->flag & PFACE_CONNECTED)
+			continue;
+
+		*stack = first->edge;
+		stack++;
+
+		while (stack != stackbase) {
+			stack--;
+			e = *stack;
+			e1 = e->next;
+			e2 = e1->next;
+
+			f = e->face;
+			f->flag |= PFACE_CONNECTED;
+
+			/* assign verts to charts so we can sort them later */
+			f->u.chart = ncharts;
+
+			if (!p_edge_connect_pair(chart, e, &stack, impl))
+				e->vert->edge = e;
+			if (!p_edge_connect_pair(chart, e1, &stack, impl))
+				e1->vert->edge = e1;
+			if (!p_edge_connect_pair(chart, e2, &stack, impl))
+				e2->vert->edge = e2;
+		}
+
+		ncharts++;
+	}
+
+	MEM_freeN(stackbase);
+
+	return ncharts;
+}
+
+static void p_split_vert(PChart *chart, PEdge *e)
+{
+	PEdge *we, *lastwe = NULL;
+	PVert *v = e->vert;
+	PBool copy = P_TRUE;
+
+	if (e->flag & PEDGE_VERTEX_SPLIT)
+		return;
+
+	/* rewind to start */
+	lastwe = e;
+	for (we = p_wheel_edge_prev(e); we && (we != e); we = p_wheel_edge_prev(we))
+		lastwe = we;
+	
+	/* go over all edges in wheel */
+	for (we = lastwe; we; we = p_wheel_edge_next(we)) {
+		if (we->flag & PEDGE_VERTEX_SPLIT)
+			break;
+
+		we->flag |= PEDGE_VERTEX_SPLIT;
+
+		if (we == v->edge) {
+			/* found it, no need to copy */
+			copy = P_FALSE;
+			phash_insert(chart->verts, (PHashLink*)v);
+		}
+	}
+
+	if (copy) {
+		/* not found, copying */
+		v = p_vert_copy(chart, v);
+		v->edge = lastwe;
+
+		we = lastwe;
+		do {
+			we->vert = v;
+			we = p_wheel_edge_next(we);
+		} while (we && (we != lastwe));
+	}
+}
+
+static PChart **p_split_charts(PHandle *handle, PChart *chart, int ncharts)
+{
+	PChart **charts = MEM_mallocN(sizeof*charts * ncharts, "PCharts"), *nchart;
+	PFace *f, *nextf;
+	int i;
+
+	for (i = 0; i < ncharts; i++)
+		charts[i] = p_chart_new(handle);
+
+	f = (PFace*)chart->faces->first;
+	while (f) {
+		PEdge *e1 = f->edge, *e2 = e1->next, *e3 = e2->next;
+		nextf = f->link.next;
+
+		nchart = charts[f->u.chart];
+
+		phash_insert(nchart->faces, (PHashLink*)f);
+		phash_insert(nchart->edges, (PHashLink*)e1);
+		phash_insert(nchart->edges, (PHashLink*)e2);
+		phash_insert(nchart->edges, (PHashLink*)e3);
+
+		p_split_vert(nchart, e1);
+		p_split_vert(nchart, e2);
+		p_split_vert(nchart, e3);
+
+		f = nextf;
+	}
+
+	return charts;
+}
+
+static void p_face_backup_uvs(PFace *f)
+{
+	PEdge *e1 = f->edge, *e2 = e1->next, *e3 = e2->next;
+
+	e1->old_uv[0] = e1->orig_uv[0];
+	e1->old_uv[1] = e1->orig_uv[1];
+	e2->old_uv[0] = e2->orig_uv[0];
+	e2->old_uv[1] = e2->orig_uv[1];
+	e3->old_uv[0] = e3->orig_uv[0];
+	e3->old_uv[1] = e3->orig_uv[1];
+}
+
+static void p_face_restore_uvs(PFace *f)
+{
+	PEdge *e1 = f->edge, *e2 = e1->next, *e3 = e2->next;
+
+	e1->orig_uv[0] = e1->old_uv[0];
+	e1->orig_uv[1] = e1->old_uv[1];
+	e2->orig_uv[0] = e2->old_uv[0];
+	e2->orig_uv[1] = e2->old_uv[1];
+	e3->orig_uv[0] = e3->old_uv[0];
+	e3->orig_uv[1] = e3->old_uv[1];
+}
+
+static PFace *p_face_add(PChart *chart, ParamKey key, ParamKey *vkeys,
+                         float *co[3], float *uv[3], int i1, int i2, int i3,
+                         ParamBool *pin, ParamBool *select)
+{
+	PFace *f;
+	PEdge *e1, *e2, *e3;
+
+	/* allocate */
+	f = (PFace*)BLI_memarena_alloc(chart->handle->arena, sizeof *f);
+
+	e1 = (PEdge*)BLI_memarena_alloc(chart->handle->arena, sizeof *e1);
+	e2 = (PEdge*)BLI_memarena_alloc(chart->handle->arena, sizeof *e2);
+	e3 = (PEdge*)BLI_memarena_alloc(chart->handle->arena, sizeof *e3);
+
+	/* set up edges */
+	f->edge = e1;
+	e1->face = e2->face = e3->face = f;
+
+	e1->next = e2;
+	e2->next = e3;
+	e3->next = e1;
+
+	if (co && uv) {
+		e1->vert = p_vert_lookup(chart, vkeys[i1], co[i1], e1);
+		e2->vert = p_vert_lookup(chart, vkeys[i2], co[i2], e2);
+		e3->vert = p_vert_lookup(chart, vkeys[i3], co[i3], e3);
+
+		e1->orig_uv = uv[i1];
+		e2->orig_uv = uv[i2];
+		e3->orig_uv = uv[i3];
+	}
+	else {
+		/* internal call to add face */
+		e1->vert = e2->vert = e3->vert = NULL;
+		e1->orig_uv = e2->orig_uv = e3->orig_uv = NULL;
+	}
+
+	if (pin) {
+		if (pin[i1]) e1->flag |= PEDGE_PIN;
+		if (pin[i2]) e2->flag |= PEDGE_PIN;
+		if (pin[i3]) e3->flag |= PEDGE_PIN;
+	}
+
+	if (select) {
+		if (select[i1]) e1->flag |= PEDGE_SELECT;
+		if (select[i2]) e2->flag |= PEDGE_SELECT;
+		if (select[i3]) e3->flag |= PEDGE_SELECT;
+	}
+
+	/* insert into hash */
+	f->link.key = key;
+	phash_insert(chart->faces, (PHashLink*)f);
+
+	e1->link.key = vkeys[i1]^vkeys[i2];
+	e2->link.key = vkeys[i2]^vkeys[i3];
+	e3->link.key = vkeys[i3]^vkeys[i1];
+
+	phash_insert(chart->edges, (PHashLink*)e1);
+	phash_insert(chart->edges, (PHashLink*)e2);
+	phash_insert(chart->edges, (PHashLink*)e3);
+
+	return f;
+}
+
+static PBool p_quad_split_direction(float **co)
+{
+    float a1, a2;
+	
+	a1 = p_vec_angle_cos(co[0], co[1], co[2]);
+	a1 += p_vec_angle_cos(co[1], co[0], co[2]);
+	a1 += p_vec_angle_cos(co[2], co[0], co[1]);
+
+	a2 = p_vec_angle_cos(co[0], co[1], co[3]);
+	a2 += p_vec_angle_cos(co[1], co[0], co[3]);
+	a2 += p_vec_angle_cos(co[3], co[0], co[1]);
+
+	return (a1 > a2);
+}
+
+static float p_edge_length(PEdge *e)
+{
+    PVert *v1 = e->vert, *v2 = e->next->vert;
+    float d[3];
+
+    d[0] = v2->co[0] - v1->co[0];
+    d[1] = v2->co[1] - v1->co[1];
+    d[2] = v2->co[2] - v1->co[2];
+
+    return sqrt(d[0]*d[0] + d[1]*d[1] + d[2]*d[2]);
+}
+
+static float p_edge_uv_length(PEdge *e)
+{
+    PVert *v1 = e->vert, *v2 = e->next->vert;
+    float d[3];
+
+    d[0] = v2->uv[0] - v1->uv[0];
+    d[1] = v2->uv[1] - v1->uv[1];
+
+    return sqrt(d[0]*d[0] + d[1]*d[1]);
+}
+
+void p_chart_uv_bbox(PChart *chart, float *minv, float *maxv)
+{
+	PVert *v;
+
+	INIT_MINMAX2(minv, maxv);
+
+	for (v=(PVert*)chart->verts->first; v; v=v->link.next) {
+		DO_MINMAX2(v->uv, minv, maxv);
+	}
+}
+
+static void p_chart_uv_scale(PChart *chart, float scale)
+{
+	PVert *v;
+
+	for (v=(PVert*)chart->verts->first; v; v=v->link.next) {
+		v->uv[0] *= scale;
+		v->uv[1] *= scale;
+	}
+}
+
+static void p_chart_uv_translate(PChart *chart, float trans[2])
+{
+	PVert *v;
+
+	for (v=(PVert*)chart->verts->first; v; v=v->link.next) {
+		v->uv[0] += trans[0];
+		v->uv[1] += trans[1];
+	}
+}
+
+static void p_chart_boundaries(PChart *chart, int *nboundaries, PEdge **outer)
+{   
+    PEdge *e, *be;
+    float len, maxlen = -1.0;
+
+	*nboundaries = 0;
+	*outer = NULL;
+
+	for (e=(PEdge*)chart->edges->first; e; e=e->link.next) {
+        if (e->pair || (e->flag & PEDGE_DONE))
+            continue;
+
+		(*nboundaries)++;
+        len = 0.0f;
+    
+		be = e;
+		do {
+            be->flag |= PEDGE_DONE;
+            len += p_edge_length(be);
+			be = be->next->vert->edge;
+        } while(be != e);
+
+        if (len > maxlen) {
+			*outer = e;
+            maxlen = len;
+        }
+    }
+
+	for (e=(PEdge*)chart->edges->first; e; e=e->link.next)
+        e->flag &= ~PEDGE_DONE;
+}
+
+static float p_edge_boundary_angle(PEdge *e)
+{
+	PEdge *we;
+	PVert *v, *v1, *v2;
+	float angle;
+	int n = 0;
+
+	v = e->vert;
+
+	/* concave angle check -- could be better */
+	angle = M_PI;
+
+	we = v->edge;
+	do {
+		v1 = we->next->vert;
+		v2 = we->next->next->vert;
+		angle -= p_vec_angle(v1->co, v->co, v2->co);
+
+		we = we->next->next->pair;
+		n++;
+	} while (we && (we != v->edge));
+
+	return angle;
+}
+
+static PEdge *p_boundary_edge_next(PEdge *e)
+{
+	return e->next->vert->edge;
+}
+
+static PEdge *p_boundary_edge_prev(PEdge *e)
+{
+    PEdge *we = e, *last;
+
+	do {
+		last = we;
+		we = p_wheel_edge_next(we);
+	} while (we && (we != e));
+
+	return last->next->next;
+}
+
+static void p_chart_fill_boundary(PChart *chart, PEdge *be, int nedges)
+{
+	PEdge *e, *e1, *e2;
+	PHashKey vkeys[3];
+	PFace *f;
+	struct PHeap *heap = pheap_new(nedges);
+	float angle;
+
+	e = be;
+	do {
+		angle = p_edge_boundary_angle(e);
+		e->u.heaplink = pheap_insert(heap, angle, e);
+
+		e = e->next->vert->edge;
+	} while(e != be);
+
+	if (nedges == 2) {
+		/* no real boundary, but an isolated seam */
+		e = be->next->vert->edge;
+		e->pair = be;
+		be->pair = e;
+
+		pheap_remove(heap, e->u.heaplink);
+		pheap_remove(heap, be->u.heaplink);
+	}
+	else {
+		while (nedges > 2) {
+			PEdge *ne, *ne1, *ne2;
+
+			e = pheap_popmin(heap);
+
+			e1 = p_boundary_edge_prev(e);
+			e2 = p_boundary_edge_next(e);
+
+			pheap_remove(heap, e1->u.heaplink);
+			pheap_remove(heap, e2->u.heaplink);
+			e->u.heaplink = e1->u.heaplink = e2->u.heaplink = NULL;
+
+			e->flag |= PEDGE_FILLED;
+			e1->flag |= PEDGE_FILLED;
+
+			vkeys[0] = e->vert->link.key;
+			vkeys[1] = e1->vert->link.key;
+			vkeys[2] = e2->vert->link.key;
+
+			f = p_face_add(chart, -1, vkeys, NULL, NULL, 0, 1, 2, NULL, NULL);
+			f->flag |= PFACE_FILLED;
+
+			ne = f->edge->next->next;
+			ne1 = f->edge;
+			ne2 = f->edge->next;
+
+			ne->flag = ne1->flag = ne2->flag = PEDGE_FILLED;
+
+			e->pair = ne;
+			ne->pair = e;
+			e1->pair = ne1;
+			ne1->pair = e1;
+
+			ne->vert = e2->vert;
+			ne1->vert = e->vert;
+			ne2->vert = e1->vert;
+
+			if (nedges == 3) {
+				e2->pair = ne2;
+				ne2->pair = e2;
+			}
+			else {
+				ne2->vert->edge = ne2;
+				
+				ne2->u.heaplink = pheap_insert(heap, p_edge_boundary_angle(ne2), ne2);
+				e2->u.heaplink = pheap_insert(heap, p_edge_boundary_angle(e2), e2);
+			}
+
+			nedges--;
+		}
+	}
+
+	pheap_delete(heap);
+}
+
+static void p_chart_fill_boundaries(PChart *chart, PEdge *outer)
+{
+	PEdge *e, *enext, *be;
+	int nedges;
+
+	for (e=(PEdge*)chart->edges->first; e; e=e->link.next) {
+		enext = e->link.next;
+
+        if (e->pair || (e->flag & PEDGE_FILLED))
+            continue;
+
+		nedges = 0;
+		be = e;
+		do {
+			be->flag |= PEDGE_FILLED;
+			be = be->next->vert->edge;
+			nedges++;
+		} while(be != e);
+
+		if (e != outer)
+			p_chart_fill_boundary(chart, e, nedges);
+    }
+}
+
+static void p_flush_uvs(PChart *chart)
+{
+	PEdge *e;
+
+	for (e=(PEdge*)chart->edges->first; e; e=e->link.next) {
+		if (e->orig_uv) {
+			e->orig_uv[0] = e->vert->uv[0];
+			e->orig_uv[1] = e->vert->uv[1];
+		}
+	}
+}
+
+/* Exported */
+
+ParamHandle *param_construct_begin()
+{
+	PHandle *handle = MEM_callocN(sizeof*handle, "PHandle");
+	handle->construction_chart = p_chart_new(handle);
+	handle->state = PHANDLE_STATE_ALLOCATED;
+	handle->arena = BLI_memarena_new((1<<16));
+	
+	return (ParamHandle*)handle;
+}
+
+void param_delete(ParamHandle *handle)
+{
+	PHandle *phandle = (PHandle*)handle;
+	int i;
+
+	param_assert((phandle->state == PHANDLE_STATE_ALLOCATED) ||
+	             (phandle->state == PHANDLE_STATE_CONSTRUCTED));
+
+	for (i = 0; i < phandle->ncharts; i++)
+		p_chart_delete(phandle->charts[i]);
+	
+	if (phandle->charts)
+		MEM_freeN(phandle->charts);
+
+	if (phandle->construction_chart)
+		p_chart_delete(phandle->construction_chart);
+
+	BLI_memarena_free(phandle->arena);
+	MEM_freeN(phandle);
+}
+
+void param_face_add(ParamHandle *handle, ParamKey key, int nverts,
+                    ParamKey *vkeys, float **co, float **uv,
+                    ParamBool *pin, ParamBool *select)
+{
+	PHandle *phandle = (PHandle*)handle;
+	PChart *chart = phandle->construction_chart;
+
+	param_assert(phash_lookup(chart->faces, key) == NULL);
+	param_assert(phandle->state == PHANDLE_STATE_ALLOCATED);
+	param_assert((nverts == 3) || (nverts == 4));
+
+	if (nverts == 4) {
+		if (p_quad_split_direction(co)) {
+			p_face_add(chart, key, vkeys, co, uv, 0, 1, 2, pin, select);
+			p_face_add(chart, key, vkeys, co, uv, 0, 2, 3, pin, select);
+		}
+		else {
+			p_face_add(chart, key, vkeys, co, uv, 0, 1, 3, pin, select);
+			p_face_add(chart, key, vkeys, co, uv, 1, 2, 3, pin, select);
+		}
+	}
+	else
+		p_face_add(chart, key, vkeys, co, uv, 0, 1, 2, pin, select);
+}
+
+void param_edge_set_seam(ParamHandle *handle, ParamKey *vkeys)
+{
+	PHandle *phandle = (PHandle*)handle;
+	PChart *chart = phandle->construction_chart;
+	PEdge *e;
+
+	param_assert(phandle->state == PHANDLE_STATE_ALLOCATED);
+
+	e = p_edge_lookup(chart, vkeys);
+	if (e)
+		e->flag |= PEDGE_SEAM;
+}
+
+void param_construct_end(ParamHandle *handle, ParamBool fill, ParamBool impl)
+{
+	PHandle *phandle = (PHandle*)handle;
+	PChart *chart = phandle->construction_chart;
+	int i, j, nboundaries = 0;
+	PEdge *outer;
+
+	param_assert(phandle->state == PHANDLE_STATE_ALLOCATED);
+
+	phandle->ncharts = p_connect_pairs(chart, impl);
+	phandle->charts = p_split_charts(phandle, chart, phandle->ncharts);
+
+	p_chart_delete(chart);
+	phandle->construction_chart = NULL;
+
+	for (i = j = 0; i < phandle->ncharts; i++) {
+		p_chart_boundaries(phandle->charts[i], &nboundaries, &outer);
+
+		if (nboundaries > 0) {
+			phandle->charts[j] = phandle->charts[i];
+			j++;
+
+			if (fill && (nboundaries > 1))
+				p_chart_fill_boundaries(phandle->charts[i], outer);
+		}
+		else
+			p_chart_delete(phandle->charts[i]);
+	}
+
+	phandle->ncharts = j;
+
+	phandle->state = PHANDLE_STATE_CONSTRUCTED;
+}
+
+/* Least Squares Conformal Maps */
+
+static void p_chart_lscm_load_solution(PChart *chart)
+{
+	PVert *pin = chart->u.lscm.singlepin, *v;
+	float translation[2];
+
+	if (pin) {
+		translation[0] = pin->uv[0] - nlGetVariable(2*pin->u.index);
+		translation[1] = pin->uv[1] - nlGetVariable(2*pin->u.index + 1);
+	}
+	else
+		translation[0] = translation[1] = 0.0f;
+
+	for (v=(PVert*)chart->verts->first; v; v=v->link.next) {
+		v->uv[0] = nlGetVariable(2*v->u.index) + translation[0];
+		v->uv[1] = nlGetVariable(2*v->u.index + 1) + translation[1];
+	}
+}
+
+static void p_chart_lscm_begin(PChart *chart)
+{
+	PVert *v, *pin1, *pin2;
+	int npins = 0, id = 0;
+
+	nlNewContext();
+	nlSolverParameteri(NL_NB_VARIABLES, 2*phash_size(chart->verts));
+	nlSolverParameteri(NL_LEAST_SQUARES, NL_TRUE);
+
+	/* give vertices matrix indices and count pins */
+	for (v=(PVert*)chart->verts->first; v; v=v->link.next) {
+		p_vert_load_pin_uvs(v);
+
+		if (v->flag & PVERT_PIN) {
+			npins++;
+			chart->u.lscm.singlepin = v;
+		}
+
+		v->u.index = id++;
+	}
+
+	if (npins <= 1) {
+		/* not enough pins, lets find some ourself */
+		p_extrema_verts(chart, &pin1, &pin2);
+
+		chart->u.lscm.pin1 = pin1;
+		chart->u.lscm.pin2 = pin2;
+	}
+	else
+		chart->u.lscm.singlepin = NULL;
+
+	chart->u.lscm.context = nlGetCurrent();
+}
+
+static PBool p_chart_lscm_solve(PChart *chart)
+{
+	PVert *v, *pin1 = chart->u.lscm.pin1, *pin2 = chart->u.lscm.pin2;
+	PFace *f;
+
+	nlMakeCurrent(chart->u.lscm.context);
+
+	nlBegin(NL_SYSTEM);
+
+	if (chart->u.lscm.pin1) {
+		nlLockVariable(2*pin1->u.index);
+		nlLockVariable(2*pin1->u.index + 1);
+		nlLockVariable(2*pin2->u.index);
+		nlLockVariable(2*pin2->u.index + 1);
+	
+		nlSetVariable(2*pin1->u.index, 0.0f);
+		nlSetVariable(2*pin1->u.index + 1, 0.5f);
+		nlSetVariable(2*pin2->u.index, 1.0f);
+		nlSetVariable(2*pin2->u.index + 1, 0.5f);
+	}
+	else {
+		/* set and lock the pins */
+		for (v=(PVert*)chart->verts->first; v; v=v->link.next) {
+			if (v->flag & PVERT_PIN) {
+				nlLockVariable(2*v->u.index);
+				nlLockVariable(2*v->u.index + 1);
+
+				nlSetVariable(2*v->u.index, v->uv[0]);
+				nlSetVariable(2*v->u.index + 1, v->uv[1]);
+			}
+		}
+	}
+
+	/* construct matrix */
+
+	nlBegin(NL_MATRIX);
+
+	for (f=(PFace*)chart->faces->first; f; f=f->link.next) {
+		PEdge *e1 = f->edge, *e2 = e1->next, *e3 = e2->next;
+		PVert *v1 = e1->vert, *v2 = e2->vert, *v3 = e3->vert;
+		float a1, a2, a3, ratio, cosine, sine;
+		float sina1, sina2, sina3, sinmax;
+
+		if (chart->u.lscm.abf_alpha) {
+			/* use abf angles if passed on */
+			a1 = *(chart->u.lscm.abf_alpha++);
+			a2 = *(chart->u.lscm.abf_alpha++);
+			a3 = *(chart->u.lscm.abf_alpha++);
+		}
+		else
+			p_face_angles(f, &a1, &a2, &a3);
+
+		sina1 = sin(a1);
+		sina2 = sin(a2);
+		sina3 = sin(a3);
+
+		sinmax = MAX3(sina1, sina2, sina3);
+
+		/* shift vertices to find most stable order */
+		#define SHIFT3(type, a, b, c) \
+			{ type tmp; tmp = a; a = c; c = b; b = tmp; }
+
+		if (sina3 != sinmax) {
+			SHIFT3(PVert*, v1, v2, v3);
+			SHIFT3(float, a1, a2, a3);
+			SHIFT3(float, sina1, sina2, sina3);
+
+			if (sina2 == sinmax) {
+				SHIFT3(PVert*, v1, v2, v3);
+				SHIFT3(float, a1, a2, a3);
+				SHIFT3(float, sina1, sina2, sina3);
+			}
+		}
+
+		/* angle based lscm formulation */
+		ratio = sina2/sina3;
+		cosine = cos(a1)*ratio;
+		sine = sina1*ratio;
+
+		nlBegin(NL_ROW);
+		nlCoefficient(2*v1->u.index,   cosine - 1.0);
+		nlCoefficient(2*v1->u.index+1, -sine);
+		nlCoefficient(2*v2->u.index,   -cosine);
+		nlCoefficient(2*v2->u.index+1, sine);
+		nlCoefficient(2*v3->u.index,   1.0);
+		nlEnd(NL_ROW);
+
+		nlBegin(NL_ROW);
+		nlCoefficient(2*v1->u.index,   sine);
+		nlCoefficient(2*v1->u.index+1, cosine - 1.0);
+		nlCoefficient(2*v2->u.index,   -sine);
+		nlCoefficient(2*v2->u.index+1, -cosine);
+		nlCoefficient(2*v3->u.index+1, 1.0);
+		nlEnd(NL_ROW);
+	}
+
+	nlEnd(NL_MATRIX);
+
+	nlEnd(NL_SYSTEM);
+
+	if (nlSolveAdvanced(NULL, NL_TRUE)) {
+		p_chart_lscm_load_solution(chart);
+		p_flush_uvs(chart);
+
+		return P_TRUE;
+	}
+
+	return P_FALSE;
+}
+
+static void p_chart_lscm_end(PChart *chart)
+{
+	if (chart->u.lscm.context)
+		nlDeleteContext(chart->u.lscm.context);
+
+	chart->u.lscm.context = NULL;
+	chart->u.lscm.singlepin = NULL;
+	chart->u.lscm.pin1 = NULL;
+	chart->u.lscm.pin2 = NULL;
+}
+
+void param_lscm_begin(ParamHandle *handle)
+{
+	PHandle *phandle = (PHandle*)handle;
+	int i;
+
+	param_assert(phandle->state == PHANDLE_STATE_CONSTRUCTED);
+	phandle->state = PHANDLE_STATE_LSCM;
+
+	for (i = 0; i < phandle->ncharts; i++)
+		p_chart_lscm_begin(phandle->charts[i]);
+}
+
+void param_lscm_solve(ParamHandle *handle)
+{
+	PHandle *phandle = (PHandle*)handle;
+	PChart *chart;
+	int i;
+	PBool result;
+
+	param_assert(phandle->state == PHANDLE_STATE_LSCM);
+
+	for (i = 0; i < phandle->ncharts; i++) {
+		chart = phandle->charts[i];
+
+		if (chart->u.lscm.context) {
+			result = p_chart_lscm_solve(chart);
+
+			if (!result || (chart->u.lscm.pin1))
+				p_chart_lscm_end(chart);
+		}
+	}
+}
+
+void param_lscm_end(ParamHandle *handle)
+{
+	PHandle *phandle = (PHandle*)handle;
+	int i;
+
+	param_assert(phandle->state == PHANDLE_STATE_LSCM);
+
+	for (i = 0; i < phandle->ncharts; i++)
+		p_chart_lscm_end(phandle->charts[i]);
+
+	phandle->state = PHANDLE_STATE_CONSTRUCTED;
+}
+
+/* Stretch */
+
+#define P_STRETCH_ITER 20
+
+static void p_stretch_pin_boundary(PChart *chart)
+{
+	PVert *v;
+
+	for(v=(PVert*)chart->verts->first; v; v=v->link.next)
+		if (v->edge->pair == NULL)
+			v->flag |= PVERT_PIN;
+		else
+			v->flag &= ~PVERT_PIN;
+}
+
+static float p_face_stretch(PFace *f)
+{
+	float T, w, tmp[3];
+	float Ps[3], Pt[3];
+	float a, c, area;
+	PEdge *e1 = f->edge, *e2 = e1->next, *e3 = e2->next;
+	PVert *v1 = e1->vert, *v2 = e2->vert, *v3 = e3->vert;
+
+	area = p_face_uv_area_signed(f);
+
+	if (area <= 0.0f) /* flipped face -> infinite stretch */
+		return 1e10f;
+	
+	if (f->flag & PFACE_FILLED)
+		return 0.0f;
+
+	w= 1.0f/(2.0f*area);
+
+	/* compute derivatives */
+	VecCopyf(Ps, v1->co);
+	VecMulf(Ps, (v2->uv[1] - v3->uv[1]));
+
+	VecCopyf(tmp, v2->co);
+	VecMulf(tmp, (v3->uv[1] - v1->uv[1]));
+	VecAddf(Ps, Ps, tmp);
+
+	VecCopyf(tmp, v3->co);
+	VecMulf(tmp, (v1->uv[1] - v2->uv[1]));
+	VecAddf(Ps, Ps, tmp);
+
+	VecMulf(Ps, w);
+
+	VecCopyf(Pt, v1->co);
+	VecMulf(Pt, (v3->uv[0] - v2->uv[0]));
+
+	VecCopyf(tmp, v2->co);
+	VecMulf(tmp, (v1->uv[0] - v3->uv[0]));
+	VecAddf(Pt, Pt, tmp);
+
+	VecCopyf(tmp, v3->co);
+	VecMulf(tmp, (v2->uv[0] - v1->uv[0]));
+	VecAddf(Pt, Pt, tmp);
+
+	VecMulf(Pt, w);
+
+	/* Sander Tensor */
+	a= Inpf(Ps, Ps);
+	c= Inpf(Pt, Pt);
+
+	T = sqrt(0.5f*(a + c)*f->u.area3d);
+
+	return T;
+}
+
+static float p_stretch_compute_vertex(PVert *v)
+{
+	PEdge *e = v->edge;
+	float sum = 0.0f;
+
+	do {
+		sum += p_face_stretch(e->face);
+		e = p_wheel_edge_next(e);
+	} while (e && e != (v->edge));
+
+	return sum;
+}
+
+static void p_chart_stretch_minimize(PChart *chart, RNG *rng)
+{
+	PVert *v;
+	PEdge *e;
+	int j, nedges;
+	float orig_stretch, low, stretch_low, high, stretch_high, mid, stretch;
+	float orig_uv[2], dir[2], random_angle, trusted_radius;
+
+	for(v=(PVert*)chart->verts->first; v; v=v->link.next) {
+		if((v->flag & PVERT_PIN) || !(v->flag & PVERT_SELECT))
+			continue;
+
+		orig_stretch = p_stretch_compute_vertex(v);
+		orig_uv[0] = v->uv[0];
+		orig_uv[1] = v->uv[1];
+
+		/* move vertex in a random direction */
+		trusted_radius = 0.0f;
+		nedges = 0;
+		e = v->edge;
+
+		do {
+			trusted_radius += p_edge_uv_length(e);
+			nedges++;
+
+			e = p_wheel_edge_next(e);
+		} while (e && e != (v->edge));
+
+		trusted_radius /= 2 * nedges;
+
+		random_angle = rng_getFloat(rng) * 2.0 * M_PI;
+		dir[0] = trusted_radius * cos(random_angle);
+		dir[1] = trusted_radius * sin(random_angle);
+
+		/* calculate old and new stretch */
+		low = 0;
+		stretch_low = orig_stretch;
+
+		Vec2Addf(v->uv, orig_uv, dir);
+		high = 1;
+		stretch = stretch_high = p_stretch_compute_vertex(v);
+
+		/* binary search for lowest stretch position */
+		for (j = 0; j < P_STRETCH_ITER; j++) {
+			mid = 0.5 * (low + high);
+			v->uv[0]= orig_uv[0] + mid*dir[0];
+			v->uv[1]= orig_uv[1] + mid*dir[1];
+			stretch = p_stretch_compute_vertex(v);
+
+			if (stretch_low < stretch_high) {
+				high = mid;
+				stretch_high = stretch;
+			}
+			else {
+				low = mid;
+				stretch_low = stretch;
+			}
+		}
+
+		/* no luck, stretch has increased, reset to old values */
+		if(stretch >= orig_stretch)
+			Vec2Copyf(v->uv, orig_uv);
+	}
+}
+
+void param_stretch_begin(ParamHandle *handle)
+{
+	PHandle *phandle = (PHandle*)handle;
+	PChart *chart;
+	PVert *v;
+	PFace *f;
+	int i;
+
+	param_assert(phandle->state == PHANDLE_STATE_CONSTRUCTED);
+	phandle->state = PHANDLE_STATE_STRETCH;
+
+	phandle->rng = rng_new(31415926);
+
+	for (i = 0; i < phandle->ncharts; i++) {
+		chart = phandle->charts[i];
+
+		for (v=(PVert*)chart->verts->first; v; v=v->link.next)
+			p_vert_load_select_uvs(v);
+
+		p_stretch_pin_boundary(chart);
+
+		for (f=(PFace*)chart->faces->first; f; f=f->link.next) {
+			p_face_backup_uvs(f);
+			f->u.area3d = p_face_area(f);
+		}
+	}
+}
+
+void param_stretch_iter(ParamHandle *handle)
+{
+	PHandle *phandle = (PHandle*)handle;
+	PChart *chart;
+	int i;
+
+	param_assert(phandle->state == PHANDLE_STATE_STRETCH);
+
+	for (i = 0; i < phandle->ncharts; i++) {
+		chart = phandle->charts[i];
+
+		p_chart_stretch_minimize(chart, phandle->rng);
+		p_flush_uvs(chart);
+	}
+}
+
+void param_stretch_end(ParamHandle *handle, ParamBool restore)
+{
+	PHandle *phandle = (PHandle*)handle;
+	PChart *chart;
+	PFace *f;
+	int i;
+
+	param_assert(phandle->state == PHANDLE_STATE_STRETCH);
+	phandle->state = PHANDLE_STATE_CONSTRUCTED;
+
+	rng_free(phandle->rng);
+	phandle->rng = NULL;
+
+	if (restore) {
+		for (i = 0; i < phandle->ncharts; i++) {
+			chart = phandle->charts[i];
+
+			for (f=(PFace*)chart->faces->first; f; f=f->link.next)
+				p_face_restore_uvs(f);
+		}
+	}
+}
+
+/* Packing */
+
+static PBool p_pack_try(PHandle *handle, float side)
+{
+	PChart *chart;
+	float packx, packy, rowh, groupw, w, h;
+	int i;
+
+	packx= packy= 0.0;
+	rowh= 0.0;
+	groupw= 1.0/sqrt(handle->ncharts);
+
+	for (i = 0; i < handle->ncharts; i++) {
+		chart = handle->charts[i];
+		w = chart->u.pack.size[0];
+		h = chart->u.pack.size[1];
+
+		if(w <= (1.0-packx)) {
+			chart->u.pack.trans[0] = packx;
+			chart->u.pack.trans[1] = packy;
+
+			packx += w;
+			rowh= MAX2(rowh, h);
+		}
+		else {
+			packy += rowh;
+			packx = w;
+			rowh = h;
+
+			chart->u.pack.trans[0] = 0.0;
+			chart->u.pack.trans[1] = packy;
+		}
+
+		if (rowh > side)
+			return P_FALSE;
+	}
+
+	return P_TRUE;
+}
+
+#define PACK_SEARCH_DEPTH 7
+
+void param_pack(ParamHandle *handle)
+{
+	PHandle *phandle = (PHandle*)handle;
+	PChart *chart;
+	float uv_area, area, trans[2], minside, maxside, totarea, side;
+	int i;
+
+	/* very simple rectangle packing */
+
+	if (phandle->ncharts == 0)
+		return;
+
+	totarea = 0.0f;
+	maxside = 0.0f;
+
+	for (i = 0; i < phandle->ncharts; i++) {
+		chart = phandle->charts[i];
+
+		p_chart_area(chart, &uv_area, &area);
+		chart->u.pack.rescale = (uv_area > 0.0f)? area/uv_area: 0.0f;
+		totarea += uv_area*chart->u.pack.rescale;
+
+		p_chart_uv_bbox(chart, trans, chart->u.pack.size);
+		trans[0] = -trans[0];
+		trans[1] = -trans[1];
+		p_chart_uv_translate(chart, trans);
+		p_chart_uv_scale(chart, chart->u.pack.rescale);
+
+		chart->u.pack.size[0] -= trans[0];
+		chart->u.pack.size[1] -= trans[1];
+		chart->u.pack.size[0] *= chart->u.pack.rescale;
+		chart->u.pack.size[1] *= chart->u.pack.rescale;
+
+		maxside = MAX3(maxside, chart->u.pack.size[0], chart->u.pack.size[1]);
+	}
+
+	return;
+
+	printf("%f\n", maxside);
+
+	/* binary search over pack region size */
+	minside = sqrt(totarea);
+	maxside = (((int)sqrt(phandle->ncharts-1))+1)*maxside;
+
+	for (i = 0; i < PACK_SEARCH_DEPTH; i++) {
+		printf("%f %f\n", minside, maxside);
+		if (p_pack_try(phandle, (minside+maxside)*0.5f))
+			minside = (minside+maxside)*0.5f;
+		else
+			maxside = (minside+maxside)*0.5f;
+	}
+
+	side = maxside + 1e-5; /* prevent floating point errors */
+	if (!p_pack_try(phandle, side))
+		printf("impossible\n");
+
+	for (i = 0; i < phandle->ncharts; i++) {
+		chart = phandle->charts[i];
+
+		p_chart_uv_scale(chart, chart->u.pack.rescale/side);
+		trans[0] = chart->u.pack.trans[0]/side;
+		trans[1] = chart->u.pack.trans[1]/side;
+		p_chart_uv_translate(chart, trans);
+	}
+}
+