d2/image_weighted_median.h

   1 // Copyright 2002, 2003, 2004 David Hilvert <dhilvert@auricle.dyndns.org>,
   2 //                                          <dhilvert@ugcs.caltech.edu>
   3
   4 /*  This file is part of the Anti-Lamenessing Engine.
   5
   6     The Anti-Lamenessing Engine is free software; you can redistribute it and/or modify
   7     it under the terms of the GNU General Public License as published by
   8     the Free Software Foundation; either version 3 of the License, or
   9     (at your option) any later version.
  10
  11     The Anti-Lamenessing Engine is distributed in the hope that it will be useful,
  12     but WITHOUT ANY WARRANTY; without even the implied warranty of
  13     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  14     GNU General Public License for more details.
  15
  16     You should have received a copy of the GNU General Public License
  17     along with the Anti-Lamenessing Engine; if not, write to the Free Software
  18     Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
  19 */
  20
  21 /*
  22  * image_weighted_median.h: Image representing a weighted median of inputs.
  23  */
  24
  25 #ifndef __image_weighted_median_h__
  26 #define __image_weighted_median_h__
  27
  28 #include "exposure/exposure.h"
  29 #include "point.h"
  30 #include "image.h"
  31
  32 class image_weighted_median : public image_weighted_avg {
  33 private:
  34
  35         /*
  36          * Array 'colors' stores image colors, sorted by intensity for each
  37          * channel at each pixel location.
  38          *
  39          * Array 'weights' stores the weights associated with each color, where
  40          * the weights are represented cumulatively, so that for weights and
  41          * intensities:
  42          *
  43          * Color:  1  2  3  6  7
  44          * Weight: 2  2  1  1  1
  45          *
  46          * The (cumulative) representation would be:
  47          *
  48          * Color:  1  2  3  6  7
  49          * Weight: 2  4  5  6  7
  50          *
  51          * XXX: This storage approach may have poor cache characteristics.
  52          * It might be better to localize elements having identical spatial
  53          * coordinates.
  54          */
  55         image **colors;
  56         image **weights;
  57         unsigned int capacity;
  58
  59 public:
  60         image_weighted_median (unsigned int dimy, unsigned int dimx, unsigned int
  61                         depth, int capacity = -1, const char *name = "anonymous")
  62                         : image_weighted_avg(dimy, dimx, depth, name) {
  63
  64                 if (capacity == -1) {
  65                         this->capacity = image_rw::count();
  66                 } else if (capacity >= 0) {
  67                         this->capacity = (unsigned int) capacity;
  68                 } else
  69                         assert(0);
  70
  71                 colors  = (image **) malloc(this->capacity * sizeof(image *));
  72                 weights = (image **) malloc(this->capacity * sizeof(image *));
  73
  74                 assert(colors);
  75                 assert(weights);
  76
  77                 if (!colors || !weights) {
  78                         fprintf(stderr, "Could not allocate memory for image data.\n");
  79                         exit(1);
  80                 }
  81
  82                 for (unsigned int f = 0; f < this->capacity; f++) {
  83                         colors[f] = new_image_ale_real(dimy, dimx, depth);
  84                         weights[f] = new_image_ale_real(dimy, dimx, depth);
  85
  86                         assert(colors[f]);
  87                         assert(weights[f]);
  88
  89                         if (!colors[f] || !weights[f]) {
  90                                 fprintf(stderr, "Could not allocate memory for image data.\n");
  91                                 exit(1);
  92                         }
  93                 }
  94         }
  95
  96         virtual ~image_weighted_median() {
  97                 for (unsigned int f = 0; f < capacity; f++) {
  98                         delete colors[f];
  99                         delete weights[f];
 100                 }
 101
 102                 free(colors);
 103                 free(weights);
 104         }
 105
 106         /*
 107          * Extend the image area to the top, bottom, left, and right,
 108          * initializing the new image areas with black pixels.  Negative values
 109          * shrink the image.
 110          */
 111         image *_extend(int top, int bottom, int left, int right) {
 112
 113                 for (unsigned int f = 0; f < capacity; f++) {
 114                         extend(&colors[f], top, bottom, left, right);
 115                         extend(&weights[f], top, bottom, left, right);
 116                 }
 117
 118                 _dimx = colors[0]->width();
 119                 _dimy = colors[0]->height();
 120                 _offset = colors[0]->offset();
 121
 122                 return NULL;
 123         }
 124
 125         int accumulate_norender(int i, int j) {
 126                 return 0;
 127         }
 128
 129         /*
 130          * Perform insertion sort on the arrays, where sort is by color.
 131          *
 132          * XXX: This does a poor job of handling multiple contributions from
 133          * the same frame, especially when the number of frames is 1.
 134          */
 135         void accumulate(int i, int j, int f, pixel new_value, pixel new_weight) {
 136                 for (unsigned int k = 0; k < 3; k++) {
 137
 138                         if (fabs(new_weight[k]) > ale_real_weight_floor)
 139                                 new_value[k] /= new_weight[k];
 140                         else
 141                                 continue;
 142
 143                         /*
 144                          * XXX: This initialization should not be necessary.
 145                          */
 146                         if (f == 0)
 147                         for (unsigned int ff = 0; ff < capacity; ff++)
 148                                 weights[ff]->set_chan(i, j, k, 0);
 149
 150                         assert (finite(new_weight[k]));
 151
 152                         if (new_weight[k] <= 0)
 153                                 continue;
 154
 155                         for (unsigned int ff = 0; ff < capacity; ff++) {
 156                                 assert (ff <= (unsigned int) f);
 157                                 if (ff == capacity - 1) {
 158                                         colors[ff]->set_chan(i, j, k, new_value[k]);
 159                                         weights[ff]->set_chan(i, j, k,
 160                                                 weights[ff]->get_chan(i, j, k) + new_weight[k]);
 161                                         break;
 162                                 }
 163                                 if ((ff == 0 && weights[ff]->get_chan(i, j, k) == 0)
 164                                  || (ff >  0 && weights[ff]->get_chan(i, j, k) == weights[ff - 1]->get_chan(i, j, k))) {
 165                                         colors[ff]->set_chan(i, j, k, new_value[k]);
 166                                         for (unsigned int fff = ff; fff < capacity; fff++)
 167                                                 weights[fff]->set_chan(i, j, k,
 168                                                         weights[fff]->get_chan(i, j, k) + new_weight[k]);
 169                                         break;
 170                                 }
 171                                 if (colors[ff]->get_chan(i, j, k) == (ale_sreal) new_value[k]) {
 172                                         for (unsigned int fff = ff; fff < capacity; fff++)
 173                                                 weights[fff]->set_chan(i, j, k,
 174                                                         weights[fff]->get_chan(i, j, k) + new_weight[k]);
 175                                         break;
 176                                 }
 177                                 if (colors[ff]->get_chan(i, j, k) > (ale_sreal) new_value[k]) {
 178                                         for (unsigned int fff = capacity - 1; fff > ff; fff--) {
 179                                                 weights[fff]->set_chan(i, j, k,  weights[fff - 1]->get_pixel(i, j)[k] + new_weight[k]);
 180                                                 colors[fff]->set_chan(i, j, k, colors[fff - 1]->get_pixel(i, j)[k]);
 181                                         }
 182                                         colors[ff]->set_chan(i, j, k, new_value[k]);
 183                                         weights[ff]->set_chan(i, j, k, new_weight[k]);
 184                                         if (ff > 0)
 185                                                 weights[ff]->set_chan(i, j, k,
 186                                                         weights[ff]->get_chan(i, j, k)
 187                                                       + weights[ff - 1]->get_chan(i, j, k));
 188
 189                                         break;
 190                                 }
 191                         }
 192                 }
 193         }
 194
 195         /*
 196          * XXX: This is inefficient in cases where only one channel is desired.
 197          */
 198         spixel get_pixel(unsigned int y, unsigned int x) const {
 199                 pixel result;
 200
 201                 for (int k = 0; k < 3; k++) {
 202                         ale_real midpoint = weights[capacity - 1]->get_chan(y, x, k) / 2;
 203
 204                         if (midpoint == 0)
 205                                 return pixel::zero();
 206
 207                         /*
 208                          * Binary search.
 209                          */
 210
 211                         unsigned int l = 0;
 212                         unsigned int h = capacity - 1;
 213                         unsigned int m = h / 2;
 214
 215                         while (h > l + 1) {
 216                                 if ((ale_real) weights[m]->get_chan(y, x, k) < midpoint)
 217                                         l = m;
 218                                 else
 219                                         h = m;
 220                                 m = (h + l) / 2;
 221                         }
 222
 223                         if ((ale_real) weights[l]->get_chan(y, x, k) < midpoint)
 224                                 l = h;
 225                         if ((ale_real) weights[l]->get_chan(y, x, k) > midpoint)
 226                                 result[k] = colors[l]->get_chan(y, x, k);
 227                         else if ((ale_real) weights[l]->get_chan(y, x, k) == midpoint)
 228                                 result[k] = (colors[l]->get_chan(y, x, k)
 229                                            + colors[l + 1]->get_chan(y, x, k)) / 2;
 230                         else
 231                                 assert(0);
 232
 233                 }
 234
 235                 return result;
 236         }
 237
 238         image *get_weights() {
 239                 return weights[capacity - 1];
 240         }
 241
 242         image *get_colors() {
 243                 return this;
 244         }
 245
 246 };
 247
 248 #endif