- added log() and exp() functions to transform2()
[imager.git] / regmach.c
1 #include "regmach.h"
2 #include <float.h>
3
4 /*#define DEBUG*/
5 #ifdef DEBUG
6 #define DBG(x) printf x
7 #else
8 #define DBG(x) 
9 #endif
10
11 static float MAX_EXP_ARG; /* = log(DBL_MAX); */
12
13
14 /* these functions currently assume RGB images - there seems to be some 
15    support for other color spaces, but I can't tell how you find what 
16    space an image is using.
17
18    HSV conversions from pages 401-403 "Procedural Elements for Computer 
19    Graphics", 1985, ISBN 0-07-053534-5.  The algorithm presents to produce
20    an HSV color calculates all components at once - I don't, so I've 
21    simiplified the algorithm to avoid unnecessary calculation (any errors 
22    (of which I had a few ;) are mine).
23 */
24
25 /* returns the value (brightness) of color from 0 to 1 */
26 static double hsv_value(i_color color) {
27   return i_max(i_max(color.rgb.r, color.rgb.g), color.rgb.b) / 255.0;
28 }
29
30 /* returns the hue (color) of color from 0 to 360 */
31 static double hsv_hue(i_color color) {
32   int val;
33   int temp;
34   temp = i_min(i_min(color.rgb.r, color.rgb.g), color.rgb.b);
35   val = i_max(color.rgb.r, i_max(color.rgb.g, color.rgb.b));
36   if (val == 0 || val==temp) {
37     return 0;
38   }
39   else {
40     double cr = (val - color.rgb.r) / (double)(val - temp);
41     double cg = (val - color.rgb.g) / (double)(val - temp);
42     double cb = (val - color.rgb.b) / (double)(val - temp);
43     double hue;
44     if (color.rgb.r == val) {
45       hue = cb-cg;
46     }
47     else if (color.rgb.g == val) {
48       hue = 2.0 + cr-cb;
49     }
50     else { /* if (blue == val) */
51       hue = 4.0 + cg - cr;
52     }
53     hue *= 60.0; /* to degrees */
54     if (hue < 0) 
55       hue += 360;
56
57     return hue;
58   }
59 }
60
61 /* return the saturation of color from 0 to 1 */
62 static double hsv_sat(i_color color) {
63   int value = i_max(i_max(color.rgb.r, color.rgb.g), color.rgb.b);
64   if (value == 0) {
65     return 0;
66   }
67   else {
68     int temp = i_min(i_max(color.rgb.r, color.rgb.g), color.rgb.b);
69     return (value - temp) / (double)value;
70   }
71 }
72
73 static i_color make_hsv(double hue, double sat, double val, int alpha) {
74   int i;
75   i_color c;
76   for( i=0; i< MAXCHANNELS; i++) c.channel[i]=0;
77   DBG(("hsv=%f %f %f\n", hue, sat, val));
78   if (sat <= 0) { /* handle -ve in case someone supplies a bad value */
79     /* should this be * 256? */
80     c.rgb.r = c.rgb.g = c.rgb.b = 255 * val;
81   }
82   else {
83     int i, m, n, k, v;
84     double f;
85     
86     if (val < 0) val = 0;
87     if (val > 1) val = 1;
88     if (sat > 1) sat = 1;
89
90     /* I want to handle -360 <= hue < 720 so that the caller can
91        fiddle with colour 
92     */
93     if (hue >= 360)
94       hue -= 360;
95     else if (hue < 0) 
96       hue += 360;
97     hue /= 60;
98     i = hue; /* floor */
99     f = hue - i;
100     val *= 255; 
101     m = val * (1.0 - sat);
102     n = val * (1.0 - sat * f);
103     k = val * (1.0 - sat * (1 - f));
104     v = val;
105     switch (i) {
106     case 0:
107       c.rgb.r = v; c.rgb.g = k; c.rgb.b = m;
108       break;
109     case 1:
110       c.rgb.r = n; c.rgb.g = v; c.rgb.b = m;
111       break;
112     case 2:
113       c.rgb.r = m; c.rgb.g = v; c.rgb.b = k;
114       break;
115     case 3:
116       c.rgb.r = m; c.rgb.g = n; c.rgb.b = v;
117       break;
118     case 4:
119       c.rgb.r = k; c.rgb.g = m; c.rgb.b = v;
120       break;
121     case 5:
122       c.rgb.r = v; c.rgb.g = m; c.rgb.b = n;
123       break;
124     }
125   }
126   c.rgba.a = alpha;
127
128   return c;
129 }
130
131 static i_color make_rgb(int r, int g, int b, int a) {
132   i_color c;
133   if (r < 0)
134     r = 0;
135   if (r > 255)
136     r = 255;
137   c.rgb.r = r;
138   if (g < 0)
139     g = 0;
140   if (g > 255)
141     g = 255;
142   c.rgb.g = g;
143   if (b < 0)
144     b = 0;
145   if (b > 255)
146     b = 255;
147   c.rgb.b = b;
148
149   c.rgba.a = a;
150
151   return c;
152 }
153
154 /* greatly simplifies the code */
155 #define nout n_regs[codes->rout]
156 #define na n_regs[codes->ra]
157 #define nb n_regs[codes->rb]
158 #define nc n_regs[codes->rc]
159 #define nd n_regs[codes->rd]
160 #define cout c_regs[codes->rout]
161 #define ca c_regs[codes->ra]
162 #define cb c_regs[codes->rb]
163 #define cc c_regs[codes->rc]
164 #define cd c_regs[codes->rd]
165
166 /* this is a pretty poor epsilon used for loosening up equality comparisons
167    It isn't currently used for inequalities 
168 */
169
170 #define n_epsilon(x, y) (abs(x)+abs(y))*0.001
171 static i_color bcol = {{ 0 }};
172
173 i_color i_rm_run(struct rm_op codes[], size_t code_count, 
174                double n_regs[],  size_t n_regs_count,
175                i_color c_regs[], size_t c_regs_count,
176                i_img *images[],  size_t image_count) {
177   double dx, dy;
178   struct rm_op *codes_base = codes;
179   size_t count_base = code_count;
180
181   DBG(("rm_run(%p, %d)\n", codes, code_count));
182   while (code_count) {
183     DBG((" rm_code %d\n", codes->code));
184     switch (codes->code) {
185     case rbc_add:
186       nout = na + nb;
187       break;
188       
189     case rbc_subtract:
190       nout = na - nb;
191       break;
192       
193     case rbc_mult:
194       nout = na * nb;
195       break;
196       
197     case rbc_div:
198       if (abs(nb) < 1e-10)
199         nout = 1e10;
200       else
201         nout = na / nb;
202       break;
203       
204     case rbc_mod:
205       if (abs(nb) > 1e-10) {
206         nout = fmod(na, nb);
207       }
208       else {
209         nout = 0; /* close enough ;) */
210       }
211       break;
212
213     case rbc_pow:
214       nout = pow(na, nb);
215       break;
216
217     case rbc_uminus:
218       nout = -na;
219
220     case rbc_multp:
221       cout = make_rgb(ca.rgb.r * nb, ca.rgb.g * nb, ca.rgb.b * nb, 255);
222       break;
223
224     case rbc_addp:
225       cout = make_rgb(ca.rgb.r + cb.rgb.r, ca.rgb.g + cb.rgb.g, 
226                       ca.rgb.b + cb.rgb.b, 255);
227       break;
228
229     case rbc_subtractp:
230       cout = make_rgb(ca.rgb.r - cb.rgb.r, ca.rgb.g - cb.rgb.g, 
231                       ca.rgb.b - cb.rgb.b, 255);
232       break;
233
234     case rbc_sin:
235       nout = sin(na);
236       break;
237
238     case rbc_cos:
239       nout = cos(na);
240       break;
241
242     case rbc_atan2:
243       nout = atan2(na, nb);
244       break;
245
246     case rbc_sqrt:
247       nout = sqrt(na);
248       break;
249
250     case rbc_distance:
251       dx = na-nc;
252       dy = nb-nd;
253       nout = sqrt(dx*dx+dy*dy);
254       break;
255
256     case rbc_getp1:
257       i_gpix(images[0], na, nb, c_regs+codes->rout);
258       if (images[0]->channels < 4) cout.rgba.a = 255;
259       break;
260
261     case rbc_getp2:
262       i_gpix(images[1], na, nb, c_regs+codes->rout);
263       if (images[1]->channels < 4) cout.rgba.a = 255;
264       break;
265
266     case rbc_getp3:
267       i_gpix(images[2], na, nb, c_regs+codes->rout);
268       if (images[2]->channels < 4) cout.rgba.a = 255;
269       break;
270
271     case rbc_value:
272       nout = hsv_value(ca);
273       break;
274
275     case rbc_hue:
276       nout = hsv_hue(ca);
277       break;
278
279     case rbc_sat:
280       nout = hsv_sat(ca);
281       break;
282       
283     case rbc_hsv:
284       cout = make_hsv(na, nb, nc, 255);
285       break;
286
287     case rbc_hsva:
288       cout = make_hsv(na, nb, nc, nd);
289       break;
290
291     case rbc_red:
292       nout = ca.rgb.r;
293       break;
294
295     case rbc_green:
296       nout = ca.rgb.g;
297       break;
298
299     case rbc_blue:
300       nout = ca.rgb.b;
301       break;
302
303     case rbc_alpha:
304       nout = ca.rgba.a;
305       break;
306
307     case rbc_rgb:
308       cout = make_rgb(na, nb, nc, 255);
309       break;
310
311     case rbc_rgba:
312       cout = make_rgb(na, nb, nc, nd);
313       break;
314
315     case rbc_int:
316       nout = (int)(na);
317       break;
318
319     case rbc_if:
320       nout = na ? nb : nc;
321       break;
322
323     case rbc_ifp:
324       cout = na ? cb : cc;
325       break;
326
327     case rbc_le:
328       nout = na <= nb + n_epsilon(na,nb);
329       break;
330
331     case rbc_lt:
332       nout = na < nb;
333       break;
334
335     case rbc_ge:
336       nout = na >= nb - n_epsilon(na,nb);
337       break;
338
339     case rbc_gt:
340       nout = na > nb;
341       break;
342
343     case rbc_eq:
344       nout = abs(na-nb) <= n_epsilon(na,nb);
345       break;
346
347     case rbc_ne:
348       nout = abs(na-nb) > n_epsilon(na,nb);
349       break;
350
351     case rbc_and:
352       nout = na && nb;
353       break;
354  
355     case rbc_or:
356       nout = na || nb;
357       break;
358
359     case rbc_not:
360       nout = !na;
361       break;
362
363     case rbc_abs:
364       nout = abs(na);
365       break;
366
367     case rbc_ret:
368       return ca;
369       break;
370
371     case rbc_jump:
372       /* yes, order is important here */
373       code_count = count_base - codes->ra;
374       codes = codes_base + codes->ra;
375       continue;
376     
377     case rbc_jumpz:
378       if (!na) {        
379         /* yes, order is important here */
380         code_count = count_base - codes->rb;
381         codes = codes_base + codes->rb;
382         continue;
383       }
384       break;
385
386     case rbc_jumpnz:
387       if (na) {
388         /* yes, order is important here */
389         code_count = count_base - codes->rb;
390         codes = codes_base + codes->rb;
391         continue;
392       }
393       break;
394
395     case rbc_set:
396       nout = na;
397       break;
398
399     case rbc_setp:
400       cout = ca;
401       break;
402
403     case rbc_log:
404       if (na > 0) {
405         nout = log(na);
406       }
407       else {
408         nout = DBL_MAX;
409       }
410       break;
411
412     case rbc_exp:
413       if (!MAX_EXP_ARG) MAX_EXP_ARG = log(DBL_MAX);
414       if (na <= MAX_EXP_ARG) {
415         nout = exp(na);
416       }
417       else {
418         nout = DBL_MAX;
419       }
420       break;
421
422     case rbc_print:
423       printf("r%d is %g\n", codes->ra, na);
424       break;
425
426     default:
427       /*croak("bad opcode"); */
428       printf("bad op %d\n", codes->code);
429       return bcol;
430     }
431     --code_count;
432     ++codes;
433   }
434   return bcol;
435   /* croak("no return opcode"); */
436 }