based on discussion with lathos on IRC
[imager.git] / filters.c
CommitLineData
02d1d628
AMH
1#include "image.h"
2#include <stdlib.h>
3#include <math.h>
4
5
6/*
7=head1 NAME
8
9filters.c - implements filters that operate on images
10
11=head1 SYNOPSIS
12
13
14 i_contrast(im, 0.8);
15 i_hardinvert(im);
16 // and more
17
18=head1 DESCRIPTION
19
20filters.c implements basic filters for Imager. These filters
21should be accessible from the filter interface as defined in
22the pod for Imager.
23
24=head1 FUNCTION REFERENCE
25
26Some of these functions are internal.
27
28=over 4
29
30=cut
31*/
32
33
34
35
36
37
38
39/*
40=item i_contrast(im, intensity)
41
42Scales the pixel values by the amount specified.
43
44 im - image object
45 intensity - scalefactor
46
47=cut
48*/
49
50void
51i_contrast(i_img *im, float intensity) {
52 int x, y;
53 unsigned char ch;
54 unsigned int new_color;
55 i_color rcolor;
56
57 mm_log((1,"i_contrast(im %p, intensity %f)\n", im, intensity));
58
59 if(intensity < 0) return;
60
61 for(y = 0; y < im->ysize; y++) for(x = 0; x < im->xsize; x++) {
62 i_gpix(im, x, y, &rcolor);
63
64 for(ch = 0; ch < im->channels; ch++) {
65 new_color = (unsigned int) rcolor.channel[ch];
66 new_color *= intensity;
67
68 if(new_color > 255) {
69 new_color = 255;
70 }
71 rcolor.channel[ch] = (unsigned char) new_color;
72 }
73 i_ppix(im, x, y, &rcolor);
74 }
75}
76
77
78/*
79=item i_hardinvert(im)
80
81Inverts the pixel values of the input image.
82
83 im - image object
84
85=cut
86*/
87
88void
89i_hardinvert(i_img *im) {
90 int x, y;
91 unsigned char ch;
92
93 i_color rcolor;
94
95 mm_log((1,"i_hardinvert(im %p)\n", im));
96
97 for(y = 0; y < im->ysize; y++) {
98 for(x = 0; x < im->xsize; x++) {
99 i_gpix(im, x, y, &rcolor);
100
101 for(ch = 0; ch < im->channels; ch++) {
102 rcolor.channel[ch] = 255 - rcolor.channel[ch];
103 }
104
105 i_ppix(im, x, y, &rcolor);
106 }
107 }
108}
109
110
111
112/*
113=item i_noise(im, amount, type)
114
115Inverts the pixel values by the amount specified.
116
117 im - image object
118 amount - deviation in pixel values
119 type - noise individual for each channel if true
120
121=cut
122*/
123
124#ifdef _MSC_VER
125/* random() is non-ASCII, even if it is better than rand() */
126#define random() rand()
127#endif
128
129void
130i_noise(i_img *im, float amount, unsigned char type) {
131 int x, y;
132 unsigned char ch;
133 int new_color;
134 float damount = amount * 2;
135 i_color rcolor;
136 int color_inc = 0;
137
138 mm_log((1,"i_noise(im %p, intensity %.2f\n", im, amount));
139
140 if(amount < 0) return;
141
142 for(y = 0; y < im->ysize; y++) for(x = 0; x < im->xsize; x++) {
143 i_gpix(im, x, y, &rcolor);
144
145 if(type == 0) {
146 color_inc = (amount - (damount * ((float)random() / RAND_MAX)));
147 }
148
149 for(ch = 0; ch < im->channels; ch++) {
150 new_color = (int) rcolor.channel[ch];
151
152 if(type != 0) {
153 new_color += (amount - (damount * ((float)random() / RAND_MAX)));
154 } else {
155 new_color += color_inc;
156 }
157
158 if(new_color < 0) {
159 new_color = 0;
160 }
161 if(new_color > 255) {
162 new_color = 255;
163 }
164
165 rcolor.channel[ch] = (unsigned char) new_color;
166 }
167
168 i_ppix(im, x, y, &rcolor);
169 }
170}
171
172
173/*
174=item i_noise(im, amount, type)
175
176Inverts the pixel values by the amount specified.
177
178 im - image object
179 amount - deviation in pixel values
180 type - noise individual for each channel if true
181
182=cut
183*/
184
185
186/*
187=item i_applyimage(im, add_im, mode)
188
189Apply's an image to another image
190
191 im - target image
192 add_im - image that is applied to target
193 mode - what method is used in applying:
194
195 0 Normal
196 1 Multiply
197 2 Screen
198 3 Overlay
199 4 Soft Light
200 5 Hard Light
201 6 Color dodge
202 7 Color Burn
203 8 Darker
204 9 Lighter
205 10 Add
206 11 Subtract
207 12 Difference
208 13 Exclusion
209
210=cut
211*/
212
213void i_applyimage(i_img *im, i_img *add_im, unsigned char mode) {
214 int x, y;
215 int mx, my;
216
217 mm_log((1, "i_applyimage(im %p, add_im %p, mode %d", im, add_im, mode));
218
219 mx = (add_im->xsize <= im->xsize) ? add_im->xsize : add_im->xsize;
220 my = (add_im->ysize <= im->ysize) ? add_im->ysize : add_im->ysize;
221
222 for(x = 0; x < mx; x++) {
223 for(y = 0; y < my; y++) {
224 }
225 }
226}
227
228
229/*
230=item i_bumpmap(im, bump, channel, light_x, light_y, st)
231
232Makes a bumpmap on image im using the bump image as the elevation map.
233
234 im - target image
235 bump - image that contains the elevation info
236 channel - to take the elevation information from
237 light_x - x coordinate of light source
238 light_y - y coordinate of light source
239 st - length of shadow
240
241=cut
242*/
243
244void
245i_bumpmap(i_img *im, i_img *bump, int channel, int light_x, int light_y, int st) {
246 int x, y, ch;
247 int mx, my;
248 i_color x1_color, y1_color, x2_color, y2_color, dst_color;
249 double nX, nY;
250 double tX, tY, tZ;
251 double aX, aY, aL;
252 double fZ;
253 unsigned char px1, px2, py1, py2;
254
255 i_img new_im;
256
257 mm_log((1, "i_bumpmap(im %p, add_im %p, channel %d, light_x %d, light_y %d, st %d)\n",
258 im, bump, channel, light_x, light_y, st));
259
260
261 if(channel >= bump->channels) {
262 mm_log((1, "i_bumpmap: channel = %d while bump image only has %d channels\n", channel, bump->channels));
263 return;
264 }
265
266 mx = (bump->xsize <= im->xsize) ? bump->xsize : im->xsize;
267 my = (bump->ysize <= im->ysize) ? bump->ysize : im->ysize;
268
269 i_img_empty_ch(&new_im, im->xsize, im->ysize, im->channels);
270
271 aX = (light_x > (mx >> 1)) ? light_x : mx - light_x;
272 aY = (light_y > (my >> 1)) ? light_y : my - light_y;
273
274 aL = sqrt((aX * aX) + (aY * aY));
275
276 for(y = 1; y < my - 1; y++) {
277 for(x = 1; x < mx - 1; x++) {
278 i_gpix(bump, x + st, y, &x1_color);
279 i_gpix(bump, x, y + st, &y1_color);
280 i_gpix(bump, x - st, y, &x2_color);
281 i_gpix(bump, x, y - st, &y2_color);
282
283 i_gpix(im, x, y, &dst_color);
284
285 px1 = x1_color.channel[channel];
286 py1 = y1_color.channel[channel];
287 px2 = x2_color.channel[channel];
288 py2 = y2_color.channel[channel];
289
290 nX = px1 - px2;
291 nY = py1 - py2;
292
293 nX += 128;
294 nY += 128;
295
296 fZ = (sqrt((nX * nX) + (nY * nY)) / aL);
297
298 tX = abs(x - light_x) / aL;
299 tY = abs(y - light_y) / aL;
300
301 tZ = 1 - (sqrt((tX * tX) + (tY * tY)) * fZ);
302
303 if(tZ < 0) tZ = 0;
304 if(tZ > 2) tZ = 2;
305
306 for(ch = 0; ch < im->channels; ch++)
307 dst_color.channel[ch] = (unsigned char) (float)(dst_color.channel[ch] * tZ);
308
309 i_ppix(&new_im, x, y, &dst_color);
310 }
311 }
312
313 i_copyto(im, &new_im, 0, 0, (int)im->xsize, (int)im->ysize, 0, 0);
314
315 i_img_exorcise(&new_im);
316}
317
318
319
320/*
321=item i_postlevels(im, levels)
322
323Quantizes Images to fewer levels.
324
325 im - target image
326 levels - number of levels
327
328=cut
329*/
330
331void
332i_postlevels(i_img *im, int levels) {
333 int x, y, ch;
334 float pv;
335 int rv;
336 float av;
337
338 i_color rcolor;
339
340 rv = (int) ((float)(256 / levels));
341 av = (float)levels;
342
343 for(y = 0; y < im->ysize; y++) for(x = 0; x < im->xsize; x++) {
344 i_gpix(im, x, y, &rcolor);
345
346 for(ch = 0; ch < im->channels; ch++) {
347 pv = (((float)rcolor.channel[ch] / 255)) * av;
348 pv = (int) ((int)pv * rv);
349
350 if(pv < 0) pv = 0;
351 else if(pv > 255) pv = 255;
352
353 rcolor.channel[ch] = (unsigned char) pv;
354 }
355 i_ppix(im, x, y, &rcolor);
356 }
357}
358
359
360/*
361=item i_mosaic(im, size)
362
363Makes an image looks like a mosaic with tilesize of size
364
365 im - target image
366 size - size of tiles
367
368=cut
369*/
370
371void
372i_mosaic(i_img *im, int size) {
373 int x, y, ch;
374 int lx, ly, z;
375 long sqrsize;
376
377 i_color rcolor;
378 long col[256];
379
380 sqrsize = size * size;
381
382 for(y = 0; y < im->ysize; y += size) for(x = 0; x < im->xsize; x += size) {
383 for(z = 0; z < 256; z++) col[z] = 0;
384
385 for(lx = 0; lx < size; lx++) {
386 for(ly = 0; ly < size; ly++) {
387 i_gpix(im, (x + lx), (y + ly), &rcolor);
388
389 for(ch = 0; ch < im->channels; ch++) {
390 col[ch] += rcolor.channel[ch];
391 }
392 }
393 }
394
395 for(ch = 0; ch < im->channels; ch++)
396 rcolor.channel[ch] = (int) ((float)col[ch] / sqrsize);
397
398
399 for(lx = 0; lx < size; lx++)
400 for(ly = 0; ly < size; ly++)
401 i_ppix(im, (x + lx), (y + ly), &rcolor);
402
403 }
404}
405
406/*
407=item saturate(in)
408
409Clamps the input value between 0 and 255. (internal)
410
411 in - input integer
412
413=cut
414*/
415
416static
417unsigned char
418saturate(int in) {
419 if (in>255) { return 255; }
420 else if (in>0) return in;
421 return 0;
422}
423
424
425/*
426=item i_watermark(im, wmark, tx, ty, pixdiff)
427
428Applies a watermark to the target image
429
430 im - target image
431 wmark - watermark image
432 tx - x coordinate of where watermark should be applied
433 ty - y coordinate of where watermark should be applied
434 pixdiff - the magnitude of the watermark, controls how visible it is
435
436=cut
437*/
438
439void
440i_watermark(i_img *im, i_img *wmark, int tx, int ty, int pixdiff) {
441 int vx, vy, ch;
442 i_color val, wval;
443
444 for(vx=0;vx<128;vx++) for(vy=0;vy<110;vy++) {
445
446 i_gpix(im, tx+vx, ty+vy,&val );
447 i_gpix(wmark, vx, vy, &wval);
448
449 for(ch=0;ch<im->channels;ch++)
450 val.channel[ch] = saturate( val.channel[ch] + (pixdiff* (wval.channel[0]-128) )/128 );
451
452 i_ppix(im,tx+vx,ty+vy,&val);
453 }
454}
455
456
457/*
458=item i_autolevels(im, lsat, usat, skew)
459
460Scales and translates each color such that it fills the range completely.
461Skew is not implemented yet - purpose is to control the color skew that can
462occur when changing the contrast.
463
464 im - target image
465 lsat - fraction of pixels that will be truncated at the lower end of the spectrum
466 usat - fraction of pixels that will be truncated at the higher end of the spectrum
467 skew - not used yet
468
469=cut
470*/
471
472void
473i_autolevels(i_img *im, float lsat, float usat, float skew) {
474 i_color val;
475 int i, x, y, rhist[256], ghist[256], bhist[256];
476 int rsum, rmin, rmax;
477 int gsum, gmin, gmax;
478 int bsum, bmin, bmax;
479 int rcl, rcu, gcl, gcu, bcl, bcu;
480
481 mm_log((1,"i_autolevels(im %p, lsat %f,usat %f,skew %f)\n", im, lsat,usat,skew));
482
483 rsum=gsum=bsum=0;
484 for(i=0;i<256;i++) rhist[i]=ghist[i]=bhist[i] = 0;
485 /* create histogram for each channel */
486 for(y = 0; y < im->ysize; y++) for(x = 0; x < im->xsize; x++) {
487 i_gpix(im, x, y, &val);
488 rhist[val.channel[0]]++;
489 ghist[val.channel[1]]++;
490 bhist[val.channel[2]]++;
491 }
492
493 for(i=0;i<256;i++) {
494 rsum+=rhist[i];
495 gsum+=ghist[i];
496 bsum+=bhist[i];
497 }
498
499 rmin = gmin = bmin = 0;
500 rmax = gmax = bmax = 255;
501
502 rcu = rcl = gcu = gcl = bcu = bcl = 0;
503
504 for(i=0; i<256; i++) {
505 rcl += rhist[i]; if ( (rcl<rsum*lsat) ) rmin=i;
506 rcu += rhist[255-i]; if ( (rcu<rsum*usat) ) rmax=255-i;
507
508 gcl += ghist[i]; if ( (gcl<gsum*lsat) ) gmin=i;
509 gcu += ghist[255-i]; if ( (gcu<gsum*usat) ) gmax=255-i;
510
511 bcl += bhist[i]; if ( (bcl<bsum*lsat) ) bmin=i;
512 bcu += bhist[255-i]; if ( (bcu<bsum*usat) ) bmax=255-i;
513 }
514
515 for(y = 0; y < im->ysize; y++) for(x = 0; x < im->xsize; x++) {
516 i_gpix(im, x, y, &val);
517 val.channel[0]=saturate((val.channel[0]-rmin)*255/(rmax-rmin));
518 val.channel[1]=saturate((val.channel[1]-gmin)*255/(gmax-gmin));
519 val.channel[2]=saturate((val.channel[2]-bmin)*255/(bmax-bmin));
520 i_ppix(im, x, y, &val);
521 }
522}
523
524/*
525=item Noise(x,y)
526
527Pseudo noise utility function used to generate perlin noise. (internal)
528
529 x - x coordinate
530 y - y coordinate
531
532=cut
533*/
534
535static
536float
537Noise(int x, int y) {
538 int n = x + y * 57;
539 n = (n<<13) ^ n;
540 return ( 1.0 - ( (n * (n * n * 15731 + 789221) + 1376312589) & 0x7fffffff) / 1073741824.0);
541}
542
543/*
544=item SmoothedNoise1(x,y)
545
546Pseudo noise utility function used to generate perlin noise. (internal)
547
548 x - x coordinate
549 y - y coordinate
550
551=cut
552*/
553
554static
555float
556SmoothedNoise1(float x, float y) {
557 float corners = ( Noise(x-1, y-1)+Noise(x+1, y-1)+Noise(x-1, y+1)+Noise(x+1, y+1) ) / 16;
558 float sides = ( Noise(x-1, y) +Noise(x+1, y) +Noise(x, y-1) +Noise(x, y+1) ) / 8;
559 float center = Noise(x, y) / 4;
560 return corners + sides + center;
561}
562
563
564/*
565=item G_Interpolate(a, b, x)
566
567Utility function used to generate perlin noise. (internal)
568
569=cut
570*/
571
572static
573float C_Interpolate(float a, float b, float x) {
574 /* float ft = x * 3.1415927; */
575 float ft = x * PI;
576 float f = (1 - cos(ft)) * .5;
577 return a*(1-f) + b*f;
578}
579
580
581/*
582=item InterpolatedNoise(x, y)
583
584Utility function used to generate perlin noise. (internal)
585
586=cut
587*/
588
589static
590float
591InterpolatedNoise(float x, float y) {
592
593 int integer_X = x;
594 float fractional_X = x - integer_X;
595 int integer_Y = y;
596 float fractional_Y = y - integer_Y;
597
598 float v1 = SmoothedNoise1(integer_X, integer_Y);
599 float v2 = SmoothedNoise1(integer_X + 1, integer_Y);
600 float v3 = SmoothedNoise1(integer_X, integer_Y + 1);
601 float v4 = SmoothedNoise1(integer_X + 1, integer_Y + 1);
602
603 float i1 = C_Interpolate(v1 , v2 , fractional_X);
604 float i2 = C_Interpolate(v3 , v4 , fractional_X);
605
606 return C_Interpolate(i1 , i2 , fractional_Y);
607}
608
609
610
611/*
612=item PerlinNoise_2D(x, y)
613
614Utility function used to generate perlin noise. (internal)
615
616=cut
617*/
618
619static
620float
621PerlinNoise_2D(float x, float y) {
622 int i,frequency;
623 float amplitude;
624 float total = 0;
625 int Number_Of_Octaves=6;
626 int n = Number_Of_Octaves - 1;
627
628 for(i=0;i<n;i++) {
629 frequency = 2*i;
630 amplitude = PI;
631 total = total + InterpolatedNoise(x * frequency, y * frequency) * amplitude;
632 }
633
634 return total;
635}
636
637
638/*
639=item i_radnoise(im, xo, yo, rscale, ascale)
640
641Perlin-like radial noise.
642
643 im - target image
644 xo - x coordinate of center
645 yo - y coordinate of center
646 rscale - radial scale
647 ascale - angular scale
648
649=cut
650*/
651
652void
653i_radnoise(i_img *im, int xo, int yo, float rscale, float ascale) {
654 int x, y, ch;
655 i_color val;
656 unsigned char v;
657 float xc, yc, r;
658 double a;
659
660 for(y = 0; y < im->ysize; y++) for(x = 0; x < im->xsize; x++) {
661 xc = (float)x-xo+0.5;
662 yc = (float)y-yo+0.5;
663 r = rscale*sqrt(xc*xc+yc*yc)+1.2;
664 a = (PI+atan2(yc,xc))*ascale;
665 v = saturate(128+100*(PerlinNoise_2D(a,r)));
666 /* v=saturate(120+12*PerlinNoise_2D(xo+(float)x/scale,yo+(float)y/scale)); Good soft marble */
667 for(ch=0; ch<im->channels; ch++) val.channel[ch]=v;
668 i_ppix(im, x, y, &val);
669 }
670}
671
672
673/*
674=item i_turbnoise(im, xo, yo, scale)
675
676Perlin-like 2d noise noise.
677
678 im - target image
679 xo - x coordinate translation
680 yo - y coordinate translation
681 scale - scale of noise
682
683=cut
684*/
685
686void
687i_turbnoise(i_img *im, float xo, float yo, float scale) {
688 int x,y,ch;
689 unsigned char v;
690 i_color val;
691
692 for(y = 0; y < im->ysize; y++) for(x = 0; x < im->xsize; x++) {
693 /* v=saturate(125*(1.0+PerlinNoise_2D(xo+(float)x/scale,yo+(float)y/scale))); */
694 v = saturate(120*(1.0+sin(xo+(float)x/scale+PerlinNoise_2D(xo+(float)x/scale,yo+(float)y/scale))));
695 for(ch=0; ch<im->channels; ch++) val.channel[ch] = v;
696 i_ppix(im, x, y, &val);
697 }
698}
699
700
701
702/*
703=item i_gradgen(im, num, xo, yo, ival, dmeasure)
704
705Gradient generating function.
706
707 im - target image
708 num - number of points given
709 xo - array of x coordinates
710 yo - array of y coordinates
711 ival - array of i_color objects
712 dmeasure - distance measure to be used.
713 0 = Euclidean
714 1 = Euclidean squared
715 2 = Manhattan distance
716
717=cut
718*/
719
720
721void
722i_gradgen(i_img *im, int num, int *xo, int *yo, i_color *ival, int dmeasure) {
723
724 i_color val;
725 int p, x, y, ch;
726 int channels = im->channels;
727 int xsize = im->xsize;
728 int ysize = im->ysize;
729
730 float *fdist;
731
732 mm_log((1,"i_gradgen(im %p, num %d, xo %p, yo %p, ival %p, dmeasure %d)\n", im, num, xo, yo, ival, dmeasure));
733
734 for(p = 0; p<num; p++) {
735 mm_log((1,"i_gradgen: (%d, %d)\n", xo[p], yo[p]));
736 ICL_info(&ival[p]);
737 }
738
739 fdist = mymalloc( sizeof(float) * num );
740
741 for(y = 0; y<ysize; y++) for(x = 0; x<xsize; x++) {
742 float cs = 0;
743 float csd = 0;
744 for(p = 0; p<num; p++) {
745 int xd = x-xo[p];
746 int yd = y-yo[p];
747 switch (dmeasure) {
748 case 0: /* euclidean */
749 fdist[p] = sqrt(xd*xd + yd*yd); /* euclidean distance */
750 break;
751 case 1: /* euclidean squared */
752 fdist[p] = xd*xd + yd*yd; /* euclidean distance */
753 break;
754 case 2: /* euclidean squared */
755 fdist[p] = max(xd*xd, yd*yd); /* manhattan distance */
756 break;
757 default:
758 m_fatal(3,"i_gradgen: Unknown distance measure\n");
759 }
760 cs += fdist[p];
761 }
762
763 csd = 1/((num-1)*cs);
764
765 for(p = 0; p<num; p++) fdist[p] = (cs-fdist[p])*csd;
766
767 for(ch = 0; ch<channels; ch++) {
768 int tres = 0;
769 for(p = 0; p<num; p++) tres += ival[p].channel[ch] * fdist[p];
770 val.channel[ch] = saturate(tres);
771 }
772 i_ppix(im, x, y, &val);
773 }
774
775}
776
02d1d628
AMH
777void
778i_nearest_color_foo(i_img *im, int num, int *xo, int *yo, i_color *ival, int dmeasure) {
779
a743c0a6 780 int p, x, y;
02d1d628
AMH
781 int xsize = im->xsize;
782 int ysize = im->ysize;
783
784 mm_log((1,"i_gradgen(im %p, num %d, xo %p, yo %p, ival %p, dmeasure %d)\n", im, num, xo, yo, ival, dmeasure));
785
786 for(p = 0; p<num; p++) {
787 mm_log((1,"i_gradgen: (%d, %d)\n", xo[p], yo[p]));
788 ICL_info(&ival[p]);
789 }
790
791 for(y = 0; y<ysize; y++) for(x = 0; x<xsize; x++) {
792 int midx = 0;
793 float mindist = 0;
794 float curdist = 0;
795
796 int xd = x-xo[0];
797 int yd = y-yo[0];
798
799 switch (dmeasure) {
800 case 0: /* euclidean */
801 mindist = sqrt(xd*xd + yd*yd); /* euclidean distance */
802 break;
803 case 1: /* euclidean squared */
804 mindist = xd*xd + yd*yd; /* euclidean distance */
805 break;
806 case 2: /* euclidean squared */
807 mindist = max(xd*xd, yd*yd); /* manhattan distance */
808 break;
809 default:
810 m_fatal(3,"i_nearest_color: Unknown distance measure\n");
811 }
812
813 for(p = 1; p<num; p++) {
814 xd = x-xo[p];
815 yd = y-yo[p];
816 switch (dmeasure) {
817 case 0: /* euclidean */
818 curdist = sqrt(xd*xd + yd*yd); /* euclidean distance */
819 break;
820 case 1: /* euclidean squared */
821 curdist = xd*xd + yd*yd; /* euclidean distance */
822 break;
823 case 2: /* euclidean squared */
824 curdist = max(xd*xd, yd*yd); /* manhattan distance */
825 break;
826 default:
827 m_fatal(3,"i_nearest_color: Unknown distance measure\n");
828 }
829 if (curdist < mindist) {
830 mindist = curdist;
831 midx = p;
832 }
833 }
834 i_ppix(im, x, y, &ival[midx]);
835 }
836}
837
02d1d628
AMH
838void
839i_nearest_color(i_img *im, int num, int *xo, int *yo, i_color *oval, int dmeasure) {
840 i_color *ival;
841 float *tval;
842 float c1, c2;
843 i_color val;
844 int p, x, y, ch;
02d1d628
AMH
845 int xsize = im->xsize;
846 int ysize = im->ysize;
847 int *cmatch;
848
d7f707d8 849 mm_log((1,"i_nearest_color(im %p, num %d, xo %p, yo %p, ival %p, dmeasure %d)\n", im, num, xo, yo, oval, dmeasure));
02d1d628
AMH
850
851 tval = mymalloc( sizeof(float)*num*im->channels );
852 ival = mymalloc( sizeof(i_color)*num );
853 cmatch = mymalloc( sizeof(int)*num );
854
855 for(p = 0; p<num; p++) {
856 for(ch = 0; ch<im->channels; ch++) tval[ p * im->channels + ch] = 0;
857 cmatch[p] = 0;
858 }
859
860
861 for(y = 0; y<ysize; y++) for(x = 0; x<xsize; x++) {
862 int midx = 0;
863 float mindist = 0;
864 float curdist = 0;
865
866 int xd = x-xo[0];
867 int yd = y-yo[0];
868
869 switch (dmeasure) {
870 case 0: /* euclidean */
871 mindist = sqrt(xd*xd + yd*yd); /* euclidean distance */
872 break;
873 case 1: /* euclidean squared */
874 mindist = xd*xd + yd*yd; /* euclidean distance */
875 break;
876 case 2: /* euclidean squared */
877 mindist = max(xd*xd, yd*yd); /* manhattan distance */
878 break;
879 default:
880 m_fatal(3,"i_nearest_color: Unknown distance measure\n");
881 }
882
883 for(p = 1; p<num; p++) {
884 xd = x-xo[p];
885 yd = y-yo[p];
886 switch (dmeasure) {
887 case 0: /* euclidean */
888 curdist = sqrt(xd*xd + yd*yd); /* euclidean distance */
889 break;
890 case 1: /* euclidean squared */
891 curdist = xd*xd + yd*yd; /* euclidean distance */
892 break;
893 case 2: /* euclidean squared */
894 curdist = max(xd*xd, yd*yd); /* manhattan distance */
895 break;
896 default:
897 m_fatal(3,"i_nearest_color: Unknown distance measure\n");
898 }
899 if (curdist < mindist) {
900 mindist = curdist;
901 midx = p;
902 }
903 }
904
905 cmatch[midx]++;
906 i_gpix(im, x, y, &val);
907 c2 = 1.0/(float)(cmatch[midx]);
908 c1 = 1.0-c2;
909
02d1d628
AMH
910 for(ch = 0; ch<im->channels; ch++)
911 tval[midx*im->channels + ch] = c1*tval[midx*im->channels + ch] + c2 * (float) val.channel[ch];
3bb1c1f3 912
02d1d628
AMH
913
914 }
915
916 for(p = 0; p<num; p++) for(ch = 0; ch<im->channels; ch++) ival[p].channel[ch] = tval[p*im->channels + ch];
917
918 i_nearest_color_foo(im, num, xo, yo, ival, dmeasure);
919}
6607600c
TC
920
921/*
922 Keep state information used by each type of fountain fill
923*/
924struct fount_state {
925 /* precalculated for the equation of the line perpendicular to the line AB */
926 double lA, lB, lC;
927 double AB;
928 double sqrtA2B2;
929 double mult;
930 double cos;
931 double sin;
932 double theta;
933 int xa, ya;
934 void *ssample_data;
935};
936
937static double linear_fount_f(double x, double y, struct fount_state *state);
938static double bilinear_fount_f(double x, double y, struct fount_state *state);
939static double radial_fount_f(double x, double y, struct fount_state *state);
940static double square_fount_f(double x, double y, struct fount_state *state);
941static double revolution_fount_f(double x, double y,
942 struct fount_state *state);
943static double conical_fount_f(double x, double y, struct fount_state *state);
944
945typedef double (*fount_func)(double, double, struct fount_state *);
946static fount_func fount_funcs[] =
947{
948 linear_fount_f,
949 bilinear_fount_f,
950 radial_fount_f,
951 square_fount_f,
952 revolution_fount_f,
953 conical_fount_f,
954};
955
956static double linear_interp(double pos, i_fountain_seg *seg);
957static double sine_interp(double pos, i_fountain_seg *seg);
958static double sphereup_interp(double pos, i_fountain_seg *seg);
959static double spheredown_interp(double pos, i_fountain_seg *seg);
960typedef double (*fount_interp)(double pos, i_fountain_seg *seg);
961static fount_interp fount_interps[] =
962{
963 linear_interp,
964 linear_interp,
965 sine_interp,
966 sphereup_interp,
967 spheredown_interp,
968};
969
970static void direct_cinterp(i_fcolor *out, double pos, i_fountain_seg *seg);
971static void hue_up_cinterp(i_fcolor *out, double pos, i_fountain_seg *seg);
972static void hue_down_cinterp(i_fcolor *out, double pos, i_fountain_seg *seg);
973typedef void (*fount_cinterp)(i_fcolor *out, double pos, i_fountain_seg *seg);
974static fount_cinterp fount_cinterps[] =
975{
976 direct_cinterp,
977 hue_up_cinterp,
978 hue_down_cinterp,
979};
980
981typedef double (*fount_repeat)(double v);
982static double fount_r_none(double v);
983static double fount_r_sawtooth(double v);
984static double fount_r_triangle(double v);
985static double fount_r_saw_both(double v);
986static double fount_r_tri_both(double v);
987static fount_repeat fount_repeats[] =
988{
989 fount_r_none,
990 fount_r_sawtooth,
991 fount_r_triangle,
992 fount_r_saw_both,
993 fount_r_tri_both,
994};
995
996static int simple_ssample(i_fcolor *out, double parm, double x, double y,
997 struct fount_state *state,
998 fount_func ffunc, fount_repeat rpfunc,
999 i_fountain_seg *segs, int count);
1000static int random_ssample(i_fcolor *out, double parm, double x, double y,
1001 struct fount_state *state,
1002 fount_func ffunc, fount_repeat rpfunc,
1003 i_fountain_seg *segs, int count);
1004static int circle_ssample(i_fcolor *out, double parm, double x, double y,
1005 struct fount_state *state,
1006 fount_func ffunc, fount_repeat rpfunc,
1007 i_fountain_seg *segs, int count);
1008typedef int (*fount_ssample)(i_fcolor *out, double parm, double x, double y,
1009 struct fount_state *state,
1010 fount_func ffunc, fount_repeat rpfunc,
1011 i_fountain_seg *segs, int count);
1012static fount_ssample fount_ssamples[] =
1013{
1014 NULL,
1015 simple_ssample,
1016 random_ssample,
1017 circle_ssample,
1018};
1019
1020static int
1021fount_getat(i_fcolor *out, double x, double y, fount_func ffunc,
1022 fount_repeat rpfunc, struct fount_state *state,
1023 i_fountain_seg *segs, int count);
1024
1025#define EPSILON (1e-6)
1026
1027/*
1028=item i_fountain(im, xa, ya, xb, yb, type, repeat, combine, super_sample, ssample_param, count, segs)
1029
1030Draws a fountain fill using A(xa, ya) and B(xb, yb) as reference points.
1031
1032I<type> controls how the reference points are used:
1033
1034=over
1035
1036=item i_ft_linear
1037
1038linear, where A is 0 and B is 1.
1039
1040=item i_ft_bilinear
1041
1042linear in both directions from A.
1043
1044=item i_ft_radial
1045
1046circular, where A is the centre of the fill, and B is a point
1047on the radius.
1048
1049=item i_ft_radial_square
1050
1051where A is the centre of the fill and B is the centre of
1052one side of the square.
1053
1054=item i_ft_revolution
1055
1056where A is the centre of the fill and B defines the 0/1.0
1057angle of the fill.
1058
1059=item i_ft_conical
1060
1061similar to i_ft_revolution, except that the revolution goes in both
1062directions
1063
1064=back
1065
1066I<repeat> can be one of:
1067
1068=over
1069
1070=item i_fr_none
1071
1072values < 0 are treated as zero, values > 1 are treated as 1.
1073
1074=item i_fr_sawtooth
1075
1076negative values are treated as 0, positive values are modulo 1.0
1077
1078=item i_fr_triangle
1079
1080negative values are treated as zero, if (int)value is odd then the value is treated as 1-(value
1081mod 1.0), otherwise the same as for sawtooth.
1082
1083=item i_fr_saw_both
1084
1085like i_fr_sawtooth, except that the sawtooth pattern repeats into
1086negative values.
1087
1088=item i_fr_tri_both
1089
1090Like i_fr_triangle, except that negative values are handled as their
1091absolute values.
1092
1093=back
1094
1095If combine is non-zero then non-opaque values are combined with the
1096underlying color.
1097
1098I<super_sample> controls super sampling, if any. At some point I'll
1099probably add a adaptive super-sampler. Current possible values are:
1100
1101=over
1102
1103=item i_fts_none
1104
1105No super-sampling is done.
1106
1107=item i_fts_grid
1108
1109A square grid of points withing the pixel are sampled.
1110
1111=item i_fts_random
1112
1113Random points within the pixel are sampled.
1114
1115=item i_fts_circle
1116
1117Points on the radius of a circle are sampled. This produces fairly
1118good results, but is fairly slow since sin() and cos() are evaluated
1119for each point.
1120
1121=back
1122
1123I<ssample_param> is intended to be roughly the number of points
1124sampled within the pixel.
1125
1126I<count> and I<segs> define the segments of the fill.
1127
1128=cut
1129
1130*/
1131
1132void
1133i_fountain(i_img *im, double xa, double ya, double xb, double yb,
1134 i_fountain_type type, i_fountain_repeat repeat,
1135 int combine, int super_sample, double ssample_param,
1136 int count, i_fountain_seg *segs) {
1137 struct fount_state state;
1138 fount_func ffunc;
1139 fount_ssample ssfunc;
1140 fount_repeat rpfunc;
1141 int x, y;
1142 i_fcolor *line = mymalloc(sizeof(i_fcolor) * im->xsize);
1143 int i, j;
1144 i_fountain_seg *my_segs = mymalloc(sizeof(i_fountain_seg) * count);
1145 int have_alpha = im->channels == 2 || im->channels == 4;
1146 int ch;
1147
1148 /* we keep a local copy that we can adjust for speed */
1149 for (i = 0; i < count; ++i) {
1150 i_fountain_seg *seg = my_segs + i;
1151
1152 *seg = segs[i];
1153 if (seg->type < 0 || type >= i_ft_end)
1154 seg->type = i_ft_linear;
1155 if (seg->color < 0 || seg->color >= i_fc_end)
1156 seg->color = i_fc_direct;
1157 if (seg->color == i_fc_hue_up || seg->color == i_fc_hue_down) {
1158 /* so we don't have to translate to HSV on each request, do it here */
1159 for (j = 0; j < 2; ++j) {
1160 i_rgb_to_hsvf(seg->c+j);
1161 }
1162 if (seg->color == i_fc_hue_up) {
1163 if (seg->c[1].channel[0] <= seg->c[0].channel[0])
1164 seg->c[1].channel[0] += 1.0;
1165 }
1166 else {
1167 if (seg->c[0].channel[0] <= seg->c[0].channel[1])
1168 seg->c[0].channel[0] += 1.0;
1169 }
1170 }
1171 /*printf("start %g mid %g end %g c0(%g,%g,%g,%g) c1(%g,%g,%g,%g) type %d color %d\n",
1172 seg->start, seg->middle, seg->end, seg->c[0].channel[0],
1173 seg->c[0].channel[1], seg->c[0].channel[2], seg->c[0].channel[3],
1174 seg->c[1].channel[0], seg->c[1].channel[1], seg->c[1].channel[2],
1175 seg->c[1].channel[3], seg->type, seg->color);*/
1176
1177 }
1178
1179 /* initialize each engine */
1180 /* these are so common ... */
1181 state.lA = xb - xa;
1182 state.lB = yb - ya;
1183 state.AB = sqrt(state.lA * state.lA + state.lB * state.lB);
1184 state.xa = xa;
1185 state.ya = ya;
1186 switch (type) {
1187 default:
1188 type = i_ft_linear; /* make the invalid value valid */
1189 case i_ft_linear:
1190 case i_ft_bilinear:
1191 state.lC = ya * ya - ya * yb + xa * xa - xa * xb;
1192 state.mult = 1;
1193 state.mult = 1/linear_fount_f(xb, yb, &state);
1194 break;
1195
1196 case i_ft_radial:
1197 state.mult = 1.0 / sqrt((double)(xb-xa)*(xb-xa)
1198 + (double)(yb-ya)*(yb-ya));
1199 break;
1200
1201 case i_ft_radial_square:
1202 state.cos = state.lA / state.AB;
1203 state.sin = state.lB / state.AB;
1204 state.mult = 1.0 / state.AB;
1205 break;
1206
1207 case i_ft_revolution:
1208 state.theta = atan2(yb-ya, xb-xa);
1209 state.mult = 1.0 / (PI * 2);
1210 break;
1211
1212 case i_ft_conical:
1213 state.theta = atan2(yb-ya, xb-xa);
1214 state.mult = 1.0 / PI;
1215 break;
1216 }
1217 ffunc = fount_funcs[type];
1218 if (super_sample < 0
1219 || super_sample >= (sizeof(fount_ssamples)/sizeof(*fount_ssamples))) {
1220 super_sample = 0;
1221 }
1222 state.ssample_data = NULL;
1223 switch (super_sample) {
1224 case i_fts_grid:
1225 ssample_param = floor(0.5 + sqrt(ssample_param));
1226 state.ssample_data = mymalloc(sizeof(i_fcolor) * ssample_param * ssample_param);
1227 break;
1228
1229 case i_fts_random:
1230 case i_fts_circle:
1231 ssample_param = floor(0.5+ssample_param);
1232 state.ssample_data = mymalloc(sizeof(i_fcolor) * ssample_param);
1233 break;
1234 }
1235 ssfunc = fount_ssamples[super_sample];
1236 if (repeat < 0 || repeat >= (sizeof(fount_repeats)/sizeof(*fount_repeats)))
1237 repeat = 0;
1238 rpfunc = fount_repeats[repeat];
1239
1240 for (y = 0; y < im->ysize; ++y) {
1241 i_glinf(im, 0, im->xsize, y, line);
1242 for (x = 0; x < im->xsize; ++x) {
1243 i_fcolor c;
1244 int got_one;
1245 double v;
1246 if (super_sample == i_fts_none)
1247 got_one = fount_getat(&c, x, y, ffunc, rpfunc, &state, my_segs, count);
1248 else
1249 got_one = ssfunc(&c, ssample_param, x, y, &state, ffunc, rpfunc,
1250 my_segs, count);
1251 if (got_one) {
1252 i_fountain_seg *seg = my_segs + i;
1253 if (combine) {
1254 for (ch = 0; ch < im->channels; ++ch) {
1255 line[x].channel[ch] = line[x].channel[ch] * (1.0 - c.channel[3])
1256 + c.channel[ch] * c.channel[3];
1257 }
1258 }
1259 else
1260 line[x] = c;
1261 }
1262 }
1263 i_plinf(im, 0, im->xsize, y, line);
1264 }
1265 myfree(line);
1266 myfree(my_segs);
1267 if (state.ssample_data)
1268 myfree(state.ssample_data);
1269}
1270
1271/*
1272=back
1273
1274=head1 INTERNAL FUNCTIONS
1275
1276=over
1277
1278=item fount_getat(out, x, y, ffunc, rpfunc, state, segs, count)
1279
1280Evaluates the fountain fill at the given point.
1281
1282This is called by both the non-super-sampling and super-sampling code.
1283
1284You might think that it would make sense to sample the fill parameter
1285instead, and combine those, but this breaks badly.
1286
1287=cut
1288*/
1289
1290static int
1291fount_getat(i_fcolor *out, double x, double y, fount_func ffunc,
1292 fount_repeat rpfunc, struct fount_state *state,
1293 i_fountain_seg *segs, int count) {
1294 double v = rpfunc(ffunc(x, y, state));
1295 int i;
1296
1297 i = 0;
1298 while (i < count && (v < segs[i].start || v > segs[i].end)) {
1299 ++i;
1300 }
1301 if (i < count) {
1302 v = (fount_interps[segs[i].type])(v, segs+i);
1303 (fount_cinterps[segs[i].color])(out, v, segs+i);
1304 return 1;
1305 }
1306 else
1307 return 0;
1308}
1309
1310/*
1311=item linear_fount_f(x, y, state)
1312
1313Calculate the fill parameter for a linear fountain fill.
1314
1315Uses the point to line distance function, with some precalculation
1316done in i_fountain().
1317
1318=cut
1319*/
1320static double
1321linear_fount_f(double x, double y, struct fount_state *state) {
1322 return (state->lA * x + state->lB * y + state->lC) / state->AB * state->mult;
1323}
1324
1325/*
1326=item bilinear_fount_f(x, y, state)
1327
1328Calculate the fill parameter for a bi-linear fountain fill.
1329
1330=cut
1331*/
1332static double
1333bilinear_fount_f(double x, double y, struct fount_state *state) {
1334 return fabs((state->lA * x + state->lB * y + state->lC) / state->AB * state->mult);
1335}
1336
1337/*
1338=item radial_fount_f(x, y, state)
1339
1340Calculate the fill parameter for a radial fountain fill.
1341
1342Simply uses the distance function.
1343
1344=cut
1345 */
1346static double
1347radial_fount_f(double x, double y, struct fount_state *state) {
1348 return sqrt((double)(state->xa-x)*(state->xa-x)
1349 + (double)(state->ya-y)*(state->ya-y)) * state->mult;
1350}
1351
1352/*
1353=item square_fount_f(x, y, state)
1354
1355Calculate the fill parameter for a square fountain fill.
1356
1357Works by rotating the reference co-ordinate around the centre of the
1358square.
1359
1360=cut
1361*/
1362static double
1363square_fount_f(double x, double y, struct fount_state *state) {
1364 int xc, yc; /* centred on A */
1365 double xt, yt; /* rotated by theta */
1366 xc = x - state->xa;
1367 yc = y - state->ya;
1368 xt = fabs(xc * state->cos + yc * state->sin);
1369 yt = fabs(-xc * state->sin + yc * state->cos);
1370 return (xt > yt ? xt : yt) * state->mult;
1371}
1372
1373/*
1374=item revolution_fount_f(x, y, state)
1375
1376Calculates the fill parameter for the revolution fountain fill.
1377
1378=cut
1379*/
1380static double
1381revolution_fount_f(double x, double y, struct fount_state *state) {
1382 double angle = atan2(y - state->ya, x - state->xa);
1383
1384 angle -= state->theta;
1385 if (angle < 0) {
1386 angle = fmod(angle+ PI * 4, PI*2);
1387 }
1388
1389 return angle * state->mult;
1390}
1391
1392/*
1393=item conical_fount_f(x, y, state)
1394
1395Calculates the fill parameter for the conical fountain fill.
1396
1397=cut
1398*/
1399static double
1400conical_fount_f(double x, double y, struct fount_state *state) {
1401 double angle = atan2(y - state->ya, x - state->xa);
1402
1403 angle -= state->theta;
1404 if (angle < -PI)
1405 angle += PI * 2;
1406 else if (angle > PI)
1407 angle -= PI * 2;
1408
1409 return fabs(angle) * state->mult;
1410}
1411
1412/*
1413=item linear_interp(pos, seg)
1414
1415Calculates linear interpolation on the fill parameter. Breaks the
1416segment into 2 regions based in the I<middle> value.
1417
1418=cut
1419*/
1420static double
1421linear_interp(double pos, i_fountain_seg *seg) {
1422 if (pos < seg->middle) {
1423 double len = seg->middle - seg->start;
1424 if (len < EPSILON)
1425 return 0.0;
1426 else
1427 return (pos - seg->start) / len / 2;
1428 }
1429 else {
1430 double len = seg->end - seg->middle;
1431 if (len < EPSILON)
1432 return 1.0;
1433 else
1434 return 0.5 + (pos - seg->middle) / len / 2;
1435 }
1436}
1437
1438/*
1439=item sine_interp(pos, seg)
1440
1441Calculates sine function interpolation on the fill parameter.
1442
1443=cut
1444*/
1445static double
1446sine_interp(double pos, i_fountain_seg *seg) {
1447 /* I wonder if there's a simple way to smooth the transition for this */
1448 double work = linear_interp(pos, seg);
1449
1450 return (1-cos(work * PI))/2;
1451}
1452
1453/*
1454=item sphereup_interp(pos, seg)
1455
1456Calculates spherical interpolation on the fill parameter, with the cusp
1457at the low-end.
1458
1459=cut
1460*/
1461static double
1462sphereup_interp(double pos, i_fountain_seg *seg) {
1463 double work = linear_interp(pos, seg);
1464
1465 return sqrt(1.0 - (1-work) * (1-work));
1466}
1467
1468/*
1469=item spheredown_interp(pos, seg)
1470
1471Calculates spherical interpolation on the fill parameter, with the cusp
1472at the high-end.
1473
1474=cut
1475*/
1476static double
1477spheredown_interp(double pos, i_fountain_seg *seg) {
1478 double work = linear_interp(pos, seg);
1479
1480 return 1-sqrt(1.0 - work * work);
1481}
1482
1483/*
1484=item direct_cinterp(out, pos, seg)
1485
1486Calculates the fountain color based on direct scaling of the channels
1487of the color channels.
1488
1489=cut
1490*/
1491static void
1492direct_cinterp(i_fcolor *out, double pos, i_fountain_seg *seg) {
1493 int ch;
1494 for (ch = 0; ch < MAXCHANNELS; ++ch) {
1495 out->channel[ch] = seg->c[0].channel[ch] * (1 - pos)
1496 + seg->c[1].channel[ch] * pos;
1497 }
1498}
1499
1500/*
1501=item hue_up_cinterp(put, pos, seg)
1502
1503Calculates the fountain color based on scaling a HSV value. The hue
1504increases as the fill parameter increases.
1505
1506=cut
1507*/
1508static void
1509hue_up_cinterp(i_fcolor *out, double pos, i_fountain_seg *seg) {
1510 int ch;
1511 for (ch = 0; ch < MAXCHANNELS; ++ch) {
1512 out->channel[ch] = seg->c[0].channel[ch] * (1 - pos)
1513 + seg->c[1].channel[ch] * pos;
1514 }
1515 i_hsv_to_rgbf(out);
1516}
1517
1518/*
1519=item hue_down_cinterp(put, pos, seg)
1520
1521Calculates the fountain color based on scaling a HSV value. The hue
1522decreases as the fill parameter increases.
1523
1524=cut
1525*/
1526static void
1527hue_down_cinterp(i_fcolor *out, double pos, i_fountain_seg *seg) {
1528 int ch;
1529 for (ch = 0; ch < MAXCHANNELS; ++ch) {
1530 out->channel[ch] = seg->c[0].channel[ch] * (1 - pos)
1531 + seg->c[1].channel[ch] * pos;
1532 }
1533 i_hsv_to_rgbf(out);
1534}
1535
1536/*
1537=item simple_ssample(out, parm, x, y, state, ffunc, rpfunc, segs, count)
1538
1539Simple grid-based super-sampling.
1540
1541=cut
1542*/
1543static int
1544simple_ssample(i_fcolor *out, double parm, double x, double y,
1545 struct fount_state *state,
1546 fount_func ffunc, fount_repeat rpfunc, i_fountain_seg *segs,
1547 int count) {
1548 i_fcolor *work = state->ssample_data;
1549 int dx, dy;
1550 int grid = parm;
1551 double base = -0.5 + 0.5 / grid;
1552 double step = 1.0 / grid;
1553 int ch, i;
1554 int samp_count = 0;
1555
1556 for (dx = 0; dx < grid; ++dx) {
1557 for (dy = 0; dy < grid; ++dy) {
1558 if (fount_getat(work+samp_count, x + base + step * dx,
1559 y + base + step * dy, ffunc, rpfunc, state,
1560 segs, count)) {
1561 ++samp_count;
1562 }
1563 }
1564 }
1565 for (ch = 0; ch < MAXCHANNELS; ++ch) {
1566 out->channel[ch] = 0;
1567 for (i = 0; i < samp_count; ++i) {
1568 out->channel[ch] += work[i].channel[ch];
1569 }
1570 /* we divide by 4 rather than samp_count since if there's only one valid
1571 sample it should be mostly transparent */
1572 out->channel[ch] /= grid * grid;
1573 }
1574 return samp_count;
1575}
1576
1577/*
1578=item random_ssample(out, parm, x, y, state, ffunc, rpfunc, segs, count)
1579
1580Random super-sampling.
1581
1582=cut
1583*/
1584static int
1585random_ssample(i_fcolor *out, double parm, double x, double y,
1586 struct fount_state *state,
1587 fount_func ffunc, fount_repeat rpfunc, i_fountain_seg *segs,
1588 int count) {
1589 i_fcolor *work = state->ssample_data;
1590 int i, ch;
1591 int maxsamples = parm;
1592 double rand_scale = 1.0 / RAND_MAX;
1593 int samp_count = 0;
1594 for (i = 0; i < maxsamples; ++i) {
1595 if (fount_getat(work+samp_count, x - 0.5 + rand() * rand_scale,
1596 y - 0.5 + rand() * rand_scale, ffunc, rpfunc, state,
1597 segs, count)) {
1598 ++samp_count;
1599 }
1600 }
1601 for (ch = 0; ch < MAXCHANNELS; ++ch) {
1602 out->channel[ch] = 0;
1603 for (i = 0; i < samp_count; ++i) {
1604 out->channel[ch] += work[i].channel[ch];
1605 }
1606 /* we divide by maxsamples rather than samp_count since if there's
1607 only one valid sample it should be mostly transparent */
1608 out->channel[ch] /= maxsamples;
1609 }
1610 return samp_count;
1611}
1612
1613/*
1614=item circle_ssample(out, parm, x, y, state, ffunc, rpfunc, segs, count)
1615
1616Super-sampling around the circumference of a circle.
1617
1618I considered saving the sin()/cos() values and transforming step-size
1619around the circle, but that's inaccurate, though it may not matter
1620much.
1621
1622=cut
1623 */
1624static int
1625circle_ssample(i_fcolor *out, double parm, double x, double y,
1626 struct fount_state *state,
1627 fount_func ffunc, fount_repeat rpfunc, i_fountain_seg *segs,
1628 int count) {
1629 i_fcolor *work = state->ssample_data;
1630 int i, ch;
1631 int maxsamples = parm;
1632 double angle = 2 * PI / maxsamples;
1633 double radius = 0.3; /* semi-random */
1634 int samp_count = 0;
1635 for (i = 0; i < maxsamples; ++i) {
1636 if (fount_getat(work+samp_count, x + radius * cos(angle * i),
1637 y + radius * sin(angle * i), ffunc, rpfunc, state,
1638 segs, count)) {
1639 ++samp_count;
1640 }
1641 }
1642 for (ch = 0; ch < MAXCHANNELS; ++ch) {
1643 out->channel[ch] = 0;
1644 for (i = 0; i < samp_count; ++i) {
1645 out->channel[ch] += work[i].channel[ch];
1646 }
1647 /* we divide by maxsamples rather than samp_count since if there's
1648 only one valid sample it should be mostly transparent */
1649 out->channel[ch] /= maxsamples;
1650 }
1651 return samp_count;
1652}
1653
1654/*
1655=item fount_r_none(v)
1656
1657Implements no repeats. Simply clamps the fill value.
1658
1659=cut
1660*/
1661static double
1662fount_r_none(double v) {
1663 return v < 0 ? 0 : v > 1 ? 1 : v;
1664}
1665
1666/*
1667=item fount_r_sawtooth(v)
1668
1669Implements sawtooth repeats. Clamps negative values and uses fmod()
1670on others.
1671
1672=cut
1673*/
1674static double
1675fount_r_sawtooth(double v) {
1676 return v < 0 ? 0 : fmod(v, 1.0);
1677}
1678
1679/*
1680=item fount_r_triangle(v)
1681
1682Implements triangle repeats. Clamps negative values, uses fmod to get
1683a range 0 through 2 and then adjusts values > 1.
1684
1685=cut
1686*/
1687static double
1688fount_r_triangle(double v) {
1689 if (v < 0)
1690 return 0;
1691 else {
1692 v = fmod(v, 2.0);
1693 return v > 1.0 ? 2.0 - v : v;
1694 }
1695}
1696
1697/*
1698=item fount_r_saw_both(v)
1699
1700Implements sawtooth repeats in the both postive and negative directions.
1701
1702Adjusts the value to be postive and then just uses fmod().
1703
1704=cut
1705*/
1706static double
1707fount_r_saw_both(double v) {
1708 if (v < 0)
1709 v += 1+(int)(-v);
1710 return fmod(v, 1.0);
1711}
1712
1713/*
1714=item fount_r_tri_both(v)
1715
1716Implements triangle repeats in the both postive and negative directions.
1717
1718Uses fmod on the absolute value, and then adjusts values > 1.
1719
1720=cut
1721*/
1722static double
1723fount_r_tri_both(double v) {
1724 v = fmod(fabs(v), 2.0);
1725 return v > 1.0 ? 2.0 - v : v;
1726}
1727
1728/*
1729=back
1730
1731=head1 AUTHOR
1732
1733Arnar M. Hrafnkelsson <addi@umich.edu>
1734
1735Tony Cook <tony@develop-help.com> (i_fountain())
1736
1737=head1 SEE ALSO
1738
1739Imager(3)
1740
1741=cut
1742*/