50#define PAR_SHAPES_T uint16_t
103 int stacks,
void* userdata);
121 float const* center,
float const* normal);
193#define PAR_PI (3.14159265359)
194#define PAR_MIN(a, b) (a > b ? b : a)
195#define PAR_MAX(a, b) (a > b ? a : b)
196#define PAR_CLAMP(v, lo, hi) PAR_MAX(lo, PAR_MIN(hi, v))
197#define PAR_SWAP(T, A, B) { T tmp = B; B = A; A = tmp; }
198#define PAR_SQR(a) ((a) * (a))
202#define PAR_MALLOC(T, N) ((T*) malloc(N * sizeof(T)))
203#define PAR_CALLOC(T, N) ((T*) calloc(N * sizeof(T), 1))
204#define PAR_REALLOC(T, BUF, N) ((T*) realloc(BUF, sizeof(T) * (N)))
205#define PAR_FREE(BUF) free(BUF)
216#ifdef PAR_SHAPES_IMPLEMENTATION
225static float par_shapes__epsilon_welded_normals = 0.001;
226static float par_shapes__epsilon_degenerate_sphere = 0.0001;
228static void par_shapes__sphere(
float const* uv,
float* xyz,
void*);
229static void par_shapes__hemisphere(
float const* uv,
float* xyz,
void*);
230static void par_shapes__plane(
float const* uv,
float* xyz,
void*);
231static void par_shapes__klein(
float const* uv,
float* xyz,
void*);
232static void par_shapes__cylinder(
float const* uv,
float* xyz,
void*);
233static void par_shapes__cone(
float const* uv,
float* xyz,
void*);
234static void par_shapes__torus(
float const* uv,
float* xyz,
void*);
235static void par_shapes__trefoil(
float const* uv,
float* xyz,
void*);
238static int par__simplex_noise(
int64_t seed,
struct osn_context** ctx);
239static void par__simplex_noise_free(
struct osn_context* ctx);
240static double par__simplex_noise2(
struct osn_context* ctx,
double x,
double y);
242static void par_shapes__copy3(
float* result,
float const* a)
249static float par_shapes__dot3(
float const* a,
float const* b)
251 return b[0] * a[0] + b[1] * a[1] + b[2] * a[2];
254static void par_shapes__transform3(
float* p,
float const* x,
float const* y,
257 float px = par_shapes__dot3(p, x);
258 float py = par_shapes__dot3(p, y);
259 float pz = par_shapes__dot3(p, z);
265static void par_shapes__cross3(
float* result,
float const* a,
float const* b)
267 float x = (a[1] * b[2]) - (a[2] * b[1]);
268 float y = (a[2] * b[0]) - (a[0] * b[2]);
269 float z = (a[0] * b[1]) - (a[1] * b[0]);
275static void par_shapes__mix3(
float* d,
float const* a,
float const* b,
float t)
277 float x = b[0] * t + a[0] * (1 - t);
278 float y = b[1] * t + a[1] * (1 - t);
279 float z = b[2] * t + a[2] * (1 - t);
285static void par_shapes__scale3(
float* result,
float a)
292static void par_shapes__normalize3(
float* v)
294 float lsqr = sqrt(v[0]*v[0] + v[1]*v[1] + v[2]*v[2]);
296 par_shapes__scale3(v, 1.0f / lsqr);
300static void par_shapes__subtract3(
float* result,
float const* a)
307static void par_shapes__add3(
float* result,
float const* a)
314static float par_shapes__sqrdist3(
float const* a,
float const* b)
316 float dx = a[0] - b[0];
317 float dy = a[1] - b[1];
318 float dz = a[2] - b[2];
319 return dx * dx + dy * dy + dz * dz;
324 const float epsilon = par_shapes__epsilon_welded_normals;
330 for (
int i = 0; i < m->
npoints; i++, pdst += 3) {
332 float const* pnormal = welded->
normals + d * 3;
333 pdst[0] = pnormal[0];
334 pdst[1] = pnormal[1];
335 pdst[2] = pnormal[2];
343 if (slices < 3 || stacks < 1) {
352 if (slices < 3 || stacks < 1) {
370 if (slices < 3 || stacks < 3) {
381 if (slices < 3 || stacks < 3) {
392 if (slices < 3 || stacks < 3) {
395 assert(radius <= 1.0 &&
"Use smaller radius to avoid self-intersection.");
396 assert(radius >= 0.1 &&
"Use larger radius to avoid self-intersection.");
397 void* userdata = (
void*) &radius;
404 if (slices < 3 || stacks < 3) {
408 par_shapes__klein, slices, stacks, 0);
410 for (
int stack = 0; stack < stacks; stack++) {
411 for (
int slice = 0; slice < slices; slice++, face += 2) {
412 if (stack < 27 * stacks / 32) {
424 if (slices < 3 || stacks < 3) {
427 assert(radius <= 3.0 &&
"Use smaller radius to avoid self-intersection.");
428 assert(radius >= 0.5 &&
"Use larger radius to avoid self-intersection.");
429 void* userdata = (
void*) &radius;
436 if (slices < 1 || stacks < 1) {
444 int slices,
int stacks,
void* userdata)
449 mesh->
npoints = (slices + 1) * (stacks + 1);
453 float* points = mesh->
points;
454 for (
int stack = 0; stack < stacks + 1; stack++) {
455 uv[0] = (float) stack / stacks;
456 for (
int slice = 0; slice < slices + 1; slice++) {
457 uv[1] = (float) slice / slices;
458 fn(uv, xyz, userdata);
468 for (
int stack = 0; stack < stacks + 1; stack++) {
469 uv[0] = (float) stack / stacks;
470 for (
int slice = 0; slice < slices + 1; slice++) {
471 uv[1] = (float) slice / slices;
482 for (
int stack = 0; stack < stacks; stack++) {
483 for (
int slice = 0; slice < slices; slice++) {
484 int next = slice + 1;
485 *face++ = v + slice + slices + 1;
488 *face++ = v + slice + slices + 1;
489 *face++ = v + next + slices + 1;
510 FILE* objfile = fopen(filename,
"wt");
511 float const* points = mesh->
points;
512 float const* tcoords = mesh->
tcoords;
513 float const* norms = mesh->
normals;
515 if (tcoords && norms) {
516 for (
int nvert = 0; nvert < mesh->
npoints; nvert++) {
517 fprintf(objfile,
"v %f %f %f\n", points[0], points[1], points[2]);
518 fprintf(objfile,
"vt %f %f\n", tcoords[0], tcoords[1]);
519 fprintf(objfile,
"vn %f %f %f\n", norms[0], norms[1], norms[2]);
524 for (
int nface = 0; nface < mesh->
ntriangles; nface++) {
525 int a = 1 + *indices++;
526 int b = 1 + *indices++;
527 int c = 1 + *indices++;
528 fprintf(objfile,
"f %d/%d/%d %d/%d/%d %d/%d/%d\n",
529 a, a, a, b, b, b, c, c, c);
532 for (
int nvert = 0; nvert < mesh->
npoints; nvert++) {
533 fprintf(objfile,
"v %f %f %f\n", points[0], points[1], points[2]);
534 fprintf(objfile,
"vn %f %f %f\n", norms[0], norms[1], norms[2]);
538 for (
int nface = 0; nface < mesh->
ntriangles; nface++) {
539 int a = 1 + *indices++;
540 int b = 1 + *indices++;
541 int c = 1 + *indices++;
542 fprintf(objfile,
"f %d//%d %d//%d %d//%d\n", a, a, b, b, c, c);
544 }
else if (tcoords) {
545 for (
int nvert = 0; nvert < mesh->
npoints; nvert++) {
546 fprintf(objfile,
"v %f %f %f\n", points[0], points[1], points[2]);
547 fprintf(objfile,
"vt %f %f\n", tcoords[0], tcoords[1]);
551 for (
int nface = 0; nface < mesh->
ntriangles; nface++) {
552 int a = 1 + *indices++;
553 int b = 1 + *indices++;
554 int c = 1 + *indices++;
555 fprintf(objfile,
"f %d/%d %d/%d %d/%d\n", a, a, b, b, c, c);
558 for (
int nvert = 0; nvert < mesh->
npoints; nvert++) {
559 fprintf(objfile,
"v %f %f %f\n", points[0], points[1], points[2]);
562 for (
int nface = 0; nface < mesh->
ntriangles; nface++) {
563 int a = 1 + *indices++;
564 int b = 1 + *indices++;
565 int c = 1 + *indices++;
566 fprintf(objfile,
"f %d %d %d\n", a, b, c);
572static void par_shapes__sphere(
float const* uv,
float* xyz,
void* userdata)
574 float phi = uv[0] *
PAR_PI;
575 float theta = uv[1] * 2 *
PAR_PI;
576 xyz[0] = cosf(theta) * sinf(phi);
577 xyz[1] = sinf(theta) * sinf(phi);
581static void par_shapes__hemisphere(
float const* uv,
float* xyz,
void* userdata)
583 float phi = uv[0] *
PAR_PI;
584 float theta = uv[1] *
PAR_PI;
585 xyz[0] = cosf(theta) * sinf(phi);
586 xyz[1] = sinf(theta) * sinf(phi);
590static void par_shapes__plane(
float const* uv,
float* xyz,
void* userdata)
597static void par_shapes__klein(
float const* uv,
float* xyz,
void* userdata)
600 float v = uv[1] * 2 *
PAR_PI;
603 xyz[0] = 3 * cosf(u) * (1 + sinf(u)) + (2 * (1 - cosf(u) / 2)) *
605 xyz[2] = -8 * sinf(u) - 2 * (1 - cosf(u) / 2) * sinf(u) * cosf(v);
607 xyz[0] = 3 * cosf(u) * (1 + sinf(u)) + (2 * (1 - cosf(u) / 2)) *
609 xyz[2] = -8 * sinf(u);
611 xyz[1] = -2 * (1 - cosf(u) / 2) * sinf(v);
614static void par_shapes__cylinder(
float const* uv,
float* xyz,
void* userdata)
616 float theta = uv[1] * 2 *
PAR_PI;
617 xyz[0] = sinf(theta);
618 xyz[1] = cosf(theta);
622static void par_shapes__cone(
float const* uv,
float* xyz,
void* userdata)
624 float r = 1.0f - uv[0];
625 float theta = uv[1] * 2 *
PAR_PI;
626 xyz[0] = r * sinf(theta);
627 xyz[1] = r * cosf(theta);
631static void par_shapes__torus(
float const* uv,
float* xyz,
void* userdata)
634 float minor = *((
float*) userdata);
635 float theta = uv[0] * 2 *
PAR_PI;
636 float phi = uv[1] * 2 *
PAR_PI;
637 float beta = major + minor * cosf(phi);
638 xyz[0] = cosf(theta) * beta;
639 xyz[1] = sinf(theta) * beta;
640 xyz[2] = sinf(phi) * minor;
643static void par_shapes__trefoil(
float const* uv,
float* xyz,
void* userdata)
645 float minor = *((
float*) userdata);
646 const float a = 0.5f;
647 const float b = 0.3f;
648 const float c = 0.5f;
649 const float d = minor * 0.1f;
650 const float u = (1 - uv[0]) * 4 *
PAR_PI;
651 const float v = uv[1] * 2 *
PAR_PI;
652 const float r = a + b * cos(1.5f * u);
653 const float x = r * cos(u);
654 const float y = r * sin(u);
655 const float z = c * sin(1.5f * u);
658 -1.5f * b * sin(1.5f * u) * cos(u) - (a + b * cos(1.5f * u)) * sin(u);
660 -1.5f * b * sin(1.5f * u) * sin(u) + (a + b * cos(1.5f * u)) * cos(u);
661 q[2] = 1.5f * c * cos(1.5f * u);
662 par_shapes__normalize3(q);
663 float qvn[3] = {q[1], -q[0], 0};
664 par_shapes__normalize3(qvn);
666 par_shapes__cross3(ww, q, qvn);
667 xyz[0] = x + d * (qvn[0] * cos(v) + ww[0] * sin(v));
668 xyz[1] = y + d * (qvn[1] * cos(v) + ww[1] * sin(v));
669 xyz[2] = z + d * ww[2] * sin(v);
673 par_shapes__epsilon_welded_normals = epsilon;
677 par_shapes__epsilon_degenerate_sphere = epsilon;
684 int vecsize =
sizeof(float) * 3;
696 int uvsize =
sizeof(float) * 2;
708 *ptriangles++ = offset + *striangles++;
709 *ptriangles++ = offset + *striangles++;
710 *ptriangles++ = offset + *striangles++;
716 float const* center,
float const* normal)
721 float* points = mesh->
points;
725 for (
int i = 0; i < slices; i++) {
726 float theta = i *
PAR_PI * 2 / slices;
727 *points++ = radius * cos(theta);
728 *points++ = radius * sin(theta);
731 float nnormal[3] = {normal[0], normal[1], normal[2]};
732 par_shapes__normalize3(nnormal);
735 for (
int i = 0; i < mesh->
npoints; i++) {
736 *norms++ = nnormal[0];
737 *norms++ = nnormal[1];
738 *norms++ = nnormal[2];
743 for (
int i = 0; i < slices; i++) {
745 *triangles++ = 1 + i;
746 *triangles++ = 1 + (i + 1) % slices;
748 float k[3] = {0, 0, -1};
750 par_shapes__cross3(axis, nnormal, k);
751 par_shapes__normalize3(axis);
764 float* points = m->
points;
765 for (
int i = 0; i < m->
npoints; i++) {
774 float s = sinf(radians);
775 float c = cosf(radians);
782 float oneMinusC = 1.0f - c;
784 (((x * x) * oneMinusC) + c),
785 ((xy * oneMinusC) + (z * s)), ((zx * oneMinusC) - (y * s))
788 ((xy * oneMinusC) - (z * s)),
789 (((y * y) * oneMinusC) + c), ((yz * oneMinusC) + (x * s))
792 ((zx * oneMinusC) + (y * s)),
793 ((yz * oneMinusC) - (x * s)), (((z * z) * oneMinusC) + c)
796 for (
int i = 0; i < mesh->
npoints; i++, p += 3) {
797 float x = col0[0] * p[0] + col1[0] * p[1] + col2[0] * p[2];
798 float y = col0[1] * p[0] + col1[1] * p[1] + col2[1] * p[2];
799 float z = col0[2] * p[0] + col1[2] * p[1] + col2[2] * p[2];
806 for (
int i = 0; i < mesh->
npoints; i++, n += 3) {
807 float x = col0[0] * n[0] + col1[0] * n[1] + col2[0] * n[2];
808 float y = col0[1] * n[0] + col1[1] * n[1] + col2[1] * n[2];
809 float z = col0[2] * n[0] + col1[2] * n[1] + col2[2] * n[2];
819 float* points = m->
points;
820 for (
int i = 0; i < m->
npoints; i++) {
826 if (n && !(x == y && y == z)) {
827 bool x_zero = x == 0;
828 bool y_zero = y == 0;
829 bool z_zero = z == 0;
830 if (!x_zero && !y_zero && !z_zero) {
835 x = x_zero && !y_zero && !z_zero;
836 y = y_zero && !x_zero && !z_zero;
837 z = z_zero && !x_zero && !y_zero;
839 for (
int i = 0; i < m->
npoints; i++, n += 3) {
843 par_shapes__normalize3(n);
856 float* points = m->
points;
857 aabb[0] = aabb[3] = points[0];
858 aabb[1] = aabb[4] = points[1];
859 aabb[2] = aabb[5] = points[2];
861 for (
int i = 1; i < m->
npoints; i++, points += 3) {
862 aabb[0] =
PAR_MIN(points[0], aabb[0]);
863 aabb[1] =
PAR_MIN(points[1], aabb[1]);
864 aabb[2] =
PAR_MIN(points[2], aabb[2]);
865 aabb[3] =
PAR_MAX(points[0], aabb[3]);
866 aabb[4] =
PAR_MAX(points[1], aabb[4]);
867 aabb[5] =
PAR_MAX(points[2], aabb[5]);
875 for (
int i = 0; i < nfaces; i++) {
883 static float verts[] = {
887 -0.724, 0.526, 0.447,
888 -0.724, -0.526, 0.447,
889 0.276, -0.851, 0.447,
890 0.724, 0.526, -0.447,
891 -0.276, 0.851, -0.447,
892 -0.894, 0.000, -0.447,
893 -0.276, -0.851, -0.447,
894 0.724, -0.526, -0.447,
920 mesh->
npoints =
sizeof(verts) /
sizeof(verts[0]) / 3;
922 memcpy(mesh->
points, verts,
sizeof(verts));
923 mesh->
ntriangles =
sizeof(faces) /
sizeof(faces[0]) / 3;
925 memcpy(mesh->
triangles, faces,
sizeof(faces));
931 static float verts[20 * 3] = {
934 -0.491, 0.357, 0.795,
935 -0.491, -0.357, 0.795,
936 0.188, -0.577, 0.795,
939 -0.795, 0.577, 0.188,
940 -0.795, -0.577, 0.188,
941 0.304, -0.934, 0.188,
942 0.795, 0.577, -0.188,
943 -0.304, 0.934, -0.188,
944 -0.982, 0.000, -0.188,
945 -0.304, -0.934, -0.188,
946 0.795, -0.577, -0.188,
947 0.491, 0.357, -0.795,
948 -0.188, 0.577, -0.795,
949 -0.607, 0.000, -0.795,
950 -0.188, -0.577, -0.795,
951 0.491, -0.357, -0.795,
967 int npentagons =
sizeof(pentagons) /
sizeof(pentagons[0]) / 5;
969 int ncorners =
sizeof(verts) /
sizeof(verts[0]) / 3;
972 memcpy(mesh->
points, verts,
sizeof(verts));
977 for (
int p = 0; p < npentagons; p++, pentagon += 5) {
978 *tris++ = pentagon[0];
979 *tris++ = pentagon[1];
980 *tris++ = pentagon[2];
981 *tris++ = pentagon[0];
982 *tris++ = pentagon[2];
983 *tris++ = pentagon[3];
984 *tris++ = pentagon[0];
985 *tris++ = pentagon[3];
986 *tris++ = pentagon[4];
993 static float verts[6 * 3] = {
997 -1.000, 0.000, 0.000,
998 0.000, -1.000, 0.000,
1011 int ntris =
sizeof(triangles) /
sizeof(triangles[0]) / 3;
1013 int ncorners =
sizeof(verts) /
sizeof(verts[0]) / 3;
1016 memcpy(mesh->
points, verts,
sizeof(verts));
1021 for (
int p = 0; p < ntris; p++) {
1022 *tris++ = *triangle++;
1023 *tris++ = *triangle++;
1024 *tris++ = *triangle++;
1031 static float verts[4 * 3] = {
1043 int ntris =
sizeof(triangles) /
sizeof(triangles[0]) / 3;
1045 int ncorners =
sizeof(verts) /
sizeof(verts[0]) / 3;
1048 memcpy(mesh->
points, verts,
sizeof(verts));
1053 for (
int p = 0; p < ntris; p++) {
1054 *tris++ = *triangle++;
1055 *tris++ = *triangle++;
1056 *tris++ = *triangle++;
1063 static float verts[8 * 3] = {
1081 int nquads =
sizeof(quads) /
sizeof(quads[0]) / 4;
1083 int ncorners =
sizeof(verts) /
sizeof(verts[0]) / 3;
1086 memcpy(mesh->
points, verts,
sizeof(verts));
1091 for (
int p = 0; p < nquads; p++, quad += 4) {
1105} par_shapes__command;
1111 par_shapes__command* commands;
1119 par_shapes__rule* rule;
1120} par_shapes__stackframe;
1122static par_shapes__rule* par_shapes__pick_rule(
const char* name,
1123 par_shapes__rule* rules,
int nrules)
1125 par_shapes__rule* rule = 0;
1127 for (
int i = 0; i < nrules; i++) {
1129 if (!strcmp(rule->name, name)) {
1130 total += rule->weight;
1133 float r = (float) rand() / RAND_MAX;
1135 for (
int i = 0; i < nrules; i++) {
1137 if (!strcmp(rule->name, name)) {
1138 t += (float) rule->weight / total;
1149 const float xaxis[] = {1, 0, 0};
1150 const float yaxis[] = {0, 1, 0};
1151 const float zaxis[] = {0, 0, 1};
1155 par_shapes__copy3(turtle->
points + 0, xaxis);
1156 par_shapes__copy3(turtle->
points + 3, yaxis);
1157 par_shapes__copy3(turtle->
points + 6, zaxis);
1165 for (
int p = 0; p < m->
npoints; p++) {
1166 float* pt = m->
points + p * 3;
1170 par_shapes__transform3(pt,
1183 int npoints = (slices + 1) * (stacks + 1);
1184 assert(scene->
npoints >= npoints &&
"Cannot connect to empty scene.");
1187 npoints = scene->
npoints + (slices + 1);
1188 float* points =
PAR_MALLOC(
float, npoints * 3);
1189 memcpy(points, scene->
points,
sizeof(
float) * scene->
npoints * 3);
1190 float* newpts = points + scene->
npoints * 3;
1191 memcpy(newpts, cylinder->
points + (slices + 1) * 3,
1192 sizeof(
float) * (slices + 1) * 3);
1197 int ntriangles = scene->
ntriangles + 2 * slices * stacks;
1201 int v = scene->
npoints - (slices + 1);
1203 for (
int stack = 0; stack < stacks; stack++) {
1204 for (
int slice = 0; slice < slices; slice++) {
1205 int next = slice + 1;
1206 *face++ = v + slice + slices + 1;
1208 *face++ = v + slice;
1209 *face++ = v + slice + slices + 1;
1210 *face++ = v + next + slices + 1;
1226 program =
PAR_MALLOC(
char, strlen(text) + 1);
1229 strcpy(program, text);
1230 char *cmd = strtok(program,
" ");
1234 char *arg = strtok(0,
" ");
1236 puts(
"lsystem error: unexpected end of program.");
1239 if (!strcmp(cmd,
"rule")) {
1244 cmd = strtok(0,
" ");
1248 par_shapes__rule* rules =
PAR_MALLOC(par_shapes__rule, nrules);
1249 par_shapes__command* commands =
PAR_MALLOC(par_shapes__command, ncommands);
1252 par_shapes__rule* current_rule = &rules[0];
1253 par_shapes__command* current_command = &commands[0];
1254 current_rule->name =
"entry";
1255 current_rule->weight = 1;
1256 current_rule->ncommands = 0;
1257 current_rule->commands = current_command;
1260 strcpy(program, text);
1261 cmd = strtok(program,
" ");
1263 char *arg = strtok(0,
" ");
1264 if (!strcmp(cmd,
"rule")) {
1268 char* dot = strchr(arg,
'.');
1270 current_rule->weight = atoi(dot + 1);
1273 current_rule->weight = 1;
1276 current_rule->name = arg;
1277 current_rule->ncommands = 0;
1278 current_rule->commands = current_command;
1280 current_rule->ncommands++;
1281 current_command->cmd = cmd;
1282 current_command->arg = arg;
1285 cmd = strtok(0,
" ");
1290 for (
int i = 0; i < nrules; i++) {
1291 par_shapes__rule rule = rules[i];
1292 printf(
"rule %s.%d\n", rule.name, rule.weight);
1293 for (
int c = 0; c < rule.ncommands; c++) {
1294 par_shapes__command cmd = rule.commands[c];
1295 printf(
"\t%s %s\n", cmd.cmd, cmd.arg);
1312 const float xaxis[] = {1, 0, 0};
1313 const float yaxis[] = {0, 1, 0};
1314 const float zaxis[] = {0, 0, 1};
1315 const float units[] = {1, 1, 1};
1318 par_shapes__stackframe* stack =
1319 PAR_CALLOC(par_shapes__stackframe, maxdepth);
1321 stack[0].orientation = turtle;
1322 stack[0].rule = &rules[0];
1323 par_shapes__copy3(stack[0].scale, units);
1324 while (stackptr >= 0) {
1325 par_shapes__stackframe* frame = &stack[stackptr];
1326 par_shapes__rule* rule = frame->rule;
1328 float* position = frame->position;
1329 float* scale = frame->scale;
1330 if (frame->pc >= rule->ncommands) {
1336 par_shapes__command* cmd = rule->commands + (frame->pc++);
1338 printf(
"%5s %5s %5s:%d %03d\n", cmd->cmd, cmd->arg, rule->name,
1339 frame->pc - 1, stackptr);
1343 if (!strcmp(cmd->cmd,
"shape")) {
1346 if (!strcmp(cmd->arg,
"connect")) {
1352 }
else if (!strcmp(cmd->cmd,
"call") && stackptr < maxdepth - 1) {
1353 rule = par_shapes__pick_rule(cmd->arg, rules, nrules);
1354 frame = &stack[++stackptr];
1358 par_shapes__copy3(frame->scale, scale);
1359 par_shapes__copy3(frame->position, position);
1362 value = atof(cmd->arg);
1363 if (!strcmp(cmd->cmd,
"rx")) {
1365 }
else if (!strcmp(cmd->cmd,
"ry")) {
1367 }
else if (!strcmp(cmd->cmd,
"rz")) {
1369 }
else if (!strcmp(cmd->cmd,
"tx")) {
1370 float vec[3] = {value, 0, 0};
1372 par_shapes__dot3(turtle->
points + 0, vec),
1373 par_shapes__dot3(turtle->
points + 3, vec),
1374 par_shapes__dot3(turtle->
points + 6, vec)
1376 par_shapes__add3(position, t);
1377 }
else if (!strcmp(cmd->cmd,
"ty")) {
1378 float vec[3] = {0, value, 0};
1380 par_shapes__dot3(turtle->
points + 0, vec),
1381 par_shapes__dot3(turtle->
points + 3, vec),
1382 par_shapes__dot3(turtle->
points + 6, vec)
1384 par_shapes__add3(position, t);
1385 }
else if (!strcmp(cmd->cmd,
"tz")) {
1386 float vec[3] = {0, 0, value};
1388 par_shapes__dot3(turtle->
points + 0, vec),
1389 par_shapes__dot3(turtle->
points + 3, vec),
1390 par_shapes__dot3(turtle->
points + 6, vec)
1392 par_shapes__add3(position, t);
1393 }
else if (!strcmp(cmd->cmd,
"sx")) {
1395 }
else if (!strcmp(cmd->cmd,
"sy")) {
1397 }
else if (!strcmp(cmd->cmd,
"sz")) {
1399 }
else if (!strcmp(cmd->cmd,
"sa")) {
1416 float* points =
PAR_MALLOC(
float, 3 * npoints);
1417 float* dst = points;
1419 for (
int i = 0; i < npoints; i++) {
1420 float const* src = mesh->
points + 3 * (*index++);
1428 if (create_indices) {
1431 for (
int i = 0; i < mesh->
ntriangles * 3; i++) {
1444 float next[3], prev[3], cp[3];
1445 for (
int f = 0; f < m->
ntriangles; f++, triangle += 3) {
1446 float const* pa = m->
points + 3 * triangle[0];
1447 float const* pb = m->
points + 3 * triangle[1];
1448 float const* pc = m->
points + 3 * triangle[2];
1449 par_shapes__copy3(next, pb);
1450 par_shapes__subtract3(next, pa);
1451 par_shapes__copy3(prev, pc);
1452 par_shapes__subtract3(prev, pa);
1453 par_shapes__cross3(cp, next, prev);
1454 par_shapes__add3(m->
normals + 3 * triangle[0], cp);
1455 par_shapes__copy3(next, pc);
1456 par_shapes__subtract3(next, pb);
1457 par_shapes__copy3(prev, pa);
1458 par_shapes__subtract3(prev, pb);
1459 par_shapes__cross3(cp, next, prev);
1460 par_shapes__add3(m->
normals + 3 * triangle[1], cp);
1461 par_shapes__copy3(next, pa);
1462 par_shapes__subtract3(next, pc);
1463 par_shapes__copy3(prev, pb);
1464 par_shapes__subtract3(prev, pc);
1465 par_shapes__cross3(cp, next, prev);
1466 par_shapes__add3(m->
normals + 3 * triangle[2], cp);
1469 for (
int p = 0; p < m->
npoints; p++, normal += 3) {
1470 par_shapes__normalize3(normal);
1478 int npoints = ntriangles * 3;
1479 float* points =
PAR_CALLOC(
float, npoints * 3);
1480 float* dpoint = points;
1481 float const* spoint = mesh->
points;
1482 for (
int t = 0; t < mesh->
ntriangles; t++, spoint += 9, dpoint += 3) {
1483 float const* a = spoint;
1484 float const* b = spoint + 3;
1485 float const* c = spoint + 6;
1486 float const* p0 = dpoint;
1487 float const* p1 = dpoint + 3;
1488 float const* p2 = dpoint + 6;
1489 par_shapes__mix3(dpoint, a, b, 0.5);
1490 par_shapes__mix3(dpoint += 3, b, c, 0.5);
1491 par_shapes__mix3(dpoint += 3, a, c, 0.5);
1492 par_shapes__add3(dpoint += 3, a);
1493 par_shapes__add3(dpoint += 3, p0);
1494 par_shapes__add3(dpoint += 3, p2);
1495 par_shapes__add3(dpoint += 3, p0);
1496 par_shapes__add3(dpoint += 3, b);
1497 par_shapes__add3(dpoint += 3, p1);
1498 par_shapes__add3(dpoint += 3, p2);
1499 par_shapes__add3(dpoint += 3, p1);
1500 par_shapes__add3(dpoint += 3, c);
1515 par_shapes__subdivide(mesh);
1517 for (
int i = 0; i < mesh->
npoints; i++) {
1518 par_shapes__normalize3(mesh->
points + i * 3);
1521 for (
int i = 0; i < mesh->
ntriangles * 3; i++) {
1534 struct osn_context* ctx;
1535 par__simplex_noise(seed, &ctx);
1536 for (
int p = 0; p < mesh->
npoints; p++) {
1537 float* pt = mesh->
points + p * 3;
1538 float a = 0.25, f = 1.0;
1539 double n = a * par__simplex_noise2(ctx, f * pt[0], f * pt[2]);
1541 n += a * par__simplex_noise2(ctx, f * pt[0], f * pt[2]);
1546 pt[1] = -pow(-pt[1], 0.5) / 2;
1549 par__simplex_noise_free(ctx);
1571 sizeof(
float) * 3 * clone->
npoints);
1576 sizeof(
float) * 2 * clone->
npoints);
1582 float const* points;
1584} par_shapes__sort_context;
1586static int par_shapes__cmp1(
const void *arg0,
const void *arg1)
1588 const int g = par_shapes__sort_context.gridsize;
1592 float const* p0 = par_shapes__sort_context.points + d0 * 3;
1593 int i0 = (int) p0[0];
1594 int j0 = (int) p0[1];
1595 int k0 = (int) p0[2];
1596 int index0 = i0 + g * j0 + g * g * k0;
1600 float const* p1 = par_shapes__sort_context.points + d1 * 3;
1601 int i1 = (int) p1[0];
1602 int j1 = (int) p1[1];
1603 int k1 = (int) p1[2];
1604 int index1 = i1 + g * j1 + g * g * k1;
1607 if (index0 < index1)
return -1;
1608 if (index0 > index1)
return 1;
1612static void par_shapes__sort_points(
par_shapes_mesh* mesh,
int gridsize,
1617 for (
int i = 0; i < mesh->
npoints; i++) {
1620 par_shapes__sort_context.gridsize = gridsize;
1621 par_shapes__sort_context.points = mesh->
points;
1627 float* dstpt = newpts;
1628 for (
int i = 0; i < mesh->
npoints; i++) {
1629 invmap[sortmap[i]] = i;
1630 float const* srcpt = mesh->
points + 3 * sortmap[i];
1631 *dstpt++ = *srcpt++;
1632 *dstpt++ = *srcpt++;
1633 *dstpt++ = *srcpt++;
1642 for (
int i = 0; i < mesh->
ntriangles * 3; i++) {
1643 *dstind++ = invmap[*srcind++];
1653static void par_shapes__weld_points(
par_shapes_mesh* mesh,
int gridsize,
1661 gridsize * gridsize * gridsize);
1662 int prev_binindex = -1;
1663 for (
int p = 0; p < mesh->
npoints; p++) {
1664 float const* pt = mesh->
points + p * 3;
1665 int i = (int) pt[0];
1666 int j = (int) pt[1];
1667 int k = (int) pt[2];
1668 int this_binindex = i + gridsize * j + gridsize * gridsize * k;
1669 if (this_binindex != prev_binindex) {
1670 bins[this_binindex] = 1 + p;
1672 prev_binindex = this_binindex;
1677 float const* pt = mesh->
points;
1679 for (
int p = 0; p < mesh->
npoints; p++, pt += 3) {
1682 if (weldmap[p] != p) {
1689 int minp[3], maxp[3];
1690 for (
int c = 0; c < 3; c++) {
1691 minp[c] = (int) (pt[c] - epsilon);
1692 maxp[c] = (int) (pt[c] + epsilon);
1694 for (
int i = minp[0]; i <= maxp[0]; i++) {
1695 for (
int j = minp[1]; j <= maxp[1]; j++) {
1696 for (
int k = minp[2]; k <= maxp[2]; k++) {
1697 int binindex = i + gridsize * j + gridsize * gridsize * k;
1701 printf(
"Epsilon value is too large.\n");
1704 nearby[nbins++] = binindex;
1711 for (
int b = 0; b < nbins; b++) {
1712 int binindex = nearby[b];
1715 assert(nindex < mesh->npoints);
1719 if (nindex != p && weldmap[nindex] == nindex) {
1720 float const* thatpt = mesh->
points + nindex * 3;
1721 float dist2 = par_shapes__sqrdist3(thatpt, pt);
1722 if (dist2 < epsilon) {
1723 weldmap[nindex] = p;
1729 if (++nindex >= mesh->
npoints) {
1734 float const* nextpt = mesh->
points + nindex * 3;
1735 int i = (int) nextpt[0];
1736 int j = (int) nextpt[1];
1737 int k = (int) nextpt[2];
1738 int nextbinindex = i + gridsize * j + gridsize * gridsize * k;
1739 if (nextbinindex != binindex) {
1748 int npoints = mesh->
npoints - nremoved;
1749 float* newpts =
PAR_MALLOC(
float, 3 * npoints);
1750 float* dst = newpts;
1753 float const* src = mesh->
points;
1755 for (
int p = 0; p < mesh->
npoints; p++, src += 3) {
1756 if (weldmap[p] == p) {
1762 *cmap++ = condensed_map[weldmap[p]];
1765 assert(ci == npoints);
1776 for (
int i = 0; i < mesh->
ntriangles; i++, tsrc += 3) {
1780 if (a != b && a != c && b != c) {
1781 assert(a < mesh->npoints);
1782 assert(b < mesh->npoints);
1783 assert(c < mesh->npoints);
1799 float maxcell = gridsize - 1;
1802 aabb[3] == aabb[0] ? 1.0f : maxcell / (aabb[3] - aabb[0]),
1803 aabb[4] == aabb[1] ? 1.0f : maxcell / (aabb[4] - aabb[1]),
1804 aabb[5] == aabb[2] ? 1.0f : maxcell / (aabb[5] - aabb[2]),
1809 par_shapes__sort_points(clone, gridsize, sortmap);
1815 for (
int i = 0; i < mesh->
npoints; i++) {
1818 par_shapes__weld_points(clone, gridsize, epsilon, weldmap);
1823 for (
int i = 0; i < mesh->
npoints; i++) {
1824 newmap[i] = weldmap[sortmap[i]];
1839#define STRETCH_CONSTANT_2D (-0.211324865405187)
1840#define SQUISH_CONSTANT_2D (0.366025403784439)
1841#define STRETCH_CONSTANT_3D (-1.0 / 6.0)
1842#define SQUISH_CONSTANT_3D (1.0 / 3.0)
1843#define STRETCH_CONSTANT_4D (-0.138196601125011)
1844#define SQUISH_CONSTANT_4D (0.309016994374947)
1846#define NORM_CONSTANT_2D (47.0)
1847#define NORM_CONSTANT_3D (103.0)
1848#define NORM_CONSTANT_4D (30.0)
1850#define DEFAULT_SEED (0LL)
1857#define ARRAYSIZE(x) (sizeof((x)) / sizeof((x)[0]))
1863static const int8_t gradients2D[] = {
1864 5, 2, 2, 5, -5, 2, -2, 5, 5, -2, 2, -5, -5, -2, -2, -5,
1873static const signed char gradients3D[] = {
1874 -11, 4, 4, -4, 11, 4, -4, 4, 11, 11, 4, 4, 4, 11, 4, 4, 4, 11, -11, -4, 4,
1875 -4, -11, 4, -4, -4, 11, 11, -4, 4, 4, -11, 4, 4, -4, 11, -11, 4, -4, -4, 11,
1876 -4, -4, 4, -11, 11, 4, -4, 4, 11, -4, 4, 4, -11, -11, -4, -4, -4, -11, -4,
1877 -4, -4, -11, 11, -4, -4, 4, -11, -4, 4, -4, -11,
1886static const signed char gradients4D[] = {
1887 3, 1, 1, 1, 1, 3, 1, 1, 1, 1, 3, 1, 1, 1, 1, 3, -3, 1, 1, 1, -1, 3, 1, 1,
1888 -1, 1, 3, 1, -1, 1, 1, 3, 3, -1, 1, 1, 1, -3, 1, 1, 1, -1, 3, 1, 1, -1, 1,
1889 3, -3, -1, 1, 1, -1, -3, 1, 1, -1, -1, 3, 1, -1, -1, 1, 3, 3, 1, -1, 1, 1,
1890 3, -1, 1, 1, 1, -3, 1, 1, 1, -1, 3, -3, 1, -1, 1, -1, 3, -1, 1, -1, 1, -3,
1891 1, -1, 1, -1, 3, 3, -1, -1, 1, 1, -3, -1, 1, 1, -1, -3, 1, 1, -1, -1, 3, -3,
1892 -1, -1, 1, -1, -3, -1, 1, -1, -1, -3, 1, -1, -1, -1, 3, 3, 1, 1, -1, 1, 3,
1893 1, -1, 1, 1, 3, -1, 1, 1, 1, -3, -3, 1, 1, -1, -1, 3, 1, -1, -1, 1, 3, -1,
1894 -1, 1, 1, -3, 3, -1, 1, -1, 1, -3, 1, -1, 1, -1, 3, -1, 1, -1, 1, -3, -3,
1895 -1, 1, -1, -1, -3, 1, -1, -1, -1, 3, -1, -1, -1, 1, -3, 3, 1, -1, -1, 1, 3,
1896 -1, -1, 1, 1, -3, -1, 1, 1, -1, -3, -3, 1, -1, -1, -1, 3, -1, -1, -1, 1, -3,
1897 -1, -1, 1, -1, -3, 3, -1, -1, -1, 1, -3, -1, -1, 1, -1, -3, -1, 1, -1, -1,
1898 -3, -3, -1, -1, -1, -1, -3, -1, -1, -1, -1, -3, -1, -1, -1, -1, -3,
1901static double extrapolate2(
1902 struct osn_context* ctx,
int xsb,
int ysb,
double dx,
double dy)
1905 int index = perm[(perm[xsb & 0xFF] + ysb) & 0xFF] & 0x0E;
1906 return gradients2D[index] * dx + gradients2D[index + 1] * dy;
1909static inline int fastFloor(
double x)
1912 return x < xi ? xi - 1 : xi;
1915static int allocate_perm(
struct osn_context* ctx,
int nperm,
int ngrad)
1924 if (!ctx->permGradIndex3D) {
1931static int par__simplex_noise(
int64_t seed,
struct osn_context** ctx)
1942 (*ctx)->perm =
NULL;
1943 (*ctx)->permGradIndex3D =
NULL;
1944 rc = allocate_perm(*ctx, 256, 256);
1949 perm = (*ctx)->perm;
1950 permGradIndex3D = (*ctx)->permGradIndex3D;
1951 for (i = 0; i < 256; i++) {
1954 seed = seed * 6364136223846793005LL + 1442695040888963407LL;
1955 seed = seed * 6364136223846793005LL + 1442695040888963407LL;
1956 seed = seed * 6364136223846793005LL + 1442695040888963407LL;
1957 for (i = 255; i >= 0; i--) {
1958 seed = seed * 6364136223846793005LL + 1442695040888963407LL;
1959 int r = (int) ((seed + 31) % (i + 1));
1962 perm[i] = source[r];
1963 permGradIndex3D[i] =
1964 (short) ((perm[i] % (ARRAYSIZE(gradients3D) / 3)) * 3);
1965 source[r] = source[i];
1970static void par__simplex_noise_free(
struct osn_context* ctx)
1978 if (ctx->permGradIndex3D) {
1980 ctx->permGradIndex3D =
NULL;
1985static double par__simplex_noise2(
struct osn_context* ctx,
double x,
double y)
1988 double stretchOffset = (x + y) * STRETCH_CONSTANT_2D;
1989 double xs = x + stretchOffset;
1990 double ys = y + stretchOffset;
1994 int xsb = fastFloor(xs);
1995 int ysb = fastFloor(ys);
1999 double squishOffset = (xsb + ysb) * SQUISH_CONSTANT_2D;
2000 double xb = xsb + squishOffset;
2001 double yb = ysb + squishOffset;
2004 double xins = xs - xsb;
2005 double yins = ys - ysb;
2008 double inSum = xins + yins;
2011 double dx0 = x - xb;
2012 double dy0 = y - yb;
2015 double dx_ext, dy_ext;
2016 int xsv_ext, ysv_ext;
2021 double dx1 = dx0 - 1 - SQUISH_CONSTANT_2D;
2022 double dy1 = dy0 - 0 - SQUISH_CONSTANT_2D;
2023 double attn1 = 2 - dx1 * dx1 - dy1 * dy1;
2026 value += attn1 * attn1 * extrapolate2(ctx, xsb + 1, ysb + 0, dx1, dy1);
2030 double dx2 = dx0 - 0 - SQUISH_CONSTANT_2D;
2031 double dy2 = dy0 - 1 - SQUISH_CONSTANT_2D;
2032 double attn2 = 2 - dx2 * dx2 - dy2 * dy2;
2035 value += attn2 * attn2 * extrapolate2(ctx, xsb + 0, ysb + 1, dx2, dy2);
2039 double zins = 1 - inSum;
2040 if (zins > xins || zins > yins) {
2055 dx_ext = dx0 - 1 - 2 * SQUISH_CONSTANT_2D;
2056 dy_ext = dy0 - 1 - 2 * SQUISH_CONSTANT_2D;
2059 double zins = 2 - inSum;
2060 if (zins < xins || zins < yins) {
2064 dx_ext = dx0 - 2 - 2 * SQUISH_CONSTANT_2D;
2065 dy_ext = dy0 + 0 - 2 * SQUISH_CONSTANT_2D;
2069 dx_ext = dx0 + 0 - 2 * SQUISH_CONSTANT_2D;
2070 dy_ext = dy0 - 2 - 2 * SQUISH_CONSTANT_2D;
2080 dx0 = dx0 - 1 - 2 * SQUISH_CONSTANT_2D;
2081 dy0 = dy0 - 1 - 2 * SQUISH_CONSTANT_2D;
2085 double attn0 = 2 - dx0 * dx0 - dy0 * dy0;
2088 value += attn0 * attn0 * extrapolate2(ctx, xsb, ysb, dx0, dy0);
2092 double attn_ext = 2 - dx_ext * dx_ext - dy_ext * dy_ext;
2094 attn_ext *= attn_ext;
2095 value += attn_ext * attn_ext *
2096 extrapolate2(ctx, xsv_ext, ysv_ext, dx_ext, dy_ext);
2099 return value / NORM_CONSTANT_2D;
2108 float next[3], prev[3], cp[3];
2109 float mincplen2 = (mintriarea * 2) * (mintriarea * 2);
2110 for (
int f = 0; f < mesh->
ntriangles; f++, src += 3) {
2111 float const* pa = mesh->
points + 3 * src[0];
2112 float const* pb = mesh->
points + 3 * src[1];
2113 float const* pc = mesh->
points + 3 * src[2];
2114 par_shapes__copy3(next, pb);
2115 par_shapes__subtract3(next, pa);
2116 par_shapes__copy3(prev, pc);
2117 par_shapes__subtract3(prev, pa);
2118 par_shapes__cross3(cp, next, prev);
2119 float cplen2 = par_shapes__dot3(cp, cp);
2120 if (cplen2 >= mincplen2) {
void par_shapes_compute_aabb(par_shapes_mesh const *mesh, float *aabb)
void par_shapes_merge_and_free(par_shapes_mesh *dst, par_shapes_mesh *src)
par_shapes_mesh * par_shapes_create_dodecahedron()
void par_shapes_invert(par_shapes_mesh *, int startface, int nfaces)
par_shapes_mesh * par_shapes_create_parametric_sphere(int slices, int stacks)
par_shapes_mesh * par_shapes_create_cube()
void par_shapes_remove_degenerate(par_shapes_mesh *, float minarea)
void par_shapes_scale(par_shapes_mesh *, float x, float y, float z)
par_shapes_mesh * par_shapes_create_cone(int slices, int stacks)
void par_shapes_free_mesh(par_shapes_mesh *)
void par_shapes__compute_welded_normals(par_shapes_mesh *m)
par_shapes_mesh * par_shapes_create_cylinder(int slices, int stacks)
par_shapes_mesh * par_shapes_create_parametric_disk(int slices, int stacks)
#define PAR_REALLOC(T, BUF, N)
par_shapes_mesh * par_shapes_create_empty()
par_shapes_mesh * par_shapes_create_torus(int slices, int stacks, float radius)
void par_shapes_unweld(par_shapes_mesh *mesh, bool create_indices)
par_shapes_mesh * par_shapes_create_plane(int slices, int stacks)
par_shapes_mesh * par_shapes_create_disk(float radius, int slices, float const *center, float const *normal)
par_shapes_mesh * par_shapes_create_trefoil_knot(int slices, int stacks, float radius)
void par_shapes_set_epsilon_welded_normals(float epsilon)
void par_shapes_set_epsilon_degenerate_sphere(float epsilon)
void par_shapes_export(par_shapes_mesh const *, char const *objfile)
par_shapes_mesh * par_shapes_create_rock(int seed, int nsubdivisions)
void par_shapes_compute_normals(par_shapes_mesh *m)
void par_shapes_translate(par_shapes_mesh *, float x, float y, float z)
par_shapes_mesh * par_shapes_create_hemisphere(int slices, int stacks)
par_shapes_mesh * par_shapes_create_tetrahedron()
par_shapes_mesh * par_shapes_create_octahedron()
par_shapes_mesh * par_shapes_clone(par_shapes_mesh const *mesh, par_shapes_mesh *target)
void par_shapes__connect(par_shapes_mesh *scene, par_shapes_mesh *cylinder, int slices)
struct par_shapes_mesh_s par_shapes_mesh
void par_shapes_merge(par_shapes_mesh *dst, par_shapes_mesh const *src)
#define PAR_SWAP(T, A, B)
par_shapes_mesh * par_shapes_create_parametric(par_shapes_fn, int slices, int stacks, void *userdata)
par_shapes_mesh * par_shapes_create_subdivided_sphere(int nsubdivisions)
par_shapes_mesh * par_shapes_create_lsystem(char const *program, int slices, int maxdepth)
void(* par_shapes_fn)(float const *, float *, void *)
par_shapes_mesh * par_shapes_create_icosahedron()
void par_shapes_rotate(par_shapes_mesh *, float radians, float const *axis)
par_shapes_mesh * par_shapes_create_klein_bottle(int slices, int stacks)
par_shapes_mesh * par_shapes_weld(par_shapes_mesh const *, float epsilon, PAR_SHAPES_T *mapping)