/** demosaic_sharpen.c * when demosaicing the unbayered image, * don't just use bilinear interpolation, but weigh * the to be guessed values according to the differences * of the known value. * Of course, smaller differences mean higher weights. * * © Kurt Garloff <garloff@suse.de>, 2002/01/15 * License: GNU GPL * * Note: Interpolation techniques more intelligent than * bilinear inerpolation have been subject to investigations. * Those two links came from Jérôme Fenal: * http://ise.stanford.edu/class/psych221/99/dixiedog/vision.htm * http://www-ise.stanford.edu/~tingchen/main.htm * * From the evaluation, it seems that we should strive to do something like * "Linear Interpolation with Laplacian 2nd order color correction terms I". * * Here; I went my own way. * The reason for this is twofold: * - Avoid all possible patent issues * (If I would be American, I would probably want to patent this algo) * - Be conservative: By using an interpolation method (with weights bound * below 1), we avoid the risk to get artefacts, like e.g. seen in the * smooth hue algorithms, as our interpolated values are always in between * those from the neighbours, just a little closer to the ones we think * fits better. * * Our algorithm: * Trying to predict the known value based on the next neighbours with * the same colour component will yield weights. To achieve this: * - Measure the difference |dv| of the known colour component * - Choose a function f(|dv|) which prefers points with smaller |dv| * - Function should be monotonic and ]0;1[ * - Sum of weights must be 1, of course * * We choose f(|dv|) = N / (ALPHA + |dv|) * and scale with the reciprocal total sum, so we fulfill the conditions. * The algorithm is integer only (unless you enable DEBUG) * and therefore suitable for FPU-less machines and/or the kernel. * * I've chosen ALPHA = 2 and the results look really good. * * ToDo: * - A similar algorithm in HSI space might be slightly better. * - Different weighing functions might do better * - Do rigorous performance analysis (quality and computation cost) * in comparison to other algos as in papers cited above. * - There's the bilinear interpolation included for reference * (debugging purposes). Use it or get rid of it for slightly better * performance. * * I've implemented this algo here as a testbed. It's only tested for * BAYER_TILE_GBRG_INTERLACED, though it's been designed to be general * for all BAYERs in gphoto2. (Hence all those tables ...) * It should be moved to the gphoto2 infrastructure to help all * cameras, not just mine. It includes the demosaicing, so it should * be merged with the gp_bayer_decode (or the bilinear demosaicing * could be removed from the latter) to avoid double work. * * History: * 2001-01-15, 0.90, KG, * working for inner points (2,2)-(width-3,height-3) * 2001-01-15, 1.00, KG, * handle boundary points * */ #include <stdlib.h> #include "demosaic_sharpen.h" /* we define bayer as * +---> x * | 0 1 * v 2 3 * y */ typedef enum { RED = 0, GREEN = 1, BLUE = 2 } col; /* Don't get confused reading this code; there's lots of * indirection through the tables to avoid branches in the code; * maybe I love long pipelines too much. * If I look at the code long enough, I get confused myself. * The boundary special cases unfortunately do introduce some extra * branching. (KG) */ /* relative postition */ typedef struct _off { signed char dx, dy; } off; typedef enum { NB_DIAG = 0, NB_TLRB, NB_LR, NB_TB, NB_TLRB2 } nb_pat; /* locations of neighbour points with the same colour */ typedef struct _neighbours { unsigned char num; off nb_pts[4]; } neighbours; /* possible locations */ static const neighbours n_pos[8] = { { /* NB_DIAG */ 4, { {-1,-1}, { 1,-1}, {-1, 1}, { 1, 1} } },{ /* NB_TLRB */ 4, { { 0,-1}, {-1, 0}, { 1, 0}, { 0, 1} } },{ /* NB_LR */ 2, { {-1, 0}, { 1, 0}, { 0, 0}, { 0, 0} } },{ /* NB_TB */ 2, { { 0,-1}, { 0, 1}, { 0, 0}, { 0, 0} } },{ /* NB_TLRB2 */ 4, { { 0,-2}, {-2, 0}, { 2, 0}, { 0, 2}, } } }; typedef struct _t_coeff { unsigned char cf[4][4]; unsigned char num; } t_coeff; typedef enum { DIAG_TO_LR = 0, DIAG_TO_TB, TLRB2_TO_DIAG, TLRB2_TO_TLRB, PATCONV_NONE } patconv; /* Transfer matrix pattern to pattern */ static const t_coeff pat_to_pat[4] = { { /* DIAG_TO_LR */ { {2, 0, 2, 0}, {0, 2, 0, 2}, {0, 0, 0, 0}, {0, 0, 0, 0} }, 2 },{ /* DIAG_TO_TB */ { { 2, 2, 0, 0}, { 0, 0, 2, 2}, { 0, 0, 0, 0}, { 0, 0, 0, 0}, }, 2 },{ /* TLRB2_TO_DIAG */ { {1, 1, 0, 0}, {1, 0, 1, 0}, {0, 1, 0, 1}, {0, 0, 1, 1} }, 4 },{ /* TLRB2_TO_TLRB (trivial) */ { {2, 0, 0, 0}, {0, 2, 0, 0}, {0, 0, 2, 0}, {0, 0, 0, 2} }, 4 } }; static const patconv pconvmap[5][5] = { { PATCONV_NONE, PATCONV_NONE, DIAG_TO_LR, DIAG_TO_TB, PATCONV_NONE }, { PATCONV_NONE, PATCONV_NONE, PATCONV_NONE, PATCONV_NONE, PATCONV_NONE }, { PATCONV_NONE, PATCONV_NONE, PATCONV_NONE, PATCONV_NONE, PATCONV_NONE }, { PATCONV_NONE, PATCONV_NONE, PATCONV_NONE, PATCONV_NONE, PATCONV_NONE }, { TLRB2_TO_DIAG, TLRB2_TO_TLRB, PATCONV_NONE, PATCONV_NONE, PATCONV_NONE } }; /* Next mapping: col of pixel (0,1,2 = RGB & * index into n_pos for own, own+1, own+2 */ typedef struct _bayer_desc { col colour; nb_pat idx_pts[3]; } bayer_desc; /* T = Bayer Tile, P = Bayer point no * T P */ static const bayer_desc bayers[4][4] = { { /* TILE_RGGB */ { RED, {NB_TLRB2, NB_TLRB, NB_DIAG} }, { GREEN, {NB_DIAG, NB_TB, NB_LR} }, { GREEN, {NB_DIAG, NB_LR, NB_TB} }, { BLUE, {NB_TLRB2, NB_DIAG, NB_TLRB} }, },{ /* TILE_GRBG */ { GREEN, {NB_DIAG, NB_TB, NB_LR} }, { RED, {NB_TLRB2, NB_TLRB, NB_DIAG} }, { BLUE, {NB_TLRB2, NB_DIAG, NB_TLRB} }, { GREEN, {NB_DIAG, NB_LR, NB_TB} }, },{ /* TILE_BGGR */ { BLUE, {NB_TLRB2, NB_DIAG, NB_TLRB} }, { GREEN, {NB_DIAG, NB_LR, NB_TB} }, { GREEN, {NB_DIAG, NB_TB, NB_LR} }, { RED, {NB_TLRB2, NB_TLRB, NB_DIAG} }, },{ /* TILE_GBRG */ { GREEN, {NB_DIAG, NB_LR, NB_TB} }, { BLUE, {NB_TLRB2, NB_DIAG, NB_TLRB} }, { RED, {NB_TLRB2, NB_TLRB, NB_DIAG} }, { GREEN, {NB_DIAG, NB_TB, NB_LR} } } }; /* Use integer arithmetic. Accuracy is 10^-6, which is good enough */ #define SHIFT 20 static inline int weight (const unsigned char dx, const int alpha) { return (1<<SHIFT)/(alpha + dx); } /* alpha controls the strength of the weighting; 1 = strongest, 64 = weak */ void demosaic_sharpen (const int width, const int height, const unsigned char * const src_region, unsigned char * const dest_region, const int alpha, const BayerTile bt) { const unsigned char* src_ptr = src_region; unsigned char* dst_ptr = dest_region; const bayer_desc *bay_des = bayers [bt & 3]; /* Don't care about interlace */ int x, y; for (y = 0; y < height; y++) { for (x = 0; x < width; x++, src_ptr += 3, dst_ptr += 3) { /* 3 2 */ /* 1 0 */ const unsigned char bayer = (1^(x&1)) + ((1^(y&1))<<1); const col colour = bay_des[bayer].colour; /* nb_pat[0] is our own pattern */ const nb_pat * const nbpts = bay_des[bayer].idx_pts; /* less strong weighting for TLRB2 pattern */ const int myalpha = (*nbpts == NB_TLRB2? (alpha << 1): alpha); const unsigned char colval = src_ptr[colour]; int weights[4]; int sum_weights = 0.0; patconv pconv; /* Calc coeffs for prediction */ int nbs; col ncol; int othcol; int i; int skno; int nsumw; int predcol = 0; /* Only for DEBUG */ /* DPRINTF("(%i,%i)(%p): bay %i, col %i, pat %i, val %i\n", x, y, src_ptr, bayer, colour, nbpts[0], colval);*/ /* Copy own colour */ dst_ptr[colour] = colval; /* Now calc weights */ for (nbs = 0; nbs < 4; nbs++) { const off offset = n_pos[nbpts[0]].nb_pts[nbs]; const int nx = x + offset.dx; const int ny = y + offset.dy; const signed long addr_off = 3 * (offset.dx + width * offset.dy); unsigned char thisval = colval; int coeff = 0; if (nx >= 0 && nx < width && ny >= 0 && ny < height) { thisval = src_ptr[addr_off+colour]; coeff = weight (abs ((int)colval - thisval), myalpha); } else if (nbpts[0] == NB_TLRB2 && x > 0 && x < width-1 && y > 0 && y < height-1) { coeff = weight (128, myalpha); /* assign some small weight */ } /*DPRINTF(" (%i,%i)(%p): val %i, diff %i, weight %i\n", nx, ny, src_ptr+addr_off, thisval, abs ((int)colval - thisval), coeff);*/ predcol += thisval * coeff; weights[nbs] = coeff; sum_weights += coeff; }; #ifdef DEBUG printf(" Coeffs:"); for (nbs = 0; nbs < 4; nbs++) printf (" %6.4f", (double)weights[nbs]/sum_weights); printf (" -> pred %i\n", predcol/sum_weights); #endif /* Now calculate other colours */ ncol = (colour+1)%3; pconv = pconvmap[nbpts[0]][nbpts[1]]; if (pconv == PATCONV_NONE) abort (); othcol = 0; predcol = 0; nsumw = 0; skno = 0; /*DPRINTF(" Col %i: pat %i pconv %i\n", ncol, nbpts[1], pconv);*/ for (nbs = 0; nbs < n_pos[nbpts[1]].num; nbs++) { off offset = n_pos[nbpts[1]].nb_pts[nbs]; const int nx = x + offset.dx; const int ny = y + offset.dy; const signed long addr_off = 3 * (offset.dx + width * offset.dy); int eff_weight = 0; unsigned char thisval; for (i = 0; i < 4; i++) eff_weight += pat_to_pat[pconv].cf[nbs][i] * weights[i]; if (nx >= 0 && nx < width && ny >= 0 && ny < height) { thisval = src_ptr[addr_off+ncol]; nsumw += eff_weight; /*DPRINTF(" (%i,%i): val %i, eff_w %6.4f\n", nx, ny, thisval, (double)(eff_weight>>1)/sum_weights);*/ othcol += thisval * eff_weight; predcol += thisval; } else { skno++; } }; dst_ptr[ncol] = othcol/nsumw; /*DPRINTF( " -> val %i (bilin: %i)\n", dst_ptr[ncol], predcol/(n_pos[nbpts[1]].num-skno));*/ /* Third colour */ ncol = (colour+2)%3; pconv = pconvmap[nbpts[0]][nbpts[2]]; if (pconv == PATCONV_NONE) abort (); othcol = 0; predcol = 0; nsumw = 0; skno = 0; /*DPRINTF(" Col %i: pat %i pconv %i\n", ncol, nbpts[2], pconv);*/ for (nbs = 0; nbs < n_pos[nbpts[2]].num; nbs++) { off offset = n_pos[nbpts[2]].nb_pts[nbs]; const int nx = x + offset.dx; const int ny = y + offset.dy; const signed long addr_off = 3 * (offset.dx + width * offset.dy); int eff_weight = 0; unsigned char thisval; for (i = 0; i < 4; i++) eff_weight += pat_to_pat[pconv].cf[nbs][i] * weights[i]; if (nx >= 0 && nx < width && ny >= 0 && ny < height) { thisval = src_ptr[addr_off+ncol]; nsumw += eff_weight; /*DPRINTF(" (%i,%i): val %i, eff_w %6.4f\n", nx, ny, thisval, (double)(eff_weight>>1)/sum_weights);*/ othcol += thisval * eff_weight; predcol += thisval; } else { skno++; } }; dst_ptr[ncol] = othcol/nsumw; /*DPRINTF( " -> val %i (bilin: %i)\n", dst_ptr[ncol], predcol/(n_pos[nbpts[1]].num-skno));*/ } } }

Generated by Doxygen 1.6.0 Back to index