commit 616515bb47336481bcbf33bb655537dacddbdaba
parent e3e1b6f775c011a498afa32e2a401cca3341c3e3
Author: Vincent Forest <vincent.forest@meso-star.com>
Date: Mon, 28 Mar 2016 16:37:04 +0200
Add the internal helical pipe mesh
Diffstat:
6 files changed, 327 insertions(+), 39 deletions(-)
diff --git a/src/schiff.c b/src/schiff.c
@@ -278,9 +278,6 @@ main(int argc, char** argv)
int err = 0;
res_T res;
- schiff_mesh_init_helicoid(NULL, NULL);
- exit(0);
-
res = schiff_args_init(&args, argc, argv);
if(res != RES_OK) goto error;
if(!args.ngeoms_dump && !args.wavelengths) goto exit;
diff --git a/src/schiff_args.c b/src/schiff_args.c
@@ -117,6 +117,12 @@ geometry_release(struct schiff_geometry* geom)
case SCHIFF_CYLINDER_AS_SPHERE:
param_release(&geom->data.cylinder.radius);
break;
+ case SCHIFF_HELICAL_PIPE:
+ param_release(&geom->data.helical_pipe.pitch);
+ param_release(&geom->data.helical_pipe.height);
+ param_release(&geom->data.helical_pipe.helicoid_radius);
+ param_release(&geom->data.helical_pipe.circle_radius);
+ break;
case SCHIFF_SPHERE:
param_release(&geom->data.sphere.radius);
break;
@@ -859,6 +865,30 @@ parse_yaml_cylinder
}
static res_T
+parse_yaml_helical_pipe
+ (const char* filename,
+ yaml_document_t* doc,
+ const yaml_node_t* node,
+ struct schiff_geometry* geom,
+ double* geom_proba)
+{
+ ASSERT(filename && doc && node && geom && geom_proba);
+ geom->type = SCHIFF_HELICAL_PIPE;
+ geom->data.helical_pipe.nslices_helicoid = 32;
+ geom->data.helical_pipe.nslices_circle = 16;
+ geom->data.helical_pipe.pitch.data.constant = 1;
+ geom->data.helical_pipe.height.data.constant = 0.6;
+ geom->data.helical_pipe.helicoid_radius.data.constant = 3;
+ geom->data.helical_pipe.circle_radius.data.constant = 1;
+ geom->data.helical_pipe.pitch.distribution = SCHIFF_PARAM_CONSTANT;
+ geom->data.helical_pipe.height.distribution = SCHIFF_PARAM_CONSTANT;
+ geom->data.helical_pipe.helicoid_radius.distribution = SCHIFF_PARAM_CONSTANT;
+ geom->data.helical_pipe.circle_radius.distribution = SCHIFF_PARAM_CONSTANT;
+ *geom_proba = 1;
+ return RES_OK;
+}
+
+static res_T
parse_yaml_sphere
(const char* filename,
yaml_document_t* doc,
@@ -1070,6 +1100,8 @@ parse_yaml_geom_distrib
res = parse_yaml_cylinder(filename, doc, val, geom, proba);
} else if(!strcmp((char*)key->data.scalar.value, "sphere")) {
res = parse_yaml_sphere(filename, doc, val, geom, proba);
+ } else if(!strcmp((char*)key->data.scalar.value, "helical_pipe")) {
+ res = parse_yaml_helical_pipe(filename, doc, val, geom, proba);
} else if(!strcmp((char*)key->data.scalar.value, "supershape")) {
res = parse_yaml_supershape(filename, doc, val, geom, proba);
} else {
diff --git a/src/schiff_geometry.c b/src/schiff_geometry.c
@@ -49,6 +49,14 @@ struct cylinder_context {
float height;
};
+/* 3D Context of an helical pipe */
+struct helical_pipe_context {
+ double helicoid_radius;
+ double circle_radius;
+ double pitch;
+ double height;
+};
+
/* 3D Context of the supershape geometry */
struct supershape_context {
double formulas[2][6];
@@ -69,6 +77,7 @@ struct mesh_context {
union {
struct ellipsoid_context ellipsoid;
struct cylinder_context cylinder;
+ struct helical_pipe_context helical_pipe;
struct sphere_context sphere;
struct supershape_context supershape;
} data;
@@ -193,6 +202,14 @@ geometry_get_indices(const unsigned itri, unsigned ids[3], void* ctx)
}
static void
+geometry_get_position(const unsigned ivert, float vertex[3], void* ctx)
+{
+ struct mesh_context* mesh_ctx = ctx;
+ ASSERT(ctx);
+ schiff_mesh_get_cartesian_position(mesh_ctx->mesh, ivert, vertex);
+}
+
+static void
ellipsoid_get_position(const unsigned ivert, float vertex[3], void* ctx)
{
struct mesh_context* mesh_ctx = ctx;
@@ -368,6 +385,52 @@ geometry_sample_cylinder
}
static res_T
+geometry_sample_helical_pipe
+ (struct ssp_rng* rng,
+ const struct schiff_geometry* geom,
+ struct schiff_mesh* mesh,
+ struct s3d_shape* shape)
+{
+ const struct schiff_helical_pipe* helical_pipe;
+ struct mesh_context mesh_ctx;
+ struct s3d_vertex_data attrib;
+ size_t nprims;
+ size_t nverts;
+ double pitch;
+ double height;
+ double hradius;
+ double cradius;
+ res_T res = RES_OK;
+ ASSERT(geom && geom->type == SCHIFF_HELICAL_PIPE);
+
+ helical_pipe = &geom->data.helical_pipe;
+
+ /* Sample the helical pipe attributes */
+ hradius = eval_param(&helical_pipe->helicoid_radius, rng);
+ cradius = eval_param(&helical_pipe->circle_radius, rng);
+ pitch = eval_param(&helical_pipe->pitch, rng);
+ height = eval_param(&helical_pipe->height, rng);
+
+ /* Setup the mesh vertices */
+ res = schiff_mesh_setup_helical_pipe(mesh, pitch, height, hradius, cradius,
+ helical_pipe->nslices_helicoid, helical_pipe->nslices_circle);
+ if(res != RES_OK) return res;
+
+ mesh_ctx.mesh = mesh;
+ mesh_ctx.type = SCHIFF_HELICAL_PIPE;
+
+ attrib.usage = S3D_POSITION;
+ attrib.type = S3D_FLOAT3;
+ attrib.get = geometry_get_position;
+
+ nverts = darray_float_size_get(&mesh->vertices.cartesian) / 3/*#coords*/;
+ nprims = darray_uint_size_get(&mesh->indices) / 3/*#indices per prim*/;
+
+ return s3d_mesh_setup_indexed_vertices(shape, (unsigned)nprims,
+ geometry_get_indices, (unsigned)nverts, &attrib, 1, &mesh_ctx);
+}
+
+static res_T
geometry_sample_sphere
(struct ssp_rng* rng,
const struct schiff_geometry* geom,
@@ -378,7 +441,7 @@ geometry_sample_sphere
struct mesh_context mesh_ctx;
struct s3d_vertex_data attrib;
size_t nverts, nprims;
- ASSERT(rng && geom && mesh && shape);
+ ASSERT(geom && mesh && shape);
ASSERT(geom->type == SCHIFF_SPHERE);
sphere = &geom->data.sphere;
@@ -409,7 +472,7 @@ geometry_sample_supershape
struct s3d_vertex_data attrib;
size_t nverts, nprims;
int iform, iattr;
- ASSERT(rng && geom && mesh && shape);
+ ASSERT(geom && mesh && shape);
ASSERT(geom->type == SCHIFF_SUPERSHAPE);
sshape = &geom->data.supershape;
@@ -441,7 +504,7 @@ geometry_sample
{
struct geometry_distribution_context* distrib = ctx;
const struct schiff_geometry* geom;
- const struct schiff_mesh* mesh;
+ struct schiff_mesh* mesh;
size_t isamp;
res_T res = RES_OK;
ASSERT(rng && shape && ctx);
@@ -459,6 +522,9 @@ geometry_sample
case SCHIFF_CYLINDER_AS_SPHERE:
res = geometry_sample_cylinder(rng, geom, mesh, shape);
break;
+ case SCHIFF_HELICAL_PIPE:
+ res = geometry_sample_helical_pipe(rng, geom, mesh, shape);
+ break;
case SCHIFF_SPHERE:
res = geometry_sample_sphere(rng, geom, mesh, shape);
break;
@@ -556,6 +622,12 @@ schiff_geometry_distribution_init
res = schiff_mesh_init_cylinder(&mem_default_allocator, &ctx->meshes[i],
geoms[i].data.cylinder.nslices);
break;
+ case SCHIFF_HELICAL_PIPE:
+ res = schiff_mesh_init_helical_pipe
+ (&mem_default_allocator, &ctx->meshes[i],
+ geoms[i].data.helical_pipe.nslices_helicoid,
+ geoms[i].data.helical_pipe.nslices_circle);
+ break;
case SCHIFF_SPHERE:
res = schiff_mesh_init_sphere(&mem_default_allocator, &ctx->meshes[i],
geoms[i].data.sphere.nslices);
diff --git a/src/schiff_geometry.h b/src/schiff_geometry.h
@@ -43,6 +43,7 @@ enum schiff_geometry_type {
SCHIFF_ELLIPSOID_AS_SPHERE,/* The ellipsoid volume is controled by a sphere */
SCHIFF_CYLINDER,
SCHIFF_CYLINDER_AS_SPHERE,/* The cylinder volume is controled by a sphere */
+ SCHIFF_HELICAL_PIPE,
SCHIFF_SPHERE,
SCHIFF_SUPERSHAPE,
SCHIFF_NONE
@@ -72,6 +73,25 @@ struct schiff_sphere {
unsigned nslices;
};
+struct schiff_helical_pipe {
+ struct schiff_param pitch; /* Elevation distance of a full revolution */
+ struct schiff_param height; /* Total heigh of the helical pipe */
+ struct schiff_param helicoid_radius; /* Radius of the helicoid */
+ struct schiff_param circle_radius; /* Radius of the meridian circle */
+
+ unsigned nslices_helicoid; /* # Discrete steps of the helicoid */
+ unsigned nslices_circle; /* # Discrete steps along 2PI */
+};
+
+#define SCHIFF_HELICAL_PIPE_DEFAULT__ \
+ {SCHIFF_PARAM_DEFAULT__, \
+ SCHIFF_PARAM_DEFAULT__, \
+ SCHIFF_PARAM_DEFAULT__, \
+ SCHIFF_PARAM_DEFAULT__, \
+ 128, 64}
+static const struct schiff_helical_pipe SCHIFF_HELICAL_PIPE_DEFAULT =
+ SCHIFF_HELICAL_PIPE_DEFAULT__;
+
#define SCHIFF_SPHERE_DEFAULT__ {SCHIFF_PARAM_DEFAULT__, 64}
static const struct schiff_sphere SCHIFF_SPHERE_DEFAULT =
SCHIFF_SPHERE_DEFAULT__;
@@ -116,6 +136,7 @@ struct schiff_geometry {
union {
struct schiff_ellipsoid ellipsoid;
struct schiff_cylinder cylinder;
+ struct schiff_helical_pipe helical_pipe;
struct schiff_sphere sphere;
struct schiff_supershape supershape;
} data;
diff --git a/src/schiff_mesh.c b/src/schiff_mesh.c
@@ -17,6 +17,7 @@
#include <rsys/float2.h>
#include <rsys/float3.h>
+#include <rsys/stretchy_array.h>
/*******************************************************************************
* Helper functions
@@ -247,7 +248,7 @@ schiff_mesh_init_cylinder
darray_uint_init(allocator, &cylinder->indices);
nverts = nsteps*2/* #contour verts */ + 2/* #polar verts */;
- ntris = nsteps*2/* #contour tris */ + 2*nsteps/* #caop tris */;
+ ntris = nsteps*2/* #contour tris */ + 2*nsteps/* #cap tris */;
res = darray_float_resize(&cylinder->vertices.cartesian, nverts*3/*#coords*/);
if(res != RES_OK) goto error;
@@ -305,44 +306,193 @@ error:
goto exit;
}
-#include <rsys/float33.h>
res_T
-schiff_mesh_init_helicoid
+schiff_mesh_init_helical_pipe
(struct mem_allocator* allocator,
- struct schiff_mesh* mesh)
+ struct schiff_mesh* helical_pipe,
+ const unsigned nsteps_helicoid,
+ const unsigned nsteps_circle)
+
+{
+ unsigned* indices;
+ unsigned* bottom_cap;
+ unsigned* top_cap;
+ size_t ihelicoid, icircle;
+ size_t nverts, ntris;
+ size_t ibottom_cap_index, itop_cap_index;
+ size_t ibottom_cap_vertex, itop_cap_vertex;
+ size_t ilast_meridian_vertex;
+ res_T res;
+ ASSERT(allocator && helical_pipe && nsteps_helicoid>=2 && nsteps_circle>=4);
+
+ helical_pipe->coordinates = SCHIFF_CARTESIAN;
+ darray_float_init(allocator, &helical_pipe->vertices.cartesian);
+ darray_uint_init(allocator, &helical_pipe->indices);
+
+ nverts = (nsteps_helicoid+1)*nsteps_circle/*#contour*/ + 2/*#cap*/;
+ ntris = (nsteps_helicoid)*nsteps_circle*2/*#contour*/ + 2*nsteps_circle/*#cap*/;
+
+ res = darray_float_resize(&helical_pipe->vertices.cartesian, nverts*3/*#coords*/);
+ if(res != RES_OK) goto error;
+ res = darray_uint_resize(&helical_pipe->indices, ntris*3/*#indices per tri*/);
+ if(res != RES_OK) goto error;
+
+ indices = darray_uint_data_get(&helical_pipe->indices);
+
+ /* Setup the indices toward the contour primitives */
+ FOR_EACH(ihelicoid, 0, nsteps_helicoid) {
+ /* Offset toward the first index of the current meridian */
+ const size_t offset = ihelicoid * nsteps_circle;
+
+ FOR_EACH(icircle, 0, nsteps_circle) {
+ size_t prim_ids[4];
+ const size_t id = icircle + offset;
+ unsigned* iprim = indices + id*3/*#ids per prim*/*2/*#prims per step*/;
+
+ /* Indices of the quad */
+ prim_ids[0] = icircle + offset;
+ prim_ids[1] = icircle + nsteps_circle + offset;
+ prim_ids[2] = (icircle + 1) % nsteps_circle + nsteps_circle + offset;
+ prim_ids[3] = (icircle + 1) % nsteps_circle + offset;
+ /* First triangle of the quad */
+ iprim[0] = (unsigned)prim_ids[0];
+ iprim[1] = (unsigned)prim_ids[3];
+ iprim[2] = (unsigned)prim_ids[1];
+ /* Second triangule of the quad */
+ iprim[3] = (unsigned)prim_ids[1];
+ iprim[4] = (unsigned)prim_ids[3];
+ iprim[5] = (unsigned)prim_ids[2];
+ }
+ }
+
+ /* Define the index of the cap vertices */
+ ibottom_cap_vertex = nverts - 2; /* Index toward the bottom cap vertex */
+ itop_cap_vertex = nverts - 1; /* Index toward the top cap vertex */
+
+ /* Index of the first vertex of the last meridian circle */
+ ilast_meridian_vertex = nsteps_helicoid * nsteps_circle;
+
+ ibottom_cap_index = /* Index of the 1st bottom cap index */
+ nsteps_helicoid
+ * nsteps_circle
+ * 2 /*#prims per step*/
+ * 3 /*#ids per prim*/;
+ itop_cap_index = /* Index of the 1st top cap index */
+ ibottom_cap_index
+ + nsteps_circle * 3/*# ids per primitive */;
+
+ bottom_cap = indices + ibottom_cap_index ;
+ top_cap = indices + itop_cap_index;
+
+ FOR_EACH(icircle, 0, nsteps_circle) {
+ const size_t id = icircle * 3/*#ids per prim*/;
+
+ bottom_cap[id + 0] = (unsigned)icircle;
+ bottom_cap[id + 1] = (unsigned)ibottom_cap_vertex;
+ bottom_cap[id + 2] = (unsigned)((icircle + 1) % nsteps_circle);
+
+ top_cap[id + 0] = (unsigned)((icircle+1)%nsteps_circle + ilast_meridian_vertex);
+ top_cap[id + 1] = (unsigned)itop_cap_vertex;
+ top_cap[id + 2] = (unsigned)(icircle + ilast_meridian_vertex);
+ }
+
+exit:
+ return res;
+error:
+ schiff_mesh_release(helical_pipe);
+ goto exit;
+}
+
+res_T
+schiff_mesh_setup_helical_pipe
+ (struct schiff_mesh* mesh,
+ const double pitch,
+ const double height,
+ const double hradius, /* Helicoid radius */
+ const double cradius, /* Circle radius */
+ const unsigned nsteps_helicoid,
+ const unsigned nsteps_circle)
{
- const double height = 0.7;
- const double height_total = 3;
- const double c = height / 2*PI;
- const double phi_max = height_total * 2*PI / height;
- const size_t circle_nsteps = 64;
- const double circle_step = 2*PI / (double)circle_nsteps;
- const double circle_radius = 0.5;
- const size_t helicoid_nsteps = 256;
- const double helicoid_step = phi_max / (double)helicoid_nsteps;
- const double helicoid_radius = 1.0;
- size_t icircle, ihelicoid;
- /*ASSERT(allocator && mesh);*/
-
- FOR_EACH(ihelicoid, 0, helicoid_nsteps) {
- const double phi = (double)ihelicoid * helicoid_step;
+ double* meridian = NULL;
+ float* vertices;
+ double c;
+ double phi_max;
+ double step_helicoid, step_circle;
+ double rcp_sqrt_hradius2_add_c2; /* Precomputed value */
+ size_t ihelicoid, icircle;
+ size_t icoord, icoord_top, icoord_bottom;
+ res_T res = RES_OK;
+
+ ASSERT(mesh && pitch > 0 && height > 0 && hradius > 0 && cradius > 0);
+ ASSERT(mesh->coordinates == SCHIFF_CARTESIAN);
+ ASSERT(nsteps_helicoid * nsteps_circle * 2 * 3 + nsteps_circle*3*2/*cap*/
+ == darray_uint_size_get(&mesh->indices));
+ ASSERT((nsteps_helicoid+1) * nsteps_circle * 3 + 2*3/*cap*/
+ == darray_float_size_get(&mesh->vertices.cartesian));
+
+ if(!sa_add(meridian, nsteps_circle * 3/*#coords per vertex*/)) {
+ res = RES_MEM_ERR;
+ goto error;
+ }
+
+ c = pitch / 2*PI;
+ phi_max = height * 2*PI / pitch;
+ step_helicoid = phi_max / (double)nsteps_helicoid;
+ step_circle = 2*PI / (double)nsteps_circle;
+ rcp_sqrt_hradius2_add_c2 = 1.0 / sqrt(hradius*hradius + c*c);
+
+ /* Compute the meridian positions */
+ FOR_EACH(icircle, 0, nsteps_circle) {
+ const double theta = (double)icircle * step_circle;
+ const double cos_theta = cos(theta);
+ const double sin_theta = sin(theta);
+ const size_t id = icircle *3;
+
+ meridian[id + 0] = cradius * cos_theta + hradius;
+ meridian[id + 1] = -cradius * c * rcp_sqrt_hradius2_add_c2 * sin_theta;
+ meridian[id + 2] = cradius*hradius * rcp_sqrt_hradius2_add_c2 * sin_theta;
+ }
+
+ vertices = darray_float_data_get(&mesh->vertices.cartesian);
+ icoord = 0;
+ icoord_bottom =
+ (nsteps_helicoid + 1)
+ * nsteps_circle
+ * 3/*#coords per vertex*/;
+ icoord_top = icoord_bottom + 3/*#coords per vertex */;
+
+ FOR_EACH(ihelicoid, 0, nsteps_helicoid + 1) {
+ const double phi = (double)ihelicoid * step_helicoid;
const double cos_phi = cos(phi);
const double sin_phi = sin(phi);
- FOR_EACH(icircle, 0, circle_nsteps) {
- double x, y, z;
- const double theta = (double)icircle * circle_step;
- const double circle_x = circle_radius*cos(theta) + helicoid_radius;
- const double circle_y = circle_radius*sin(theta) + helicoid_radius;
-
- x = circle_x * cos_phi;
- y = circle_x * sin_phi;
- z = circle_y + c*phi;
- /* TODO rotate the x, y and z coordinates with respect to the helocoid
- * slope */
+ FOR_EACH(icircle, 0, nsteps_circle) {
+ const double* pos = meridian + icircle*3;
+ vertices[icoord + 0] = (float)(pos[0]*cos_phi - pos[1]*sin_phi);
+ vertices[icoord + 1] = (float)(pos[0]*sin_phi + pos[1]*cos_phi);
+ vertices[icoord + 2] = (float)(pos[2] + c * phi);
+ icoord += 3;
+ }
+
+ /* Top cap vertex */
+ if(ihelicoid == 0) {
+ vertices[icoord_bottom + 0] = (float)(hradius * cos_phi);
+ vertices[icoord_bottom + 1] = (float)(hradius * sin_phi);
+ vertices[icoord_bottom + 2] = (float)(c*phi);
+
+ /* Bottom cap vertex */
+ } else if(ihelicoid == nsteps_helicoid) {
+ vertices[icoord_top + 0] = (float)(hradius * cos_phi);
+ vertices[icoord_top + 1] = (float)(hradius * sin_phi);
+ vertices[icoord_top + 2] = (float)(c * phi);
}
}
- return RES_OK;
+
+exit:
+ if(meridian) sa_release(meridian);
+ return res;
+error:
+ goto exit;
}
void
diff --git a/src/schiff_mesh.h b/src/schiff_mesh.h
@@ -64,10 +64,26 @@ schiff_mesh_init_cylinder
struct schiff_mesh* mesh,
const unsigned nslices);
+/* Initialise the indices of the helical pipe. The vertices are still not
+ * defined. Use the schiff_mesh_setup_helical_pipe to define these data */
extern LOCAL_SYM res_T
-schiff_mesh_init_helicoid
+schiff_mesh_init_helical_pipe
(struct mem_allocator* allocator,
- struct schiff_mesh* mesh);
+ struct schiff_mesh* mesh,
+ const unsigned nsteps_helicoid,
+ const unsigned nsteps_circle);
+
+extern LOCAL_SYM res_T
+schiff_mesh_setup_helical_pipe
+ (struct schiff_mesh* mesh,
+ const double pitch,
+ const double height,
+ const double radius_helicoid,
+ const double radius_circle,
+ /* The following attribs are assumed to be equal to attributes used by the
+ * init function */
+ const unsigned nsteps_helicoid,
+ const unsigned nsteps_circle);
extern LOCAL_SYM void
schiff_mesh_release
