started repo again without all the data files gunking the works

1 files changed, 284 insertions, 0 deletions

diff --git a/src/fortune/main.c b/src/fortune/main.c
new file mode 100644
index 0000000..7571945
--- /dev/null
+++ b/src/fortune/main.c
@@ -0,0 +1,284 @@
+#
+#include <stdio.h>
+#include <stdint.h>
+#include "defs.h"
+#include <math.h>
+struct Site *nextone();
+
+int triangulate, sorted, plot, debug;
+double xmin, xmax, ymin, ymax, deltax, deltay;
+int vert_count, edge_count, line_count, dual_count, site_count;
+double *site_list, *vert_list, *line_list;
+unsigned int *edge_list, *dual_list;
+
+/* sort sites on y, then x, coord */
+int scomp(s1, s2) struct Point *s1, *s2;
+{
+	if (s1->y < s2->y)
+		return (-1);
+	if (s1->y > s2->y)
+		return (1);
+	if (s1->x < s2->x)
+		return (-1);
+	if (s1->x > s2->x)
+		return (1);
+	return (0);
+}
+
+/* read all sites, sort, and compute xmin, xmax, ymin, ymax */
+void readsites(double *lattice) {
+
+	sites = (struct Site *)myalloc(nsites * sizeof *sites);
+	for (unsigned int j = 0; j < nsites; j++) {
+		sites[j].coord.x = lattice[2 * j];
+		sites[j].coord.y = lattice[2 * j + 1];
+		sites[j].sitenbr = j;
+		sites[j].refcnt = 0;
+	}
+	qsort(sites, nsites, sizeof(struct Site), scomp);
+	xmin = sites[0].coord.x;
+	xmax = sites[0].coord.x;
+	for (unsigned int i = 1; i < nsites; i += 1) {
+		if (sites[i].coord.x < xmin)
+			xmin = sites[i].coord.x;
+		if (sites[i].coord.x > xmax)
+			xmax = sites[i].coord.x;
+	};
+	ymin = sites[0].coord.y;
+	ymax = sites[nsites - 1].coord.y;
+}
+
+unsigned int delete_duplicates(unsigned int ne, unsigned int *etv) {
+	unsigned int ndup = 0;
+	bool *duplicates = (bool *)calloc(ne, sizeof(bool));
+	for (unsigned int i = 0; i < ne; i++) {
+		unsigned int v1 = etv[2 * i];
+		unsigned int v2 = etv[2 * i + 1];
+		for (unsigned int j = 0; j < i; j++) {
+			unsigned int w1 = etv[2 * j];
+			unsigned int w2 = etv[2 * j + 1];
+			if (w1 == v1 && w2 == v2) {
+				duplicates[j] = true;
+				ndup++;
+				break;
+			}
+		}
+	}
+	unsigned int newsize = (int)ne - (int)ndup;
+	unsigned int count = 0;
+	for (unsigned int i = 0; i < ne; i++) {
+		if (!duplicates[i]) {
+			etv[2 * count] = etv[2 * i];
+			etv[2 * count + 1] = etv[2 * i + 1];
+			count++;
+		}
+	}
+	free(duplicates);
+	return newsize;
+}
+
+intptr_t *run_voronoi(unsigned int num, double *lattice, double xmin, double xmax, double ymin, double ymax)
+{
+	struct Site *(*next)();
+
+	sorted = 0;
+	triangulate = 0;
+	plot = 0;
+	debug = 0;
+
+	alloclist = (char **)malloc(9 * num * sizeof(char *));
+	allocnum = 0;
+
+	freeinit(&sfl, sizeof *sites);
+
+	unsigned int eff_num = 9 * num;
+	double *eff_lattice = (double *)malloc(2 * eff_num * sizeof(double));
+
+	for (unsigned int i = 0; i < num; i++) {
+		// original sites - our baby boys
+		eff_lattice[2*i] = lattice[2*i];
+		eff_lattice[2*i+1] = lattice[2*i+1];
+
+		// sites to the right
+		eff_lattice[2*(num+i)+1] = lattice[2*i+1] - xmin + xmax;
+		eff_lattice[2*(num+i)] = lattice[2*i];
+
+		// sites to the left
+		eff_lattice[2*(2*num+i)+1] = lattice[2*i+1] - xmax + xmin;
+		eff_lattice[2*(2*num+i)] = lattice[2*i];
+
+		// sites to the top
+		eff_lattice[2*(3*num+i)+1] = lattice[2*i+1];
+		eff_lattice[2*(3*num+i)] = lattice[2*i] - ymin + ymax;
+
+		// sites to the bottom
+		eff_lattice[2*(4*num+i)+1] = lattice[2*i+1];
+		eff_lattice[2*(4*num+i)] = lattice[2*i] - ymax + ymin;
+
+		// sites to the upper right
+		eff_lattice[2*(5*num+i)+1] = lattice[2*i+1] - xmin + xmax;
+		eff_lattice[2*(5*num+i)] = lattice[2*i] - ymin + ymax;
+
+		// sites to the upper left
+		eff_lattice[2*(6*num+i)+1] = lattice[2*i+1] - xmax + xmin;
+		eff_lattice[2*(6*num+i)] = lattice[2*i] - ymin + ymax;
+
+		// sites to the lower left
+		eff_lattice[2*(7*num+i)+1] = lattice[2*i+1] - xmax + xmin;
+		eff_lattice[2*(7*num+i)] = lattice[2*i] + ymin - ymax;
+
+		// sites to the lower right
+		eff_lattice[2*(8*num+i)+1] = lattice[2*i+1] + xmax - xmin;
+		eff_lattice[2*(8*num+i)] = lattice[2*i] + ymin - ymax;
+	}
+
+	nsites = eff_num;
+	readsites(eff_lattice);
+	free(eff_lattice);
+	next = nextone;
+
+	siteidx = 0;
+	geominit();
+
+	site_list = malloc(2 * nsites * sizeof(double));
+	vert_list = NULL;
+	line_list = NULL;
+	edge_list = NULL;
+	dual_list = NULL;
+	site_count = 0;
+	vert_count = 0;
+	edge_count = 0;
+	line_count = 0;
+	dual_count = 0;
+
+	voronoi(triangulate, next);
+
+	unsigned int real_vert_count = vert_count;
+	unsigned int real_edge_count = edge_count;
+	unsigned int real_dual_count = dual_count;
+	double *real_vert_list = vert_list;
+	unsigned int *real_edge_list = edge_list;
+	unsigned int *real_dual_list = dual_list;
+
+	{
+		real_vert_count = 0;
+		real_edge_count = 0;
+		real_dual_count = 0;
+		real_vert_list = (double *)malloc(2 * vert_count * sizeof(double));
+		real_edge_list = (unsigned int *)malloc(2 * edge_count * sizeof(unsigned int));
+		real_dual_list = (unsigned int *)malloc(3 * dual_count * sizeof(unsigned int));
+		unsigned int *triangle_analyzed = (unsigned int *)malloc(4 * dual_count * sizeof(unsigned int));
+		unsigned int q = 0;
+		int *new_vert_ind = (int *)malloc(vert_count * sizeof(int));
+
+		for (unsigned int i = 0; i < vert_count; i++) {
+			unsigned int t1, t2, t3;
+			t1 = dual_list[3*i]; t2 = dual_list[3*i+1]; t3 = dual_list[3*i+2];
+			if (t1 >= num && t2 >= num && t3 >= num) {
+				new_vert_ind[i] = -1;
+			}
+			else if (t1 < num && t2 < num && t3 < num) {
+				real_vert_list[2*real_vert_count] = vert_list[2*i];
+				real_vert_list[2*real_vert_count+1] = vert_list[2*i+1];
+				real_dual_list[3*real_vert_count] = t1;
+				real_dual_list[3*real_vert_count+1] = t2;
+				real_dual_list[3*real_vert_count+2] = t3;
+				new_vert_ind[i] = real_vert_count;
+				real_vert_count++;
+			}
+			else {
+				unsigned int tt1, tt2, tt3;
+				tt1 = t1 % num; tt2 = t2 % num; tt3 = t3 % num;
+				unsigned int s1, s2, s3;
+				s1 = tt1 < tt2 ? (tt1 < tt3 ? tt1 : tt3) : (tt2 < tt3 ? tt2 : tt3);
+				s3 = tt1 > tt2 ? (tt1 > tt3 ? tt1 : tt3) : (tt2 > tt3 ? tt2 : tt3);
+				s2 = tt1 < tt2 ? (tt1 > tt3 ? tt1 : (tt2 < tt3 ? tt2 : tt3)) : (tt1 < tt3 ? tt1 : (tt2 < tt3 ? tt3: tt2));
+
+				bool matched = false;
+				unsigned int ass_vert;
+
+				for (unsigned int j = 0; j < q; j++) {
+					if (s1 == triangle_analyzed[4*j] && s2 == triangle_analyzed[4*j+1] && s3 == triangle_analyzed[4*j+2]) {
+						matched = true;
+						ass_vert = triangle_analyzed[4*j+3];
+						break;
+					}
+				}
+
+				if (matched) {
+					new_vert_ind[i] = ass_vert;
+				} else {
+					triangle_analyzed[4*q] = s1;
+					triangle_analyzed[4*q+1] = s2;
+					triangle_analyzed[4*q+2] = s3;
+					triangle_analyzed[4*q+3] = real_vert_count;
+					q++;
+					real_vert_list[2*real_vert_count+1] = fmod(2*xmax+vert_list[2*i+1], xmax);
+					real_vert_list[2*real_vert_count] = fmod(2*ymax+vert_list[2*i], ymax);
+					real_dual_list[3*real_vert_count] = t1 % num;
+					real_dual_list[3*real_vert_count+1] = t2 % num;
+					real_dual_list[3*real_vert_count+2] = t3 % num;
+					new_vert_ind[i] = real_vert_count;
+					real_vert_count++;
+				}
+			}
+		}
+		for (unsigned int i = 0; i < edge_count; i++) {
+			unsigned int v1, v2;
+			v1 = edge_list[2 * i];
+			v2 = edge_list[2 * i+1];
+			if (v1 != UINT_MAX && v2 != UINT_MAX) {
+				if (new_vert_ind[v1] >= 0 && new_vert_ind[v2] >=0) {
+					unsigned int nv1 = new_vert_ind[v1];
+					unsigned int nv2 = new_vert_ind[v2];
+					real_edge_list[2*real_edge_count] = nv1 < nv2 ? nv1 : nv2;
+					real_edge_list[2*real_edge_count+1] = nv1 < nv2 ? nv2 : nv1;
+					real_edge_count++;
+				}
+			}
+		}
+		unsigned int new_edge_count = delete_duplicates(real_edge_count, real_edge_list);
+		real_edge_count = new_edge_count;
+		free(triangle_analyzed);
+		free(new_vert_ind);
+		free(dual_list);
+		free(edge_list);
+		free(vert_list);
+	}
+
+
+	intptr_t *output = (intptr_t *)malloc(6 * sizeof(intptr_t));
+	output[0] = (intptr_t)malloc(4 * sizeof(unsigned int));
+	((unsigned int *)output[0])[0] = real_vert_count;
+	((unsigned int *)output[0])[1] = real_edge_count;
+	((unsigned int *)output[0])[2] = line_count;
+	((unsigned int *)output[0])[3] = real_vert_count;
+	output[1] = (intptr_t)site_list;
+	output[2] = (intptr_t)real_vert_list;
+	output[3] = (intptr_t)real_edge_list;
+	output[4] = (intptr_t)line_list;
+	output[5] = (intptr_t)real_dual_list;
+
+	free(ELhash);
+	free(PQhash);
+	free(sites);
+	for (unsigned int i = 0; i < allocnum; i++) {
+		free(alloclist[i]);
+	}
+	free(alloclist);
+
+	return output;
+}
+
+/* return a single in-storage site */
+struct Site *nextone() {
+	for (; siteidx < nsites; siteidx += 1) {
+		if (siteidx == 0 || sites[siteidx].coord.x != sites[siteidx - 1].coord.x ||
+				sites[siteidx].coord.y != sites[siteidx - 1].coord.y) {
+			siteidx += 1;
+			return (&sites[siteidx - 1]);
+		};
+	};
+	return ((struct Site *)NULL);
+}
+