diff options
Diffstat (limited to 'apps/plugins/fractals/mandelbrot_set.c')
-rw-r--r-- | apps/plugins/fractals/mandelbrot_set.c | 392 |
1 files changed, 392 insertions, 0 deletions
diff --git a/apps/plugins/fractals/mandelbrot_set.c b/apps/plugins/fractals/mandelbrot_set.c new file mode 100644 index 0000000000..ccc65b8e91 --- /dev/null +++ b/apps/plugins/fractals/mandelbrot_set.c | |||
@@ -0,0 +1,392 @@ | |||
1 | /*************************************************************************** | ||
2 | * __________ __ ___. | ||
3 | * Open \______ \ ____ ____ | | _\_ |__ _______ ___ | ||
4 | * Source | _// _ \_/ ___\| |/ /| __ \ / _ \ \/ / | ||
5 | * Jukebox | | ( <_> ) \___| < | \_\ ( <_> > < < | ||
6 | * Firmware |____|_ /\____/ \___ >__|_ \|___ /\____/__/\_ \ | ||
7 | * \/ \/ \/ \/ \/ | ||
8 | * $Id: mandelbrot.c 24010 2009-12-15 20:51:41Z tomers $ | ||
9 | * | ||
10 | * Copyright (C) 2004 Matthias Wientapper | ||
11 | * Heavily extended 2005 Jens Arnold | ||
12 | * | ||
13 | * | ||
14 | * This program is free software; you can redistribute it and/or | ||
15 | * modify it under the terms of the GNU General Public License | ||
16 | * as published by the Free Software Foundation; either version 2 | ||
17 | * of the License, or (at your option) any later version. | ||
18 | * | ||
19 | * This software is distributed on an "AS IS" basis, WITHOUT WARRANTY OF ANY | ||
20 | * KIND, either express or implied. | ||
21 | * | ||
22 | ****************************************************************************/ | ||
23 | #include "mandelbrot_set.h" | ||
24 | |||
25 | #define BUTTON_YIELD_TIMEOUT (HZ / 4) | ||
26 | |||
27 | /* 8 entries cyclical, last entry is black (convergence) */ | ||
28 | #ifdef HAVE_LCD_COLOR | ||
29 | static const fb_data color[9] = { | ||
30 | LCD_RGBPACK(255, 0, 159), LCD_RGBPACK(159, 0, 255), LCD_RGBPACK(0, 0, 255), | ||
31 | LCD_RGBPACK(0, 159, 255), LCD_RGBPACK(0, 255, 128), LCD_RGBPACK(128, 255, 0), | ||
32 | LCD_RGBPACK(255, 191, 0), LCD_RGBPACK(255, 0, 0), LCD_RGBPACK(0, 0, 0) | ||
33 | }; | ||
34 | #else /* greyscale */ | ||
35 | static const unsigned char color[9] = { | ||
36 | 255, 223, 191, 159, 128, 96, 64, 32, 0 | ||
37 | }; | ||
38 | #endif | ||
39 | |||
40 | #if CONFIG_LCD == LCD_SSD1815 | ||
41 | /* Recorder, Ondio: pixel_height == 1.25 * pixel_width */ | ||
42 | #define MB_HEIGHT (LCD_HEIGHT*5/4) | ||
43 | #else | ||
44 | /* square pixels */ | ||
45 | #define MB_HEIGHT LCD_HEIGHT | ||
46 | #endif | ||
47 | |||
48 | #define MB_XOFS (-0x03000000L) /* -0.75 (s5.26) */ | ||
49 | #if 3000*MB_HEIGHT/LCD_WIDTH >= 2400 /* width is limiting factor */ | ||
50 | #define MB_XFAC (0x06000000LL) /* 1.5 (s5.26) */ | ||
51 | #define MB_YFAC (MB_XFAC*MB_HEIGHT/LCD_WIDTH) | ||
52 | #else /* height is limiting factor */ | ||
53 | #define MB_YFAC (0x04cccccdLL) /* 1.2 (s5.26) */ | ||
54 | #define MB_XFAC (MB_YFAC*LCD_WIDTH/MB_HEIGHT) | ||
55 | #endif | ||
56 | |||
57 | #if (LCD_DEPTH < 8) | ||
58 | #else | ||
59 | #define UPDATE_FREQ (HZ/50) | ||
60 | #endif | ||
61 | |||
62 | /* Fixed point format s5.26: sign, 5 bits integer part, 26 bits fractional part */ | ||
63 | struct mandelbrot_ctx | ||
64 | { | ||
65 | struct fractal_ops *ops; | ||
66 | long x_min; | ||
67 | long x_max; | ||
68 | long x_step; | ||
69 | long x_delta; | ||
70 | long y_min; | ||
71 | long y_max; | ||
72 | long y_step; | ||
73 | long y_delta; | ||
74 | int step_log2; | ||
75 | unsigned max_iter; | ||
76 | } ctx; | ||
77 | |||
78 | static void mandelbrot_init(void); | ||
79 | |||
80 | static int mandelbrot_calc_low_prec(struct fractal_rect *rect, | ||
81 | int (*button_yield_cb)(void *), void *button_yield_ctx); | ||
82 | |||
83 | static int mandelbrot_calc_high_prec(struct fractal_rect *rect, | ||
84 | int (*button_yield_cb)(void *), void *button_yield_ctx); | ||
85 | |||
86 | static void mandelbrot_move(int dx, int dy); | ||
87 | |||
88 | static void mandelbrot_zoom(int factor); | ||
89 | |||
90 | static int mandelbrot_precision(int d); | ||
91 | |||
92 | struct fractal_ops mandelbrot_ops = | ||
93 | { | ||
94 | .init = mandelbrot_init, | ||
95 | .calc = NULL, | ||
96 | .move = mandelbrot_move, | ||
97 | .zoom = mandelbrot_zoom, | ||
98 | .precision = mandelbrot_precision, | ||
99 | }; | ||
100 | |||
101 | static int ilog2_fp(long value) /* calculate integer log2(value_fp_6.26) */ | ||
102 | { | ||
103 | int i = 0; | ||
104 | |||
105 | if (value <= 0) | ||
106 | { | ||
107 | return -32767; | ||
108 | } | ||
109 | else if (value > (1L << 26)) | ||
110 | { | ||
111 | while (value >= (2L << 26)) | ||
112 | { | ||
113 | value >>= 1; | ||
114 | i++; | ||
115 | } | ||
116 | } | ||
117 | else | ||
118 | { | ||
119 | while (value < (1L<<26)) | ||
120 | { | ||
121 | value <<= 1; | ||
122 | i--; | ||
123 | } | ||
124 | } | ||
125 | return i; | ||
126 | } | ||
127 | |||
128 | static void recalc_parameters(void) | ||
129 | { | ||
130 | ctx.x_step = (ctx.x_max - ctx.x_min) / LCD_WIDTH; | ||
131 | ctx.x_delta = (ctx.x_step * LCD_WIDTH) / 8; | ||
132 | ctx.y_step = (ctx.y_max - ctx.y_min) / LCD_HEIGHT; | ||
133 | ctx.y_delta = (ctx.y_step * LCD_HEIGHT) / 8; | ||
134 | ctx.step_log2 = ilog2_fp(MIN(ctx.x_step, ctx.y_step)); | ||
135 | ctx.max_iter = MAX(15, -15 * ctx.step_log2 - 45); | ||
136 | |||
137 | ctx.ops->calc = (ctx.step_log2 <= -10) ? | ||
138 | mandelbrot_calc_high_prec : mandelbrot_calc_low_prec; | ||
139 | } | ||
140 | |||
141 | static void mandelbrot_init(void) | ||
142 | { | ||
143 | ctx.ops = &mandelbrot_ops; | ||
144 | |||
145 | ctx.x_min = MB_XOFS - MB_XFAC; | ||
146 | ctx.x_max = MB_XOFS + MB_XFAC; | ||
147 | ctx.y_min = -MB_YFAC; | ||
148 | ctx.y_max = MB_YFAC; | ||
149 | |||
150 | recalc_parameters(); | ||
151 | } | ||
152 | |||
153 | static int mandelbrot_calc_low_prec(struct fractal_rect *rect, | ||
154 | int (*button_yield_cb)(void *), void *button_yield_ctx) | ||
155 | { | ||
156 | #ifndef USEGSLIB | ||
157 | long next_update = *rb->current_tick; | ||
158 | int last_px = rect->px_min; | ||
159 | #endif | ||
160 | unsigned n_iter; | ||
161 | long a32, b32; | ||
162 | short x, x2, y, y2, a, b; | ||
163 | int p_x, p_y; | ||
164 | unsigned long last_yield = *rb->current_tick; | ||
165 | unsigned long last_button_yield = *rb->current_tick; | ||
166 | |||
167 | a32 = ctx.x_min + ctx.x_step * rect->px_min; | ||
168 | |||
169 | for (p_x = rect->px_min; p_x < rect->px_max; p_x++) | ||
170 | { | ||
171 | a = a32 >> 16; | ||
172 | |||
173 | b32 = ctx.y_min + ctx.y_step * (LCD_HEIGHT - rect->py_max); | ||
174 | |||
175 | for (p_y = rect->py_max - 1; p_y >= rect->py_min; p_y--) | ||
176 | { | ||
177 | b = b32 >> 16; | ||
178 | x = a; | ||
179 | y = b; | ||
180 | n_iter = 0; | ||
181 | |||
182 | while (++n_iter <= ctx.max_iter) | ||
183 | { | ||
184 | x2 = MULS16_ASR10(x, x); | ||
185 | y2 = MULS16_ASR10(y, y); | ||
186 | |||
187 | if (x2 + y2 > (4<<10)) break; | ||
188 | |||
189 | y = 2 * MULS16_ASR10(x, y) + b; | ||
190 | x = x2 - y2 + a; | ||
191 | } | ||
192 | |||
193 | if (n_iter > ctx.max_iter) | ||
194 | imgbuffer[p_y] = color[8]; | ||
195 | else | ||
196 | imgbuffer[p_y] = color[n_iter & 7]; | ||
197 | |||
198 | /* be nice to other threads: | ||
199 | * if at least one tick has passed, yield */ | ||
200 | if (TIME_AFTER(*rb->current_tick, last_yield)) | ||
201 | { | ||
202 | rb->yield(); | ||
203 | last_yield = *rb->current_tick; | ||
204 | } | ||
205 | |||
206 | if (TIME_AFTER(*rb->current_tick, last_button_yield)) | ||
207 | { | ||
208 | if (button_yield_cb(button_yield_ctx)) | ||
209 | { | ||
210 | #ifndef USEGSLIB | ||
211 | /* update screen part that was changed since last yield */ | ||
212 | rb->lcd_update_rect(last_px, rect->py_min, | ||
213 | p_x - last_px + 1, rect->py_max - rect->py_min); | ||
214 | #endif | ||
215 | rect->px_min = (p_x == 0) ? 0 : p_x - 1; | ||
216 | |||
217 | return 1; | ||
218 | } | ||
219 | |||
220 | last_button_yield = *rb->current_tick + BUTTON_YIELD_TIMEOUT; | ||
221 | } | ||
222 | |||
223 | b32 += ctx.y_step; | ||
224 | } | ||
225 | #ifdef USEGSLIB | ||
226 | grey_ub_gray_bitmap_part(imgbuffer, 0, rect->py_min, 1, | ||
227 | p_x, rect->py_min, 1, rect->py_max - rect->py_min); | ||
228 | #else | ||
229 | rb->lcd_bitmap_part(imgbuffer, 0, rect->py_min, 1, | ||
230 | p_x, rect->py_min, 1, rect->py_max - rect->py_min); | ||
231 | |||
232 | if ((p_x == rect->px_max - 1) || | ||
233 | TIME_AFTER(*rb->current_tick, next_update)) | ||
234 | { | ||
235 | next_update = *rb->current_tick + UPDATE_FREQ; | ||
236 | |||
237 | /* update screen part that was changed since last yield */ | ||
238 | rb->lcd_update_rect(last_px, rect->py_min, | ||
239 | p_x - last_px + 1, rect->py_max - rect->py_min); | ||
240 | last_px = p_x; | ||
241 | } | ||
242 | #endif | ||
243 | |||
244 | a32 += ctx.x_step; | ||
245 | } | ||
246 | |||
247 | rect->valid = 0; | ||
248 | |||
249 | return 0; | ||
250 | } | ||
251 | |||
252 | static int mandelbrot_calc_high_prec(struct fractal_rect *rect, | ||
253 | int (*button_yield_cb)(void *), void *button_yield_ctx) | ||
254 | { | ||
255 | #ifndef USEGSLIB | ||
256 | long next_update = *rb->current_tick; | ||
257 | int last_px = rect->px_min; | ||
258 | #endif | ||
259 | unsigned n_iter; | ||
260 | long x, x2, y, y2, a, b; | ||
261 | int p_x, p_y; | ||
262 | unsigned long last_yield = *rb->current_tick; | ||
263 | unsigned long last_button_yield = *rb->current_tick; | ||
264 | |||
265 | MULS32_INIT(); | ||
266 | |||
267 | a = ctx.x_min + ctx.x_step * rect->px_min; | ||
268 | |||
269 | for (p_x = rect->px_min; p_x < rect->px_max; p_x++) | ||
270 | { | ||
271 | b = ctx.y_min + ctx.y_step * (LCD_HEIGHT - rect->py_max); | ||
272 | |||
273 | for (p_y = rect->py_max - 1; p_y >= rect->py_min; p_y--) | ||
274 | { | ||
275 | x = a; | ||
276 | y = b; | ||
277 | n_iter = 0; | ||
278 | |||
279 | while (++n_iter <= ctx.max_iter) | ||
280 | { | ||
281 | x2 = MULS32_ASR26(x, x); | ||
282 | y2 = MULS32_ASR26(y, y); | ||
283 | |||
284 | if (x2 + y2 > (4L<<26)) break; | ||
285 | |||
286 | y = 2 * MULS32_ASR26(x, y) + b; | ||
287 | x = x2 - y2 + a; | ||
288 | } | ||
289 | |||
290 | if (n_iter > ctx.max_iter) | ||
291 | imgbuffer[p_y] = color[8]; | ||
292 | else | ||
293 | imgbuffer[p_y] = color[n_iter & 7]; | ||
294 | |||
295 | /* be nice to other threads: | ||
296 | * if at least one tick has passed, yield */ | ||
297 | if (TIME_AFTER(*rb->current_tick, last_yield)) | ||
298 | { | ||
299 | rb->yield(); | ||
300 | last_yield = *rb->current_tick; | ||
301 | } | ||
302 | |||
303 | if (TIME_AFTER(*rb->current_tick, last_button_yield)) | ||
304 | { | ||
305 | if (button_yield_cb(button_yield_ctx)) | ||
306 | { | ||
307 | #ifndef USEGSLIB | ||
308 | /* update screen part that was changed since last yield */ | ||
309 | rb->lcd_update_rect(last_px, rect->py_min, | ||
310 | p_x - last_px + 1, rect->py_max - rect->py_min); | ||
311 | #endif | ||
312 | rect->px_min = (p_x == 0) ? 0 : p_x - 1; | ||
313 | |||
314 | return 1; | ||
315 | } | ||
316 | |||
317 | last_button_yield = *rb->current_tick + BUTTON_YIELD_TIMEOUT; | ||
318 | } | ||
319 | |||
320 | b += ctx.y_step; | ||
321 | } | ||
322 | #ifdef USEGSLIB | ||
323 | grey_ub_gray_bitmap_part(imgbuffer, 0, rect->py_min, 1, | ||
324 | p_x, rect->py_min, 1, rect->py_max - rect->py_min); | ||
325 | #else | ||
326 | rb->lcd_bitmap_part(imgbuffer, 0, rect->py_min, 1, | ||
327 | p_x, rect->py_min, 1, rect->py_max - rect->py_min); | ||
328 | |||
329 | if ((p_x == rect->px_max - 1) || | ||
330 | TIME_AFTER(*rb->current_tick, next_update)) | ||
331 | { | ||
332 | next_update = *rb->current_tick + UPDATE_FREQ; | ||
333 | |||
334 | /* update screen part that was changed since last yield */ | ||
335 | rb->lcd_update_rect(last_px, rect->py_min, | ||
336 | p_x - last_px + 1, rect->py_max - rect->py_min); | ||
337 | last_px = p_x; | ||
338 | } | ||
339 | #endif | ||
340 | a += ctx.x_step; | ||
341 | } | ||
342 | |||
343 | rect->valid = 0; | ||
344 | |||
345 | return 0; | ||
346 | } | ||
347 | |||
348 | static void mandelbrot_move(int dx, int dy) | ||
349 | { | ||
350 | long d_x = (long)dx * ctx.x_delta; | ||
351 | long d_y = (long)dy * ctx.y_delta; | ||
352 | |||
353 | ctx.x_min += d_x; | ||
354 | ctx.x_max += d_x; | ||
355 | ctx.y_min += d_y; | ||
356 | ctx.y_max += d_y; | ||
357 | } | ||
358 | |||
359 | static void mandelbrot_zoom(int factor) | ||
360 | { | ||
361 | long factor_x = (long)factor * ctx.x_delta; | ||
362 | long factor_y = (long)factor * ctx.y_delta; | ||
363 | |||
364 | ctx.x_min += factor_x; | ||
365 | ctx.x_max -= factor_x; | ||
366 | ctx.y_min += factor_y; | ||
367 | ctx.y_max -= factor_y; | ||
368 | |||
369 | recalc_parameters(); | ||
370 | } | ||
371 | |||
372 | static int mandelbrot_precision(int d) | ||
373 | { | ||
374 | int changed = 0; | ||
375 | |||
376 | /* Precision increase */ | ||
377 | for (; d > 0; d--) | ||
378 | { | ||
379 | ctx.max_iter += ctx.max_iter / 2; | ||
380 | changed = 1; | ||
381 | } | ||
382 | |||
383 | /* Precision decrease */ | ||
384 | for (; d < 0 && ctx.max_iter >= 15; d++) | ||
385 | { | ||
386 | ctx.max_iter -= ctx.max_iter / 3; | ||
387 | changed = 1; | ||
388 | } | ||
389 | |||
390 | return changed; | ||
391 | } | ||
392 | |||