1 /* GIMP - The GNU Image Manipulation Program
2 * Copyright (C) 1995 Spencer Kimball and Peter Mattis
3 *
4 * IfsCompose is a interface for creating IFS fractals by
5 * direct manipulation.
6 * Copyright (C) 1997 Owen Taylor
7 *
8 * it under the terms of the GNU General Public License as published by
9 * the Free Software Foundation; either version 3 of the License, or
10 * (at your option) any later version.
11 *
12 * This program is distributed in the hope that it will be useful,
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 * GNU General Public License for more details.
16 *
17 * You should have received a copy of the GNU General Public License
18 * along with this program. If not, see <https://www.gnu.org/licenses/>.
19 */
20
21 #include "config.h"
22
23 #include <stdlib.h>
24 #include <string.h>
25
26 #include <gdk/gdk.h>
27
28 #include <libgimp/gimp.h>
29
30 #include "ifs-compose.h"
31
32
33 typedef struct
34 {
35 GdkPoint point;
36 gdouble angle;
37 } SortPoint;
38
39
40 /* local functions */
41 static void aff_element_compute_click_boundary (AffElement *elem,
42 gint num_elements,
43 gdouble *points_x,
44 gdouble *points_y);
45 static guchar * create_brush (IfsComposeVals *ifsvals,
46 gint *brush_size,
47 gdouble *brush_offset);
48
49
50 void
aff2_translate(Aff2 * naff,gdouble x,gdouble y)51 aff2_translate (Aff2 *naff,
52 gdouble x,
53 gdouble y)
54 {
55 naff->a11 = 1.0;
56 naff->a12 = 0;
57 naff->a21 = 0;
58 naff->a22 = 1.0;
59 naff->b1 = x;
60 naff->b2 = y;
61 }
62
63 void
aff2_rotate(Aff2 * naff,gdouble theta)64 aff2_rotate (Aff2 *naff,
65 gdouble theta)
66 {
67 naff->a11 = cos(theta);
68 naff->a12 = sin(theta);
69 naff->a21 = -naff->a12;
70 naff->a22 = naff->a11;
71 naff->b1 = 0;
72 naff->b2 = 0;
73 }
74
75 void
aff2_scale(Aff2 * naff,gdouble s,gboolean flip)76 aff2_scale (Aff2 *naff,
77 gdouble s,
78 gboolean flip)
79 {
80 if (flip)
81 naff->a11 = -s;
82 else
83 naff->a11 = s;
84
85 naff->a12 = 0;
86 naff->a21 = 0;
87 naff->a22 = s;
88 naff->b1 = 0;
89 naff->b2 = 0;
90 }
91
92 /* Create a unitary transform with given x-y asymmetry and shear */
93 void
aff2_distort(Aff2 * naff,gdouble asym,gdouble shear)94 aff2_distort (Aff2 *naff,
95 gdouble asym,
96 gdouble shear)
97 {
98 naff->a11 = asym;
99 naff->a22 = 1/asym;
100 naff->a12 = shear;
101 naff->a21 = 0;
102 naff->b1 = 0;
103 naff->b2 = 0;
104 }
105
106 /* Find a pure stretch in some direction that brings xo,yo to xn,yn */
107 void
aff2_compute_stretch(Aff2 * naff,gdouble xo,gdouble yo,gdouble xn,gdouble yn)108 aff2_compute_stretch (Aff2 *naff,
109 gdouble xo,
110 gdouble yo,
111 gdouble xn,
112 gdouble yn)
113 {
114 gdouble denom = xo*xn + yo*yn;
115
116 if (denom == 0.0) /* singular */
117 {
118 naff->a11 = 1.0;
119 naff->a12 = 0.0;
120 naff->a21 = 0.0;
121 naff->a22 = 1.0;
122 }
123 else
124 {
125 naff->a11 = (SQR(xn) + SQR(yo)) / denom;
126 naff->a22 = (SQR(xo) + SQR(yn)) / denom;
127 naff->a12 = naff->a21 = (xn * yn - xo * yo) / denom;
128 }
129
130 naff->b1 = 0.0;
131 naff->b2 = 0.0;
132 }
133
134 void
aff2_compose(Aff2 * naff,Aff2 * aff1,Aff2 * aff2)135 aff2_compose (Aff2 *naff,
136 Aff2 *aff1,
137 Aff2 *aff2)
138 {
139 naff->a11 = aff1->a11 * aff2->a11 + aff1->a12 * aff2->a21;
140 naff->a12 = aff1->a11 * aff2->a12 + aff1->a12 * aff2->a22;
141 naff->b1 = aff1->a11 * aff2->b1 + aff1->a12 * aff2->b2 + aff1->b1;
142 naff->a21 = aff1->a21 * aff2->a11 + aff1->a22 * aff2->a21;
143 naff->a22 = aff1->a21 * aff2->a12 + aff1->a22 * aff2->a22;
144 naff->b2 = aff1->a21 * aff2->b1 + aff1->a22 * aff2->b2 + aff1->b2;
145 }
146
147 /* Returns the identity matrix if the original matrix was singular */
148 void
aff2_invert(Aff2 * naff,Aff2 * aff)149 aff2_invert (Aff2 *naff,
150 Aff2 *aff)
151 {
152 gdouble det = aff->a11 * aff->a22 - aff->a12 * aff->a21;
153
154 if (det==0)
155 {
156 aff2_scale (naff, 1.0, 0);
157 }
158 else
159 {
160 naff->a11 = aff->a22 / det;
161 naff->a22 = aff->a11 / det;
162 naff->a21 = - aff->a21 / det;
163 naff->a12 = - aff->a12 / det;
164 naff->b1 = - naff->a11 * aff->b1 - naff->a12 * aff->b2;
165 naff->b2 = - naff->a21 * aff->b1 - naff->a22 * aff->b2;
166 }
167 }
168
169 void
aff2_apply(Aff2 * aff,gdouble x,gdouble y,gdouble * xf,gdouble * yf)170 aff2_apply (Aff2 *aff,
171 gdouble x,
172 gdouble y,
173 gdouble *xf,
174 gdouble *yf)
175 {
176 gdouble xt = aff->a11 * x + aff->a12 * y + aff->b1;
177 gdouble yt = aff->a21 * x + aff->a22 * y + aff->b2;
178
179 *xf = xt;
180 *yf = yt;
181 }
182
183 /* Find the fixed point of an affine transformation
184 (Will return garbage for pure translations) */
185
186 void
aff2_fixed_point(Aff2 * aff,gdouble * xf,gdouble * yf)187 aff2_fixed_point (Aff2 *aff,
188 gdouble *xf,
189 gdouble *yf)
190 {
191 Aff2 t1,t2;
192
193 t1.a11 = 1 - aff->a11;
194 t1.a22 = 1 - aff->a22;
195 t1.a12 = -aff->a12;
196 t1.a21 = -aff->a21;
197 t1.b1 = 0;
198 t1.b2 = 0;
199
200 aff2_invert (&t2, &t1);
201 aff2_apply (&t2, aff->b1, aff->b2, xf, yf);
202 }
203
204 void
aff3_apply(Aff3 * t,gdouble x,gdouble y,gdouble z,gdouble * xf,gdouble * yf,gdouble * zf)205 aff3_apply (Aff3 *t,
206 gdouble x,
207 gdouble y,
208 gdouble z,
209 gdouble *xf,
210 gdouble *yf,
211 gdouble *zf)
212 {
213 gdouble xt = (t->vals[0][0] * x +
214 t->vals[0][1] * y +
215 t->vals[0][2] * z + t->vals[0][3]);
216 gdouble yt = (t->vals[1][0] * x +
217 t->vals[1][1] * y +
218 t->vals[1][2] * z + t->vals[1][3]);
219 gdouble zt = (t->vals[2][0] * x +
220 t->vals[2][1] * y +
221 t->vals[2][2] * z + t->vals[2][3]);
222
223 *xf = xt;
224 *yf = yt;
225 *zf = zt;
226 }
227
228 static int
ipolygon_sort_func(const void * a,const void * b)229 ipolygon_sort_func (const void *a,
230 const void *b)
231 {
232 if (((SortPoint *)a)->angle < ((SortPoint *)b)->angle)
233 return -1;
234 else if (((SortPoint *)a)->angle > ((SortPoint *)b)->angle)
235 return 1;
236 else
237 return 0;
238 }
239
240 /* Return a newly-allocated polygon which is the convex hull
241 of the given polygon.
242
243 Uses the Graham scan. see
244 http://www.cs.curtin.edu.au/units/cg201/notes/node77.html
245
246 for a description
247 */
248
249 IPolygon *
ipolygon_convex_hull(IPolygon * poly)250 ipolygon_convex_hull (IPolygon *poly)
251 {
252 gint num_new = poly->npoints;
253 GdkPoint *new_points = g_new (GdkPoint, num_new);
254 SortPoint *sort_points = g_new (SortPoint, num_new);
255 IPolygon *new_poly = g_new (IPolygon, 1);
256
257 gint i, j;
258 gint x1, x2, y1, y2;
259 gint lowest;
260 GdkPoint lowest_pt;
261
262 new_poly->points = new_points;
263 if (num_new <= 3)
264 {
265 memcpy (new_points, poly->points, num_new * sizeof (GdkPoint));
266 new_poly->npoints = num_new;
267 g_free (sort_points);
268 return new_poly;
269 }
270
271 /* scan for the lowest point */
272 lowest_pt = poly->points[0];
273 lowest = 0;
274
275 for (i = 1; i < num_new; i++)
276 if (poly->points[i].y < lowest_pt.y)
277 {
278 lowest_pt = poly->points[i];
279 lowest = i;
280 }
281
282 /* sort by angle from lowest point */
283
284 for (i = 0, j = 0; i < num_new; i++, j++)
285 {
286 if (i==lowest)
287 j--;
288 else
289 {
290 gdouble dy = poly->points[i].y - lowest_pt.y;
291 gdouble dx = poly->points[i].x - lowest_pt.x;
292
293 if (dy == 0 && dx == 0)
294 {
295 j--;
296 num_new--;
297 continue;
298 }
299 sort_points[j].point = poly->points[i];
300 sort_points[j].angle = atan2 (dy, dx);
301 }
302 }
303
304 qsort (sort_points, num_new - 1, sizeof (SortPoint), ipolygon_sort_func);
305
306 /* now ensure that all turns as we trace the perimeter are
307 counter-clockwise */
308
309 new_points[0] = lowest_pt;
310 new_points[1] = sort_points[0].point;
311 x1 = new_points[1].x - new_points[0].x;
312 y1 = new_points[1].y - new_points[0].y;
313
314 for (i = 1, j = 2; j < num_new; i++, j++)
315 {
316 x2 = sort_points[i].point.x - new_points[j - 1].x;
317 y2 = sort_points[i].point.y - new_points[j - 1].y;
318
319 if (x2 == 0 && y2 == 0)
320 {
321 num_new--;
322 j--;
323 continue;
324 }
325
326 while (x1 * y2 - x2 * y1 < 0) /* clockwise rotation */
327 {
328 num_new--;
329 j--;
330 x1 = new_points[j - 1].x - new_points[j - 2].x;
331 y1 = new_points[j - 1].y - new_points[j - 2].y;
332 x2 = sort_points[i].point.x - new_points[j - 1].x;
333 y2 = sort_points[i].point.y - new_points[j - 1].y;
334 }
335 new_points[j] = sort_points[i].point;
336 x1 = x2;
337 y1 = y2;
338 }
339
340 g_free (sort_points);
341
342 new_poly->npoints = num_new;
343
344 return new_poly;
345 }
346
347 /* Determines whether a specified point is in the given polygon.
348 Based on
349
350 inpoly.c by Bob Stein and Craig Yap.
351
352 (Linux Journal, Issue 35 (March 1997), p 68)
353 */
354
355 gint
ipolygon_contains(IPolygon * poly,gint xt,gint yt)356 ipolygon_contains (IPolygon *poly,
357 gint xt,
358 gint yt)
359 {
360 gint xnew, ynew;
361 gint xold, yold;
362 gint x1,y1;
363 gint x2,y2;
364
365 gint i;
366 gint inside = 0;
367
368 if (poly->npoints < 3)
369 return 0;
370
371 xold=poly->points[poly->npoints - 1].x;
372 yold=poly->points[poly->npoints - 1].y;
373 for (i = 0; i < poly->npoints; i++)
374 {
375 xnew = poly->points[i].x;
376 ynew = poly->points[i].y;
377 if (xnew > xold)
378 {
379 x1 = xold;
380 x2 = xnew;
381 y1 = yold;
382 y2 = ynew;
383 }
384 else
385 {
386 x1 = xnew;
387 x2 = xold;
388 y1 = ynew;
389 y2 = yold;
390 }
391 if ((xnew < xt) == (xt <= xold) &&
392 (yt - y1)*(x2 - x1) < (y2 - y1)*(xt - x1))
393 inside = !inside;
394 xold = xnew;
395 yold = ynew;
396 }
397 return inside;
398 }
399
400 void
aff_element_compute_color_trans(AffElement * elem)401 aff_element_compute_color_trans (AffElement *elem)
402 {
403 int i, j;
404
405 if (elem->v.simple_color)
406 {
407 gdouble mag2;
408
409 mag2 = SQR (elem->v.target_color.r);
410 mag2 += SQR (elem->v.target_color.g);
411 mag2 += SQR (elem->v.target_color.b);
412
413 /* For mag2 == 0, the transformation blows up in general
414 but is well defined for hue_scale == value_scale, so
415 we assume that special case. */
416 if (mag2 == 0)
417 for (i = 0; i < 3; i++)
418 {
419 for (j = 0; j < 4; j++)
420 elem->color_trans.vals[i][j] = 0.0;
421
422 elem->color_trans.vals[i][i] = elem->v.hue_scale;
423 }
424 else
425 {
426 /* red */
427 for (j = 0; j < 3; j++)
428 {
429 elem->color_trans.vals[0][j] = elem->v.target_color.r
430 / mag2 * (elem->v.value_scale - elem->v.hue_scale);
431 }
432
433 /* green */
434 for (j = 0; j < 3; j++)
435 {
436 elem->color_trans.vals[1][j] = elem->v.target_color.g
437 / mag2 * (elem->v.value_scale - elem->v.hue_scale);
438 }
439
440 /* blue */
441 for (j = 0; j < 3; j++)
442 {
443 elem->color_trans.vals[2][j] = elem->v.target_color.g
444 / mag2 * (elem->v.value_scale - elem->v.hue_scale);
445 }
446
447 elem->color_trans.vals[0][0] += elem->v.hue_scale;
448 elem->color_trans.vals[1][1] += elem->v.hue_scale;
449 elem->color_trans.vals[2][2] += elem->v.hue_scale;
450
451 elem->color_trans.vals[0][3] =
452 (1 - elem->v.value_scale) * elem->v.target_color.r;
453 elem->color_trans.vals[1][3] =
454 (1 - elem->v.value_scale) * elem->v.target_color.g;
455 elem->color_trans.vals[2][3] =
456 (1 - elem->v.value_scale) * elem->v.target_color.b;
457
458 }
459
460
461 aff3_apply (&elem->color_trans, 1.0, 0.0, 0.0,
462 &elem->v.red_color.r,
463 &elem->v.red_color.g,
464 &elem->v.red_color.b);
465 aff3_apply (&elem->color_trans, 0.0, 1.0, 0.0,
466 &elem->v.green_color.r,
467 &elem->v.green_color.g,
468 &elem->v.green_color.b);
469 aff3_apply (&elem->color_trans, 0.0, 0.0, 1.0,
470 &elem->v.blue_color.r,
471 &elem->v.blue_color.g,
472 &elem->v.blue_color.b);
473 aff3_apply (&elem->color_trans, 0.0, 0.0, 0.0,
474 &elem->v.black_color.r,
475 &elem->v.black_color.g,
476 &elem->v.black_color.b);
477 }
478 else
479 {
480 elem->color_trans.vals[0][0] =
481 elem->v.red_color.r - elem->v.black_color.r;
482 elem->color_trans.vals[1][0] =
483 elem->v.red_color.g - elem->v.black_color.g;
484 elem->color_trans.vals[2][0] =
485 elem->v.red_color.b - elem->v.black_color.b;
486
487 elem->color_trans.vals[0][1] =
488 elem->v.green_color.r - elem->v.black_color.r;
489 elem->color_trans.vals[1][1] =
490 elem->v.green_color.g - elem->v.black_color.g;
491 elem->color_trans.vals[2][1] =
492 elem->v.green_color.b - elem->v.black_color.b;
493
494 elem->color_trans.vals[0][2] =
495 elem->v.blue_color.r - elem->v.black_color.r;
496 elem->color_trans.vals[1][2] =
497 elem->v.blue_color.g - elem->v.black_color.g;
498 elem->color_trans.vals[2][2] =
499 elem->v.blue_color.b - elem->v.black_color.b;
500
501 elem->color_trans.vals[0][3] = elem->v.black_color.r;
502 elem->color_trans.vals[1][3] = elem->v.black_color.g;
503 elem->color_trans.vals[2][3] = elem->v.black_color.b;
504 }
505 }
506
507 void
aff_element_compute_trans(AffElement * elem,gdouble width,gdouble height,gdouble center_x,gdouble center_y)508 aff_element_compute_trans (AffElement *elem,
509 gdouble width,
510 gdouble height,
511 gdouble center_x,
512 gdouble center_y)
513 {
514 Aff2 t1, t2, t3;
515
516 /* create the rotation, scaling and shearing part of the transform */
517 aff2_distort (&t1, elem->v.asym, elem->v.shear);
518 aff2_scale (&t2, elem->v.scale, elem->v.flip);
519 aff2_compose (&t3, &t2, &t1);
520 aff2_rotate (&t2, elem->v.theta);
521 aff2_compose (&t1, &t2, &t3);
522
523 /* now create the translational part */
524 aff2_translate (&t2, -center_x*width, -center_y*width);
525 aff2_compose (&t3, &t1, &t2);
526 aff2_translate (&t2, elem->v.x*width, elem->v.y*width);
527 aff2_compose (&elem->trans, &t2, &t3);
528 }
529
530 void
aff_element_decompose_trans(AffElement * elem,Aff2 * aff,gdouble width,gdouble height,gdouble center_x,gdouble center_y)531 aff_element_decompose_trans (AffElement *elem,
532 Aff2 *aff,
533 gdouble width,
534 gdouble height,
535 gdouble center_x,
536 gdouble center_y)
537 {
538 Aff2 t1, t2;
539 gdouble det, scale, sign;
540
541 /* pull of the translational parts */
542 aff2_translate (&t1,center_x * width, center_y * width);
543 aff2_compose (&t2, aff, &t1);
544
545 elem->v.x = t2.b1 / width;
546 elem->v.y = t2.b2 / width;
547
548 det = t2.a11 * t2.a22 - t2.a12 * t2.a21;
549
550 if (det == 0.0)
551 {
552 elem->v.scale = 0.0;
553 elem->v.theta = 0.0;
554 elem->v.asym = 1.0;
555 elem->v.shear = 0.0;
556 elem->v.flip = 0;
557 }
558 else
559 {
560 if (det >= 0)
561 {
562 scale = elem->v.scale = sqrt (det);
563 sign = 1;
564 elem->v.flip = 0;
565 }
566 else
567 {
568 scale = elem->v.scale = sqrt (-det);
569 sign = -1;
570 elem->v.flip = 1;
571 }
572
573 elem->v.theta = atan2 (-t2.a21, t2.a11);
574
575 if (cos (elem->v.theta) == 0.0)
576 {
577 elem->v.asym = - t2.a21 / scale / sin (elem->v.theta);
578 elem->v.shear = - sign * t2.a22 / scale / sin (elem->v.theta);
579 }
580 else
581 {
582 elem->v.asym = sign * t2.a11 / scale / cos (elem->v.theta);
583 elem->v.shear = sign *
584 (t2.a12/scale - sin (elem->v.theta)/elem->v.asym)
585 / cos (elem->v.theta);
586 }
587 }
588 }
589
590 static void
aff_element_compute_click_boundary(AffElement * elem,int num_elements,gdouble * points_x,gdouble * points_y)591 aff_element_compute_click_boundary (AffElement *elem,
592 int num_elements,
593 gdouble *points_x,
594 gdouble *points_y)
595 {
596 gint i;
597 gdouble xtot = 0;
598 gdouble ytot = 0;
599 gdouble xc, yc;
600 gdouble theta;
601 gdouble sth, cth; /* sin(theta), cos(theta) */
602 gdouble axis1, axis2;
603 gdouble axis1max, axis2max, axis1min, axis2min;
604
605 /* compute the center of mass of the points */
606 for (i = 0; i < num_elements; i++)
607 {
608 xtot += points_x[i];
609 ytot += points_y[i];
610 }
611 xc = xtot / num_elements;
612 yc = ytot / num_elements;
613
614 /* compute the sum of the (x+iy)^2, and take half the the resulting
615 angle (xtot+iytot = A*exp(2i*theta)), to get an average direction */
616
617 xtot = 0;
618 ytot = 0;
619 for (i = 0; i < num_elements; i++)
620 {
621 xtot += SQR (points_x[i] - xc) - SQR (points_y[i] - yc);
622 ytot += 2 * (points_x[i] - xc) * (points_y[i] - yc);
623 }
624 theta = 0.5 * atan2 (ytot, xtot);
625 sth = sin (theta);
626 cth = cos (theta);
627
628 /* compute the minimum rectangle at angle theta that bounds the points,
629 1/2 side lengths left in axis1, axis2, center in xc, yc */
630
631 axis1max = axis1min = 0.0;
632 axis2max = axis2min = 0.0;
633 for (i = 0; i < num_elements; i++)
634 {
635 gdouble proj1 = (points_x[i] - xc) * cth + (points_y[i] - yc) * sth;
636 gdouble proj2 = -(points_x[i] - xc) * sth + (points_y[i] - yc) * cth;
637 if (proj1 < axis1min)
638 axis1min = proj1;
639 if (proj1 > axis1max)
640 axis1max = proj1;
641 if (proj2 < axis2min)
642 axis2min = proj2;
643 if (proj2 > axis2max)
644 axis2max = proj2;
645 }
646 axis1 = 0.5 * (axis1max - axis1min);
647 axis2 = 0.5 * (axis2max - axis2min);
648 xc += 0.5 * ((axis1max + axis1min) * cth - (axis2max + axis2min) * sth);
649 yc += 0.5 * ((axis1max + axis1min) * sth + (axis2max + axis2min) * cth);
650
651 /* if the the rectangle is less than 10 pixels in any dimension,
652 make it click_boundary, otherwise set click_boundary = draw_boundary */
653
654 if (axis1 < 8.0 || axis2 < 8.0)
655 {
656 GdkPoint *points = g_new (GdkPoint, 4);
657
658 elem->click_boundary = g_new (IPolygon, 1);
659 elem->click_boundary->points = points;
660 elem->click_boundary->npoints = 4;
661
662 if (axis1 < 8.0) axis1 = 8.0;
663 if (axis2 < 8.0) axis2 = 8.0;
664
665 points[0].x = xc + axis1 * cth - axis2 * sth;
666 points[0].y = yc + axis1 * sth + axis2 * cth;
667 points[1].x = xc - axis1 * cth - axis2 * sth;
668 points[1].y = yc - axis1 * sth + axis2 * cth;
669 points[2].x = xc - axis1 * cth + axis2 * sth;
670 points[2].y = yc - axis1 * sth - axis2 * cth;
671 points[3].x = xc + axis1 * cth + axis2 * sth;
672 points[3].y = yc + axis1 * sth - axis2 * cth;
673 }
674 else
675 elem->click_boundary = elem->draw_boundary;
676 }
677
678 void
aff_element_compute_boundary(AffElement * elem,gint width,gint height,AffElement ** elements,gint num_elements)679 aff_element_compute_boundary (AffElement *elem,
680 gint width,
681 gint height,
682 AffElement **elements,
683 gint num_elements)
684 {
685 gint i;
686 IPolygon tmp_poly;
687 gdouble *points_x;
688 gdouble *points_y;
689
690 if (elem->click_boundary && elem->click_boundary != elem->draw_boundary)
691 g_free (elem->click_boundary);
692 if (elem->draw_boundary)
693 g_free (elem->draw_boundary);
694
695 tmp_poly.npoints = num_elements;
696 tmp_poly.points = g_new (GdkPoint, num_elements);
697 points_x = g_new (gdouble, num_elements);
698 points_y = g_new (gdouble, num_elements);
699
700 for (i = 0; i < num_elements; i++)
701 {
702 aff2_apply (&elem->trans,
703 elements[i]->v.x * width, elements[i]->v.y * width,
704 &points_x[i],&points_y[i]);
705 tmp_poly.points[i].x = (gint)points_x[i];
706 tmp_poly.points[i].y = (gint)points_y[i];
707 }
708
709 elem->draw_boundary = ipolygon_convex_hull (&tmp_poly);
710 aff_element_compute_click_boundary (elem, num_elements, points_x, points_y);
711
712 g_free (tmp_poly.points);
713 }
714
715 void
aff_element_draw(AffElement * elem,gboolean selected,gint width,gint height,cairo_t * cr,GdkColor * color,PangoLayout * layout)716 aff_element_draw (AffElement *elem,
717 gboolean selected,
718 gint width,
719 gint height,
720 cairo_t *cr,
721 GdkColor *color,
722 PangoLayout *layout)
723 {
724 PangoRectangle rect;
725 gint i;
726
727 pango_layout_set_text (layout, elem->name, -1);
728 pango_layout_get_pixel_extents (layout, NULL, &rect);
729
730 gdk_cairo_set_source_color (cr, color);
731
732 cairo_move_to (cr,
733 elem->v.x * width - rect.width / 2,
734 elem->v.y * width + rect.height / 2);
735 pango_cairo_show_layout (cr, layout);
736 cairo_fill (cr);
737
738 cairo_set_line_width (cr, 1.0);
739
740 if (elem->click_boundary != elem->draw_boundary)
741 {
742 cairo_move_to (cr,
743 elem->click_boundary->points[0].x,
744 elem->click_boundary->points[0].y);
745
746 for (i = 1; i < elem->click_boundary->npoints; i++)
747 cairo_line_to (cr,
748 elem->click_boundary->points[i].x,
749 elem->click_boundary->points[i].y);
750
751 cairo_close_path (cr);
752
753 cairo_stroke (cr);
754 }
755
756 if (selected)
757 cairo_set_line_width (cr, 3.0);
758
759 cairo_move_to (cr,
760 elem->draw_boundary->points[0].x,
761 elem->draw_boundary->points[0].y);
762
763 for (i = 1; i < elem->draw_boundary->npoints; i++)
764 cairo_line_to (cr,
765 elem->draw_boundary->points[i].x,
766 elem->draw_boundary->points[i].y);
767
768 cairo_close_path (cr);
769
770 cairo_stroke (cr);
771 }
772
773 AffElement *
aff_element_new(gdouble x,gdouble y,GimpRGB * color,gint count)774 aff_element_new (gdouble x,
775 gdouble y,
776 GimpRGB *color,
777 gint count)
778 {
779 AffElement *elem = g_new (AffElement, 1);
780 gchar buffer[16];
781
782 elem->v.x = x;
783 elem->v.y = y;
784 elem->v.theta = 0.0;
785 elem->v.scale = 0.5;
786 elem->v.asym = 1.0;
787 elem->v.shear = 0.0;
788 elem->v.flip = 0;
789
790 elem->v.red_color = *color;
791 elem->v.blue_color = *color;
792 elem->v.green_color = *color;
793 elem->v.black_color = *color;
794
795 elem->v.target_color = *color;
796 elem->v.hue_scale = 0.5;
797 elem->v.value_scale = 0.5;
798
799 elem->v.simple_color = TRUE;
800
801 elem->draw_boundary = NULL;
802 elem->click_boundary = NULL;
803
804 aff_element_compute_color_trans (elem);
805
806 elem->v.prob = 1.0;
807
808 sprintf (buffer,"%d", count);
809 elem->name = g_strdup (buffer);
810
811 return elem;
812 }
813
814 void
aff_element_free(AffElement * elem)815 aff_element_free (AffElement *elem)
816 {
817 if (elem->click_boundary != elem->draw_boundary)
818 g_free (elem->click_boundary);
819
820 g_free (elem->draw_boundary);
821 g_free (elem);
822 }
823
824 #ifdef DEBUG_BRUSH
825 static brush_chars[] = {' ',':','*','@'};
826 #endif
827
828 static guchar *
create_brush(IfsComposeVals * ifsvals,gint * brush_size,gdouble * brush_offset)829 create_brush (IfsComposeVals *ifsvals,
830 gint *brush_size,
831 gdouble *brush_offset)
832 {
833 gint i, j;
834 gint ii, jj;
835 guchar *brush;
836 #ifdef DEBUG_BRUSH
837 gdouble totpix = 0.0;
838 #endif
839
840 gdouble radius = ifsvals->radius * ifsvals->subdivide;
841
842 *brush_size = ceil (2 * radius);
843 *brush_offset = 0.5 * (*brush_size - 1);
844
845 brush = g_new (guchar, SQR (*brush_size));
846
847 for (i = 0; i < *brush_size; i++)
848 {
849 for (j = 0; j < *brush_size; j++)
850 {
851 gdouble pixel = 0.0;
852 gdouble d = sqrt (SQR (i - *brush_offset) +
853 SQR (j - *brush_offset));
854
855 if (d - 0.5 * G_SQRT2 > radius)
856 pixel = 0.0;
857 else if (d + 0.5 * G_SQRT2 < radius)
858 pixel = 1.0;
859 else
860 for (ii = 0; ii < 10; ii++)
861 for (jj = 0; jj < 10; jj++)
862 {
863 d = sqrt (SQR (i - *brush_offset + ii * 0.1 - 0.45) +
864 SQR (j - *brush_offset + jj * 0.1 - 0.45));
865 pixel += (d < radius) / 100.0;
866 }
867
868 brush[i * *brush_size + j] = 255.999 * pixel;
869
870 #ifdef DEBUG_BRUSH
871 putchar(brush_chars[(gint)(pixel * 3.999)]);
872 totpix += pixel;
873 #endif /* DEBUG_BRUSH */
874 }
875 #ifdef DEBUG_BRUSH
876 putchar('\n');
877 #endif /* DEBUG_BRUSH */
878 }
879 #ifdef DEBUG_BRUSH
880 printf ("Brush total / area = %f\n", totpix / SQR (ifsvals->subdivide));
881 #endif /* DEBUG_BRUSH */
882 return brush;
883 }
884
885 void
ifs_render(AffElement ** elements,gint num_elements,gint width,gint height,gint nsteps,IfsComposeVals * vals,gint band_y,gint band_height,guchar * data,guchar * mask,guchar * nhits,gboolean preview)886 ifs_render (AffElement **elements,
887 gint num_elements,
888 gint width,
889 gint height,
890 gint nsteps,
891 IfsComposeVals *vals,
892 gint band_y,
893 gint band_height,
894 guchar *data,
895 guchar *mask,
896 guchar *nhits,
897 gboolean preview)
898 {
899 gint i, k, n;
900 gdouble x, y;
901 gdouble r, g, b;
902 gint ri, gi, bi;
903 guint32 p0, psum;
904 gdouble pt;
905 guchar *ptr;
906 guint32 *prob;
907 gdouble *fprob;
908 gint subdivide;
909 guchar *brush = NULL;
910 gint brush_size = 1;
911 gdouble brush_offset = 0.0;
912
913 if (preview)
914 subdivide = 1;
915 else
916 subdivide = vals->subdivide;
917
918 /* compute the probabilities and transforms */
919 fprob = g_new (gdouble, num_elements);
920 prob = g_new (guint32, num_elements);
921 pt = 0.0;
922
923 for (i = 0; i < num_elements; i++)
924 {
925 aff_element_compute_trans(elements[i],
926 width * subdivide,
927 height * subdivide,
928 vals->center_x,
929 vals->center_y);
930 fprob[i] = fabs(
931 elements[i]->trans.a11 * elements[i]->trans.a22
932 - elements[i]->trans.a12 * elements[i]->trans.a21);
933
934 /* As a heuristic, if the determinant is really small, it's
935 probably a line element, so increase the probability so
936 it gets rendered */
937
938 /* FIXME: figure out what 0.01 really should be */
939 if (fprob[i] < 0.01)
940 fprob[i] = 0.01;
941
942 fprob[i] *= elements[i]->v.prob;
943
944 pt += fprob[i];
945 }
946
947 psum = 0;
948 for (i = 0; i < num_elements; i++)
949 {
950 psum += (guint32) -1 * (fprob[i] / pt);
951 prob[i] = psum;
952 }
953
954 prob[i - 1] = (guint32) -1; /* make sure we don't get bitten by roundoff */
955
956 /* create the brush */
957 if (!preview)
958 brush = create_brush (vals, &brush_size, &brush_offset);
959
960 x = y = 0;
961 r = g = b = 0;
962
963 /* n is used to limit the number of progress updates */
964 n = nsteps / 32;
965
966 /* now run the iteration */
967 for (i = 0; i < nsteps; i++)
968 {
969 if (!preview && ((i % n) == 0))
970 gimp_progress_update ((gdouble) i / (gdouble) nsteps);
971
972 p0 = g_random_int ();
973 k = 0;
974
975 while (p0 > prob[k])
976 k++;
977
978 aff2_apply (&elements[k]->trans, x, y, &x, &y);
979 aff3_apply (&elements[k]->color_trans, r, g, b, &r, &g, &b);
980
981 if (i < 50)
982 continue;
983
984 ri = (gint) (255.0 * r + 0.5);
985 gi = (gint) (255.0 * g + 0.5);
986 bi = (gint) (255.0 * b + 0.5);
987
988 if ((ri < 0) || (ri > 255) ||
989 (gi < 0) || (gi > 255) ||
990 (bi < 0) || (bi > 255))
991 continue;
992
993 if (preview)
994 {
995 if ((x < width) && (y < (band_y + band_height)) &&
996 (x >= 0) && (y >= band_y))
997 {
998 ptr = data + 3 * (((gint) (y - band_y)) * width + (gint) x);
999
1000 *ptr++ = ri;
1001 *ptr++ = gi;
1002 *ptr = bi;
1003 }
1004 }
1005 else
1006 {
1007 if ((x < width * subdivide) && (y < height * subdivide) &&
1008 (x >= 0) && (y >= 0))
1009 {
1010 gint ii;
1011 gint jj;
1012 gint jj0 = floor (y - brush_offset - band_y * subdivide);
1013 gint ii0 = floor (x - brush_offset);
1014 gint jjmin = 0;
1015 gint iimin = 0;
1016 gint jjmax;
1017 gint iimax;
1018
1019 if (ii0 < 0)
1020 iimin = - ii0;
1021 else
1022 iimin = 0;
1023
1024 if (jj0 < 0)
1025 jjmin = - jj0;
1026 else
1027 jjmin = 0;
1028
1029 if (jj0 + brush_size >= subdivide * band_height)
1030 jjmax = subdivide * band_height - jj0;
1031 else
1032 jjmax = brush_size;
1033
1034 if (ii0 + brush_size >= subdivide * width)
1035 iimax = subdivide * width - ii0;
1036 else
1037 iimax = brush_size;
1038
1039 for (jj = jjmin; jj < jjmax; jj++)
1040 for (ii = iimin; ii < iimax; ii++)
1041 {
1042 guint m_old;
1043 guint m_new;
1044 guint m_pix;
1045 guint n_hits;
1046 guint old_scale;
1047 guint pix_scale;
1048 gint index = (jj0 + jj) * width * subdivide + ii0 + ii;
1049
1050 n_hits = nhits[index];
1051 if (n_hits == 255)
1052 continue;
1053
1054 m_pix = brush[jj * brush_size + ii];
1055 if (!m_pix)
1056 continue;
1057
1058 nhits[index] = ++n_hits;
1059 m_old = mask[index];
1060 m_new = m_old + m_pix - m_old * m_pix / 255;
1061 mask[index] = m_new;
1062
1063 /* relative probability that old colored pixel is on top */
1064 old_scale = m_old * (255 * n_hits - m_pix);
1065
1066 /* relative probability that new colored pixel is on top */
1067 pix_scale = m_pix * ((255 - m_old) * n_hits + m_old);
1068
1069 ptr = data + 3 * index;
1070
1071 *ptr = ((old_scale * (*ptr) + pix_scale * ri) /
1072 (old_scale + pix_scale));
1073 ptr++;
1074
1075 *ptr = ((old_scale * (*ptr) + pix_scale * gi) /
1076 (old_scale + pix_scale));
1077 ptr++;
1078
1079 *ptr = ((old_scale * (*ptr) + pix_scale * bi) /
1080 (old_scale + pix_scale));
1081 }
1082 }
1083 }
1084 } /* main iteration */
1085
1086 if (!preview )
1087 gimp_progress_update (1.0);
1088
1089 g_free (brush);
1090 g_free (prob);
1091 g_free (fprob);
1092 }
1093