Commit 3909dfd7 authored by Sander Mathijs van Veen's avatar Sander Mathijs van Veen

Merge branch 'master' of ssh://vo20.nl/git/uva

parents ef1707f2 9499f863
CC=gcc
WARNING_FLAGS=-Wall -Wextra -Werror-implicit-function-declaration -Wshadow -Wstrict-prototypes -pedantic-errors
CFLAGS=-g -O2 -std=c99 $(WARNING_FLAGS)
LDFLAGS=-g -lGL -lglut -lGLU
.c.o:
$(CC) -c $(CFLAGS) $<
all: main
main: main.o marching_tetrahedra.o volume.o
$(CC) $(LDFLAGS) -o main main.o marching_tetrahedra.o volume.o
clean:
rm -f *.o
rm -f main
volume.o : v3math.h volume.h
volume.o : volume.h volume.c
v3math.o : v3math.h
marching_tetrahedra.o : marching_tetrahedra.h marching_tetrahedra.c
marching_tetrahedra.o : v3math.h volume.h marching_tetrahedra.h
main.o : marching_tetrahedra.h v3math.h volume.h main.c
This diff is collapsed.
/* Computer Graphics, Assignment, Volume rendering with cubes/points/isosurface
*
* Student name ....
* Student email ...
* Collegekaart ....
* Date ............
* Comments ........
*
* (always fill in these fields before submitting!!)
*/
#include <stdlib.h>
#include <stdio.h>
#include <assert.h>
#include "marching_tetrahedra.h"
/* Compute a linearly interpolated position where an isosurface cuts
an edge between two vertices (p1 and p2), each with their own
scalar value (v1 and v2) */
static vec3
interpolate_points(unsigned char isovalue, vec3 p1, vec3 p2,
unsigned char v1, unsigned char v2)
{
/* Initially, simply return the midpoint between p1 and p2.
So no real interpolation is done yet */
return v3_add(v3_multiply(p1, 0.5), v3_multiply(p2, 0.5));
}
/* Using the given iso-value generate triangles for the tetrahedron
defined by corner vertices v0, v1, v2, v3 of cell c.
Store the resulting triangles in the "triangles" array.
Return the number of triangles created (either 0, 1, or 2).
Note: the output array "triangles" should have space for at least
2 triangles.
*/
static int
generate_tetrahedron_triangles(triangle *triangles, unsigned char isovalue,
cell c, int v0, int v1, int v2, int v3)
{
unsigned char *v = c.value;
vec3 *p = triangles[0].p;
int bitsystem = 0;
if(c.value[v0] > isovalue) bitsystem += 1;
if(c.value[v1] > isovalue) bitsystem += 2;
if(c.value[v2] > isovalue) bitsystem += 4;
if(c.value[v3] > isovalue) bitsystem += 8;
switch( bitsystem )
{
case 1: case 14: // 0001
p[0] = interpolate_points(isovalue, c.p[v0], c.p[v1], v[v0], v[v1]);
p[1] = interpolate_points(isovalue, c.p[v0], c.p[v2], v[v0], v[v2]);
p[2] = interpolate_points(isovalue, c.p[v0], c.p[v3], v[v0], v[v3]);
return 1;
case 2: case 13: // 0010
p[0] = interpolate_points(isovalue, c.p[v1], c.p[v0], v[v1], v[v0]);
p[1] = interpolate_points(isovalue, c.p[v1], c.p[v3], v[v1], v[v3]);
p[2] = interpolate_points(isovalue, c.p[v1], c.p[v2], v[v1], v[v2]);
return 1;
case 4: case 11: // 0100
p[0] = interpolate_points(isovalue, c.p[v2], c.p[v0], v[v2], v[v0]);
p[1] = interpolate_points(isovalue, c.p[v2], c.p[v1], v[v2], v[v1]);
p[2] = interpolate_points(isovalue, c.p[v2], c.p[v3], v[v2], v[v3]);
return 1;
case 8: case 7: // 1000
p[0] = interpolate_points(isovalue, c.p[v3], c.p[v0], v[v3], v[v0]);
p[1] = interpolate_points(isovalue, c.p[v3], c.p[v2], v[v3], v[v2]);
p[2] = interpolate_points(isovalue, c.p[v3], c.p[v1], v[v3], v[v1]);
return 1;
case 5: case 10: // 0101
p[0] = interpolate_points(isovalue, c.p[v0], c.p[v1], v[v0], v[v1]);
p[1] = interpolate_points(isovalue, c.p[v2], c.p[v3], v[v2], v[v3]);
p[2] = interpolate_points(isovalue, c.p[v0], c.p[v3], v[v0], v[v3]);
p = triangles[1].p;
p[0] = interpolate_points(isovalue, c.p[v0], c.p[v1], v[v0], v[v1]);
p[1] = interpolate_points(isovalue, c.p[v1], c.p[v2], v[v1], v[v2]);
p[2] = interpolate_points(isovalue, c.p[v2], c.p[v3], v[v2], v[v3]);
return 2;
case 3: case 12: // 0011
p[0] = interpolate_points(isovalue, c.p[v0], c.p[v3], v[v0], v[v3]);
p[1] = interpolate_points(isovalue, c.p[v0], c.p[v2], v[v0], v[v2]);
p[2] = interpolate_points(isovalue, c.p[v1], c.p[v3], v[v1], v[v3]);
p = triangles[1].p;
p[0] = interpolate_points(isovalue, c.p[v1], c.p[v3], v[v1], v[v3]);
p[2] = interpolate_points(isovalue, c.p[v1], c.p[v2], v[v1], v[v2]);
p[1] = interpolate_points(isovalue, c.p[v0], c.p[v2], v[v0], v[v2]);
return 2;
case 6: case 9: // 0110
p[0] = interpolate_points(isovalue, c.p[v0], c.p[v1], v[v0], v[v1]);
p[2] = interpolate_points(isovalue, c.p[v1], c.p[v3], v[v1], v[v3]);
p[1] = interpolate_points(isovalue, c.p[v2], c.p[v3], v[v2], v[v3]);
p = triangles[1].p;
p[0] = interpolate_points(isovalue, c.p[v0], c.p[v1], v[v0], v[v1]);
p[1] = interpolate_points(isovalue, c.p[v0], c.p[v2], v[v0], v[v2]);
p[2] = interpolate_points(isovalue, c.p[v2], c.p[v3], v[v2], v[v3]);
return 2;
default: // 0000
return 0;
}
}
/* Generate triangles for a single cell for the given iso-value. This function
should produce at most 6 * 2 triangles (for which the "triangles" array
should have enough space).
Use calls to generate_tetrahedron_triangles().
Return the number of triangles produced
*/
int
generate_cell_triangles(triangle *triangles, cell c, unsigned char isovalue)
{
const int points[24] = {0,1,3,7,0,2,6,7,0,1,5,7,0,2,3,7,0,4,5,7,0,4,6,7};
int i, tri_count = 0;
for( i = 0; i < 21; i += 4 )
{
tri_count += generate_tetrahedron_triangles(triangles + tri_count,
isovalue, c, points[i], points[i+1], points[i+2], points[i+3]);
}
return tri_count;
}
#ifndef MARCHING_TETRAHEDRA_H
#define MARCHING_TETRAHEDRA_H
#include "v3math.h"
#include "volume.h"
typedef struct
{
vec3 p[3];
vec3 n[3];
}
triangle;
int generate_cell_triangles(triangle *triangles, cell c, unsigned char isovalue);
#endif
#ifndef V3MATH_H
#define V3MATH_H
#include <math.h>
typedef struct
{
float x, y, z;
}
vec3;
// Create a new 3-vector of floats
static inline vec3
v3_create(float x, float y, float z)
{
vec3 res;
res.x = x;
res.y = y;
res.z = z;
return res;
}
// Return -a
static inline vec3
v3_negate(vec3 a)
{
vec3 res;
res.x = - a.x;
res.y = - a.y;
res.z = - a.z;
return res;
}
// Return a + b
static inline vec3
v3_add(vec3 a,vec3 b)
{
vec3 res;
res.x = a.x+b.x;
res.y = a.y+b.y;
res.z = a.z+b.z;
return res;
}
// Return a - b
static inline vec3
v3_subtract(vec3 a, vec3 b)
{
vec3 res;
res.x = a.x-b.x;
res.y = a.y-b.y;
res.z = a.z-b.z;
return res;
}
// Return a / |a|
static inline vec3
v3_normalize(vec3 a)
{
vec3 res;
double l = sqrt(a.x*a.x + a.y*a.y + a.z*a.z);
res = a;
res.x /= l;
res.y /= l;
res.z /= l;
return res;
}
// Return a ^ b
static inline vec3
v3_crossprod(vec3 a, vec3 b)
{
vec3 res;
res.x = a.y*b.z - a.z*b.y;
res.y = a.z*b.x - a.x*b.z;
res.z = a.x*b.y - a.y*b.x;
return res;
}
// Return a * b
static inline float
v3_dotprod(vec3 a, vec3 b)
{
return a.x*b.x + a.y*b.y + a.z*b.z;
}
// Return a*s
static inline vec3
v3_multiply(vec3 a, float s)
{
vec3 res;
res.x = a.x*s;
res.y = a.y*s;
res.z = a.z*s;
return res;
}
// Return |a|
static inline float
v3_length(vec3 a)
{
return sqrt(a.x*a.x + a.y*a.y + a.z*a.z);
}
// Return the i-th component of a, i.e. for i=0
// this returns a.x
static inline float
v3_component(vec3 a, int i)
{
if (i == 0)
return a.x;
else if (i == 1)
return a.y;
else
return a.z;
}
// Set the i-th component of a to v
static inline vec3
v3_set_component(vec3 a, int i, float v)
{
vec3 res = a;
if (i == 0)
res.x = v;
else if (i == 1)
res.y = v;
else
res.z = v;
return res;
}
#endif
/* Computer Graphics, Assignment, Volume rendering with cubes/points/isosurface
*
* Student name ....
* Student email ...
* Collegekaart ....
* Date ............
* Comments ........
*
* (always fill in these fields before submitting!!)
*/
#include <stdio.h>
#include <stdlib.h>
#include "volume.h"
/* The voxels of the volume dataset, stored as a one-dimensional array */
unsigned char *volume;
/* The dimensions of the volume dataset */
int nx, ny, nz;
/* The size of a voxel */
float sizex, sizey, sizez;
/* Utility function to convert the index of a voxel
into an index in the volume array above */
int
voxel2idx(int i, int j, int k)
{
return (k*ny + j)*nx + i;
}
/* Extract a cell from the volume, so that datapoint 0 of the
cell corresponds to voxel (i, j, k), datapoint 1 to voxel (i+1, j, k),
etc. See Figure 3 in the assignment. */
cell
get_cell(int i, int j, int k)
{
cell c;
c.p[0] = v3_create(i, j, k);
c.p[1] = v3_create(i+1, j, k);
c.p[2] = v3_create(i, j+1, k);
c.p[3] = v3_create(i+1, j+1, k);
c.p[4] = v3_create(i, j, k+1);
c.p[5] = v3_create(i+1, j, k+1);
c.p[6] = v3_create(i, j+1, k+1);
c.p[7] = v3_create(i+1, j+1, k+1);
c.value[0] = volume[voxel2idx(i, j, k)];
c.value[1] = volume[voxel2idx(i+1, j, k)];
c.value[2] = volume[voxel2idx(i, j+1, k)];
c.value[3] = volume[voxel2idx(i+1, j+1, k)];
c.value[4] = volume[voxel2idx(i, j, k+1)];
c.value[5] = volume[voxel2idx(i+1, j, k+1)];
c.value[6] = volume[voxel2idx(i, j+1, k+1)];
c.value[7] = volume[voxel2idx(i+1, j+1, k+1)];
return c;
}
/* Utility function to read a volume dataset from a VTK file.
This will store the data in the "volume" array and update the dimension
and size values. */
void
read_volume(const char *fname)
{
FILE *f;
char s[256];
int nvoxels;
printf("Reading %s\n", fname);
f = fopen(fname, "rb");
if (!f)
{
fprintf(stderr, "read_volume(): Could not open file '%s' for reading!\n", fname);
exit(-1);
}
// header line
fgets(s, 255, f);
// comment line
fgets(s, 255, f);
// BINARY
fgets(s, 255, f);
// DATASET STRUCTURED_POINTS
fgets(s, 255, f);
// DIMENSIONS %d %d %d
fscanf(f, "%s %d %d %d\n", s, &nx, &ny, &nz);
printf("%d x %d x %d voxels\n", nx, ny, nz);
// ASPECT_RATIO/SPACING %f %f %f
fscanf(f, "%s %f %f %f\n", s, &sizex, &sizey, &sizez);
printf("voxel sizes: %.3f, %.3f, %.3f\n", sizex, sizey, sizez);
// ORIGIN ...
fgets(s, 255, f);
// POINT_DATA ...
fgets(s, 255, f);
// SCALARS ...
fgets(s, 255, f);
// LOOKUP_TABLE ...
fgets(s, 255, f);
// allocate memory to hold the volume data and read it from file
nvoxels = nx * ny * nz;
volume = (unsigned char*)malloc(nvoxels);
if (fread(volume, 1, nvoxels, f) < (size_t)nvoxels)
{
printf("WARNING: not all data could be read!\n");
}
fclose(f);
}
#ifndef VOLUME_H
#define VOLUME_H
#include "v3math.h"
/* A cell in the volume dataset, consisting of 8 neighbouring datapoints;
p[] contains the corner positions, value[] the associated scalar values */
typedef struct
{
vec3 p[8];
vec3 n[8]; // Note: we experiment with normals in the solution
unsigned char value[8];
}
cell;
/* The data points in the volume dataset, stored as a one-dimensional array */
extern unsigned char *volume;
/* The dimensions of the volume dataset in number of voxels in each
dimension*/
extern int nx, ny, nz;
/* The size of a voxel for each dimension */
extern float sizex, sizey, sizez;
/* Utility function to convert the index of a datapoint
into an index in the volume array above */
int voxel2idx(int i, int j, int k);
/* Extract a cell from the volume, so that datapoint 0 of the
cell corresponds to voxel (i, j, k) */
cell get_cell(int i, int j, int k);
/* Utility function to read a volume dataset from a file.
This will store the data in the "volume" array and update the dimension
and size values. */
void read_volume(const char *fname);
#endif
PROG = ass4
PROG = philosophers
CC = gcc
CFLAGS = -std=c99 -pedantic -Wall -Wextra -D_GNU_SOURCE # -D_POSIX_SOURCE
LFLAGS = -lpthread
OFILES = main.o
CFLAGS = -ansi -Wall -Wextra -O3 -pthread
OFILES = host.o
RM = rm -f
ifdef DEBUG
CFLAGS +=-ggdb
else
CFLAGS +=-O3
endif
$(PROG): $(OFILES)
$(CC) $(CFLAGS) $(LFLAGS) -o $@ $(OFILES)
$(CC) -o $(PROG) $(OFILES) $(CFLAGS)
%.o: %.c
$(CC) $(CFLAGS) -o $@ -c $<
$(CC) -c $< $(CFLAGS)
clean:
$(RM) $(OFILES) $(PROG)
$(RM) *.o $(PROG)
Markdown is supported
0%
or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment