10 return rand() / (RAND_MAX + 1.0);
13 typedef void (*mapfunc_t) (double x_in, double y_in, double *x_out, double *y_out, double *z_out);
14 typedef void (*colorfunc_t) (double x, double y, double z, double *r, double *g, double *b);
16 void color_test(double x, double y, double z, double *r, double *g, double *b)
18 // put in a nice function here
24 double mandelbrot_iter(double zx, double zy, double cx, double cy, int maxiter)
33 for(i = 1; i < maxiter; ++i)
36 zx = zx * zx - zy * zy + cx;
37 zy = 2 * tmp * zy + cy;
39 f = zx * zx + zy * zy;
48 // f: the greater, the more in 0 direction
49 // the smaller, the more in 1 direction
51 // the greater, the more in 0 direction
52 return i + 1 / (f - 4 + 1); // f = 16: + 0, f = 4: + 1
54 // i.e. 0 for i=1, 1 for i=maxiter
57 double mandelbrot_range(double zx, double zy, double cx, double cy, int maxiter, double offset)
59 double i = mandelbrot_iter(zx, zy, cx, cy, maxiter);
60 // map linearily 1/(offset + iter) so that:
64 // f(iter) = A/(offset + iter) + B
66 // f(0) = A/offset + B = 0
67 // f(maxiter) = A/(offset + maxiter) + B = 1
69 // 1/(1/(offset + maxiter) - 1/offset) = A
70 // B = 1 + offset / maxiter
71 // A = -offset (1 + offset / maxiter)
73 // f(iter) = -offset (1 + offset / maxiter) / (offset + iter) + 1 + offset / maxiter
74 // = -offset (1 + offset / maxiter) / (offset + iter) + 1 + offset / maxiter
75 // = iter (offset + maxiter) / maxiter (offset + iter)
76 return (i * (offset + maxiter)) / ((i + offset) * maxiter);
79 double color_mandelbrot_parms[13];
80 double mandelbrot_miniter = -1;
83 double iter_mandelbrot_raw(double x, double y, double z)
85 z -= color_mandelbrot_parms[6];
90 return mandelbrot_range(color_mandelbrot_parms[4], color_mandelbrot_parms[5], color_mandelbrot_parms[0] + x * color_mandelbrot_parms[2], color_mandelbrot_parms[1] + y * color_mandelbrot_parms[3], MAXITER, color_mandelbrot_parms[9]);
95 void iter_mandelbrot_raw_initialize_min()
97 if(mandelbrot_miniter >= 0)
99 // randomly sample 256 points
101 // set that as miniter
104 mandelbrot_miniter = MAXITER;
105 for(i = 0; i < 8192; ++i)
110 double f = sqrt(x*x + y*y + z*z);
114 double a = (z - color_mandelbrot_parms[6]) / (color_mandelbrot_parms[7] - color_mandelbrot_parms[6]);
115 a = (a - color_mandelbrot_parms[8]) / (1 - color_mandelbrot_parms[8]);
118 double iterations = iter_mandelbrot_raw(x, y, z);
121 if(iterations < mandelbrot_miniter)
122 mandelbrot_miniter = iterations;
126 void color_mandelbrot(double x, double y, double z, double *r, double *g, double *b)
128 iter_mandelbrot_raw_initialize_min();
130 double iterations = iter_mandelbrot_raw(x, y, z);
131 //printf("iter = %f\n", iterations);
132 double a = (z - color_mandelbrot_parms[6]) / (color_mandelbrot_parms[7] - color_mandelbrot_parms[6]);
133 a = (a - color_mandelbrot_parms[8]) / (1 - color_mandelbrot_parms[8]);
138 iterations = iterations * a + mandelbrot_miniter * (1-a);
139 *r = pow(iterations, color_mandelbrot_parms[10]);
140 *g = pow(iterations, color_mandelbrot_parms[11]);
141 *b = pow(iterations, color_mandelbrot_parms[12]);
144 void map_back(double x_in, double y_in, double *x_out, double *y_out, double *z_out)
146 *x_out = 2 * x_in - 1;
148 *z_out = 1 - 2 * y_in;
151 void map_right(double x_in, double y_in, double *x_out, double *y_out, double *z_out)
154 *y_out = 1 - 2 * x_in;
155 *z_out = 1 - 2 * y_in;
158 void map_front(double x_in, double y_in, double *x_out, double *y_out, double *z_out)
160 *x_out = 1 - 2 * x_in;
162 *z_out = 1 - 2 * y_in;
165 void map_left(double x_in, double y_in, double *x_out, double *y_out, double *z_out)
168 *y_out = 2 * x_in - 1;
169 *z_out = 1 - 2 * y_in;
172 void map_up(double x_in, double y_in, double *x_out, double *y_out, double *z_out)
174 *x_out = 2 * y_in - 1;
175 *y_out = 1 - 2 * x_in;
179 void map_down(double x_in, double y_in, double *x_out, double *y_out, double *z_out)
181 *x_out = 1 - 2 * y_in;
182 *y_out = 1 - 2 * x_in;
186 void writepic(colorfunc_t f, mapfunc_t m, const char *fn, int width, int height)
191 FILE *file = fopen(fn, "wb");
193 err(1, "fopen %s", fn);
195 memset(tga, 0, sizeof(tga));
196 tga[2] = 2; // uncompressed type
197 tga[12] = (width >> 0) & 0xFF;
198 tga[13] = (width >> 8) & 0xFF;
199 tga[14] = (height >> 0) & 0xFF;
200 tga[15] = (height >> 8) & 0xFF;
201 tga[16] = 24; // pixel size
203 fwrite(&tga, sizeof(tga), 1, file);
204 for(y = height-1; y >= 0; --y)
205 for(x = 0; x < width; ++x)
210 double xxx, yyy, zzz;
212 xx = (x + 0.5) / width;
213 yy = (y + 0.5) / height;
214 m(xx, yy, &xxx, &yyy, &zzz);
215 r = sqrt(xxx*xxx + yyy*yyy + zzz*zzz);
219 f(xxx, yyy, zzz, &rr, &gg, &bb);
220 rgb[2] = floor(rnd() + rr * 255);
221 rgb[1] = floor(rnd() + gg * 255);
222 rgb[0] = floor(rnd() + bb * 255);
223 fwrite(rgb, sizeof(rgb), 1, file);
229 void map_all(const char *fn, colorfunc_t f, int width, int height)
232 snprintf(buf, sizeof(buf), "%s_bk.tga", fn); buf[sizeof(buf) - 1] = 0; writepic(f, map_back, buf, width, height);
233 snprintf(buf, sizeof(buf), "%s_ft.tga", fn); buf[sizeof(buf) - 1] = 0; writepic(f, map_front, buf, width, height);
234 snprintf(buf, sizeof(buf), "%s_rt.tga", fn); buf[sizeof(buf) - 1] = 0; writepic(f, map_right, buf, width, height);
235 snprintf(buf, sizeof(buf), "%s_lf.tga", fn); buf[sizeof(buf) - 1] = 0; writepic(f, map_left, buf, width, height);
236 snprintf(buf, sizeof(buf), "%s_up.tga", fn); buf[sizeof(buf) - 1] = 0; writepic(f, map_up, buf, width, height);
237 snprintf(buf, sizeof(buf), "%s_dn.tga", fn); buf[sizeof(buf) - 1] = 0; writepic(f, map_down, buf, width, height);
240 int main(int argc, char **argv)
244 errx(1, "usage: %s filename res func parms...", *argv);
245 int res = atoi(argv[2]);
246 if(!strcmp(argv[3], "mandel"))
248 f = color_mandelbrot;
249 color_mandelbrot_parms[0] = argc<= 4 ? -0.740 : atof(argv[4]); // shift xy
250 color_mandelbrot_parms[1] = argc<= 5 ? -0.314 : atof(argv[5]);
251 color_mandelbrot_parms[2] = argc<= 6 ? -0.003 : atof(argv[6]); // mul xy
252 color_mandelbrot_parms[3] = argc<= 7 ? -0.003 : atof(argv[7]);
253 color_mandelbrot_parms[4] = argc<= 8 ? 0.420 : atof(argv[8]); // shift z
254 color_mandelbrot_parms[5] = argc<= 9 ? 0.000 : atof(argv[9]);
255 color_mandelbrot_parms[6] = argc<=10 ? -0.8 : atof(argv[10]); // horizon
256 color_mandelbrot_parms[7] = argc<=11 ? -0.7 : atof(argv[11]);
257 color_mandelbrot_parms[8] = argc<=12 ? 0.5 : atof(argv[12]);
258 color_mandelbrot_parms[9] = argc<=13 ? 400 : atof(argv[13]); // coloring
259 color_mandelbrot_parms[10] = argc<=14 ? 0.6 : atof(argv[14]);
260 color_mandelbrot_parms[11] = argc<=15 ? 0.5 : atof(argv[15]);
261 color_mandelbrot_parms[12] = argc<=16 ? 0.2 : atof(argv[16]);
267 map_all(argv[1], color_mandelbrot, res, res);