voxfract

view src/main.c @ 0:8171de5a000b

initial commit
author John Tsiombikas <nuclear@member.fsf.org>
date Wed, 14 Jun 2017 18:03:57 +0300
parents
children d0297d001505
line source
1 #include <stdio.h>
2 #include <stdlib.h>
3 #include <math.h>
4 #include <assert.h>
5 #include <GL/glut.h>
6 #include <metasurf.h>
8 struct quat {
9 float x, y, z, w;
10 };
12 int init(void);
13 void cleanup(void);
14 void display(void);
15 void reshape(int x, int y);
16 void keypress(unsigned char key, int x, int y);
17 void mouse(int bn, int st, int x, int y);
18 void motion(int x, int y);
20 float eval(struct metasurface *ms, float x, float y, float z);
21 void vertex(struct metasurface *ms, float x, float y, float z);
23 int julia(struct quat *qres, struct quat *qprime, struct quat *q, struct quat *seed, int max_iter);
24 float julia_dist(struct quat *z, struct quat *seed, int max_iter);
25 void julia_grad(float *grad, float dist, struct quat *q, struct quat *seed, int max_iter);
27 void show_waitscr(void);
29 float cam_theta, cam_phi = 25, cam_dist = 5;
30 int grid_size = 350;
31 float grid_scale = 1.7;
32 struct quat seed = {0.4, 0.0, 0.0, -0.8};
33 int max_iter = 9;
35 struct metasurface *msurf;
36 int dlist;
37 FILE *fp;
39 int main(int argc, char **argv)
40 {
41 glutInit(&argc, argv);
42 glutInitWindowSize(800, 600);
43 glutInitDisplayMode(GLUT_RGB | GLUT_DEPTH | GLUT_DOUBLE);
44 glutCreateWindow("voxel fractals");
46 glutDisplayFunc(display);
47 glutReshapeFunc(reshape);
48 glutKeyboardFunc(keypress);
49 glutMouseFunc(mouse);
50 glutMotionFunc(motion);
52 if(init() == -1) {
53 return 1;
54 }
55 atexit(cleanup);
57 glutMainLoop();
58 return 0;
59 }
61 int init(void)
62 {
63 glEnable(GL_DEPTH_TEST);
64 glEnable(GL_CULL_FACE);
65 glEnable(GL_LIGHTING);
66 glEnable(GL_LIGHT0);
67 /*glShadeModel(GL_FLAT);*/
69 fp = fopen("julia.obj", "wb");
71 if(!(msurf = msurf_create())) {
72 return -1;
73 }
74 msurf_eval_func(msurf, eval);
75 msurf_vertex_func(msurf, vertex);
76 msurf_set_resolution(msurf, grid_size, grid_size, grid_size);
77 msurf_set_threshold(msurf, 0.0);
78 msurf_set_inside(msurf, MSURF_LESS);
80 show_waitscr();
82 dlist = glGenLists(1);
83 glNewList(dlist, GL_COMPILE);
85 glBegin(GL_TRIANGLES);
86 msurf_polygonize(msurf);
87 glEnd();
89 glEndList();
91 if(fp) {
92 fclose(fp);
93 }
95 glClearColor(0.05, 0.05, 0.05, 1);
97 return 0;
98 }
100 void cleanup(void)
101 {
102 msurf_free(msurf);
103 }
105 void display(void)
106 {
107 glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT);
109 glMatrixMode(GL_MODELVIEW);
110 glLoadIdentity();
111 glTranslatef(0, 0, -cam_dist);
112 glRotatef(cam_phi, 1, 0, 0);
113 glRotatef(cam_theta, 0, 1, 0);
115 glCallList(dlist);
117 assert(glGetError() == GL_NO_ERROR);
118 glutSwapBuffers();
119 }
121 void reshape(int x, int y)
122 {
123 glViewport(0, 0, x, y);
124 glMatrixMode(GL_PROJECTION);
125 glLoadIdentity();
126 gluPerspective(50.0, (float)x / (float)y, 0.5, 500.0);
127 }
129 void keypress(unsigned char key, int x, int y)
130 {
131 switch(key) {
132 case 27:
133 exit(0);
134 }
135 }
137 int prev_x, prev_y;
138 int bnstate[8];
140 void mouse(int bn, int st, int x, int y)
141 {
142 bnstate[bn - GLUT_LEFT_BUTTON] = st == GLUT_DOWN;
143 prev_x = x;
144 prev_y = y;
145 }
147 void motion(int x, int y)
148 {
149 int dx = x - prev_x;
150 int dy = y - prev_y;
151 prev_x = x;
152 prev_y = y;
154 if(!dx && !dy) return;
156 if(bnstate[0]) {
157 cam_theta += dx * 0.5;
158 cam_phi += dy * 0.5;
160 if(cam_phi < -90) cam_phi = -90;
161 if(cam_phi > 90) cam_phi = 90;
162 glutPostRedisplay();
163 }
164 if(bnstate[2]) {
165 cam_dist += dy * 0.1;
167 if(cam_dist < 0.0) cam_dist = 0.0;
168 glutPostRedisplay();
169 }
170 }
172 float eval(struct metasurface *ms, float x, float y, float z)
173 {
174 struct quat q;
176 q.w = grid_scale * x;
177 q.x = grid_scale * y;
178 q.y = grid_scale * z;
179 q.z = 0.0f;
181 return julia_dist(&q, &seed, max_iter);
182 }
184 void vertex(struct metasurface *ms, float x, float y, float z)
185 {
186 struct quat q;
187 float dist;
188 float norm[3];
189 float norm_len, nscale = 1.0;
191 q.w = grid_scale * x;
192 q.x = grid_scale * y;
193 q.y = grid_scale * z;
194 q.z = 0.0f;
196 dist = julia_dist(&q, &seed, max_iter);
198 julia_grad(norm, dist, &q, &seed, max_iter);
199 norm_len = sqrt(norm[0] * norm[0] + norm[1] * norm[1] + norm[2] * norm[2]);
201 if(norm_len != 0.0f) {
202 nscale = 1.0f / norm_len;
203 }
205 glNormal3f(norm[0] * nscale, norm[1] * nscale, norm[2] * nscale);
206 glVertex3f(x, y, z);
208 if(fp) {
209 static int nverts;
211 fprintf(fp, "v %f %f %f\n", x, y, z);
213 if((++nverts) % 3 == 0) {
214 fprintf(fp, "f %d %d %d\n", nverts - 2, nverts - 1, nverts);
215 }
216 }
217 }
219 void quat_mul(struct quat *res, struct quat *a, struct quat *b)
220 {
221 float w = a->w * b->w - (a->x * b->x + a->y * b->y + a->z * b->z);
222 float x = a->w * b->x + b->w * a->x + (a->y * b->z - a->z * b->y);
223 float y = a->w * b->y + b->w * a->y + (a->z * b->x - a->x * b->z);
224 res->z = a->w * b->z + b->w * a->z + (a->x * b->y - a->y * b->x);
225 res->x = x;
226 res->y = y;
227 res->w = w;
228 }
230 void quat_sq(struct quat *q)
231 {
232 float w = q->w * q->w - (q->x * q->x + q->y * q->y + q->z * q->z);
233 float x = 2.0 * q->w * q->x;
234 float y = 2.0 * q->w * q->y;
235 q->z = 2.0 * q->w * q->z;
236 q->x = x;
237 q->y = y;
238 q->w = w;
239 }
241 float quat_lensq(struct quat *q)
242 {
243 return q->x * q->x + q->y * q->y + q->z * q->z + q->w * q->w;
244 }
246 int julia(struct quat *qres, struct quat *qprime, struct quat *q, struct quat *seed, int max_iter)
247 {
248 int i;
250 *qres = *q;
251 qprime->x = qprime->y = qprime->z = 0.0f;
252 qprime->w = 1.0f;
254 for(i=0; i<max_iter; i++) {
255 quat_mul(qprime, qres, qprime);
256 qprime->x *= 2.0;
257 qprime->y *= 2.0;
258 qprime->z *= 2.0;
259 qprime->w *= 2.0;
261 quat_sq(qres);
262 qres->x += seed->x;
263 qres->y += seed->y;
264 qres->z += seed->z;
265 qres->w += seed->w;
267 if(quat_lensq(qres) > 8.0) {
268 return 0;
269 }
270 }
272 return 1;
273 }
275 float julia_dist(struct quat *z, struct quat *seed, int max_iter)
276 {
277 struct quat qprime, q;
279 /* calc julia at z */
280 /*int inside = */julia(&q, &qprime, z, seed, max_iter);
282 float lenq = sqrt(quat_lensq(&q));
283 float lenqp = sqrt(quat_lensq(&qprime));
284 return 0.5 * lenq * log(lenq) / lenqp;
285 }
287 #define OFFS 1e-4
288 void julia_grad(float *grad, float dist, struct quat *q, struct quat *seed, int max_iter)
289 {
290 struct quat qnext = *q;
291 struct quat qprev = *q;
293 qnext.w += OFFS;
294 qprev.w -= OFFS;
295 grad[0] = julia_dist(&qnext, seed, max_iter) - julia_dist(&qprev, seed, max_iter);
296 qnext.w = qprev.w = q->w;
298 qnext.x += OFFS;
299 qprev.x -= OFFS;
300 grad[1] = julia_dist(&qnext, seed, max_iter) - julia_dist(&qprev, seed, max_iter);
301 qnext.x = qprev.x = q->x;
303 qnext.y += OFFS;
304 qprev.y -= OFFS;
305 grad[2] = julia_dist(&qnext, seed, max_iter) - julia_dist(&qprev, seed, max_iter);
306 }
308 void show_waitscr(void)
309 {
310 const char *text = "Please wait, generating fractal...";
311 glClear(GL_COLOR_BUFFER_BIT);
313 glPushAttrib(GL_ENABLE_BIT);
314 glDisable(GL_LIGHTING);
315 glDisable(GL_DEPTH_TEST);
317 glMatrixMode(GL_PROJECTION);
318 glPushMatrix();
319 glLoadIdentity();
320 glTranslatef(-0.75, 0, 0);
321 glScalef(0.00075, 0.00075, 0.00075);
322 glMatrixMode(GL_MODELVIEW);
324 glEnable(GL_LINE_SMOOTH);
325 glEnable(GL_BLEND);
326 glBlendFunc(GL_SRC_ALPHA, GL_ONE_MINUS_SRC_ALPHA);
328 glLineWidth(1.5);
330 glColor3f(1, 1, 1);
331 while(*text) {
332 glutStrokeCharacter(GLUT_STROKE_ROMAN, *text++);
333 }
335 glMatrixMode(GL_PROJECTION);
336 glPopMatrix();
338 glPopAttrib();
340 glutSwapBuffers();
341 }