#include #include #include #include #include "ImagePalc.h" /*** * Functions and structures relating to color reduction. * Uses ImageMagick's "octree" algorithm: * http://www.cs.uni.edu/Help/ImageMagick/www/quantize.html */ struct FlatTree { int ncolors; /* values stored in and below this root */ int nhere; /* values stored *in* this root */ int Sr; /* Red sum *in* this root */ int Sg; /* Green sum *in* this root */ int Sb; /* Blue sum *in* this root */ int Serror; /* Sum of squares of error in *and below* this root */ int nchild; /* Number of non-NULL child entries */ struct FlatTree *child[8]; struct FlatTree *parent; }; /* The number of colors left in the image "at this moment" */ static int TotalNumberOfFlatNodes; static int dist3d(int,int,int, int,int,int); static struct FlatTree *newFlatTree(struct FlatTree *parent); static int addtoFlatTree(struct FlatTree **p, unsigned char color[3], int maxdepth); static void pruneFlatTree(struct FlatTree **p); static void lookupFlatTree(struct FlatTree **p, unsigned char (*color)[3]); static void clearFlatTree(struct FlatTree **p); static struct FlatTree **findMinError(struct FlatTree **root); static void palettizeFlatTree(struct FlatTree *p, unsigned char (*palette)[3], int *pal_idx); static int createFlatTreePalette(struct FlatTree **cubes, unsigned char (*palette)[3], int npalette); static int dist3d(int x1, int y1, int z1, int x2, int y2, int z2) { x1 -= x2; y1 -= y2; z1 -= z2; return (int) sqrt(x1*x1 + y1*y1 + z1*z1); } static struct FlatTree *newFlatTree(struct FlatTree *parent) { struct FlatTree *p = malloc(sizeof *p); if (p == NULL) return NULL; p->ncolors = 0; p->nhere = 0; p->Sr = p->Sg = p->Sb = 0; p->Serror = 0; p->nchild = 0; p->child[0] = NULL; p->child[1] = NULL; p->child[2] = NULL; p->child[3] = NULL; p->child[4] = NULL; p->child[5] = NULL; p->child[6] = NULL; p->child[7] = NULL; p->parent = parent; return p; } static int addtoFlatTree(struct FlatTree **p, unsigned char color[3], int maxdepth) { struct FlatTree *cur; int i, c; if (*p == NULL) *p = newFlatTree(NULL); if (*p == NULL) return -1; if (maxdepth > 8) maxdepth = 8; cur = *p; for (i=0; i < maxdepth; ++i) { ++cur->ncolors; cur->Serror += dist3d(color[0]%(1<<(8-i)), color[1]%(1<<(8-i)), color[2]%(1<<(8-i)), color[0]&(1<<(7-i)), color[1]&(1<<(7-i)), color[2]&(1<<(7-i))); c = ((color[0]>>(7-i)) & 1) << 0 | ((color[1]>>(7-i)) & 1) << 1 | ((color[2]>>(7-i)) & 1) << 2; if (cur->child[c] == NULL) { cur->child[c] = newFlatTree(cur); if (cur->child[c] == NULL) return -1; ++cur->nchild; } cur = cur->child[c]; } ++cur->ncolors; if (cur->nhere == 0) ++TotalNumberOfFlatNodes; ++cur->nhere; cur->Sr += color[0]; cur->Sg += color[1]; cur->Sb += color[2]; cur->Serror += dist3d(color[0]%(1<<(8-i)), color[1]%(1<<(8-i)), color[2]%(1<<(8-i)), color[0]&(1<<(7-i)), color[1]&(1<<(7-i)), color[2]&(1<<(7-i))); return 0; } static void pruneFlatTree(struct FlatTree **p) { /* Find the node rooted at p with the lowest error */ struct FlatTree **minError = p; if (!p || !*p) return; /* While that node has children, repeat */ while ((*minError)->nchild) { findMinError(NULL); minError = findMinError(minError); } /* Remove that leaf node from the tree by propagating upwards */ if ((*minError)->parent == NULL) { /* Just delete the node. */ if ((*minError)->nhere) --TotalNumberOfFlatNodes; clearFlatTree(minError); } else { /* Propagate the statistics upwards */ (*minError)->parent->Sr += (*minError)->Sr; (*minError)->parent->Sg += (*minError)->Sg; (*minError)->parent->Sb += (*minError)->Sb; if ((*minError)->nhere) { --TotalNumberOfFlatNodes; if ((*minError)->parent->nhere == 0) ++TotalNumberOfFlatNodes; } (*minError)->parent->nhere += (*minError)->nhere; /* Remove one child of this one's parent. */ --(*minError)->parent->nchild; clearFlatTree(minError); } } static struct FlatTree **findMinError(struct FlatTree **root) { static struct FlatTree **absMinimum = NULL; int i; if (root == NULL) { /* Flush */ absMinimum = NULL; return NULL; } if ((*root)->nchild == 0) { if (absMinimum == NULL || (*absMinimum)->Serror > (*root)->Serror) absMinimum = root; return absMinimum; } for (i=0; i < 8; ++i) { if ((*root)->child[i] != NULL) { (void)findMinError(&(*root)->child[i]); } } return absMinimum; } static void lookupFlatTree(struct FlatTree **p, unsigned char (*color)[3]) { struct FlatTree *cur; int i, c; if (*p == NULL) { (*color)[0] = 0; (*color)[1] = 0; (*color)[2] = 0; return; } cur = *p; for (i=0; i < 8; ++i) { c = (((*color)[0]>>(7-i)) & 1) << 0 | (((*color)[1]>>(7-i)) & 1) << 1 | (((*color)[2]>>(7-i)) & 1) << 2; if (cur->child[c] == NULL) { if (cur->nhere) { (*color)[0] = cur->Sr / cur->nhere; (*color)[1] = cur->Sg / cur->nhere; (*color)[2] = cur->Sb / cur->nhere; return; } else { for (c=0; c < 8; ++c) if (cur->child[c] != NULL) break; if (c == 8) { fprintf(stderr, "Yipes!\n"); exit(1); } } } cur = cur->child[c]; } (*color)[0] = cur->Sr / cur->nhere; (*color)[1] = cur->Sg / cur->nhere; (*color)[2] = cur->Sb / cur->nhere; } /* Record all the colors in this tree. */ static void palettizeFlatTree(struct FlatTree *p, unsigned char (*palette)[3], int *pal_idx) { int i; if (p == NULL) return; if (p->nhere) { palette[*pal_idx][0] = p->Sr / p->nhere; palette[*pal_idx][1] = p->Sg / p->nhere; palette[*pal_idx][2] = p->Sb / p->nhere; ++*pal_idx; } for (i=0; i < 8; ++i) { palettizeFlatTree(p->child[i], palette, pal_idx); } } static void clearFlatTree(struct FlatTree **p) { int i; if (!p || !*p) return; for (i=0; i < 8; ++i) { if ((*p)->child[i] != NULL) clearFlatTree(&(*p)->child[i]); } *p = NULL; return; } static int createFlatTreePalette(struct FlatTree **cubes, unsigned char (*palette)[3], int npalette) { int i; TotalNumberOfFlatNodes = 0; for (i=0; i < npalette; ++i) { if (addtoFlatTree(cubes, palette[i], 8) == -1) { clearFlatTree(cubes); return -3; } } return 0; } /*** * Functions and structures relating to image dithering. * Uses a slight modification of Victor Ostromoukhov's algorithm, * version 1 (August 2001), as described in the article * "A Simple and Efficient Error-Diffusion Algorithm" (SIGGRAPH'01) * http://www.iro.umontreal.ca/~ostrom/ */ /*------ variable coefficients, packed for readability ------*/ static unsigned vc_tab[256] = { 0x0400d000, 0x0400d000, 0x07015000, 0x02007000, 0x03008000, 0x1302f003, 0x09017003, 0x0600f003, 0x09016006, 0x1302b00f, 0x03007003, 0xea1f50e0, 0x750f9074, 0x4e0a5050, 0x3a07b03e, 0xea1e9100, 0x2705102c, 0xea1e3110, 0x1d03c023, 0x1a035020, 0x750ed094, 0xea1d7130, 0x01003002, 0xe71cb130, 0x13026019, 0xe71c5128, 0x730e1092, 0x4d095060, 0x3906f047, 0x2103f028, 0x2604902e, 0xe71b3110, 0x3906c043, 0x0700d008, 0x730d5082, 0xe71a7100, 0x02005003, 0x9a1190ad, 0x4d08d059, 0x9a11b0b7, 0x2604702f, 0x9a11d0c1, 0x0700d009, 0x1602901d, 0x1302401a, 0x9a1210d5, 0x4d09106d, 0x9a1230df, 0x26049039, 0x9a1250e9, 0x0b015011, 0x9a1270f3, 0x1302501f, 0x0e01b017, 0x4d095081, 0x9a12b107, 0x2604b043, 0x1602b027, 0x4d09708b, 0x9a12f11b, 0x13026024, 0x9a131125, 0x4d099095, 0x9a13312f, 0x00001001, 0x34065069, 0x1a031035, 0x3405f06b, 0x0d01701b, 0x3405906d, 0x1a02b037, 0x3405306f, 0x03005007, 0x610ac0b5, 0x3006104c, 0x20048029, 0x3007702f, 0x01004001, 0x01004001, 0x01004001, 0x01004001, 0x01004001, 0x01004001, 0x01004001, 0x01004001, 0x01004001, 0x19041012, 0x2505f01d, 0x4b0b903e, 0x0c01e00b, 0x0f02300e, 0x25055025, 0x1903701a, 0x25050029, 0x4b09b056, 0x02005003, 0x02005003, 0x02005003, 0x02005003, 0x02005003, 0x02005003, 0x02005003, 0x02005003, 0x02005003, 0x02005003, 0x02005003, 0x02005003, 0x02005003, 0x961310b0, 0x4b09b056, 0x32069038, 0x25050029, 0x1e041020, 0x1903701a, 0x9614f098, 0x25055025, 0x32073030, 0x0f02300e, 0x96163088, 0x0c01e00b, 0x9616d080, 0x4b0b903e, 0x0a019008, 0x2505f01d, 0x96181070, 0x19041012, 0x9618b068, 0x01004001, 0x01004001, 0x9618b068, 0x19041012, 0x96181070, 0x2505f01d, 0x0a019008, 0x4b0b903e, 0x9616d080, 0x0c01e00b, 0x96163088, 0x0f02300e, 0x32073030, 0x25055025, 0x9614f098, 0x1903701a, 0x1e041020, 0x25050029, 0x32069038, 0x4b09b056, 0x961310b0, 0x02005003, 0x02005003, 0x02005003, 0x02005003, 0x02005003, 0x02005003, 0x02005003, 0x02005003, 0x02005003, 0x02005003, 0x02005003, 0x02005003, 0x02005003, 0x4b09b056, 0x25050029, 0x1903701a, 0x25055025, 0x0f02300e, 0x0c01e00b, 0x4b0b903e, 0x2505f01d, 0x19041012, 0x01004001, 0x01004001, 0x01004001, 0x01004001, 0x01004001, 0x01004001, 0x01004001, 0x01004001, 0x01004001, 0x3007702f, 0x20048029, 0x3006104c, 0x610ac0b5, 0x03005007, 0x3405306f, 0x1a02b037, 0x3405906d, 0x0d01701b, 0x3405f06b, 0x1a031035, 0x34065069, 0x00001001, 0x9a13312f, 0x4d099095, 0x9a131125, 0x13026024, 0x9a12f11b, 0x4d09708b, 0x1602b027, 0x2604b043, 0x9a12b107, 0x4d095081, 0x0e01b017, 0x1302501f, 0x9a1270f3, 0x0b015011, 0x9a1250e9, 0x26049039, 0x9a1230df, 0x4d09106d, 0x9a1210d5, 0x1302401a, 0x1602901d, 0x0700d009, 0x9a11d0c1, 0x2604702f, 0x9a11b0b7, 0x4d08d059, 0x9a1190ad, 0x02005003, 0xe71a7100, 0x730d5082, 0x0700d008, 0x3906c043, 0xe71b3110, 0x2604902e, 0x2103f028, 0x3906f047, 0x4d095060, 0x730e1092, 0xe71c5128, 0x13026019, 0xe71cb130, 0x01003002, 0xea1d7130, 0x750ed094, 0x1a035020, 0x1d03c023, 0xea1e3110, 0x2705102c, 0xea1e9100, 0x3a07b03e, 0x4e0a5050, 0x750f9074, 0xea1f50e0, 0x03007003, 0x1302b00f, 0x09016006, 0x0600f003, 0x09017003, 0x1302f003, 0x03008000, 0x02007000, 0x07015000, 0x0400d000, 0x0400d000, }; static void shift_carry_buffers(double **carry_line_0, double **carry_line_1, int w) { int i; double *tmp = *carry_line_0; *carry_line_0 = *carry_line_1; *carry_line_1 = tmp; for (i=0; i < w; ++i) (*carry_line_1)[i] = 0.0; } static void shift_carry_buffers_RGB(double (**carry_line_0)[3], double (**carry_line_1)[3], int w) { int i; double (*tmp)[3] = *carry_line_0; *carry_line_0 = *carry_line_1; *carry_line_1 = tmp; for (i=0; i < w; ++i) { (*carry_line_1)[i][0] = 0.0; (*carry_line_1)[i][1] = 0.0; (*carry_line_1)[i][2] = 0.0; } } static void distribute_error(double *cl0, double *cl1, int x, double diff, int dir, int input_level) { unsigned long coefs = vc_tab[input_level]; unsigned sum = (((coefs >> 24) & 0xFF) << 2) + 2; unsigned term_r = ((coefs >> 12) & 0xFFF); unsigned term_dl = (coefs & 0xFFF); unsigned term_d = sum - (term_r + term_dl); if (term_d > sum) { sum = term_r + term_dl; term_d = 0; } cl0[x+dir] += term_r * diff/sum; cl1[x-dir] += term_dl * diff/sum; cl1[x] += term_d * diff/sum; } static void distribute_error_RGB(double (*cl0)[3], double (*cl1)[3], int x, double diff[3], int dir, unsigned char input[3]) { int i; for (i=0; i < 3; ++i) { unsigned long coefs = vc_tab[(int)input[i]]; unsigned sum = (((coefs >> 24) & 0xFF) << 2) + 2; unsigned term_r = ((coefs >> 12) & 0xFFF); unsigned term_dl = (coefs & 0xFFF); unsigned term_d = sum - (term_r + term_dl); if (term_d > sum) { sum = term_r + term_dl; term_d = 0; } cl0[x+dir][i] += term_r * diff[i]/sum; cl1[x-dir][i] += term_dl * diff[i]/sum; cl1[x][i] += term_d * diff[i]/sum; } } static void carryFlatTree(struct FlatTree **p, double dcorr[3], int (*color)[3]) { struct FlatTree *cur; int i, c; int corr[3]; /* Convert from floating point to integer so we can use the bit patterns */ corr[0] = dcorr[0]+0.5; corr[0] = (corr[0]>255)? 255: (corr[0]<0)? 0: corr[0]; corr[1] = dcorr[1]+0.5; corr[1] = (corr[1]>255)? 255: (corr[1]<0)? 0: corr[1]; corr[2] = dcorr[2]+0.5; corr[2] = (corr[2]>255)? 255: (corr[2]<0)? 0: corr[2]; if (*p == NULL) { memset(*color, 0x00, sizeof *color); return; } cur = *p; for (i=0; i < 8; ++i) { c = ((corr[0]>>(7-i)) & 1) << 0 | ((corr[1]>>(7-i)) & 1) << 1 | ((corr[2]>>(7-i)) & 1) << 2; if (cur->child[c] == NULL) { if (cur->nhere == 0) { int j; for (j=0; j < 8; ++j) if (cur->child[j]) break; cur = cur->child[j]; continue; } (*color)[0] = cur->Sr / cur->nhere; (*color)[1] = cur->Sg / cur->nhere; (*color)[2] = cur->Sb / cur->nhere; return; } cur = cur->child[c]; } (*color)[0] = cur->Sr / cur->nhere; (*color)[1] = cur->Sg / cur->nhere; (*color)[2] = cur->Sb / cur->nhere; } static int lookupPalette(unsigned char color[3], unsigned char (*pal)[3], int npal) { int i; for (i=0; i < npal; ++i) { if (memcmp(color, pal[i], 3) == 0) return i; } return -1; } /************************************************************************* * PUBLIC FUNCTIONS */ /* Create a palette with a maximum of nColors colors for this image. */ int CreateRGBPalette(unsigned char (*data)[3], int w, int h, unsigned char (**palette)[3], int nColors) { struct FlatTree *cubes = NULL; int picSize = w*h; int i; int maxdepth; int pal_idx; *palette = NULL; if (nColors == 0) { maxdepth = 8; } else { for (maxdepth = 0; (1 << (1 << maxdepth)) < nColors; ++maxdepth) /* continue */ ; maxdepth += 1; } TotalNumberOfFlatNodes = 0; for (i=0; i < picSize; ++i) { if (addtoFlatTree(&cubes, data[i], maxdepth) == -1) { clearFlatTree(&cubes); return -3; } } if (nColors) { while (TotalNumberOfFlatNodes > nColors) { pruneFlatTree(&cubes); } } *palette = malloc(TotalNumberOfFlatNodes * sizeof **palette); if (*palette == NULL) { clearFlatTree(&cubes); return -3; } pal_idx = 0; palettizeFlatTree(cubes, *palette, &pal_idx); clearFlatTree(&cubes); return pal_idx; } /* 'Posterize' the image using the specified palette. */ int RGB2Pal(unsigned char (*rgb)[3], unsigned char **pal, int w, int h, unsigned char (*palette)[3], int nColors) { struct FlatTree *cubes = NULL; int picSize = w*h; int i; int TempPaletteFlag = 0; *pal = NULL; if (nColors > 256 || (palette != NULL && nColors == 0)) return -1; /* Create the palette if necessary. */ if (palette == NULL) { int rc = CreateRGBPalette(rgb, w, h, &palette, nColors); if (rc < 0) return rc; nColors = rc; TempPaletteFlag = 1; } i = createFlatTreePalette(&cubes, palette, nColors); if (i != 0) { if (TempPaletteFlag) free(palette); return i; } *pal = malloc(picSize * sizeof **pal); if (*pal == NULL) { if (TempPaletteFlag) free(palette); return -3; } for (i=0; i < picSize; ++i) { unsigned char temp[3]; memcpy(temp, rgb[i], 3); lookupFlatTree(&cubes, &temp); (*pal)[i] = lookupPalette(temp, palette, nColors); } clearFlatTree(&cubes); if (TempPaletteFlag) free(palette); return 0; } int DitherGray2BW(unsigned char *input_image, int xdim_in, int ydim_in, unsigned char **output_image, int xdim_out, int ydim_out) { int OriginalImage = 0; int x, y; int *x_tab; int *y_tab; int xstart, xstop, xstep; double x_scale_factor; double y_scale_factor; double *carry_line_0; double *carry_line_1; double threshold = 127.5; /* allocate output image */ if (xdim_in==xdim_out && ydim_in==ydim_out && *output_image==input_image) { /* Use the original image as the output image! */ OriginalImage = 1; } else { /* Create a new working area. */ *output_image = malloc(xdim_out*ydim_out * sizeof **output_image); } if (*output_image == NULL) return -3; /* allocate carry_line_0 and carry_line_1 */ carry_line_0 = calloc(xdim_out+2, sizeof *carry_line_0); carry_line_1 = calloc(xdim_out+2, sizeof *carry_line_1); /* allocate and init x_tab and y_tab */ x_tab = calloc(xdim_out, sizeof *x_tab); y_tab = calloc(ydim_out, sizeof *y_tab); if (carry_line_0==NULL || carry_line_1==NULL || x_tab==NULL || y_tab==NULL) { if (!OriginalImage) { free(*output_image); *output_image = NULL; } free(carry_line_0); free(carry_line_1); free(x_tab); free(y_tab); return -3; } ++carry_line_0; ++carry_line_1; x_scale_factor = (double)xdim_out / xdim_in; y_scale_factor = (double)ydim_out / ydim_in; for (x = 0; x < xdim_out; ++x) x_tab[x] = (int)(x / x_scale_factor); for (y = 0; y < ydim_out; y++) y_tab[y] = (int)(y / y_scale_factor); for (y = 0; y < ydim_out; ++y) { int ysrc = y_tab[y]; if (y % 2) { xstart = xdim_out-1; xstop = -1; xstep = -1; } else { xstart = 0; xstop = xdim_out; xstep = 1; } for (x = xstart; x != xstop; x += xstep) { int xsrc = x_tab[x]; int input = input_image[ysrc*xdim_in+xsrc]; double corrected_level = input + carry_line_0[x]; int intensity = (corrected_level <= threshold)? 0: 255; double diff = corrected_level - intensity; distribute_error(carry_line_0, carry_line_1, x, diff, xstep, input); if (input == 0 || intensity == 0) (*output_image)[y*xdim_out+x] = 0; else (*output_image)[y*xdim_out+x] = 255; } shift_carry_buffers(&carry_line_0, &carry_line_1, xdim_out); } free(carry_line_0 - 1); free(carry_line_1 - 1); free(x_tab); free(y_tab); return 0; } /* * Dither a 24-bit image onto the provided palette. * If 'palette' is NULL, then CreateRGBPalette(data, w, h, [temp], nColors) * before dithering. */ int DitherRGB(unsigned char (*data)[3], int w, int h, unsigned char (*palette)[3], int nColors) { struct FlatTree *cubes = NULL; double (*carry_line_0)[3]; double (*carry_line_1)[3]; int x, y; int xstart, xstop, xstep; int TempPaletteFlag = 0; /* allocate carry_line_0 and carry_line_1 */ carry_line_0 = calloc(w+2, sizeof *carry_line_0); carry_line_1 = calloc(w+2, sizeof *carry_line_1); if (carry_line_0==NULL || carry_line_1==NULL) { free(carry_line_0); free(carry_line_1); return -3; } /* Create the palette if necessary. */ if (palette == NULL) { int rc = CreateRGBPalette(data, w, h, &palette, nColors); if (rc < 0) { free(carry_line_0); free(carry_line_1); return rc; } nColors = rc; TempPaletteFlag = 1; } if (createFlatTreePalette(&cubes, palette, nColors)) { free(carry_line_0); free(carry_line_1); if (TempPaletteFlag) free(palette); return -3; } ++carry_line_0; ++carry_line_1; for (y = 0; y < h; ++y) { if (y % 2) { xstart = w-1; xstop = -1; xstep = -1; } else { xstart = 0; xstop = w; xstep = 1; } for (x = xstart; x != xstop; x += xstep) { double corrected[3]; int intensity[3]; double diff[3]; corrected[0] = data[y*w+x][0] + carry_line_0[x][0]; corrected[1] = data[y*w+x][1] + carry_line_0[x][1]; corrected[2] = data[y*w+x][2] + carry_line_0[x][2]; carryFlatTree(&cubes, corrected, &intensity); diff[0] = corrected[0] - intensity[0]; diff[1] = corrected[1] - intensity[1]; diff[2] = corrected[2] - intensity[2]; distribute_error_RGB(carry_line_0, carry_line_1, x, diff, xstep, data[y*w+x]); data[y*w+x][0] = intensity[0]; data[y*w+x][1] = intensity[1]; data[y*w+x][2] = intensity[2]; } shift_carry_buffers_RGB(&carry_line_0, &carry_line_1, w); } free(carry_line_0 - 1); free(carry_line_1 - 1); clearFlatTree(&cubes); if (TempPaletteFlag) free(palette); return 0; } /* Convert palettized to RGB data. */ int Pal2RGB(unsigned char *pal, unsigned char (**rgb)[3], int w, int h, unsigned char (*palette)[3], int nColors) { int i; int picSize = w*h; /* Used only for bounds checking and consistency */ (void)nColors; *rgb = malloc(picSize * sizeof **rgb); if (*rgb == NULL) return -3; for (i=0; i < picSize; ++i) { memcpy((*rgb)[i], palette[pal[i]], 3); } return 0; }