From 77440b774d5ab3b1235ace39320b601168e06617 Mon Sep 17 00:00:00 2001 From: Alessandro Ranellucci Date: Sun, 23 Jun 2013 18:18:38 +0200 Subject: [PATCH] Include admesh code for STL repair --- xs/MANIFEST | 7 + xs/MANIFEST.bak | 4 +- xs/src/admesh/connect.c | 1123 +++++++++++++++++++++++++++++++++++++++ xs/src/admesh/normals.c | 391 ++++++++++++++ xs/src/admesh/shared.c | 243 +++++++++ xs/src/admesh/stl.h | 173 ++++++ xs/src/admesh/stl_io.c | 467 ++++++++++++++++ xs/src/admesh/stlinit.c | 360 +++++++++++++ xs/src/admesh/util.c | 372 +++++++++++++ xs/src/myinit.h | 2 + xs/xsp/TriangleMesh.xsp | 23 + 11 files changed, 3164 insertions(+), 1 deletion(-) create mode 100644 xs/src/admesh/connect.c create mode 100644 xs/src/admesh/normals.c create mode 100644 xs/src/admesh/shared.c create mode 100644 xs/src/admesh/stl.h create mode 100644 xs/src/admesh/stl_io.c create mode 100644 xs/src/admesh/stlinit.c create mode 100644 xs/src/admesh/util.c diff --git a/xs/MANIFEST b/xs/MANIFEST index 43e2324d1..22cd010db 100644 --- a/xs/MANIFEST +++ b/xs/MANIFEST @@ -1,6 +1,13 @@ Build.PL lib/Slic3r/XS.pm MANIFEST This list of files +src/admesh/connect.c +src/admesh/normals.c +src/admesh/shared.c +src/admesh/stl.h +src/admesh/stl_io.c +src/admesh/stlinit.c +src/admesh/util.c src/myinit.h src/ppport.h t/01_trianglemesh.t diff --git a/xs/MANIFEST.bak b/xs/MANIFEST.bak index a473c6c89..43e2324d1 100644 --- a/xs/MANIFEST.bak +++ b/xs/MANIFEST.bak @@ -1,11 +1,13 @@ Build.PL -lib/Slic3r/TriangleMesh/XS.pm lib/Slic3r/XS.pm MANIFEST This list of files src/myinit.h src/ppport.h t/01_trianglemesh.t +t/02_object.t +xsp/my.map xsp/mytype.map +xsp/Object.xsp xsp/TriangleMesh.xsp xsp/typemap.xspt xsp/XS.xsp diff --git a/xs/src/admesh/connect.c b/xs/src/admesh/connect.c new file mode 100644 index 000000000..1a455a628 --- /dev/null +++ b/xs/src/admesh/connect.c @@ -0,0 +1,1123 @@ +/* ADMesh -- process triangulated solid meshes + * Copyright (C) 1995 Anthony D. Martin + * + * 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 + * the Free Software Foundation; either version 2, or (at your option) + * any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program; if not, write to the Free Software + * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. + * + * Questions, comments, suggestions, etc to + */ + +#include +#include +#include +#include + +#include "stl.h" + + +static void stl_match_neighbors_exact(stl_file *stl, + stl_hash_edge *edge_a, stl_hash_edge *edge_b); +static void stl_match_neighbors_nearby(stl_file *stl, + stl_hash_edge *edge_a, stl_hash_edge *edge_b); +static void stl_record_neighbors(stl_file *stl, + stl_hash_edge *edge_a, stl_hash_edge *edge_b); +static void stl_initialize_facet_check_exact(stl_file *stl); +static void stl_initialize_facet_check_nearby(stl_file *stl); +static void stl_load_edge_exact(stl_file *stl, stl_hash_edge *edge, + stl_vertex *a, stl_vertex *b); +static int stl_load_edge_nearby(stl_file *stl, stl_hash_edge *edge, + stl_vertex *a, stl_vertex *b, float tolerance); +static void insert_hash_edge(stl_file *stl, stl_hash_edge edge, + void (*match_neighbors)(stl_file *stl, + stl_hash_edge *edge_a, stl_hash_edge *edge_b)); +static int stl_get_hash_for_edge(int M, stl_hash_edge *edge); +static int stl_compare_function(stl_hash_edge *edge_a, stl_hash_edge *edge_b); +static void stl_free_edges(stl_file *stl); +static void stl_remove_facet(stl_file *stl, int facet_number); +static void stl_change_vertices(stl_file *stl, int facet_num, int vnot, + stl_vertex new_vertex); +static void stl_which_vertices_to_change(stl_file *stl, stl_hash_edge *edge_a, + stl_hash_edge *edge_b, int *facet1, int *vertex1, + int *facet2, int *vertex2, + stl_vertex *new_vertex1, stl_vertex *new_vertex2); +static void stl_remove_degenerate(stl_file *stl, int facet); +static void stl_add_facet(stl_file *stl, stl_facet *new_facet); +extern int stl_check_normal_vector(stl_file *stl, + int facet_num, int normal_fix_flag); +static void stl_update_connects_remove_1(stl_file *stl, int facet_num); + + +void +stl_check_facets_exact(stl_file *stl) +{ +/* This function builds the neighbors list. No modifications are made + * to any of the facets. The edges are said to match only if all six + * floats of the first edge matches all six floats of the second edge. + */ + + stl_hash_edge edge; + stl_facet facet; + int i; + int j; + + stl->stats.connected_edges = 0; + stl->stats.connected_facets_1_edge = 0; + stl->stats.connected_facets_2_edge = 0; + stl->stats.connected_facets_3_edge = 0; + + stl_initialize_facet_check_exact(stl); + + for(i = 0; i < stl->stats.number_of_facets; i++) + { + facet = stl->facet_start[i]; + + if( !memcmp(&facet.vertex[0], &facet.vertex[1], + sizeof(stl_vertex)) + || !memcmp(&facet.vertex[1], &facet.vertex[2], + sizeof(stl_vertex)) + || !memcmp(&facet.vertex[0], &facet.vertex[2], + sizeof(stl_vertex))) + { + stl->stats.degenerate_facets += 1; + stl_remove_facet(stl, i); + i--; + continue; + + } + for(j = 0; j < 3; j++) + { + edge.facet_number = i; + edge.which_edge = j; + stl_load_edge_exact(stl, &edge, &facet.vertex[j], + &facet.vertex[(j + 1) % 3]); + + insert_hash_edge(stl, edge, stl_match_neighbors_exact); + } + } + stl_free_edges(stl); +} + +static void +stl_load_edge_exact(stl_file *stl, stl_hash_edge *edge, + stl_vertex *a, stl_vertex *b) +{ + + float diff_x; + float diff_y; + float diff_z; + float max_diff; + + + diff_x = ABS(a->x - b->x); + diff_y = ABS(a->y - b->y); + diff_z = ABS(a->z - b->z); + max_diff = STL_MAX(diff_x, diff_y); + max_diff = STL_MAX(diff_z, max_diff); + stl->stats.shortest_edge = STL_MIN(max_diff, stl->stats.shortest_edge); + + if(diff_x == max_diff) + { + if(a->x > b->x) + { + memcpy(&edge->key[0], a, sizeof(stl_vertex)); + memcpy(&edge->key[3], b, sizeof(stl_vertex)); + } + else + { + memcpy(&edge->key[0], b, sizeof(stl_vertex)); + memcpy(&edge->key[3], a, sizeof(stl_vertex)); + edge->which_edge += 3; /* this edge is loaded backwards */ + } + } + else if(diff_y == max_diff) + { + if(a->y > b->y) + { + memcpy(&edge->key[0], a, sizeof(stl_vertex)); + memcpy(&edge->key[3], b, sizeof(stl_vertex)); + } + else + { + memcpy(&edge->key[0], b, sizeof(stl_vertex)); + memcpy(&edge->key[3], a, sizeof(stl_vertex)); + edge->which_edge += 3; /* this edge is loaded backwards */ + } + } + else + { + if(a->z > b->z) + { + memcpy(&edge->key[0], a, sizeof(stl_vertex)); + memcpy(&edge->key[3], b, sizeof(stl_vertex)); + } + else + { + memcpy(&edge->key[0], b, sizeof(stl_vertex)); + memcpy(&edge->key[3], a, sizeof(stl_vertex)); + edge->which_edge += 3; /* this edge is loaded backwards */ + } + } +} + +static void +stl_initialize_facet_check_exact(stl_file *stl) +{ + int i; + + stl->stats.malloced = 0; + stl->stats.freed = 0; + stl->stats.collisions = 0; + + + stl->M = 81397; + + for(i = 0; i < stl->stats.number_of_facets ; i++) + { + /* initialize neighbors list to -1 to mark unconnected edges */ + stl->neighbors_start[i].neighbor[0] = -1; + stl->neighbors_start[i].neighbor[1] = -1; + stl->neighbors_start[i].neighbor[2] = -1; + } + + stl->heads = (stl_hash_edge**)calloc(stl->M, sizeof(*stl->heads)); + if(stl->heads == NULL) perror("stl_initialize_facet_check_exact"); + + stl->tail = (stl_hash_edge*)malloc(sizeof(stl_hash_edge)); + if(stl->tail == NULL) perror("stl_initialize_facet_check_exact"); + + stl->tail->next = stl->tail; + + for(i = 0; i < stl->M; i++) + { + stl->heads[i] = stl->tail; + } +} + +static void +insert_hash_edge(stl_file *stl, stl_hash_edge edge, + void (*match_neighbors)(stl_file *stl, + stl_hash_edge *edge_a, stl_hash_edge *edge_b)) +{ + stl_hash_edge *link; + stl_hash_edge *new_edge; + stl_hash_edge *temp; + int chain_number; + + chain_number = stl_get_hash_for_edge(stl->M, &edge); + + link = stl->heads[chain_number]; + + if(link == stl->tail) + { + /* This list doesn't have any edges currently in it. Add this one. */ + new_edge = (stl_hash_edge*)malloc(sizeof(stl_hash_edge)); + if(new_edge == NULL) perror("insert_hash_edge"); + stl->stats.malloced++; + *new_edge = edge; + new_edge->next = stl->tail; + stl->heads[chain_number] = new_edge; + return; + } + else if(!stl_compare_function(&edge, link)) + { + /* This is a match. Record result in neighbors list. */ + match_neighbors(stl, &edge, link); + /* Delete the matched edge from the list. */ + stl->heads[chain_number] = link->next; + free(link); + stl->stats.freed++; + return; + } + else + { /* Continue through the rest of the list */ + for(;;) + { + if(link->next == stl->tail) + { + /* This is the last item in the list. Insert a new edge. */ + new_edge = (stl_hash_edge*)malloc(sizeof(stl_hash_edge)); + if(new_edge == NULL) perror("insert_hash_edge"); + stl->stats.malloced++; + *new_edge = edge; + new_edge->next = stl->tail; + link->next = new_edge; + stl->stats.collisions++; + return; + } + else if(!stl_compare_function(&edge, link->next)) + { + /* This is a match. Record result in neighbors list. */ + match_neighbors(stl, &edge, link->next); + + /* Delete the matched edge from the list. */ + temp = link->next; + link->next = link->next->next; + free(temp); + stl->stats.freed++; + return; + } + else + { + /* This is not a match. Go to the next link */ + link = link->next; + stl->stats.collisions++; + } + } + } +} + + +static int +stl_get_hash_for_edge(int M, stl_hash_edge *edge) +{ + return ((edge->key[0] / 23 + edge->key[1] / 19 + edge->key[2] / 17 + + edge->key[3] /13 + edge->key[4] / 11 + edge->key[5] / 7 ) % M); +} + +static int +stl_compare_function(stl_hash_edge *edge_a, stl_hash_edge *edge_b) +{ + if(edge_a->facet_number == edge_b->facet_number) + { + return 1; /* Don't match edges of the same facet */ + } + else + { + return memcmp(edge_a, edge_b, SIZEOF_EDGE_SORT); + } +} + +void +stl_check_facets_nearby(stl_file *stl, float tolerance) +{ + stl_hash_edge edge[3]; + stl_facet facet; + int i; + int j; + + + if( (stl->stats.connected_facets_1_edge == stl->stats.number_of_facets) + && (stl->stats.connected_facets_2_edge == stl->stats.number_of_facets) + && (stl->stats.connected_facets_3_edge == stl->stats.number_of_facets)) + { + /* No need to check any further. All facets are connected */ + return; + } + + stl_initialize_facet_check_nearby(stl); + + for(i = 0; i < stl->stats.number_of_facets; i++) + { + facet = stl->facet_start[i]; + for(j = 0; j < 3; j++) + { + if(stl->neighbors_start[i].neighbor[j] == -1) + { + edge[j].facet_number = i; + edge[j].which_edge = j; + if(stl_load_edge_nearby(stl, &edge[j], &facet.vertex[j], + &facet.vertex[(j + 1) % 3], + tolerance)) + { + /* only insert edges that have different keys */ + insert_hash_edge(stl, edge[j], stl_match_neighbors_nearby); + } + } + } + } + + stl_free_edges(stl); +} + +static int +stl_load_edge_nearby(stl_file *stl, stl_hash_edge *edge, + stl_vertex *a, stl_vertex *b, float tolerance) +{ + float diff_x; + float diff_y; + float diff_z; + float max_diff; + unsigned vertex1[3]; + unsigned vertex2[3]; + + + diff_x = ABS(a->x - b->x); + diff_y = ABS(a->y - b->y); + diff_z = ABS(a->z - b->z); + max_diff = STL_MAX(diff_x, diff_y); + max_diff = STL_MAX(diff_z, max_diff); + + vertex1[0] = (unsigned)((a->x - stl->stats.min.x) / tolerance); + vertex1[1] = (unsigned)((a->y - stl->stats.min.y) / tolerance); + vertex1[2] = (unsigned)((a->z - stl->stats.min.z) / tolerance); + vertex2[0] = (unsigned)((b->x - stl->stats.min.x) / tolerance); + vertex2[1] = (unsigned)((b->y - stl->stats.min.y) / tolerance); + vertex2[2] = (unsigned)((b->z - stl->stats.min.z) / tolerance); + + if( (vertex1[0] == vertex2[0]) + && (vertex1[1] == vertex2[1]) + && (vertex1[2] == vertex2[2])) + { + /* Both vertices hash to the same value */ + return 0; + } + + if(diff_x == max_diff) + { + if(a->x > b->x) + { + memcpy(&edge->key[0], vertex1, sizeof(stl_vertex)); + memcpy(&edge->key[3], vertex2, sizeof(stl_vertex)); + } + else + { + memcpy(&edge->key[0], vertex2, sizeof(stl_vertex)); + memcpy(&edge->key[3], vertex1, sizeof(stl_vertex)); + edge->which_edge += 3; /* this edge is loaded backwards */ + } + } + else if(diff_y == max_diff) + { + if(a->y > b->y) + { + memcpy(&edge->key[0], vertex1, sizeof(stl_vertex)); + memcpy(&edge->key[3], vertex2, sizeof(stl_vertex)); + } + else + { + memcpy(&edge->key[0], vertex2, sizeof(stl_vertex)); + memcpy(&edge->key[3], vertex1, sizeof(stl_vertex)); + edge->which_edge += 3; /* this edge is loaded backwards */ + } + } + else + { + if(a->z > b->z) + { + memcpy(&edge->key[0], vertex1, sizeof(stl_vertex)); + memcpy(&edge->key[3], vertex2, sizeof(stl_vertex)); + } + else + { + memcpy(&edge->key[0], vertex2, sizeof(stl_vertex)); + memcpy(&edge->key[3], vertex1, sizeof(stl_vertex)); + edge->which_edge += 3; /* this edge is loaded backwards */ + } + } + return 1; +} + +static void +stl_free_edges(stl_file *stl) +{ + int i; + stl_hash_edge *temp; + + if(stl->stats.malloced != stl->stats.freed) + { + for(i = 0; i < stl->M; i++) + { + for(temp = stl->heads[i]; stl->heads[i] != stl->tail; + temp = stl->heads[i]) + { + stl->heads[i] = stl->heads[i]->next; + free(temp); + stl->stats.freed++; + } + } + } + free(stl->heads); + free(stl->tail); +} + +static void +stl_initialize_facet_check_nearby(stl_file *stl) +{ + int i; + + stl->stats.malloced = 0; + stl->stats.freed = 0; + stl->stats.collisions = 0; + + /* tolerance = STL_MAX(stl->stats.shortest_edge, tolerance);*/ + /* tolerance = STL_MAX((stl->stats.bounding_diameter / 500000.0), tolerance);*/ + /* tolerance *= 0.5;*/ + + stl->M = 81397; + + stl->heads = (stl_hash_edge**)calloc(stl->M, sizeof(*stl->heads)); + if(stl->heads == NULL) perror("stl_initialize_facet_check_nearby"); + + stl->tail = (stl_hash_edge*)malloc(sizeof(stl_hash_edge)); + if(stl->tail == NULL) perror("stl_initialize_facet_check_nearby"); + + stl->tail->next = stl->tail; + + for(i = 0; i < stl->M; i++) + { + stl->heads[i] = stl->tail; + } +} + + + +static void +stl_record_neighbors(stl_file *stl, + stl_hash_edge *edge_a, stl_hash_edge *edge_b) +{ + int i; + int j; + + /* Facet a's neighbor is facet b */ + stl->neighbors_start[edge_a->facet_number].neighbor[edge_a->which_edge % 3] = + edge_b->facet_number; /* sets the .neighbor part */ + + stl->neighbors_start[edge_a->facet_number]. + which_vertex_not[edge_a->which_edge % 3] = + (edge_b->which_edge + 2) % 3; /* sets the .which_vertex_not part */ + + /* Facet b's neighbor is facet a */ + stl->neighbors_start[edge_b->facet_number].neighbor[edge_b->which_edge % 3] = + edge_a->facet_number; /* sets the .neighbor part */ + + stl->neighbors_start[edge_b->facet_number]. + which_vertex_not[edge_b->which_edge % 3] = + (edge_a->which_edge + 2) % 3; /* sets the .which_vertex_not part */ + + if( ((edge_a->which_edge < 3) && (edge_b->which_edge < 3)) + || ((edge_a->which_edge > 2) && (edge_b->which_edge > 2))) + { + /* these facets are oriented in opposite directions. */ + /* their normals are probably messed up. */ + stl->neighbors_start[edge_a->facet_number]. + which_vertex_not[edge_a->which_edge % 3] += 3; + stl->neighbors_start[edge_b->facet_number]. + which_vertex_not[edge_b->which_edge % 3] += 3; + } + + + /* Count successful connects */ + /* Total connects */ + stl->stats.connected_edges += 2; + /* Count individual connects */ + i = ((stl->neighbors_start[edge_a->facet_number].neighbor[0] == -1) + + (stl->neighbors_start[edge_a->facet_number].neighbor[1] == -1) + + (stl->neighbors_start[edge_a->facet_number].neighbor[2] == -1)); + j = ((stl->neighbors_start[edge_b->facet_number].neighbor[0] == -1) + + (stl->neighbors_start[edge_b->facet_number].neighbor[1] == -1) + + (stl->neighbors_start[edge_b->facet_number].neighbor[2] == -1)); + if(i == 2) + { + stl->stats.connected_facets_1_edge +=1; + } + else if(i == 1) + { + stl->stats.connected_facets_2_edge +=1; + } + else + { + stl->stats.connected_facets_3_edge +=1; + } + if(j == 2) + { + stl->stats.connected_facets_1_edge +=1; + } + else if(j == 1) + { + stl->stats.connected_facets_2_edge +=1; + } + else + { + stl->stats.connected_facets_3_edge +=1; + } +} + +static void +stl_match_neighbors_exact(stl_file *stl, + stl_hash_edge *edge_a, stl_hash_edge *edge_b) +{ + stl_record_neighbors(stl, edge_a, edge_b); +} + +static void +stl_match_neighbors_nearby(stl_file *stl, + stl_hash_edge *edge_a, stl_hash_edge *edge_b) +{ + int facet1; + int facet2; + int vertex1; + int vertex2; + int vnot1; + int vnot2; + stl_vertex new_vertex1; + stl_vertex new_vertex2; + + stl_record_neighbors(stl, edge_a, edge_b); + stl_which_vertices_to_change(stl, edge_a, edge_b, &facet1, &vertex1, + &facet2, &vertex2, &new_vertex1, &new_vertex2); + if(facet1 != -1) + { + if(facet1 == edge_a->facet_number) + { + vnot1 = (edge_a->which_edge + 2) % 3; + } + else + { + vnot1 = (edge_b->which_edge + 2) % 3; + } + if(((vnot1 + 2) % 3) == vertex1) + { + vnot1 += 3; + } + stl_change_vertices(stl, facet1, vnot1, new_vertex1); + } + if(facet2 != -1) + { + if(facet2 == edge_a->facet_number) + { + vnot2 = (edge_a->which_edge + 2) % 3; + } + else + { + vnot2 = (edge_b->which_edge + 2) % 3; + } + if(((vnot2 + 2) % 3) == vertex2) + { + vnot2 += 3; + } + stl_change_vertices(stl, facet2, vnot2, new_vertex2); + } + stl->stats.edges_fixed += 2; +} + + +static void +stl_change_vertices(stl_file *stl, int facet_num, int vnot, + stl_vertex new_vertex) +{ + int first_facet; + int direction; + int next_edge; + int pivot_vertex; + + first_facet = facet_num; + direction = 0; + + for(;;) + { + if(vnot > 2) + { + if(direction == 0) + { + pivot_vertex = (vnot + 2) % 3; + next_edge = pivot_vertex; + direction = 1; + } + else + { + pivot_vertex = (vnot + 1) % 3; + next_edge = vnot % 3; + direction = 0; + } + } + else + { + if(direction == 0) + { + pivot_vertex = (vnot + 1) % 3; + next_edge = vnot; + } + else + { + pivot_vertex = (vnot + 2) % 3; + next_edge = pivot_vertex; + } + } + stl->facet_start[facet_num].vertex[pivot_vertex] = new_vertex; + vnot = stl->neighbors_start[facet_num].which_vertex_not[next_edge]; + facet_num = stl->neighbors_start[facet_num].neighbor[next_edge]; + + if(facet_num == -1) + { + break; + } + + if(facet_num == first_facet) + { + /* back to the beginning */ + printf("\ +Back to the first facet changing vertices: probably a mobius part.\n\ +Try using a smaller tolerance or don't do a nearby check\n"); + exit(1); + break; + } + } +} + + +static void +stl_which_vertices_to_change(stl_file *stl, stl_hash_edge *edge_a, + stl_hash_edge *edge_b, int *facet1, int *vertex1, + int *facet2, int *vertex2, + stl_vertex *new_vertex1, stl_vertex *new_vertex2) +{ + int v1a; /* pair 1, facet a */ + int v1b; /* pair 1, facet b */ + int v2a; /* pair 2, facet a */ + int v2b; /* pair 2, facet b */ + + + /* Find first pair */ + if(edge_a->which_edge < 3) + { + v1a = edge_a->which_edge; + v2a = (edge_a->which_edge + 1) % 3; + } + else + { + v2a = edge_a->which_edge % 3; + v1a = (edge_a->which_edge + 1) % 3; + } + if(edge_b->which_edge < 3) + { + v1b = edge_b->which_edge; + v2b = (edge_b->which_edge + 1) % 3; + } + else + { + v2b = edge_b->which_edge % 3; + v1b = (edge_b->which_edge + 1) % 3; + } + + /* Of the first pair, which vertex, if any, should be changed */ + if(!memcmp(&stl->facet_start[edge_a->facet_number].vertex[v1a], + &stl->facet_start[edge_b->facet_number].vertex[v1b], + sizeof(stl_vertex))) + { + /* These facets are already equal. No need to change. */ + *facet1 = -1; + } + else + { + if( (stl->neighbors_start[edge_a->facet_number].neighbor[v1a] == -1) + && (stl->neighbors_start[edge_a->facet_number]. + neighbor[(v1a + 2) % 3] == -1)) + { + /* This vertex has no neighbors. This is a good one to change */ + *facet1 = edge_a->facet_number; + *vertex1 = v1a; + *new_vertex1 = stl->facet_start[edge_b->facet_number].vertex[v1b]; + } + else + { + *facet1 = edge_b->facet_number; + *vertex1 = v1b; + *new_vertex1 = stl->facet_start[edge_a->facet_number].vertex[v1a]; + } + } + + /* Of the second pair, which vertex, if any, should be changed */ + if(!memcmp(&stl->facet_start[edge_a->facet_number].vertex[v2a], + &stl->facet_start[edge_b->facet_number].vertex[v2b], + sizeof(stl_vertex))) + { + /* These facets are already equal. No need to change. */ + *facet2 = -1; + } + else + { + if( (stl->neighbors_start[edge_a->facet_number].neighbor[v2a] == -1) + && (stl->neighbors_start[edge_a->facet_number]. + neighbor[(v2a + 2) % 3] == -1)) + { + /* This vertex has no neighbors. This is a good one to change */ + *facet2 = edge_a->facet_number; + *vertex2 = v2a; + *new_vertex2 = stl->facet_start[edge_b->facet_number].vertex[v2b]; + } + else + { + *facet2 = edge_b->facet_number; + *vertex2 = v2b; + *new_vertex2 = stl->facet_start[edge_a->facet_number].vertex[v2a]; + } + } +} + +static void +stl_remove_facet(stl_file *stl, int facet_number) +{ + int neighbor[3]; + int vnot[3]; + int i; + int j; + + stl->stats.facets_removed += 1; + /* Update list of connected edges */ + j = ((stl->neighbors_start[facet_number].neighbor[0] == -1) + + (stl->neighbors_start[facet_number].neighbor[1] == -1) + + (stl->neighbors_start[facet_number].neighbor[2] == -1)); + if(j == 2) + { + stl->stats.connected_facets_1_edge -= 1; + } + else if(j == 1) + { + stl->stats.connected_facets_2_edge -= 1; + stl->stats.connected_facets_1_edge -= 1; + } + else if(j == 0) + { + stl->stats.connected_facets_3_edge -= 1; + stl->stats.connected_facets_2_edge -= 1; + stl->stats.connected_facets_1_edge -= 1; + } + + stl->facet_start[facet_number] = + stl->facet_start[stl->stats.number_of_facets - 1]; + /* I could reallocate at this point, but it is not really necessary. */ + stl->neighbors_start[facet_number] = + stl->neighbors_start[stl->stats.number_of_facets - 1]; + stl->stats.number_of_facets -= 1; + + for(i = 0; i < 3; i++) + { + neighbor[i] = stl->neighbors_start[facet_number].neighbor[i]; + vnot[i] = stl->neighbors_start[facet_number].which_vertex_not[i]; + } + + for(i = 0; i < 3; i++) + { + if(neighbor[i] != -1) + { + if(stl->neighbors_start[neighbor[i]].neighbor[(vnot[i] + 1)% 3] != + stl->stats.number_of_facets) + { + printf("\ +in stl_remove_facet: neighbor = %d numfacets = %d this is wrong\n", + stl->neighbors_start[neighbor[i]].neighbor[(vnot[i] + 1)% 3], + stl->stats.number_of_facets); + exit(1); + } + stl->neighbors_start[neighbor[i]].neighbor[(vnot[i] + 1)% 3] + = facet_number; + } + } +} + +void +stl_remove_unconnected_facets(stl_file *stl) +{ + /* A couple of things need to be done here. One is to remove any */ + /* completely unconnected facets (0 edges connected) since these are */ + /* useless and could be completely wrong. The second thing that needs to */ + /* be done is to remove any degenerate facets that were created during */ + /* stl_check_facets_nearby(). */ + + int i; + + /* remove degenerate facets */ + for(i = 0; i < stl->stats.number_of_facets; i++) + { + if( !memcmp(&stl->facet_start[i].vertex[0], + &stl->facet_start[i].vertex[1], sizeof(stl_vertex)) + || !memcmp(&stl->facet_start[i].vertex[1], + &stl->facet_start[i].vertex[2], sizeof(stl_vertex)) + || !memcmp(&stl->facet_start[i].vertex[0], + &stl->facet_start[i].vertex[2], sizeof(stl_vertex))) + { + stl_remove_degenerate(stl, i); + i--; + } + } + + if(stl->stats.connected_facets_1_edge < stl->stats.number_of_facets) + { + /* remove completely unconnected facets */ + for(i = 0; i < stl->stats.number_of_facets; i++) + { + if( (stl->neighbors_start[i].neighbor[0] == -1) + && (stl->neighbors_start[i].neighbor[1] == -1) + && (stl->neighbors_start[i].neighbor[2] == -1)) + { + /* This facet is completely unconnected. Remove it. */ + stl_remove_facet(stl, i); + i--; + } + } + } +} + +static void +stl_remove_degenerate(stl_file *stl, int facet) +{ + int edge1; + int edge2; + int edge3; + int neighbor1; + int neighbor2; + int neighbor3; + int vnot1; + int vnot2; + int vnot3; + + if( !memcmp(&stl->facet_start[facet].vertex[0], + &stl->facet_start[facet].vertex[1], sizeof(stl_vertex)) + && !memcmp(&stl->facet_start[facet].vertex[1], + &stl->facet_start[facet].vertex[2], sizeof(stl_vertex))) + { + /* all 3 vertices are equal. Just remove the facet. I don't think*/ + /* this is really possible, but just in case... */ + printf("removing a facet in stl_remove_degenerate\n"); + + stl_remove_facet(stl, facet); + return; + } + + if(!memcmp(&stl->facet_start[facet].vertex[0], + &stl->facet_start[facet].vertex[1], sizeof(stl_vertex))) + { + edge1 = 1; + edge2 = 2; + edge3 = 0; + } + else if(!memcmp(&stl->facet_start[facet].vertex[1], + &stl->facet_start[facet].vertex[2], sizeof(stl_vertex))) + { + edge1 = 0; + edge2 = 2; + edge3 = 1; + } + else if(!memcmp(&stl->facet_start[facet].vertex[2], + &stl->facet_start[facet].vertex[0], sizeof(stl_vertex))) + { + edge1 = 0; + edge2 = 1; + edge3 = 2; + } + else + { + /* No degenerate. Function shouldn't have been called. */ + return; + } + neighbor1 = stl->neighbors_start[facet].neighbor[edge1]; + neighbor2 = stl->neighbors_start[facet].neighbor[edge2]; + + if(neighbor1 == -1) + { + stl_update_connects_remove_1(stl, neighbor2); + } + if(neighbor2 == -1) + { + stl_update_connects_remove_1(stl, neighbor1); + } + + + neighbor3 = stl->neighbors_start[facet].neighbor[edge3]; + vnot1 = stl->neighbors_start[facet].which_vertex_not[edge1]; + vnot2 = stl->neighbors_start[facet].which_vertex_not[edge2]; + vnot3 = stl->neighbors_start[facet].which_vertex_not[edge3]; + + stl->neighbors_start[neighbor1].neighbor[(vnot1 + 1) % 3] = neighbor2; + stl->neighbors_start[neighbor2].neighbor[(vnot2 + 1) % 3] = neighbor1; + stl->neighbors_start[neighbor1].which_vertex_not[(vnot1 + 1) % 3] = vnot2; + stl->neighbors_start[neighbor2].which_vertex_not[(vnot2 + 1) % 3] = vnot1; + + stl_remove_facet(stl, facet); + + if(neighbor3 != -1) + { + stl_update_connects_remove_1(stl, neighbor3); + stl->neighbors_start[neighbor3].neighbor[(vnot3 + 1) % 3] = -1; + } +} + +void +stl_update_connects_remove_1(stl_file *stl, int facet_num) +{ + int j; + + /* Update list of connected edges */ + j = ((stl->neighbors_start[facet_num].neighbor[0] == -1) + + (stl->neighbors_start[facet_num].neighbor[1] == -1) + + (stl->neighbors_start[facet_num].neighbor[2] == -1)); + if(j == 0) /* Facet has 3 neighbors */ + { + stl->stats.connected_facets_3_edge -= 1; + } + else if(j == 1) /* Facet has 2 neighbors */ + { + stl->stats.connected_facets_2_edge -= 1; + } + else if(j == 2) /* Facet has 1 neighbor */ + { + stl->stats.connected_facets_1_edge -= 1; + } +} + +void +stl_fill_holes(stl_file *stl) +{ + stl_facet facet; + stl_facet new_facet; + int neighbors_initial[3]; + stl_hash_edge edge; + int first_facet; + int direction; + int facet_num; + int vnot; + int next_edge; + int pivot_vertex; + int next_facet; + int i; + int j; + int k; + + /* Insert all unconnected edges into hash list */ + stl_initialize_facet_check_nearby(stl); + for(i = 0; i < stl->stats.number_of_facets; i++) + { + facet = stl->facet_start[i]; + for(j = 0; j < 3; j++) + { + if(stl->neighbors_start[i].neighbor[j] != -1) continue; + edge.facet_number = i; + edge.which_edge = j; + stl_load_edge_exact(stl, &edge, &facet.vertex[j], + &facet.vertex[(j + 1) % 3]); + + insert_hash_edge(stl, edge, stl_match_neighbors_exact); + } + } + + for(i = 0; i < stl->stats.number_of_facets; i++) + { + facet = stl->facet_start[i]; + neighbors_initial[0] = stl->neighbors_start[i].neighbor[0]; + neighbors_initial[1] = stl->neighbors_start[i].neighbor[1]; + neighbors_initial[2] = stl->neighbors_start[i].neighbor[2]; + first_facet = i; + for(j = 0; j < 3; j++) + { + if(stl->neighbors_start[i].neighbor[j] != -1) continue; + + new_facet.vertex[0] = facet.vertex[j]; + new_facet.vertex[1] = facet.vertex[(j + 1) % 3]; + if(neighbors_initial[(j + 2) % 3] == -1) + { + direction = 1; + } + else + { + direction = 0; + } + + facet_num = i; + vnot = (j + 2) % 3; + + for(;;) + { + if(vnot > 2) + { + if(direction == 0) + { + pivot_vertex = (vnot + 2) % 3; + next_edge = pivot_vertex; + direction = 1; + } + else + { + pivot_vertex = (vnot + 1) % 3; + next_edge = vnot % 3; + direction = 0; + } + } + else + { + if(direction == 0) + { + pivot_vertex = (vnot + 1) % 3; + next_edge = vnot; + } + else + { + pivot_vertex = (vnot + 2) % 3; + next_edge = pivot_vertex; + } + } + next_facet = stl->neighbors_start[facet_num].neighbor[next_edge]; + + if(next_facet == -1) + { + new_facet.vertex[2] = stl->facet_start[facet_num]. + vertex[vnot % 3]; + stl_add_facet(stl, &new_facet); + for(k = 0; k < 3; k++) + { + edge.facet_number = stl->stats.number_of_facets - 1; + edge.which_edge = k; + stl_load_edge_exact(stl, &edge, &new_facet.vertex[k], + &new_facet.vertex[(k + 1) % 3]); + + insert_hash_edge(stl, edge, stl_match_neighbors_exact); + } + break; + } + else + { + vnot = stl->neighbors_start[facet_num]. + which_vertex_not[next_edge]; + facet_num = next_facet; + } + + if(facet_num == first_facet) + { + /* back to the beginning */ + printf("\ +Back to the first facet filling holes: probably a mobius part.\n\ +Try using a smaller tolerance or don't do a nearby check\n"); + exit(1); + break; + } + } + } + } +} + +static void +stl_add_facet(stl_file *stl, stl_facet *new_facet) +{ + stl->stats.facets_added += 1; + if(stl->stats.facets_malloced < stl->stats.number_of_facets + 1) + { + stl->facet_start = (stl_facet*)realloc(stl->facet_start, + (sizeof(stl_facet) * (stl->stats.facets_malloced + 256))); + if(stl->facet_start == NULL) perror("stl_add_facet"); + stl->neighbors_start = (stl_neighbors*)realloc(stl->neighbors_start, + (sizeof(stl_neighbors) * (stl->stats.facets_malloced + 256))); + if(stl->neighbors_start == NULL) perror("stl_add_facet"); + stl->stats.facets_malloced += 256; + } + stl->facet_start[stl->stats.number_of_facets] = *new_facet; + + /* note that the normal vector is not set here, just initialized to 0 */ + stl->facet_start[stl->stats.number_of_facets].normal.x = 0.0; + stl->facet_start[stl->stats.number_of_facets].normal.y = 0.0; + stl->facet_start[stl->stats.number_of_facets].normal.z = 0.0; + + stl->neighbors_start[stl->stats.number_of_facets].neighbor[0] = -1; + stl->neighbors_start[stl->stats.number_of_facets].neighbor[1] = -1; + stl->neighbors_start[stl->stats.number_of_facets].neighbor[2] = -1; + stl->stats.number_of_facets += 1; +} diff --git a/xs/src/admesh/normals.c b/xs/src/admesh/normals.c new file mode 100644 index 000000000..baaab36eb --- /dev/null +++ b/xs/src/admesh/normals.c @@ -0,0 +1,391 @@ +/* ADMesh -- process triangulated solid meshes + * Copyright (C) 1995 Anthony D. Martin + * + * 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 + * the Free Software Foundation; either version 2, or (at your option) + * any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program; if not, write to the Free Software + * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. + * + * Questions, comments, suggestions, etc to + */ + +#include +#include +#include +#include + +#include "stl.h" + +static void stl_reverse_facet(stl_file *stl, int facet_num); +/* static float stl_calculate_area(stl_facet *facet); */ +static void stl_reverse_vector(float v[]); +int stl_check_normal_vector(stl_file *stl, int facet_num, int normal_fix_flag); + +static void +stl_reverse_facet(stl_file *stl, int facet_num) +{ + stl_vertex tmp_vertex; + /* int tmp_neighbor;*/ + int neighbor[3]; + int vnot[3]; + + stl->stats.facets_reversed += 1; + + neighbor[0] = stl->neighbors_start[facet_num].neighbor[0]; + neighbor[1] = stl->neighbors_start[facet_num].neighbor[1]; + neighbor[2] = stl->neighbors_start[facet_num].neighbor[2]; + vnot[0] = stl->neighbors_start[facet_num].which_vertex_not[0]; + vnot[1] = stl->neighbors_start[facet_num].which_vertex_not[1]; + vnot[2] = stl->neighbors_start[facet_num].which_vertex_not[2]; + + /* reverse the facet */ + tmp_vertex = stl->facet_start[facet_num].vertex[0]; + stl->facet_start[facet_num].vertex[0] = + stl->facet_start[facet_num].vertex[1]; + stl->facet_start[facet_num].vertex[1] = tmp_vertex; + + /* fix the vnots of the neighboring facets */ + if(neighbor[0] != -1) + stl->neighbors_start[neighbor[0]].which_vertex_not[(vnot[0] + 1) % 3] = + (stl->neighbors_start[neighbor[0]]. + which_vertex_not[(vnot[0] + 1) % 3] + 3) % 6; + if(neighbor[1] != -1) + stl->neighbors_start[neighbor[1]].which_vertex_not[(vnot[1] + 1) % 3] = + (stl->neighbors_start[neighbor[1]]. + which_vertex_not[(vnot[1] + 1) % 3] + 4) % 6; + if(neighbor[2] != -1) + stl->neighbors_start[neighbor[2]].which_vertex_not[(vnot[2] + 1) % 3] = + (stl->neighbors_start[neighbor[2]]. + which_vertex_not[(vnot[2] + 1) % 3] + 2) % 6; + + /* swap the neighbors of the facet that is being reversed */ + stl->neighbors_start[facet_num].neighbor[1] = neighbor[2]; + stl->neighbors_start[facet_num].neighbor[2] = neighbor[1]; + + /* swap the vnots of the facet that is being reversed */ + stl->neighbors_start[facet_num].which_vertex_not[1] = vnot[2]; + stl->neighbors_start[facet_num].which_vertex_not[2] = vnot[1]; + + /* reverse the values of the vnots of the facet that is being reversed */ + stl->neighbors_start[facet_num].which_vertex_not[0] = + (stl->neighbors_start[facet_num].which_vertex_not[0] + 3) % 6; + stl->neighbors_start[facet_num].which_vertex_not[1] = + (stl->neighbors_start[facet_num].which_vertex_not[1] + 3) % 6; + stl->neighbors_start[facet_num].which_vertex_not[2] = + (stl->neighbors_start[facet_num].which_vertex_not[2] + 3) % 6; +} + +void +stl_fix_normal_directions(stl_file *stl) +{ + char *norm_sw; + /* int edge_num;*/ + /* int vnot;*/ + int checked = 0; + int facet_num; + /* int next_facet;*/ + int i; + int j; + int checked_before = 0; + struct stl_normal + { + int facet_num; + struct stl_normal *next; + }; + struct stl_normal *head; + struct stl_normal *tail; + struct stl_normal *newn; + struct stl_normal *temp; + + + /* Initialize linked list. */ + head = (stl_normal*)malloc(sizeof(struct stl_normal)); + if(head == NULL) perror("stl_fix_normal_directions"); + tail = (stl_normal*)malloc(sizeof(struct stl_normal)); + if(tail == NULL) perror("stl_fix_normal_directions"); + head->next = tail; + tail->next = tail; + + /* Initialize list that keeps track of already fixed facets. */ + norm_sw = (char*)calloc(stl->stats.number_of_facets, sizeof(char)); + if(norm_sw == NULL) perror("stl_fix_normal_directions"); + + + facet_num = 0; + if(stl_check_normal_vector(stl, 0, 0) == 2) + stl_reverse_facet(stl, 0); + + norm_sw[facet_num] = 1; + /* edge_num = 0; + vnot = stl->neighbors_start[0].which_vertex_not[0]; + */ + checked++; + + for(;;) + { + /* Add neighbors_to_list. */ + for(j = 0; j < 3; j++) + { + /* Reverse the neighboring facets if necessary. */ + if(stl->neighbors_start[facet_num].which_vertex_not[j] > 2) + { + if(stl->neighbors_start[facet_num].neighbor[j] != -1) + { + stl_reverse_facet + (stl, stl->neighbors_start[facet_num].neighbor[j]); + } + } + if(stl->neighbors_start[facet_num].neighbor[j] != -1) + { + if(norm_sw[stl->neighbors_start[facet_num].neighbor[j]] != 1) + { + /* Add node to beginning of list. */ + newn = (stl_normal*)malloc(sizeof(struct stl_normal)); + if(newn == NULL) perror("stl_fix_normal_directions"); + newn->facet_num = stl->neighbors_start[facet_num].neighbor[j]; + newn->next = head->next; + head->next = newn; + } + } + } + /* Get next facet to fix from top of list. */ + if(head->next != tail) + { + facet_num = head->next->facet_num; + if(norm_sw[facet_num] != 1) /* If facet is in list mutiple times */ + { + norm_sw[facet_num] = 1; /* Record this one as being fixed. */ + checked++; + } + temp = head->next; /* Delete this facet from the list. */ + head->next = head->next->next; + free(temp); + } + else + { + /* All of the facets in this part have been fixed. */ + stl->stats.number_of_parts += 1; + /* There are (checked-checked_before) facets */ + /* in part stl->stats.number_of_parts */ + checked_before = checked; + if(checked == stl->stats.number_of_facets) + { + /* All of the facets have been checked. Bail out. */ + break; + } + else + { + /* There is another part here. Find it and continue. */ + for(i = 0; i < stl->stats.number_of_facets; i++) + { + if(norm_sw[i] == 0) + { /* This is the first facet of the next part. */ + facet_num = i; + if(stl_check_normal_vector(stl, i, 0) == 2) + { + stl_reverse_facet(stl, i); + } + + norm_sw[facet_num] = 1; + checked++; + break; + } + } + } + } + } + free(head); + free(tail); + free(norm_sw); +} + +int +stl_check_normal_vector(stl_file *stl, int facet_num, int normal_fix_flag) +{ + /* Returns 0 if the normal is within tolerance */ + /* Returns 1 if the normal is not within tolerance, but direction is OK */ + /* Returns 2 if the normal is not within tolerance and backwards */ + /* Returns 4 if the status is unknown. */ + + float normal[3]; + float test_norm[3]; + stl_facet *facet; + + facet = &stl->facet_start[facet_num]; + + stl_calculate_normal(normal, facet); + stl_normalize_vector(normal); + + if( (ABS(normal[0] - facet->normal.x) < 0.001) + && (ABS(normal[1] - facet->normal.y) < 0.001) + && (ABS(normal[2] - facet->normal.z) < 0.001)) + { + /* It is not really necessary to change the values here */ + /* but just for consistency, I will. */ + facet->normal.x = normal[0]; + facet->normal.y = normal[1]; + facet->normal.z = normal[2]; + return 0; + } + + test_norm[0] = facet->normal.x; + test_norm[1] = facet->normal.y; + test_norm[2] = facet->normal.z; + + stl_normalize_vector(test_norm); + if( (ABS(normal[0] - test_norm[0]) < 0.001) + && (ABS(normal[1] - test_norm[1]) < 0.001) + && (ABS(normal[2] - test_norm[2]) < 0.001)) + { + if(normal_fix_flag) + { + facet->normal.x = normal[0]; + facet->normal.y = normal[1]; + facet->normal.z = normal[2]; + stl->stats.normals_fixed += 1; + } + return 1; + } + + stl_reverse_vector(test_norm); + if( (ABS(normal[0] - test_norm[0]) < 0.001) + && (ABS(normal[1] - test_norm[1]) < 0.001) + && (ABS(normal[2] - test_norm[2]) < 0.001)) + { + /* Facet is backwards. */ + if(normal_fix_flag) + { + facet->normal.x = normal[0]; + facet->normal.y = normal[1]; + facet->normal.z = normal[2]; + stl->stats.normals_fixed += 1; + } + return 2; + } + if(normal_fix_flag) + { + facet->normal.x = normal[0]; + facet->normal.y = normal[1]; + facet->normal.z = normal[2]; + stl->stats.normals_fixed += 1; + } + return 4; +} + +static void +stl_reverse_vector(float v[]) +{ + v[0] *= -1; + v[1] *= -1; + v[2] *= -1; +} + + +void +stl_calculate_normal(float normal[], stl_facet *facet) +{ + float v1[3]; + float v2[3]; + + v1[0] = facet->vertex[1].x - facet->vertex[0].x; + v1[1] = facet->vertex[1].y - facet->vertex[0].y; + v1[2] = facet->vertex[1].z - facet->vertex[0].z; + v2[0] = facet->vertex[2].x - facet->vertex[0].x; + v2[1] = facet->vertex[2].y - facet->vertex[0].y; + v2[2] = facet->vertex[2].z - facet->vertex[0].z; + + normal[0] = (float)((double)v1[1] * (double)v2[2]) - ((double)v1[2] * (double)v2[1]); + normal[1] = (float)((double)v1[2] * (double)v2[0]) - ((double)v1[0] * (double)v2[2]); + normal[2] = (float)((double)v1[0] * (double)v2[1]) - ((double)v1[1] * (double)v2[0]); +} + +/* +static float +stl_calculate_area(stl_facet *facet) +{ + float cross[3][3]; + float sum[3]; + float normal[3]; + float area; + int i; + + for(i = 0; i < 3; i++) + { + cross[i][0] = ((facet->vertex[i].y * facet->vertex[(i + 1) % 3].z) - + (facet->vertex[i].z * facet->vertex[(i + 1) % 3].y)); + cross[i][1] = ((facet->vertex[i].z * facet->vertex[(i + 1) % 3].x) - + (facet->vertex[i].x * facet->vertex[(i + 1) % 3].z)); + cross[i][2] = ((facet->vertex[i].x * facet->vertex[(i + 1) % 3].y) - + (facet->vertex[i].y * facet->vertex[(i + 1) % 3].x)); + } + + sum[0] = cross[0][0] + cross[1][0] + cross[2][0]; + sum[1] = cross[0][1] + cross[1][1] + cross[2][1]; + sum[2] = cross[0][2] + cross[1][2] + cross[2][2]; + + stl_calculate_normal(normal, facet); + stl_normalize_vector(normal); + area = 0.5 * (normal[0] * sum[0] + normal[1] * sum[1] + + normal[2] * sum[2]); + return ABS(area); +} +*/ + +void stl_normalize_vector(float v[]) +{ + double length; + double factor; + float min_normal_length; + + length = sqrt((double)v[0] * (double)v[0] + (double)v[1] * (double)v[1] + (double)v[2] * (double)v[2]); + min_normal_length = 0.000000000001; + if(length < min_normal_length) + { + v[0] = 1.0; + v[1] = 0.0; + v[2] = 0.0; + return; + } + factor = 1.0 / length; + v[0] *= factor; + v[1] *= factor; + v[2] *= factor; +} + +void +stl_fix_normal_values(stl_file *stl) +{ + int i; + + for(i = 0; i < stl->stats.number_of_facets; i++) + { + stl_check_normal_vector(stl, i, 1); + } +} + +void +stl_reverse_all_facets(stl_file *stl) +{ + int i; + float normal[3]; + + for(i = 0; i < stl->stats.number_of_facets; i++) + { + stl_reverse_facet(stl, i); + stl_calculate_normal(normal, &stl->facet_start[i]); + stl_normalize_vector(normal); + stl->facet_start[i].normal.x = normal[0]; + stl->facet_start[i].normal.y = normal[1]; + stl->facet_start[i].normal.z = normal[2]; + } +} + diff --git a/xs/src/admesh/shared.c b/xs/src/admesh/shared.c new file mode 100644 index 000000000..681ec0ff9 --- /dev/null +++ b/xs/src/admesh/shared.c @@ -0,0 +1,243 @@ +/* ADMesh -- process triangulated solid meshes + * Copyright (C) 1995, 1996 Anthony D. Martin + * + * 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 + * the Free Software Foundation; either version 2, or (at your option) + * any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program; if not, write to the Free Software + * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. + * + * Questions, comments, suggestions, etc to + */ + +#include +#include + +#include "stl.h" + + +void +stl_generate_shared_vertices(stl_file *stl) +{ + int i; + int j; + int first_facet; + int direction; + int facet_num; + int vnot; + int next_edge; + int pivot_vertex; + int next_facet; + int reversed; + + stl->v_indices = (v_indices_struct*) + calloc(stl->stats.number_of_facets, sizeof(v_indices_struct)); + if(stl->v_indices == NULL) perror("stl_generate_shared_vertices"); + stl->v_shared = (stl_vertex*) + calloc((stl->stats.number_of_facets / 2), sizeof(stl_vertex)); + if(stl->v_shared == NULL) perror("stl_generate_shared_vertices"); + stl->stats.shared_malloced = stl->stats.number_of_facets / 2; + stl->stats.shared_vertices = 0; + + for(i = 0; i < stl->stats.number_of_facets; i++) + { + stl->v_indices[i].vertex[0] = -1; + stl->v_indices[i].vertex[1] = -1; + stl->v_indices[i].vertex[2] = -1; + } + + + for(i = 0; i < stl->stats.number_of_facets; i++) + { + first_facet = i; + for(j = 0; j < 3; j++) + { + if(stl->v_indices[i].vertex[j] != -1) + { + continue; + } + if(stl->stats.shared_vertices == stl->stats.shared_malloced) + { + stl->stats.shared_malloced += 1024; + stl->v_shared = (stl_vertex*)realloc(stl->v_shared, + stl->stats.shared_malloced * sizeof(stl_vertex)); + if(stl->v_shared == NULL) perror("stl_generate_shared_vertices"); + } + + stl->v_shared[stl->stats.shared_vertices] = + stl->facet_start[i].vertex[j]; + + direction = 0; + reversed = 0; + facet_num = i; + vnot = (j + 2) % 3; + + for(;;) + { + if(vnot > 2) + { + if(direction == 0) + { + pivot_vertex = (vnot + 2) % 3; + next_edge = pivot_vertex; + direction = 1; + } + else + { + pivot_vertex = (vnot + 1) % 3; + next_edge = vnot % 3; + direction = 0; + } + } + else + { + if(direction == 0) + { + pivot_vertex = (vnot + 1) % 3; + next_edge = vnot; + } + else + { + pivot_vertex = (vnot + 2) % 3; + next_edge = pivot_vertex; + } + } + stl->v_indices[facet_num].vertex[pivot_vertex] = + stl->stats.shared_vertices; + + next_facet = stl->neighbors_start[facet_num].neighbor[next_edge]; + if(next_facet == -1) + { + if(reversed) + { + break; + } + else + { + direction = 1; + vnot = (j + 1) % 3; + reversed = 1; + facet_num = first_facet; + } + } + else if(next_facet != first_facet) + { + vnot = stl->neighbors_start[facet_num]. + which_vertex_not[next_edge]; + facet_num = next_facet; + } + else + { + break; + } + } + stl->stats.shared_vertices += 1; + } + } +} + +void +stl_write_off(stl_file *stl, char *file) +{ + int i; + FILE *fp; + char *error_msg; + + + /* Open the file */ + fp = fopen(file, "w"); + if(fp == NULL) + { + error_msg = (char*) + malloc(81 + strlen(file)); /* Allow 80 chars+file size for message */ + sprintf(error_msg, "stl_write_ascii: Couldn't open %s for writing", + file); + perror(error_msg); + free(error_msg); + exit(1); + } + + fprintf(fp, "OFF\n"); + fprintf(fp, "%d %d 0\n", + stl->stats.shared_vertices, stl->stats.number_of_facets); + + for(i = 0; i < stl->stats.shared_vertices; i++) + { + fprintf(fp, "\t%f %f %f\n", + stl->v_shared[i].x, stl->v_shared[i].y, stl->v_shared[i].z); + } + for(i = 0; i < stl->stats.number_of_facets; i++) + { + fprintf(fp, "\t3 %d %d %d\n", stl->v_indices[i].vertex[0], + stl->v_indices[i].vertex[1], stl->v_indices[i].vertex[2]); + } + fclose(fp); +} + +void +stl_write_vrml(stl_file *stl, char *file) +{ + int i; + FILE *fp; + char *error_msg; + + + /* Open the file */ + fp = fopen(file, "w"); + if(fp == NULL) + { + error_msg = (char*) + malloc(81 + strlen(file)); /* Allow 80 chars+file size for message */ + sprintf(error_msg, "stl_write_ascii: Couldn't open %s for writing", + file); + perror(error_msg); + free(error_msg); + exit(1); + } + + fprintf(fp, "#VRML V1.0 ascii\n\n"); + fprintf(fp, "Separator {\n"); + fprintf(fp, "\tDEF STLShape ShapeHints {\n"); + fprintf(fp, "\t\tvertexOrdering COUNTERCLOCKWISE\n"); + fprintf(fp, "\t\tfaceType CONVEX\n"); + fprintf(fp, "\t\tshapeType SOLID\n"); + fprintf(fp, "\t\tcreaseAngle 0.0\n"); + fprintf(fp, "\t}\n"); + fprintf(fp, "\tDEF STLModel Separator {\n"); + fprintf(fp, "\t\tDEF STLColor Material {\n"); + fprintf(fp, "\t\t\temissiveColor 0.700000 0.700000 0.000000\n"); + fprintf(fp, "\t\t}\n"); + fprintf(fp, "\t\tDEF STLVertices Coordinate3 {\n"); + fprintf(fp, "\t\t\tpoint [\n"); + + for(i = 0; i < (stl->stats.shared_vertices - 1); i++) + { + fprintf(fp, "\t\t\t\t%f %f %f,\n", + stl->v_shared[i].x, stl->v_shared[i].y, stl->v_shared[i].z); + } + fprintf(fp, "\t\t\t\t%f %f %f]\n", + stl->v_shared[i].x, stl->v_shared[i].y, stl->v_shared[i].z); + fprintf(fp, "\t\t}\n"); + fprintf(fp, "\t\tDEF STLTriangles IndexedFaceSet {\n"); + fprintf(fp, "\t\t\tcoordIndex [\n"); + + for(i = 0; i < (stl->stats.number_of_facets - 1); i++) + { + fprintf(fp, "\t\t\t\t%d, %d, %d, -1,\n", stl->v_indices[i].vertex[0], + stl->v_indices[i].vertex[1], stl->v_indices[i].vertex[2]); + } + fprintf(fp, "\t\t\t\t%d, %d, %d, -1]\n", stl->v_indices[i].vertex[0], + stl->v_indices[i].vertex[1], stl->v_indices[i].vertex[2]); + fprintf(fp, "\t\t}\n"); + fprintf(fp, "\t}\n"); + fprintf(fp, "}\n"); + fclose(fp); +} diff --git a/xs/src/admesh/stl.h b/xs/src/admesh/stl.h new file mode 100644 index 000000000..370a95c26 --- /dev/null +++ b/xs/src/admesh/stl.h @@ -0,0 +1,173 @@ +/* ADMesh -- process triangulated solid meshes + * Copyright (C) 1995, 1996 Anthony D. Martin + * + * 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 + * the Free Software Foundation; either version 2, or (at your option) + * any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program; if not, write to the Free Software + * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. + * + * Questions, comments, suggestions, etc to + */ + +#include + +#define STL_MAX(A,B) ((A)>(B)? (A):(B)) +#define STL_MIN(A,B) ((A)<(B)? (A):(B)) +#define ABS(X) ((X) < 0 ? -(X) : (X)) + +#define LABEL_SIZE 80 +#define NUM_FACET_SIZE 4 +#define HEADER_SIZE 84 +#define STL_MIN_FILE_SIZE 284 +#define ASCII_LINES_PER_FACET 7 +#define SIZEOF_EDGE_SORT 24 + +typedef struct +{ + float x; + float y; + float z; +}stl_vertex; + +typedef struct +{ + float x; + float y; + float z; +}stl_normal; + +typedef char stl_extra[2]; + +typedef struct +{ + stl_normal normal; + stl_vertex vertex[3]; + stl_extra extra; +}stl_facet; +#define SIZEOF_STL_FACET 50 + +typedef enum {binary, ascii} stl_type; + +typedef struct +{ + stl_vertex p1; + stl_vertex p2; + int facet_number; +}stl_edge; + +typedef struct stl_hash_edge +{ + unsigned key[6]; + int facet_number; + int which_edge; + struct stl_hash_edge *next; +}stl_hash_edge; + +typedef struct +{ + int neighbor[3]; + char which_vertex_not[3]; +}stl_neighbors; + +typedef struct +{ + int vertex[3]; +}v_indices_struct; + +typedef struct +{ + char header[81]; + stl_type type; + int number_of_facets; + stl_vertex max; + stl_vertex min; + stl_vertex size; + float bounding_diameter; + float shortest_edge; + float volume; + unsigned number_of_blocks; + int connected_edges; + int connected_facets_1_edge; + int connected_facets_2_edge; + int connected_facets_3_edge; + int facets_w_1_bad_edge; + int facets_w_2_bad_edge; + int facets_w_3_bad_edge; + int original_num_facets; + int edges_fixed; + int degenerate_facets; + int facets_removed; + int facets_added; + int facets_reversed; + int backwards_edges; + int normals_fixed; + int number_of_parts; + int malloced; + int freed; + int facets_malloced; + int collisions; + int shared_vertices; + int shared_malloced; +}stl_stats; + +typedef struct +{ + FILE *fp; + stl_facet *facet_start; + stl_edge *edge_start; + stl_hash_edge **heads; + stl_hash_edge *tail; + int M; + stl_neighbors *neighbors_start; + v_indices_struct *v_indices; + stl_vertex *v_shared; + stl_stats stats; +}stl_file; + + +extern void stl_open(stl_file *stl, char *file); +extern void stl_close(stl_file *stl); +extern void stl_stats_out(stl_file *stl, FILE *file, char *input_file); +extern void stl_print_edges(stl_file *stl, FILE *file); +extern void stl_print_neighbors(stl_file *stl, char *file); +extern void stl_write_ascii(stl_file *stl, char *file, char *label); +extern void stl_write_binary(stl_file *stl, char *file, char *label); +extern void stl_check_facets_exact(stl_file *stl); +extern void stl_check_facets_nearby(stl_file *stl, float tolerance); +extern void stl_remove_unconnected_facets(stl_file *stl); +extern void stl_write_vertex(stl_file *stl, int facet, int vertex); +extern void stl_write_facet(stl_file *stl, char *label, int facet); +extern void stl_write_edge(stl_file *stl, char *label, stl_hash_edge edge); +extern void stl_write_neighbor(stl_file *stl, int facet); +extern void stl_write_quad_object(stl_file *stl, char *file); +extern void stl_verify_neighbors(stl_file *stl); +extern void stl_fill_holes(stl_file *stl); +extern void stl_fix_normal_directions(stl_file *stl); +extern void stl_fix_normal_values(stl_file *stl); +extern void stl_reverse_all_facets(stl_file *stl); +extern void stl_translate(stl_file *stl, float x, float y, float z); +extern void stl_scale(stl_file *stl, float factor); +extern void stl_rotate_x(stl_file *stl, float angle); +extern void stl_rotate_y(stl_file *stl, float angle); +extern void stl_rotate_z(stl_file *stl, float angle); +extern void stl_mirror_xy(stl_file *stl); +extern void stl_mirror_yz(stl_file *stl); +extern void stl_mirror_xz(stl_file *stl); +extern void stl_open_merge(stl_file *stl, char *file); +extern void stl_generate_shared_vertices(stl_file *stl); +extern void stl_write_off(stl_file *stl, char *file); +extern void stl_write_dxf(stl_file *stl, char *file, char *label); +extern void stl_write_vrml(stl_file *stl, char *file); +extern void stl_calculate_normal(float normal[], stl_facet *facet); +extern void stl_normalize_vector(float v[]); +extern void stl_calculate_volume(stl_file *stl); + diff --git a/xs/src/admesh/stl_io.c b/xs/src/admesh/stl_io.c new file mode 100644 index 000000000..1732bc4ae --- /dev/null +++ b/xs/src/admesh/stl_io.c @@ -0,0 +1,467 @@ +/* ADMesh -- process triangulated solid meshes + * Copyright (C) 1995, 1996 Anthony D. Martin + * + * 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 + * the Free Software Foundation; either version 2, or (at your option) + * any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program; if not, write to the Free Software + * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. + * + * Questions, comments, suggestions, etc to + */ + +#include +#include +#include "stl.h" + +#if !defined(SEEK_SET) +#define SEEK_SET 0 +#define SEEK_CUR 1 +#define SEEK_END 2 +#endif + +static void stl_put_little_int(FILE *fp, int value); +static void stl_put_little_float(FILE *fp, float value_in); + +void +stl_print_edges(stl_file *stl, FILE *file) +{ + int i; + int edges_allocated; + + edges_allocated = stl->stats.number_of_facets * 3; + for(i = 0; i < edges_allocated; i++) + { + fprintf(file, "%d, %f, %f, %f, %f, %f, %f\n", + stl->edge_start[i].facet_number, + stl->edge_start[i].p1.x, stl->edge_start[i].p1.y, + stl->edge_start[i].p1.z, stl->edge_start[i].p2.x, + stl->edge_start[i].p2.y, stl->edge_start[i].p2.z); + } +} + + +void +stl_stats_out(stl_file *stl, FILE *file, char *input_file) +{ + fprintf(file, "\n\ +================= Results produced by ADMesh version 0.95 ================\n"); + fprintf(file, "\ +Input file : %s\n", input_file); + if(stl->stats.type == binary) + { + fprintf(file, "\ +File type : Binary STL file\n"); + } + else + { + fprintf(file, "\ +File type : ASCII STL file\n"); + } + fprintf(file, "\ +Header : %s\n", stl->stats.header); + fprintf(file, "============== Size ==============\n"); + fprintf(file, "Min X = % f, Max X = % f\n", + stl->stats.min.x, stl->stats.max.x); + fprintf(file, "Min Y = % f, Max Y = % f\n", + stl->stats.min.y, stl->stats.max.y); + fprintf(file, "Min Z = % f, Max Z = % f\n", + stl->stats.min.z, stl->stats.max.z); + + fprintf(file, "\ +========= Facet Status ========== Original ============ Final ====\n"); + fprintf(file, "\ +Number of facets : %5d %5d\n", + stl->stats.original_num_facets, stl->stats.number_of_facets); + fprintf(file, "\ +Facets with 1 disconnected edge : %5d %5d\n", + stl->stats.facets_w_1_bad_edge, stl->stats.connected_facets_2_edge - + stl->stats.connected_facets_3_edge); + fprintf(file, "\ +Facets with 2 disconnected edges : %5d %5d\n", + stl->stats.facets_w_2_bad_edge, stl->stats.connected_facets_1_edge - + stl->stats.connected_facets_2_edge); + fprintf(file, "\ +Facets with 3 disconnected edges : %5d %5d\n", + stl->stats.facets_w_3_bad_edge, stl->stats.number_of_facets - + stl->stats.connected_facets_1_edge); + fprintf(file, "\ +Total disconnected facets : %5d %5d\n", + stl->stats.facets_w_1_bad_edge + stl->stats.facets_w_2_bad_edge + + stl->stats.facets_w_3_bad_edge, stl->stats.number_of_facets - + stl->stats.connected_facets_3_edge); + + fprintf(file, +"=== Processing Statistics === ===== Other Statistics =====\n"); + fprintf(file, "\ +Number of parts : %5d Volume : % f\n", + stl->stats.number_of_parts, stl->stats.volume); + fprintf(file, "\ +Degenerate facets : %5d\n", stl->stats.degenerate_facets); + fprintf(file, "\ +Edges fixed : %5d\n", stl->stats.edges_fixed); + fprintf(file, "\ +Facets removed : %5d\n", stl->stats.facets_removed); + fprintf(file, "\ +Facets added : %5d\n", stl->stats.facets_added); + fprintf(file, "\ +Facets reversed : %5d\n", stl->stats.facets_reversed); + fprintf(file, "\ +Backwards edges : %5d\n", stl->stats.backwards_edges); + fprintf(file, "\ +Normals fixed : %5d\n", stl->stats.normals_fixed); +} + +void +stl_write_ascii(stl_file *stl, char *file, char *label) +{ + int i; + FILE *fp; + char *error_msg; + + + /* Open the file */ + fp = fopen(file, "w"); + if(fp == NULL) + { + error_msg = (char*) + malloc(81 + strlen(file)); /* Allow 80 chars+file size for message */ + sprintf(error_msg, "stl_write_ascii: Couldn't open %s for writing", + file); + perror(error_msg); + free(error_msg); + exit(1); + } + + fprintf(fp, "solid %s\n", label); + + for(i = 0; i < stl->stats.number_of_facets; i++) + { + fprintf(fp, " facet normal % .8E % .8E % .8E\n", + stl->facet_start[i].normal.x, stl->facet_start[i].normal.y, + stl->facet_start[i].normal.z); + fprintf(fp, " outer loop\n"); + fprintf(fp, " vertex % .8E % .8E % .8E\n", + stl->facet_start[i].vertex[0].x, stl->facet_start[i].vertex[0].y, + stl->facet_start[i].vertex[0].z); + fprintf(fp, " vertex % .8E % .8E % .8E\n", + stl->facet_start[i].vertex[1].x, stl->facet_start[i].vertex[1].y, + stl->facet_start[i].vertex[1].z); + fprintf(fp, " vertex % .8E % .8E % .8E\n", + stl->facet_start[i].vertex[2].x, stl->facet_start[i].vertex[2].y, + stl->facet_start[i].vertex[2].z); + fprintf(fp, " endloop\n"); + fprintf(fp, " endfacet\n"); + } + + fprintf(fp, "endsolid %s\n", label); + + fclose(fp); +} + +void +stl_print_neighbors(stl_file *stl, char *file) +{ + int i; + FILE *fp; + char *error_msg; + + /* Open the file */ + fp = fopen(file, "w"); + if(fp == NULL) + { + error_msg = (char*) + malloc(81 + strlen(file)); /* Allow 80 chars+file size for message */ + sprintf(error_msg, "stl_print_neighbors: Couldn't open %s for writing", + file); + perror(error_msg); + free(error_msg); + exit(1); + } + + for(i = 0; i < stl->stats.number_of_facets; i++) + { + fprintf(fp, "%d, %d,%d, %d,%d, %d,%d\n", + i, + stl->neighbors_start[i].neighbor[0], + (int)stl->neighbors_start[i].which_vertex_not[0], + stl->neighbors_start[i].neighbor[1], + (int)stl->neighbors_start[i].which_vertex_not[1], + stl->neighbors_start[i].neighbor[2], + (int)stl->neighbors_start[i].which_vertex_not[2]); + } +} + +static void +stl_put_little_int(FILE *fp, int value_in) +{ + int new_value; + union + { + int int_value; + char char_value[4]; + } value; + + value.int_value = value_in; + + new_value = value.char_value[0] & 0xFF; + new_value |= (value.char_value[1] & 0xFF) << 0x08; + new_value |= (value.char_value[2] & 0xFF) << 0x10; + new_value |= (value.char_value[3] & 0xFF) << 0x18; + fwrite(&new_value, sizeof(int), 1, fp); +} + +static void +stl_put_little_float(FILE *fp, float value_in) +{ + int new_value; + union + { + float float_value; + char char_value[4]; + } value; + + value.float_value = value_in; + + new_value = value.char_value[0] & 0xFF; + new_value |= (value.char_value[1] & 0xFF) << 0x08; + new_value |= (value.char_value[2] & 0xFF) << 0x10; + new_value |= (value.char_value[3] & 0xFF) << 0x18; + fwrite(&new_value, sizeof(int), 1, fp); +} + + +void +stl_write_binary(stl_file *stl, char *file, char *label) +{ + FILE *fp; + int i; + char *error_msg; + + + /* Open the file */ + fp = fopen(file, "w"); + if(fp == NULL) + { + error_msg = (char*) + malloc(81 + strlen(file)); /* Allow 80 chars+file size for message */ + sprintf(error_msg, "stl_write_binary: Couldn't open %s for writing", + file); + perror(error_msg); + free(error_msg); + exit(1); + } + + fprintf(fp, "%s", label); + for(i = strlen(label); i < LABEL_SIZE; i++) putc(0, fp); + + fseek(fp, LABEL_SIZE, SEEK_SET); + + stl_put_little_int(fp, stl->stats.number_of_facets); + + for(i = 0; i < stl->stats.number_of_facets; i++) + { + stl_put_little_float(fp, stl->facet_start[i].normal.x); + stl_put_little_float(fp, stl->facet_start[i].normal.y); + stl_put_little_float(fp, stl->facet_start[i].normal.z); + stl_put_little_float(fp, stl->facet_start[i].vertex[0].x); + stl_put_little_float(fp, stl->facet_start[i].vertex[0].y); + stl_put_little_float(fp, stl->facet_start[i].vertex[0].z); + stl_put_little_float(fp, stl->facet_start[i].vertex[1].x); + stl_put_little_float(fp, stl->facet_start[i].vertex[1].y); + stl_put_little_float(fp, stl->facet_start[i].vertex[1].z); + stl_put_little_float(fp, stl->facet_start[i].vertex[2].x); + stl_put_little_float(fp, stl->facet_start[i].vertex[2].y); + stl_put_little_float(fp, stl->facet_start[i].vertex[2].z); + fputc(stl->facet_start[i].extra[0], fp); + fputc(stl->facet_start[i].extra[1], fp); + } + + fclose(fp); +} + +void +stl_write_vertex(stl_file *stl, int facet, int vertex) +{ + printf(" vertex %d/%d % .8E % .8E % .8E\n", vertex, facet, + stl->facet_start[facet].vertex[vertex].x, + stl->facet_start[facet].vertex[vertex].y, + stl->facet_start[facet].vertex[vertex].z); +} + +void +stl_write_facet(stl_file *stl, char *label, int facet) +{ + printf("facet (%d)/ %s\n", facet, label); + stl_write_vertex(stl, facet, 0); + stl_write_vertex(stl, facet, 1); + stl_write_vertex(stl, facet, 2); +} + +void +stl_write_edge(stl_file *stl, char *label, stl_hash_edge edge) +{ + printf("edge (%d)/(%d) %s\n", edge.facet_number, edge.which_edge, label); + if(edge.which_edge < 3) + { + stl_write_vertex(stl, edge.facet_number, edge.which_edge % 3); + stl_write_vertex(stl, edge.facet_number, (edge.which_edge + 1) % 3); + } + else + { + stl_write_vertex(stl, edge.facet_number, (edge.which_edge + 1) % 3); + stl_write_vertex(stl, edge.facet_number, edge.which_edge % 3); + } +} + +void +stl_write_neighbor(stl_file *stl, int facet) +{ + printf("Neighbors %d: %d, %d, %d ; %d, %d, %d\n", facet, + stl->neighbors_start[facet].neighbor[0], + stl->neighbors_start[facet].neighbor[1], + stl->neighbors_start[facet].neighbor[2], + stl->neighbors_start[facet].which_vertex_not[0], + stl->neighbors_start[facet].which_vertex_not[1], + stl->neighbors_start[facet].which_vertex_not[2]); +} + +void +stl_write_quad_object(stl_file *stl, char *file) +{ + FILE *fp; + int i; + int j; + char *error_msg; + stl_vertex connect_color; + stl_vertex uncon_1_color; + stl_vertex uncon_2_color; + stl_vertex uncon_3_color; + stl_vertex color; + + /* Open the file */ + fp = fopen(file, "w"); + if(fp == NULL) + { + error_msg = (char*) + malloc(81 + strlen(file)); /* Allow 80 chars+file size for message */ + sprintf(error_msg, "stl_write_quad_object: Couldn't open %s for writing", + file); + perror(error_msg); + free(error_msg); + exit(1); + } + + connect_color.x = 0.0; + connect_color.y = 0.0; + connect_color.z = 1.0; + uncon_1_color.x = 0.0; + uncon_1_color.y = 1.0; + uncon_1_color.z = 0.0; + uncon_2_color.x = 1.0; + uncon_2_color.y = 1.0; + uncon_2_color.z = 1.0; + uncon_3_color.x = 1.0; + uncon_3_color.y = 0.0; + uncon_3_color.z = 0.0; + + fprintf(fp, "CQUAD\n"); + for(i = 0; i < stl->stats.number_of_facets; i++) + { + j = ((stl->neighbors_start[i].neighbor[0] == -1) + + (stl->neighbors_start[i].neighbor[1] == -1) + + (stl->neighbors_start[i].neighbor[2] == -1)); + if(j == 0) + { + color = connect_color; + } + else if(j == 1) + { + color = uncon_1_color; + } + else if(j == 2) + { + color = uncon_2_color; + } + else + { + color = uncon_3_color; + } + fprintf(fp, "%f %f %f %1.1f %1.1f %1.1f 1\n", + stl->facet_start[i].vertex[0].x, + stl->facet_start[i].vertex[0].y, + stl->facet_start[i].vertex[0].z, color.x, color.y, color.z); + fprintf(fp, "%f %f %f %1.1f %1.1f %1.1f 1\n", + stl->facet_start[i].vertex[1].x, + stl->facet_start[i].vertex[1].y, + stl->facet_start[i].vertex[1].z, color.x, color.y, color.z); + fprintf(fp, "%f %f %f %1.1f %1.1f %1.1f 1\n", + stl->facet_start[i].vertex[2].x, + stl->facet_start[i].vertex[2].y, + stl->facet_start[i].vertex[2].z, color.x, color.y, color.z); + fprintf(fp, "%f %f %f %1.1f %1.1f %1.1f 1\n", + stl->facet_start[i].vertex[2].x, + stl->facet_start[i].vertex[2].y, + stl->facet_start[i].vertex[2].z, color.x, color.y, color.z); + } + fclose(fp); +} + +void +stl_write_dxf(stl_file *stl, char *file, char *label) +{ + int i; + FILE *fp; + char *error_msg; + + + /* Open the file */ + fp = fopen(file, "w"); + if(fp == NULL) + { + error_msg = (char*) + malloc(81 + strlen(file)); /* Allow 80 chars+file size for message */ + sprintf(error_msg, "stl_write_ascii: Couldn't open %s for writing", + file); + perror(error_msg); + free(error_msg); + exit(1); + } + + fprintf(fp, "999\n%s\n", label); + fprintf(fp, "0\nSECTION\n2\nHEADER\n0\nENDSEC\n"); + fprintf(fp, "0\nSECTION\n2\nTABLES\n0\nTABLE\n2\nLAYER\n70\n1\n\ +0\nLAYER\n2\n0\n70\n0\n62\n7\n6\nCONTINUOUS\n0\nENDTAB\n0\nENDSEC\n"); + fprintf(fp, "0\nSECTION\n2\nBLOCKS\n0\nENDSEC\n"); + + fprintf(fp, "0\nSECTION\n2\nENTITIES\n"); + + for(i = 0; i < stl->stats.number_of_facets; i++) + { + fprintf(fp, "0\n3DFACE\n8\n0\n"); + fprintf(fp, "10\n%f\n20\n%f\n30\n%f\n", + stl->facet_start[i].vertex[0].x, stl->facet_start[i].vertex[0].y, + stl->facet_start[i].vertex[0].z); + fprintf(fp, "11\n%f\n21\n%f\n31\n%f\n", + stl->facet_start[i].vertex[1].x, stl->facet_start[i].vertex[1].y, + stl->facet_start[i].vertex[1].z); + fprintf(fp, "12\n%f\n22\n%f\n32\n%f\n", + stl->facet_start[i].vertex[2].x, stl->facet_start[i].vertex[2].y, + stl->facet_start[i].vertex[2].z); + fprintf(fp, "13\n%f\n23\n%f\n33\n%f\n", + stl->facet_start[i].vertex[2].x, stl->facet_start[i].vertex[2].y, + stl->facet_start[i].vertex[2].z); + } + + fprintf(fp, "0\nENDSEC\n0\nEOF\n"); + + fclose(fp); +} diff --git a/xs/src/admesh/stlinit.c b/xs/src/admesh/stlinit.c new file mode 100644 index 000000000..cd6a8e074 --- /dev/null +++ b/xs/src/admesh/stlinit.c @@ -0,0 +1,360 @@ +/* ADMesh -- process triangulated solid meshes + * Copyright (C) 1995, 1996 Anthony D. Martin + * + * 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 + * the Free Software Foundation; either version 2, or (at your option) + * any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program; if not, write to the Free Software + * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. + * + * Questions, comments, suggestions, etc to + */ + +#include +#include +#include +#include + +#include "stl.h" + +#if !defined(SEEK_SET) +#define SEEK_SET 0 +#define SEEK_CUR 1 +#define SEEK_END 2 +#endif + +static void stl_initialize(stl_file *stl, char *file); +static void stl_allocate(stl_file *stl); +static void stl_read(stl_file *stl, int first_facet, int first); +static void stl_reallocate(stl_file *stl); +static int stl_get_little_int(FILE *fp); +static float stl_get_little_float(FILE *fp); + +void +stl_open(stl_file *stl, char *file) +{ + stl_initialize(stl, file); + stl_allocate(stl); + stl_read(stl, 0, 1); + fclose(stl->fp); +} + +static int +stl_get_little_int(FILE *fp) +{ + int value; + value = fgetc(fp) & 0xFF; + value |= (fgetc(fp) & 0xFF) << 0x08; + value |= (fgetc(fp) & 0xFF) << 0x10; + value |= (fgetc(fp) & 0xFF) << 0x18; + return(value); +} + +static float +stl_get_little_float(FILE *fp) +{ + union + { + int int_value; + float float_value; + } value; + + value.int_value = fgetc(fp) & 0xFF; + value.int_value |= (fgetc(fp) & 0xFF) << 0x08; + value.int_value |= (fgetc(fp) & 0xFF) << 0x10; + value.int_value |= (fgetc(fp) & 0xFF) << 0x18; + return(value.float_value); +} + + +static void +stl_initialize(stl_file *stl, char *file) +{ + long file_size; + int header_num_facets; + int num_facets; + int i, j; + unsigned char chtest[128]; + int num_lines = 1; + char *error_msg; + + stl->stats.degenerate_facets = 0; + stl->stats.edges_fixed = 0; + stl->stats.facets_added = 0; + stl->stats.facets_removed = 0; + stl->stats.facets_reversed = 0; + stl->stats.normals_fixed = 0; + stl->stats.number_of_parts = 0; + stl->stats.original_num_facets = 0; + stl->stats.number_of_facets = 0; + stl->stats.volume = -1.0; + + stl->neighbors_start = NULL; + stl->facet_start = NULL; + stl->v_indices = NULL; + stl->v_shared = NULL; + + + /* Open the file */ + stl->fp = fopen(file, "r"); + if(stl->fp == NULL) + { + error_msg = (char*) + malloc(81 + strlen(file)); /* Allow 80 chars+file size for message */ + sprintf(error_msg, "stl_initialize: Couldn't open %s for reading", + file); + perror(error_msg); + free(error_msg); + exit(1); + } + /* Find size of file */ + fseek(stl->fp, 0, SEEK_END); + file_size = ftell(stl->fp); + + /* Check for binary or ASCII file */ + fseek(stl->fp, HEADER_SIZE, SEEK_SET); + fread(chtest, sizeof(chtest), 1, stl->fp); + stl->stats.type = ascii; + for(i = 0; i < sizeof(chtest); i++) + { + if(chtest[i] > 127) + { + stl->stats.type = binary; + break; + } + } + rewind(stl->fp); + + /* Get the header and the number of facets in the .STL file */ + /* If the .STL file is binary, then do the following */ + if(stl->stats.type == binary) + { + /* Test if the STL file has the right size */ + if(((file_size - HEADER_SIZE) % SIZEOF_STL_FACET != 0) + || (file_size < STL_MIN_FILE_SIZE)) + { + fprintf(stderr, "The file %s has the wrong size.\n", file); + exit(1); + } + num_facets = (file_size - HEADER_SIZE) / SIZEOF_STL_FACET; + + /* Read the header */ + fread(stl->stats.header, LABEL_SIZE, 1, stl->fp); + stl->stats.header[80] = '\0'; + + /* Read the int following the header. This should contain # of facets */ + header_num_facets = stl_get_little_int(stl->fp); + if(num_facets != header_num_facets) + { + fprintf(stderr, + "Warning: File size doesn't match number of facets in the header\n"); + } + } + /* Otherwise, if the .STL file is ASCII, then do the following */ + else + { + /* Find the number of facets */ + j = 0; + for(i = 0; i < file_size ; i++) + { + j++; + if(getc(stl->fp) == '\n') + { + if(j > 4) /* don't count short lines */ + { + num_lines++; + } + j = 0; + } + } + rewind(stl->fp); + + /* Get the header */ + for(i = 0; + (i < 80) && (stl->stats.header[i] = getc(stl->fp)) != '\n'; i++); + stl->stats.header[i] = '\0'; /* Lose the '\n' */ + stl->stats.header[80] = '\0'; + + num_facets = num_lines / ASCII_LINES_PER_FACET; + } + stl->stats.number_of_facets += num_facets; + stl->stats.original_num_facets = stl->stats.number_of_facets; +} + +static void +stl_allocate(stl_file *stl) +{ + /* Allocate memory for the entire .STL file */ + stl->facet_start = (stl_facet*)calloc(stl->stats.number_of_facets, + sizeof(stl_facet)); + if(stl->facet_start == NULL) perror("stl_initialize"); + stl->stats.facets_malloced = stl->stats.number_of_facets; + + /* Allocate memory for the neighbors list */ + stl->neighbors_start = (stl_neighbors*) + calloc(stl->stats.number_of_facets, sizeof(stl_neighbors)); + if(stl->facet_start == NULL) perror("stl_initialize"); +} + +void +stl_open_merge(stl_file *stl, char *file) +{ + int first_facet; + + first_facet = stl->stats.number_of_facets; + stl_initialize(stl, file); + stl_reallocate(stl); + stl_read(stl, first_facet, 0); +} + +static void +stl_reallocate(stl_file *stl) +{ + /* Reallocate more memory for the .STL file(s) */ + stl->facet_start = (stl_facet*)realloc(stl->facet_start, stl->stats.number_of_facets * + sizeof(stl_facet)); + if(stl->facet_start == NULL) perror("stl_initialize"); + stl->stats.facets_malloced = stl->stats.number_of_facets; + + /* Reallocate more memory for the neighbors list */ + stl->neighbors_start = (stl_neighbors*) + realloc(stl->neighbors_start, stl->stats.number_of_facets * + sizeof(stl_neighbors)); + if(stl->facet_start == NULL) perror("stl_initialize"); +} + +static void +stl_read(stl_file *stl, int first_facet, int first) +{ + stl_facet facet; + int i; + float diff_x; + float diff_y; + float diff_z; + float max_diff; + + + if(stl->stats.type == binary) + { + fseek(stl->fp, HEADER_SIZE, SEEK_SET); + } + else + { + rewind(stl->fp); + /* Skip the first line of the file */ + while(getc(stl->fp) != '\n'); + } + + for(i = first_facet; i < stl->stats.number_of_facets; i++) + { + if(stl->stats.type == binary) + /* Read a single facet from a binary .STL file */ + { + facet.normal.x = stl_get_little_float(stl->fp); + facet.normal.y = stl_get_little_float(stl->fp); + facet.normal.z = stl_get_little_float(stl->fp); + facet.vertex[0].x = stl_get_little_float(stl->fp); + facet.vertex[0].y = stl_get_little_float(stl->fp); + facet.vertex[0].z = stl_get_little_float(stl->fp); + facet.vertex[1].x = stl_get_little_float(stl->fp); + facet.vertex[1].y = stl_get_little_float(stl->fp); + facet.vertex[1].z = stl_get_little_float(stl->fp); + facet.vertex[2].x = stl_get_little_float(stl->fp); + facet.vertex[2].y = stl_get_little_float(stl->fp); + facet.vertex[2].z = stl_get_little_float(stl->fp); + facet.extra[0] = fgetc(stl->fp); + facet.extra[1] = fgetc(stl->fp); + } + else + /* Read a single facet from an ASCII .STL file */ + { + fscanf(stl->fp, "%*s %*s %f %f %f\n", &facet.normal.x, + &facet.normal.y, &facet.normal.z); + fscanf(stl->fp, "%*s %*s"); + fscanf(stl->fp, "%*s %f %f %f\n", &facet.vertex[0].x, + &facet.vertex[0].y, &facet.vertex[0].z); + fscanf(stl->fp, "%*s %f %f %f\n", &facet.vertex[1].x, + &facet.vertex[1].y, &facet.vertex[1].z); + fscanf(stl->fp, "%*s %f %f %f\n", &facet.vertex[2].x, + &facet.vertex[2].y, &facet.vertex[2].z); + fscanf(stl->fp, "%*s"); + fscanf(stl->fp, "%*s"); + } + /* Write the facet into memory. */ + stl->facet_start[i] = facet; + + /* while we are going through all of the facets, let's find the */ + /* maximum and minimum values for x, y, and z */ + + /* Initialize the max and min values the first time through*/ + if(first) + { + stl->stats.max.x = facet.vertex[0].x; + stl->stats.min.x = facet.vertex[0].x; + stl->stats.max.y = facet.vertex[0].y; + stl->stats.min.y = facet.vertex[0].y; + stl->stats.max.z = facet.vertex[0].z; + stl->stats.min.z = facet.vertex[0].z; + + diff_x = ABS(facet.vertex[0].x - facet.vertex[1].x); + diff_y = ABS(facet.vertex[0].y - facet.vertex[1].y); + diff_z = ABS(facet.vertex[0].z - facet.vertex[1].z); + max_diff = STL_MAX(diff_x, diff_y); + max_diff = STL_MAX(diff_z, max_diff); + stl->stats.shortest_edge = max_diff; + + first = 0; + } + /* now find the max and min values */ + stl->stats.max.x = STL_MAX(stl->stats.max.x, facet.vertex[0].x); + stl->stats.min.x = STL_MIN(stl->stats.min.x, facet.vertex[0].x); + stl->stats.max.y = STL_MAX(stl->stats.max.y, facet.vertex[0].y); + stl->stats.min.y = STL_MIN(stl->stats.min.y, facet.vertex[0].y); + stl->stats.max.z = STL_MAX(stl->stats.max.z, facet.vertex[0].z); + stl->stats.min.z = STL_MIN(stl->stats.min.z, facet.vertex[0].z); + + stl->stats.max.x = STL_MAX(stl->stats.max.x, facet.vertex[1].x); + stl->stats.min.x = STL_MIN(stl->stats.min.x, facet.vertex[1].x); + stl->stats.max.y = STL_MAX(stl->stats.max.y, facet.vertex[1].y); + stl->stats.min.y = STL_MIN(stl->stats.min.y, facet.vertex[1].y); + stl->stats.max.z = STL_MAX(stl->stats.max.z, facet.vertex[1].z); + stl->stats.min.z = STL_MIN(stl->stats.min.z, facet.vertex[1].z); + + stl->stats.max.x = STL_MAX(stl->stats.max.x, facet.vertex[2].x); + stl->stats.min.x = STL_MIN(stl->stats.min.x, facet.vertex[2].x); + stl->stats.max.y = STL_MAX(stl->stats.max.y, facet.vertex[2].y); + stl->stats.min.y = STL_MIN(stl->stats.min.y, facet.vertex[2].y); + stl->stats.max.z = STL_MAX(stl->stats.max.z, facet.vertex[2].z); + stl->stats.min.z = STL_MIN(stl->stats.min.z, facet.vertex[2].z); + } + stl->stats.size.x = stl->stats.max.x - stl->stats.min.x; + stl->stats.size.y = stl->stats.max.y - stl->stats.min.y; + stl->stats.size.z = stl->stats.max.z - stl->stats.min.z; + stl->stats.bounding_diameter = + sqrt(stl->stats.size.x * stl->stats.size.x + + stl->stats.size.y * stl->stats.size.y + + stl->stats.size.z * stl->stats.size.z); +} + + +void +stl_close(stl_file *stl) +{ + if(stl->neighbors_start != NULL) + free(stl->neighbors_start); + if(stl->facet_start != NULL) + free(stl->facet_start); + if(stl->v_indices != NULL) + free(stl->v_indices); + if(stl->v_shared != NULL) + free(stl->v_shared); +} + diff --git a/xs/src/admesh/util.c b/xs/src/admesh/util.c new file mode 100644 index 000000000..6b951e190 --- /dev/null +++ b/xs/src/admesh/util.c @@ -0,0 +1,372 @@ +/* ADMesh -- process triangulated solid meshes + * Copyright (C) 1995, 1996 Anthony D. Martin + * + * 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 + * the Free Software Foundation; either version 2, or (at your option) + * any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program; if not, write to the Free Software + * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. + * + * Questions, comments, suggestions, etc to + */ + +#include +#include +#include +#include + +#include "stl.h" + +static void stl_rotate(float *x, float *y, float angle); +static void stl_get_size(stl_file *stl); +static float get_area(stl_facet *facet); +static float get_volume(stl_file *stl); + + +void +stl_verify_neighbors(stl_file *stl) +{ + int i; + int j; + stl_edge edge_a; + stl_edge edge_b; + int neighbor; + int vnot; + + stl->stats.backwards_edges = 0; + + for(i = 0; i < stl->stats.number_of_facets; i++) + { + for(j = 0; j < 3; j++) + { + edge_a.p1 = stl->facet_start[i].vertex[j]; + edge_a.p2 = stl->facet_start[i].vertex[(j + 1) % 3]; + neighbor = stl->neighbors_start[i].neighbor[j]; + vnot = stl->neighbors_start[i].which_vertex_not[j]; + + if(neighbor == -1) + continue; /* this edge has no neighbor... Continue. */ + if(vnot < 3) + { + edge_b.p1 = stl->facet_start[neighbor].vertex[(vnot + 2) % 3]; + edge_b.p2 = stl->facet_start[neighbor].vertex[(vnot + 1) % 3]; + } + else + { + stl->stats.backwards_edges += 1; + edge_b.p1 = stl->facet_start[neighbor].vertex[(vnot + 1) % 3]; + edge_b.p2 = stl->facet_start[neighbor].vertex[(vnot + 2) % 3]; + } + if(memcmp(&edge_a, &edge_b, SIZEOF_EDGE_SORT) != 0) + { + /* These edges should match but they don't. Print results. */ + printf("edge %d of facet %d doesn't match edge %d of facet %d\n", + j, i, vnot + 1, neighbor); + stl_write_facet(stl, (char*)"first facet", i); + stl_write_facet(stl, (char*)"second facet", neighbor); + } + } + } +} + +void +stl_translate(stl_file *stl, float x, float y, float z) +{ + int i; + int j; + + for(i = 0; i < stl->stats.number_of_facets; i++) + { + for(j = 0; j < 3; j++) + { + stl->facet_start[i].vertex[j].x -= (stl->stats.min.x - x); + stl->facet_start[i].vertex[j].y -= (stl->stats.min.y - y); + stl->facet_start[i].vertex[j].z -= (stl->stats.min.z - z); + } + } + stl->stats.max.x -= (stl->stats.min.x - x); + stl->stats.max.y -= (stl->stats.min.y - y); + stl->stats.max.z -= (stl->stats.min.z - z); + stl->stats.min.x = x; + stl->stats.min.y = y; + stl->stats.min.z = z; +} + +void +stl_scale(stl_file *stl, float factor) +{ + int i; + int j; + + stl->stats.min.x *= factor; + stl->stats.min.y *= factor; + stl->stats.min.z *= factor; + stl->stats.max.x *= factor; + stl->stats.max.y *= factor; + stl->stats.max.z *= factor; + + for(i = 0; i < stl->stats.number_of_facets; i++) + { + for(j = 0; j < 3; j++) + { + stl->facet_start[i].vertex[j].x *= factor; + stl->facet_start[i].vertex[j].y *= factor; + stl->facet_start[i].vertex[j].z *= factor; + } + } +} + +static void calculate_normals(stl_file *stl) +{ + long i; + float normal[3]; + + for(i = 0; i < stl->stats.number_of_facets; i++){ + stl_calculate_normal(normal, &stl->facet_start[i]); + stl_normalize_vector(normal); + stl->facet_start[i].normal.x = normal[0]; + stl->facet_start[i].normal.y = normal[1]; + stl->facet_start[i].normal.z = normal[2]; + } +} + +void +stl_rotate_x(stl_file *stl, float angle) +{ + int i; + int j; + + for(i = 0; i < stl->stats.number_of_facets; i++) + { + for(j = 0; j < 3; j++) + { + stl_rotate(&stl->facet_start[i].vertex[j].y, + &stl->facet_start[i].vertex[j].z, angle); + } + } + stl_get_size(stl); + calculate_normals(stl); +} + +void +stl_rotate_y(stl_file *stl, float angle) +{ + int i; + int j; + + for(i = 0; i < stl->stats.number_of_facets; i++) + { + for(j = 0; j < 3; j++) + { + stl_rotate(&stl->facet_start[i].vertex[j].z, + &stl->facet_start[i].vertex[j].x, angle); + } + } + stl_get_size(stl); + calculate_normals(stl); +} + +void +stl_rotate_z(stl_file *stl, float angle) +{ + int i; + int j; + + for(i = 0; i < stl->stats.number_of_facets; i++) + { + for(j = 0; j < 3; j++) + { + stl_rotate(&stl->facet_start[i].vertex[j].x, + &stl->facet_start[i].vertex[j].y, angle); + } + } + stl_get_size(stl); + calculate_normals(stl); +} + + + +static void +stl_rotate(float *x, float *y, float angle) +{ + double r; + double theta; + double radian_angle; + + radian_angle = (angle / 180.0) * M_PI; + + r = sqrt((*x * *x) + (*y * *y)); + theta = atan2(*y, *x); + *x = r * cos(theta + radian_angle); + *y = r * sin(theta + radian_angle); +} + +static void +stl_get_size(stl_file *stl) +{ + int i; + int j; + + stl->stats.min.x = stl->facet_start[0].vertex[0].x; + stl->stats.min.y = stl->facet_start[0].vertex[0].y; + stl->stats.min.z = stl->facet_start[0].vertex[0].z; + stl->stats.max.x = stl->facet_start[0].vertex[0].x; + stl->stats.max.y = stl->facet_start[0].vertex[0].y; + stl->stats.max.z = stl->facet_start[0].vertex[0].z; + + for(i = 0; i < stl->stats.number_of_facets; i++) + { + for(j = 0; j < 3; j++) + { + stl->stats.min.x = STL_MIN(stl->stats.min.x, + stl->facet_start[i].vertex[j].x); + stl->stats.min.y = STL_MIN(stl->stats.min.y, + stl->facet_start[i].vertex[j].y); + stl->stats.min.z = STL_MIN(stl->stats.min.z, + stl->facet_start[i].vertex[j].z); + stl->stats.max.x = STL_MAX(stl->stats.max.x, + stl->facet_start[i].vertex[j].x); + stl->stats.max.y = STL_MAX(stl->stats.max.y, + stl->facet_start[i].vertex[j].y); + stl->stats.max.z = STL_MAX(stl->stats.max.z, + stl->facet_start[i].vertex[j].z); + } + } +} + +void +stl_mirror_xy(stl_file *stl) +{ + int i; + int j; + float temp_size; + + for(i = 0; i < stl->stats.number_of_facets; i++) + { + for(j = 0; j < 3; j++) + { + stl->facet_start[i].vertex[j].z *= -1.0; + } + } + temp_size = stl->stats.min.z; + stl->stats.min.z = stl->stats.max.z; + stl->stats.max.z = temp_size; + stl->stats.min.z *= -1.0; + stl->stats.max.z *= -1.0; +} + +void +stl_mirror_yz(stl_file *stl) +{ + int i; + int j; + float temp_size; + + for(i = 0; i < stl->stats.number_of_facets; i++) + { + for(j = 0; j < 3; j++) + { + stl->facet_start[i].vertex[j].x *= -1.0; + } + } + temp_size = stl->stats.min.x; + stl->stats.min.x = stl->stats.max.x; + stl->stats.max.x = temp_size; + stl->stats.min.x *= -1.0; + stl->stats.max.x *= -1.0; +} + +void +stl_mirror_xz(stl_file *stl) +{ + int i; + int j; + float temp_size; + + for(i = 0; i < stl->stats.number_of_facets; i++) + { + for(j = 0; j < 3; j++) + { + stl->facet_start[i].vertex[j].y *= -1.0; + } + } + temp_size = stl->stats.min.y; + stl->stats.min.y = stl->stats.max.y; + stl->stats.max.y = temp_size; + stl->stats.min.y *= -1.0; + stl->stats.max.y *= -1.0; +} + +static float get_volume(stl_file *stl) +{ + long i; + stl_vertex p0; + stl_vertex p; + stl_normal n; + float height; + float area; + float volume = 0.0; + + /* Choose a point, any point as the reference */ + p0.x = stl->facet_start[0].vertex[0].x; + p0.y = stl->facet_start[0].vertex[0].y; + p0.z = stl->facet_start[0].vertex[0].z; + + for(i = 0; i < stl->stats.number_of_facets; i++){ + p.x = stl->facet_start[i].vertex[0].x - p0.x; + p.y = stl->facet_start[i].vertex[0].y - p0.y; + p.z = stl->facet_start[i].vertex[0].z - p0.z; + /* Do dot product to get distance from point to plane */ + n = stl->facet_start[i].normal; + height = (n.x * p.x) + (n.y * p.y) + (n.z * p.z); + area = get_area(&stl->facet_start[i]); + volume += (area * height) / 3.0; + } + return volume; +} + +void stl_calculate_volume(stl_file *stl) +{ + stl->stats.volume = get_volume(stl); + if(stl->stats.volume < 0.0){ + stl_reverse_all_facets(stl); + stl->stats.volume = -stl->stats.volume; + } +} + +static float get_area(stl_facet *facet) +{ + float cross[3][3]; + float sum[3]; + float n[3]; + float area; + int i; + + for(i = 0; i < 3; i++){ + cross[i][0]=((facet->vertex[i].y * facet->vertex[(i + 1) % 3].z) - + (facet->vertex[i].z * facet->vertex[(i + 1) % 3].y)); + cross[i][1]=((facet->vertex[i].z * facet->vertex[(i + 1) % 3].x) - + (facet->vertex[i].x * facet->vertex[(i + 1) % 3].z)); + cross[i][2]=((facet->vertex[i].x * facet->vertex[(i + 1) % 3].y) - + (facet->vertex[i].y * facet->vertex[(i + 1) % 3].x)); + } + + sum[0] = cross[0][0] + cross[1][0] + cross[2][0]; + sum[1] = cross[0][1] + cross[1][1] + cross[2][1]; + sum[2] = cross[0][2] + cross[1][2] + cross[2][2]; + + /* This should already be done. But just in case, let's do it again */ + stl_calculate_normal(n, facet); + stl_normalize_vector(n); + + area = 0.5 * (n[0] * sum[0] + n[1] * sum[1] + n[2] * sum[2]); + return area; +} diff --git a/xs/src/myinit.h b/xs/src/myinit.h index f2104209e..570908808 100644 --- a/xs/src/myinit.h +++ b/xs/src/myinit.h @@ -18,4 +18,6 @@ ZTable::ZTable(std::vector* ztable) : { } +#include + #endif diff --git a/xs/xsp/TriangleMesh.xsp b/xs/xsp/TriangleMesh.xsp index be6b31df6..7ba9724f1 100644 --- a/xs/xsp/TriangleMesh.xsp +++ b/xs/xsp/TriangleMesh.xsp @@ -4,6 +4,8 @@ %{ PROTOTYPES: DISABLE +#include + std::string hello_world() CODE: @@ -11,4 +13,25 @@ hello_world() OUTPUT: RETVAL +float +stl_volume(filename) + char* filename; + CODE: + stl_file stl_in; + stl_open(&stl_in, filename); + + stl_check_facets_exact(&stl_in); + stl_in.stats.facets_w_1_bad_edge = (stl_in.stats.connected_facets_2_edge - stl_in.stats.connected_facets_3_edge); + stl_in.stats.facets_w_2_bad_edge = (stl_in.stats.connected_facets_1_edge - stl_in.stats.connected_facets_2_edge); + stl_in.stats.facets_w_3_bad_edge = (stl_in.stats.number_of_facets - stl_in.stats.connected_facets_1_edge); + + stl_fix_normal_directions(&stl_in); + stl_fix_normal_values(&stl_in); + + stl_calculate_volume(&stl_in); + RETVAL = stl_in.stats.volume; + stl_close(&stl_in); + OUTPUT: + RETVAL + %}