1 //-----------------------------------------------------------------------
2 // Color Quantization Demo
3 //
4 // Author: Roman Podobedov
5 // Email: romka@ut.ee
6 // Romka Graphics: www.ut.ee/~romka
7 //
8 // Also in this file implemented Wu's Color Quantizer algorithm (v. 2)
9 // For details see Graphics Gems vol. II, pp. 126-133
10 //
11 // Wu's Color Quantizer Algorithm:
12 // Author: Xiaolin Wu
13 // Dept. of Computer Science
14 // Univ. of Western Ontario
15 // London, Ontario N6A 5B7
16 // wu@csd.uwo.ca
17 // http://www.csd.uwo.ca/faculty/wu/
18 //-----------------------------------------------------------------------
19
20 //-----------------------------------------------------------------------------
21 //
22 // ImageLib Sources
23 // by Denton Woods
24 // Last modified: 01/04/2009
25 //
26 // Filename: src-IL/src/il_quantizer.c
27 //
28 // Description: Heavily modified by Denton Woods.
29 //
30 // 20040223 XIX : Modified so it works better with color requests < 256
31 // pallete always has memory space for 256 entries
32 // used so we can quant down to 255 colors then add a transparent color in there.
33 //
34 //-----------------------------------------------------------------------------
35
36 #include "il_internal.h"
37
38 #define MAXCOLOR 256
39 #define RED 2
40 #define GREEN 1
41 #define BLUE 0
42
43 typedef struct Box
44 {
45 ILint r0; // min value, exclusive
46 ILint r1; // max value, inclusive
47 ILint g0;
48 ILint g1;
49 ILint b0;
50 ILint b1;
51 ILint vol;
52 } Box;
53
54 /* Histogram is in elements 1..HISTSIZE along each axis,
55 * element 0 is for base or marginal value
56 * NB: these must start out 0!
57 */
58
59 ILfloat gm2[33][33][33];
60 ILint wt[33][33][33], mr[33][33][33], mg[33][33][33], mb[33][33][33];
61 ILuint size; //image size
62 ILint K; //colour look-up table size
63 ILushort *Qadd;
64
65
66 ILint WindW, WindH, WindD;
67 ILint i;
68 ILubyte *buffer;
69 static ILint Width, Height, Depth, Comp;
70 /*ILint TotalColors;
71 ILint a, b;
72 ILubyte *buf1, *buf2;*/
73
n2(ILint s)74 ILuint n2(ILint s)
75 {
76 ILint i;
77 ILint res;
78
79 res = 1;
80 for (i = 0; i < s; i++) {
81 res = res*2;
82 }
83
84 return res;
85 }
86
87
88 // Build 3-D color histogram of counts, r/g/b, c^2
Hist3d(ILubyte * Ir,ILubyte * Ig,ILubyte * Ib,ILint * vwt,ILint * vmr,ILint * vmg,ILint * vmb,ILfloat * m2)89 ILboolean Hist3d(ILubyte *Ir, ILubyte *Ig, ILubyte *Ib, ILint *vwt, ILint *vmr, ILint *vmg, ILint *vmb, ILfloat *m2)
90 {
91 ILint ind, r, g, b;
92 ILint inr, ing, inb, table[2560];
93 ILuint i;
94
95 for (i = 0; i < 256; i++)
96 {
97 table[i] = i * i;
98 }
99 Qadd = (ILushort*)ialloc(sizeof(ILushort) * size);
100 if (Qadd == NULL)
101 {
102 return IL_FALSE;
103 }
104
105 imemclear(Qadd, sizeof(ILushort) * size);
106
107 for (i = 0; i < size; i++)
108 {
109 r = Ir[i]; g = Ig[i]; b = Ib[i];
110 inr = (r>>3) + 1;
111 ing = (g>>3) + 1;
112 inb = (b>>3) + 1;
113 Qadd[i] = ind = (inr<<10)+(inr<<6)+inr+(ing<<5)+ing+inb;
114 //[inr][ing][inb]
115 vwt[ind]++;
116 vmr[ind] += r;
117 vmg[ind] += g;
118 vmb[ind] += b;
119 m2[ind] += (ILfloat)(table[r]+table[g]+table[b]);
120 }
121 return IL_TRUE;
122 }
123
124 /* At conclusion of the histogram step, we can interpret
125 * wt[r][g][b] = sum over voxel of P(c)
126 * mr[r][g][b] = sum over voxel of r*P(c) , similarly for mg, mb
127 * m2[r][g][b] = sum over voxel of c^2*P(c)
128 * Actually each of these should be divided by 'size' to give the usual
129 * interpretation of P() as ranging from 0 to 1, but we needn't do that here.
130 */
131
132 /* We now convert histogram into moments so that we can rapidly calculate
133 * the sums of the above quantities over any desired Box.
134 */
135
136
137 // Compute cumulative moments
M3d(ILint * vwt,ILint * vmr,ILint * vmg,ILint * vmb,ILfloat * m2)138 void M3d(ILint *vwt, ILint *vmr, ILint *vmg, ILint *vmb, ILfloat *m2)
139 {
140 ILushort ind1, ind2;
141 ILubyte i, r, g, b;
142 ILint line, line_r, line_g, line_b, area[33], area_r[33], area_g[33], area_b[33];
143 ILfloat line2, area2[33];
144
145 for (r = 1; r <= 32; r++) {
146 for (i = 0; i <= 32; i++) {
147 area2[i] = 0.0f;
148 area[i]=area_r[i]=area_g[i]=area_b[i]=0;
149 }
150 for (g = 1; g <= 32; g++) {
151 line2 = 0.0f;
152 line = line_r = line_g = line_b = 0;
153 for (b = 1; b <= 32; b++) {
154 ind1 = (r<<10) + (r<<6) + r + (g<<5) + g + b; // [r][g][b]
155 line += vwt[ind1];
156 line_r += vmr[ind1];
157 line_g += vmg[ind1];
158 line_b += vmb[ind1];
159 line2 += m2[ind1];
160 area[b] += line;
161 area_r[b] += line_r;
162 area_g[b] += line_g;
163 area_b[b] += line_b;
164 area2[b] += line2;
165 ind2 = ind1 - 1089; // [r-1][g][b]
166 vwt[ind1] = vwt[ind2] + area[b];
167 vmr[ind1] = vmr[ind2] + area_r[b];
168 vmg[ind1] = vmg[ind2] + area_g[b];
169 vmb[ind1] = vmb[ind2] + area_b[b];
170 m2[ind1] = m2[ind2] + area2[b];
171 }
172 }
173 }
174
175 return;
176 }
177
178
179 // Compute sum over a Box of any given statistic
Vol(Box * cube,ILint mmt[33][33][33])180 ILint Vol(Box *cube, ILint mmt[33][33][33])
181 {
182 return( mmt[cube->r1][cube->g1][cube->b1]
183 -mmt[cube->r1][cube->g1][cube->b0]
184 -mmt[cube->r1][cube->g0][cube->b1]
185 +mmt[cube->r1][cube->g0][cube->b0]
186 -mmt[cube->r0][cube->g1][cube->b1]
187 +mmt[cube->r0][cube->g1][cube->b0]
188 +mmt[cube->r0][cube->g0][cube->b1]
189 -mmt[cube->r0][cube->g0][cube->b0] );
190 }
191
192 /* The next two routines allow a slightly more efficient calculation
193 * of Vol() for a proposed subBox of a given Box. The sum of Top()
194 * and Bottom() is the Vol() of a subBox split in the given direction
195 * and with the specified new upper bound.
196 */
197
198 // Compute part of Vol(cube, mmt) that doesn't depend on r1, g1, or b1
199 // (depending on dir)
Bottom(Box * cube,ILubyte dir,ILint mmt[33][33][33])200 ILint Bottom(Box *cube, ILubyte dir, ILint mmt[33][33][33])
201 {
202 switch(dir)
203 {
204 case RED:
205 return( -mmt[cube->r0][cube->g1][cube->b1]
206 +mmt[cube->r0][cube->g1][cube->b0]
207 +mmt[cube->r0][cube->g0][cube->b1]
208 -mmt[cube->r0][cube->g0][cube->b0] );
209 break;
210 case GREEN:
211 return( -mmt[cube->r1][cube->g0][cube->b1]
212 +mmt[cube->r1][cube->g0][cube->b0]
213 +mmt[cube->r0][cube->g0][cube->b1]
214 -mmt[cube->r0][cube->g0][cube->b0] );
215 break;
216 case BLUE:
217 return( -mmt[cube->r1][cube->g1][cube->b0]
218 +mmt[cube->r1][cube->g0][cube->b0]
219 +mmt[cube->r0][cube->g1][cube->b0]
220 -mmt[cube->r0][cube->g0][cube->b0] );
221 break;
222 }
223 return 0;
224 }
225
226
227 // Compute remainder of Vol(cube, mmt), substituting pos for
228 // r1, g1, or b1 (depending on dir)
Top(Box * cube,ILubyte dir,ILint pos,ILint mmt[33][33][33])229 ILint Top(Box *cube, ILubyte dir, ILint pos, ILint mmt[33][33][33])
230 {
231 switch (dir)
232 {
233 case RED:
234 return( mmt[pos][cube->g1][cube->b1]
235 -mmt[pos][cube->g1][cube->b0]
236 -mmt[pos][cube->g0][cube->b1]
237 +mmt[pos][cube->g0][cube->b0] );
238 break;
239 case GREEN:
240 return( mmt[cube->r1][pos][cube->b1]
241 -mmt[cube->r1][pos][cube->b0]
242 -mmt[cube->r0][pos][cube->b1]
243 +mmt[cube->r0][pos][cube->b0] );
244 break;
245 case BLUE:
246 return( mmt[cube->r1][cube->g1][pos]
247 -mmt[cube->r1][cube->g0][pos]
248 -mmt[cube->r0][cube->g1][pos]
249 +mmt[cube->r0][cube->g0][pos] );
250 break;
251 }
252
253 return 0;
254 }
255
256
257 // Compute the weighted variance of a Box
258 // NB: as with the raw statistics, this is really the variance * size
Var(Box * cube)259 ILfloat Var(Box *cube)
260 {
261 ILfloat dr, dg, db, xx;
262
263 dr = (ILfloat)Vol(cube, mr);
264 dg = (ILfloat)Vol(cube, mg);
265 db = (ILfloat)Vol(cube, mb);
266 xx = gm2[cube->r1][cube->g1][cube->b1]
267 -gm2[cube->r1][cube->g1][cube->b0]
268 -gm2[cube->r1][cube->g0][cube->b1]
269 +gm2[cube->r1][cube->g0][cube->b0]
270 -gm2[cube->r0][cube->g1][cube->b1]
271 +gm2[cube->r0][cube->g1][cube->b0]
272 +gm2[cube->r0][cube->g0][cube->b1]
273 -gm2[cube->r0][cube->g0][cube->b0];
274
275 return xx - (dr*dr+dg*dg+db*db) / (ILfloat)Vol(cube, wt);
276 }
277
278 /* We want to minimize the sum of the variances of two subBoxes.
279 * The sum(c^2) terms can be ignored since their sum over both subBoxes
280 * is the same (the sum for the whole Box) no matter where we split.
281 * The remaining terms have a minus sign in the variance formula,
282 * so we drop the minus sign and MAXIMIZE the sum of the two terms.
283 */
284
Maximize(Box * cube,ILubyte dir,ILint first,ILint last,ILint * cut,ILint whole_r,ILint whole_g,ILint whole_b,ILint whole_w)285 ILfloat Maximize(Box *cube, ILubyte dir, ILint first, ILint last, ILint *cut,
286 ILint whole_r, ILint whole_g, ILint whole_b, ILint whole_w)
287 {
288 ILint half_r, half_g, half_b, half_w;
289 ILint base_r, base_g, base_b, base_w;
290 ILint i;
291 ILfloat temp, max;
292
293 base_r = Bottom(cube, dir, mr);
294 base_g = Bottom(cube, dir, mg);
295 base_b = Bottom(cube, dir, mb);
296 base_w = Bottom(cube, dir, wt);
297 max = 0.0;
298 *cut = -1;
299
300 for (i = first; i < last; ++i) {
301 half_r = base_r + Top(cube, dir, i, mr);
302 half_g = base_g + Top(cube, dir, i, mg);
303 half_b = base_b + Top(cube, dir, i, mb);
304 half_w = base_w + Top(cube, dir, i, wt);
305 // Now half_x is sum over lower half of Box, if split at i
306 if (half_w == 0) { // subBox could be empty of pixels!
307 continue; // never split into an empty Box
308 }
309 else {
310 temp = ((ILfloat)half_r*half_r + (ILfloat)half_g * half_g +
311 (ILfloat)half_b*half_b) / half_w;
312 }
313
314 half_r = whole_r - half_r;
315 half_g = whole_g - half_g;
316 half_b = whole_b - half_b;
317 half_w = whole_w - half_w;
318 if (half_w == 0) { // subBox could be empty of pixels!
319 continue; // never split into an empty Box
320 }
321 else {
322 temp += ((ILfloat)half_r*half_r + (ILfloat)half_g * half_g +
323 (ILfloat)half_b*half_b) / half_w;
324 }
325
326 if (temp > max) {
327 max = temp;
328 *cut = i;
329 }
330 }
331
332 return max;
333 }
334
335
Cut(Box * set1,Box * set2)336 ILint Cut(Box *set1, Box *set2)
337 {
338 ILubyte dir;
339 ILint cutr, cutg, cutb;
340 ILfloat maxr, maxg, maxb;
341 ILint whole_r, whole_g, whole_b, whole_w;
342
343 whole_r = Vol(set1, mr);
344 whole_g = Vol(set1, mg);
345 whole_b = Vol(set1, mb);
346 whole_w = Vol(set1, wt);
347
348 maxr = Maximize(set1, RED, set1->r0+1, set1->r1, &cutr, whole_r, whole_g, whole_b, whole_w);
349 maxg = Maximize(set1, GREEN, set1->g0+1, set1->g1, &cutg, whole_r, whole_g, whole_b, whole_w);
350 maxb = Maximize(set1, BLUE, set1->b0+1, set1->b1, &cutb, whole_r, whole_g, whole_b, whole_w);
351
352 if ((maxr >= maxg) && (maxr >= maxb)) {
353 dir = RED;
354 if (cutr < 0)
355 return 0; // can't split the Box
356 }
357 else if ((maxg >= maxr) && (maxg >= maxb))
358 dir = GREEN;
359 else
360 dir = BLUE;
361
362 set2->r1 = set1->r1;
363 set2->g1 = set1->g1;
364 set2->b1 = set1->b1;
365
366 switch (dir)
367 {
368 case RED:
369 set2->r0 = set1->r1 = cutr;
370 set2->g0 = set1->g0;
371 set2->b0 = set1->b0;
372 break;
373 case GREEN:
374 set2->g0 = set1->g1 = cutg;
375 set2->r0 = set1->r0;
376 set2->b0 = set1->b0;
377 break;
378 case BLUE:
379 set2->b0 = set1->b1 = cutb;
380 set2->r0 = set1->r0;
381 set2->g0 = set1->g0;
382 break;
383 }
384
385 set1->vol = (set1->r1-set1->r0) * (set1->g1-set1->g0) * (set1->b1-set1->b0);
386 set2->vol = (set2->r1-set2->r0) * (set2->g1-set2->g0) * (set2->b1-set2->b0);
387
388 return 1;
389 }
390
391
Mark(struct Box * cube,int label,unsigned char * tag)392 void Mark(struct Box *cube, int label, unsigned char *tag)
393 {
394 ILint r, g, b;
395
396 for (r = cube->r0 + 1; r <= cube->r1; r++) {
397 for (g = cube->g0 + 1; g <= cube->g1; g++) {
398 for (b = cube->b0 + 1; b <= cube->b1; b++) {
399 tag[(r<<10) + (r<<6) + r + (g<<5) + g + b] = label;
400 }
401 }
402 }
403 return;
404 }
405
406
iQuantizeImage(ILimage * Image,ILuint NumCols)407 ILimage *iQuantizeImage(ILimage *Image, ILuint NumCols)
408 {
409 Box cube[MAXCOLOR];
410 ILubyte *tag = NULL;
411 ILubyte lut_r[MAXCOLOR], lut_g[MAXCOLOR], lut_b[MAXCOLOR];
412 ILint next;
413 ILint weight;
414 ILuint k;
415 ILfloat vv[MAXCOLOR], temp;
416 //ILint color_num;
417 ILubyte *NewData = NULL, *Palette = NULL;
418 ILimage *TempImage = NULL, *NewImage = NULL;
419 ILubyte *Ir = NULL, *Ig = NULL, *Ib = NULL;
420
421 ILint num_alloced_colors; // number of colors we allocated space for in palette, as NumCols but will not be less than 256
422
423 num_alloced_colors=NumCols;
424 if(num_alloced_colors<256) { num_alloced_colors=256; }
425
426
427 NewImage = iCurImage;
428 iCurImage = Image;
429 TempImage = iConvertImage(iCurImage, IL_RGB, IL_UNSIGNED_BYTE);
430 iCurImage = NewImage;
431
432
433
434 if (TempImage == NULL)
435 return NULL;
436
437 buffer = Image->Data;
438 WindW = Width = Image->Width;
439 WindH = Height = Image->Height;
440 WindD = Depth = Image->Depth;
441 Comp = Image->Bpp;
442 Qadd = NULL;
443
444 //color_num = ImagePrecalculate(Image);
445
446 NewData = (ILubyte*)ialloc(Image->Width * Image->Height * Image->Depth);
447 Palette = (ILubyte*)ialloc(3 * num_alloced_colors);
448 if (!NewData || !Palette) {
449 ifree(NewData);
450 ifree(Palette);
451 return NULL;
452 }
453
454 Ir = (ILubyte*)ialloc(Width * Height * Depth);
455 Ig = (ILubyte*)ialloc(Width * Height * Depth);
456 Ib = (ILubyte*)ialloc(Width * Height * Depth);
457 if (!Ir || !Ig || !Ib) {
458 ifree(Ir);
459 ifree(Ig);
460 ifree(Ib);
461 ifree(NewData);
462 ifree(Palette);
463 return NULL;
464 }
465
466 size = Width * Height * Depth;
467
468 #ifdef ALTIVEC_GCC
469 register ILuint v_size = size>>4;
470 register ILuint pos = 0;
471 v_size = v_size /3;
472 register vector unsigned char d0,d1,d2;
473 register vector unsigned char red[3],blu[3],green[3];
474
475 register union{
476 vector unsigned char vec;
477 vector unsigned int load;
478 } mask_1, mask_2, mask_3;
479
480 mask_1.load = (vector unsigned int){0xFF0000FF,0x0000FF00,0x00FF0000,0xFF0000FF};
481 mask_2.load = (vector unsigned int){0x00FF0000,0xFF0000FF,0x0000FF00,0x00FF0000};
482 mask_2.load = (vector unsigned int){0x0000FF00,0x00FF0000,0xFF0000FF,0x0000FF00};
483
484 while( v_size >= 0 ) {
485 d0 = vec_ld(pos,TempImage->Data);
486 d1 = vec_ld(pos+16,TempImage->Data);
487 d2 = vec_ld(pos+32,TempImage->Data);
488
489 red[0] = vec_and(d0,mask_1.vec);
490 green[0] = vec_and(d0,mask_2.vec);
491 blu[0] = vec_and(d0,mask_3.vec);
492
493 red[1] = vec_and(d1,mask_3.vec);
494 green[1] = vec_and(d1,mask_1.vec);
495 blu[1] = vec_and(d1,mask_2.vec);
496
497 red[2] = vec_and(d2,mask_2.vec);
498 green[2] = vec_and(d2,mask_3.vec);
499 blu[2] = vec_and(d2,mask_1.vec);
500
501 vec_st(red[0],pos,Ir);
502 vec_st(red[1],pos+16,Ir);
503 vec_st(red[2],pos+32,Ir);
504
505 vec_st(blu[0],pos,Ib);
506 vec_st(blu[1],pos+16,Ib);
507 vec_st(blu[2],pos+32,Ib);
508
509 vec_st(green[0],pos,Ig);
510 vec_st(green[1],pos+16,Ig);
511 vec_st(green[2],pos+32,Ig);
512
513 pos += 48;
514 }
515 size -= pos;
516 #endif
517
518 for (k = 0; k < size; k++) {
519 Ir[k] = TempImage->Data[k * 3];
520 Ig[k] = TempImage->Data[k * 3 + 1];
521 Ib[k] = TempImage->Data[k * 3 + 2];
522 }
523
524 #ifdef ALTIVEC_GCC
525 size = Width * Height * Depth;
526 #endif
527
528 // Set new colors number
529 K = NumCols;
530
531 if (K <= 256) {
532 // Begin Wu's color quantization algorithm
533
534 // May have "leftovers" from a previous run.
535
536 imemclear(gm2, 33 * 33 * 33 * sizeof(ILfloat));
537 imemclear(wt, 33 * 33 * 33 * sizeof(ILint));
538 imemclear(mr, 33 * 33 * 33 * sizeof(ILint));
539 imemclear(mg, 33 * 33 * 33 * sizeof(ILint));
540 imemclear(mb, 33 * 33 * 33 * sizeof(ILint));
541
542 if (!Hist3d(Ir, Ig, Ib, (ILint*)wt, (ILint*)mr, (ILint*)mg, (ILint*)mb, (ILfloat*)gm2))
543 goto error_label;
544
545 M3d((ILint*)wt, (ILint*)mr, (ILint*)mg, (ILint*)mb, (ILfloat*)gm2);
546
547 cube[0].r0 = cube[0].g0 = cube[0].b0 = 0;
548 cube[0].r1 = cube[0].g1 = cube[0].b1 = 32;
549 next = 0;
550 for (i = 1; i < K; ++i) {
551 if (Cut(&cube[next], &cube[i])) { // volume test ensures we won't try to cut one-cell Box */
552 vv[next] = (cube[next].vol>1) ? Var(&cube[next]) : 0.0f;
553 vv[i] = (cube[i].vol>1) ? Var(&cube[i]) : 0.0f;
554 }
555 else {
556 vv[next] = 0.0; // don't try to split this Box again
557 i--; // didn't create Box i
558 }
559 next = 0;
560 temp = vv[0];
561 for (k = 1; (ILint)k <= i; ++k) {
562 if (vv[k] > temp) {
563 temp = vv[k]; next = k;
564 }
565 }
566
567 if (temp <= 0.0) {
568 K = i+1;
569 // Only got K Boxes
570 break;
571 }
572 }
573
574 tag = (ILubyte*)ialloc(33 * 33 * 33 * sizeof(ILubyte));
575 if (tag == NULL)
576 goto error_label;
577 for (k = 0; (ILint)k < K; k++) {
578 Mark(&cube[k], k, tag);
579 weight = Vol(&cube[k], wt);
580 if (weight) {
581 lut_r[k] = (ILubyte)(Vol(&cube[k], mr) / weight);
582 lut_g[k] = (ILubyte)(Vol(&cube[k], mg) / weight);
583 lut_b[k] = (ILubyte)(Vol(&cube[k], mb) / weight);
584 }
585 else {
586 // Bogus Box
587 lut_r[k] = lut_g[k] = lut_b[k] = 0;
588 }
589 }
590
591 for (i = 0; i < (ILint)size; i++) {
592 NewData[i] = tag[Qadd[i]];
593 }
594 ifree(tag);
595 ifree(Qadd);
596
597 for (k = 0; k < NumCols; k++) {
598 Palette[k * 3] = lut_b[k];
599 Palette[k * 3 + 1] = lut_g[k];
600 Palette[k * 3 + 2] = lut_r[k];
601 }
602 }
603 else { // If colors more than 256
604 // Begin Octree quantization
605 //Quant_Octree(Image->Width, Image->Height, Image->Bpp, buffer2, NewData, Palette, K);
606 ilSetError(IL_INTERNAL_ERROR); // Can't get much more specific than this.
607 goto error_label;
608 }
609
610 ifree(Ig);
611 ifree(Ib);
612 ifree(Ir);
613 ilCloseImage(TempImage);
614
615 NewImage = (ILimage*)icalloc(sizeof(ILimage), 1);
616 if (NewImage == NULL) {
617 return NULL;
618 }
619 ilCopyImageAttr(NewImage, Image);
620 NewImage->Bpp = 1;
621 NewImage->Bps = Image->Width;
622 NewImage->SizeOfPlane = NewImage->Bps * Image->Height;
623 NewImage->SizeOfData = NewImage->SizeOfPlane;
624 NewImage->Format = IL_COLOUR_INDEX;
625 NewImage->Type = IL_UNSIGNED_BYTE;
626
627 NewImage->Pal.Palette = Palette;
628 NewImage->Pal.PalSize = 256 * 3;
629 NewImage->Pal.PalType = IL_PAL_BGR24;
630 NewImage->Data = NewData;
631
632 return NewImage;
633
634 error_label:
635 ifree(NewData);
636 ifree(Palette);
637 ifree(Ig);
638 ifree(Ib);
639 ifree(Ir);
640 ifree(tag);
641 ifree(Qadd);
642 return NULL;
643 }
644