rich_image.cc

   1
   2 ///////////////////////////////////////////////////////////////////////////
   3 // This program is free software: you can redistribute it and/or modify  //
   4 // it under the terms of the version 3 of the GNU General Public License //
   5 // as published by the Free Software Foundation.                         //
   6 //                                                                       //
   7 // This program is distributed in the hope that it will be useful, but   //
   8 // WITHOUT ANY WARRANTY; without even the implied warranty of            //
   9 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU      //
  10 // General Public License for more details.                              //
  11 //                                                                       //
  12 // You should have received a copy of the GNU General Public License     //
  13 // along with this program. If not, see <http://www.gnu.org/licenses/>.  //
  14 //                                                                       //
  15 // Written by Francois Fleuret                                           //
  16 // (C) Idiap Research Institute                                          //
  17 //                                                                       //
  18 // Contact <francois.fleuret@idiap.ch> for comments & bug reports        //
  19 ///////////////////////////////////////////////////////////////////////////
  20
  21 #include <string.h>
  22
  23 #include "rich_image.h"
  24
  25 #define PIXEL_DELTA(a, b) ((a) >= (b) ? (a) - (b) : (b) - (a))
  26
  27 static inline bool edge(                   unsigned char v0, unsigned char v1,
  28                          unsigned char v2, unsigned char v3, unsigned char v4, unsigned char v5,
  29                                            unsigned char v6, unsigned char v7) {
  30   unsigned char g = PIXEL_DELTA(v3, v4);
  31
  32   return
  33     g > 8 &&
  34     g > PIXEL_DELTA(v0, v3) &&
  35     g > PIXEL_DELTA(v1, v4) &&
  36     g > PIXEL_DELTA(v2, v3) &&
  37     g > PIXEL_DELTA(v4, v5) &&
  38     g > PIXEL_DELTA(v3, v6) &&
  39     g > PIXEL_DELTA(v4, v7);
  40 }
  41
  42 void RichImage::free() {
  43   delete[] _edge_map;
  44   _edge_map = 0;
  45   delete[] _line_edge_map;
  46   _line_edge_map = 0;
  47   delete[] _tag_edge_map;
  48   _tag_edge_map = 0;
  49   delete[] _scale_edge_map;
  50   _scale_edge_map = 0;
  51   delete[] _width_at_scale;
  52   _width_at_scale = 0;
  53   delete[] _height_at_scale;
  54   _height_at_scale = 0;
  55 }
  56
  57 void RichImage::compute_one_scale_edge_maps(int width, int height,
  58                                             unsigned int ***scale_edge_map,
  59                                             unsigned char *pixel_map,
  60                                             unsigned int *sum_pixel_map,
  61                                             unsigned int *sum_sq_pixel_map) {
  62
  63   unsigned char d00, d01, d02, d03, d04, d05, d06, d07;
  64   unsigned char d08, d09, d10, d11, d12, d13, d14, d15;
  65   unsigned char *local_pixel_map;
  66   unsigned int *local_sum_pixel_map, *local_sum_sq_pixel_map;
  67
  68   // Compute the integral images
  69
  70   {
  71     unsigned int *sp = sum_pixel_map, *ssp = sum_sq_pixel_map;
  72     unsigned char *p = pixel_map;
  73     for(int x = 0; x < width; x++) {
  74       *sp++ = 0;
  75       *ssp++ = 0;
  76       p++;
  77     }
  78
  79     for(int y = 1; y < height; y++) {
  80       *sp++ = 0;
  81       *ssp++ = 0;
  82       p++;
  83       for(int x = 1; x < width; x++) {
  84         *sp = *(sp - width) + *(sp - 1) - *(sp - width - 1) + ((unsigned int) (*p));
  85         *ssp = *(ssp - width) + *(ssp - 1) - *(ssp - width - 1) + sq((unsigned int) (*p));
  86         sp++;
  87         ssp++;
  88         p++;
  89       }
  90     }
  91   }
  92
  93   const int var_square_size = 16;
  94
  95   int k00 = - 2 + width * (- 2);
  96   int k01 = - 1 + width * (- 2);
  97   int k02 = + 0 + width * (- 2);
  98   int k03 = + 1 + width * (- 2);
  99   int k04 = - 2 + width * (- 1);
 100   int k05 = - 1 + width * (- 1);
 101   int k06 = + 0 + width * (- 1);
 102   int k07 = + 1 + width * (- 1);
 103   int k08 = - 2 + width * (+ 0);
 104   int k09 = - 1 + width * (+ 0);
 105   int k10 = + 0 + width * (+ 0);
 106   int k11 = + 1 + width * (+ 0);
 107   int k12 = - 2 + width * (+ 1);
 108   int k13 = - 1 + width * (+ 1);
 109   int k14 = + 0 + width * (+ 1);
 110   int k15 = + 1 + width * (+ 1);
 111
 112   unsigned int *edge_map0 = scale_edge_map[0][0];
 113   unsigned int *edge_map1 = scale_edge_map[1][0];
 114   unsigned int *edge_map2 = scale_edge_map[2][0];
 115   unsigned int *edge_map3 = scale_edge_map[3][0];
 116   unsigned int *edge_map4 = scale_edge_map[4][0];
 117   unsigned int *edge_map5 = scale_edge_map[5][0];
 118   unsigned int *edge_map6 = scale_edge_map[6][0];
 119   unsigned int *edge_map7 = scale_edge_map[7][0];
 120   unsigned int *variance_map = scale_edge_map[8][0];
 121
 122   int a = -1 - width, b = - width, c = -1, d = 0;
 123
 124   local_pixel_map = pixel_map;
 125   local_sum_pixel_map = sum_pixel_map;
 126   local_sum_sq_pixel_map = sum_sq_pixel_map;
 127
 128   const int gray_bin_width = 256 / nb_gray_tags;
 129
 130   for(int y = 0; y < height; y++) {
 131
 132     for(int x = 0; x < width; x++) {
 133
 134       if(x == 0 || y == 0) {
 135         for(int e = 0; e < _nb_tags; e++)
 136           scale_edge_map[e][0][d] = 0;
 137       } else {
 138         for(int e = 0; e < _nb_tags; e++)
 139           scale_edge_map[e][0][d] =
 140             scale_edge_map[e][0][b] +
 141             scale_edge_map[e][0][c] -
 142             scale_edge_map[e][0][a];
 143       }
 144
 145       scale_edge_map[first_gray_tag +
 146                      (local_pixel_map[0] / gray_bin_width)][0][d]++;
 147
 148       if(x - int(var_square_size/2) >= 0 &&
 149          x + int(var_square_size/2) < width &&
 150          y - int(var_square_size/2) >= 0 &&
 151          y + int(var_square_size/2) < height) {
 152
 153         int s =
 154           + local_sum_pixel_map[ - var_square_size/2 + width * ( - var_square_size / 2)]
 155           + local_sum_pixel_map[ + var_square_size/2 + width * ( + var_square_size / 2)]
 156           - local_sum_pixel_map[ - var_square_size/2 + width * ( + var_square_size / 2)]
 157           - local_sum_pixel_map[ + var_square_size/2 + width * ( - var_square_size / 2)];
 158
 159         int s_sq =
 160           + local_sum_sq_pixel_map[ - var_square_size/2 + width * ( - var_square_size / 2)]
 161           + local_sum_sq_pixel_map[ + var_square_size/2 + width * ( + var_square_size / 2)]
 162           - local_sum_sq_pixel_map[ - var_square_size/2 + width * ( + var_square_size / 2)]
 163           - local_sum_sq_pixel_map[ + var_square_size/2 + width * ( - var_square_size / 2)];
 164
 165         if(sq(var_square_size) * s_sq - sq(s) >=
 166            100 * sq(var_square_size) * (sq(var_square_size) - 1)) {
 167
 168           d00 = local_pixel_map[k00];
 169           d01 = local_pixel_map[k01];
 170           d02 = local_pixel_map[k02];
 171           d03 = local_pixel_map[k03];
 172
 173           d04 = local_pixel_map[k04];
 174           d05 = local_pixel_map[k05];
 175           d06 = local_pixel_map[k06];
 176           d07 = local_pixel_map[k07];
 177
 178           d08 = local_pixel_map[k08];
 179           d09 = local_pixel_map[k09];
 180           d10 = local_pixel_map[k10];
 181           d11 = local_pixel_map[k11];
 182
 183           d12 = local_pixel_map[k12];
 184           d13 = local_pixel_map[k13];
 185           d14 = local_pixel_map[k14];
 186           d15 = local_pixel_map[k15];
 187
 188           /*
 189
 190               #0         #1         #2         #3
 191
 192             XXXXXX     .XXXXX     ...XXX     .....X
 193             XXXXXX     ..XXXX     ...XXX     ....XX
 194             ......     ...XXX     ...XXX     ...XXX
 195             ......     ....XX     ...XXX     ..XXXX
 196
 197               #4         #5         #6         #7
 198
 199             ......     X.....     XXX...     XXXXX.
 200             ......     XX....     XXX...     XXXX..
 201             XXXXXX     XXX...     XXX...     XXX...
 202             XXXXXX     XXXX..     XXX...     XX....
 203
 204           */
 205
 206           if(edge(d04, d08, d01, d05, d09, d13, d06, d10)) {
 207             if(d05 < d09) edge_map0[d]++;
 208             else          edge_map4[d]++;
 209           }
 210
 211           if(edge(d02, d07, d00, d05, d10, d15, d08, d13)) {
 212             if(d05 < d10) edge_map7[d]++;
 213             else          edge_map3[d]++;
 214           }
 215
 216           if(edge(d01, d02, d04, d05, d06, d07, d09, d10)) {
 217             if(d05 < d06) edge_map6[d]++;
 218             else          edge_map2[d]++;
 219           }
 220
 221           if(edge(d01, d04, d03, d06, d09, d12, d11, d14)) {
 222             if(d06 < d09) edge_map1[d]++;
 223             else          edge_map5[d]++;
 224           }
 225
 226           variance_map[d]++;
 227         }
 228       }
 229
 230       a++; b++; c++; d++;
 231       local_pixel_map++;
 232       local_sum_pixel_map++;
 233       local_sum_sq_pixel_map++;
 234     }
 235   }
 236 }
 237
 238 void RichImage::compute_rich_structure() {
 239   free();
 240
 241   unsigned char *pixel_maps;
 242   unsigned char **line_pixel_maps;
 243   unsigned char ***scale_pixel_maps;
 244
 245   _nb_scales =
 246     int(global.nb_scales_per_power_of_two
 247         * log(scalar_t(min(_width, _height))/scalar_t(_image_size_min)) / log(2)) + 1;
 248
 249   _width_at_scale = new int[_nb_scales];
 250   _height_at_scale = new int[_nb_scales];
 251   int total_surface = 0, total_lines = 0;
 252
 253   // Compute the sizes of the image at various scales
 254
 255   scalar_t rho = exp(log(0.5) / scalar_t(global.nb_scales_per_power_of_two));
 256   scalar_t factor = 1.0;
 257
 258   for(int s = 0; s < _nb_scales; s++) {
 259     _width_at_scale[s] = int(scalar_t(_width) * factor);
 260     _height_at_scale[s] = int(scalar_t(_height) * factor);
 261     total_surface += _width_at_scale[s] * _height_at_scale[s];
 262     total_lines += _height_at_scale[s];
 263     factor *= rho;
 264   }
 265
 266   // Allocate the memory for the various scales and set up the pointer
 267   // arrays to speed-up accesses
 268
 269   pixel_maps = new unsigned char[total_surface];
 270   memset(pixel_maps, 0, sizeof(unsigned char) * total_surface);
 271   line_pixel_maps = new unsigned char *[total_lines];
 272   scale_pixel_maps = new unsigned char **[_nb_scales];
 273
 274   int sum_surfaces, sum_heights;
 275
 276   sum_surfaces = 0;
 277   sum_heights = 0;
 278
 279   for(int s = 0; s < _nb_scales; s++) {
 280     scale_pixel_maps[s] = line_pixel_maps + sum_heights;
 281     for(int y = 0; y < _height_at_scale[s]; y++) {
 282       scale_pixel_maps[s][y] = pixel_maps + sum_surfaces;
 283       sum_surfaces += _width_at_scale[s];
 284     }
 285     sum_heights += _height_at_scale[s];
 286   }
 287
 288   _edge_map = new unsigned int[total_surface * _nb_tags];
 289   memset(_edge_map, 0, sizeof(unsigned int) * total_surface * _nb_tags);
 290   _line_edge_map = new unsigned int *[total_lines * _nb_tags];
 291   _tag_edge_map = new unsigned int **[_nb_tags * _nb_scales];
 292   _scale_edge_map = new unsigned int ***[_nb_scales];
 293
 294   sum_surfaces = 0;
 295   sum_heights = 0;
 296
 297   for(int s = 0; s < _nb_scales; s++) {
 298     _scale_edge_map[s] = _tag_edge_map + s * _nb_tags;
 299     for(int t = 0; t < _nb_tags; t++) {
 300       _scale_edge_map[s][t] = _line_edge_map + sum_heights;
 301       for(int y = 0; y < _height_at_scale[s]; y++) {
 302         _scale_edge_map[s][t][y] = _edge_map + sum_surfaces;
 303         sum_surfaces += _width_at_scale[s];
 304       }
 305       sum_heights += _height_at_scale[s];
 306     }
 307   }
 308
 309   // Compute the scaled images per se
 310
 311   for(int s = 0; s < _nb_scales; s++) {
 312
 313     if(s < global.nb_scales_per_power_of_two) {
 314
 315       if(s == 0) {
 316
 317         // The first image is not rescaled
 318
 319         for(int y = 0; y < _height_at_scale[s]; y++)
 320           for(int x = 0; x < _width_at_scale[s]; x++)
 321             scale_pixel_maps[s][y][x] = _content[x + _width * y];
 322
 323       } else {
 324
 325         // The nb_scales_per_power_of_two - 1 first images are
 326         // generated with a 'complex' scaling, and then we just
 327         // downscale by a factor of two
 328
 329         const int ld = 3;
 330
 331         int delta_column[_width_at_scale[s] << ld];
 332         for(int xx = 0; xx < _width_at_scale[s] << ld; xx++)
 333           delta_column[xx] = (xx * _width_at_scale[0]) / (_width_at_scale[s] << ld);
 334
 335         unsigned char *line[_height_at_scale[s] << ld];
 336         for(int yy = 0; yy < _height_at_scale[s] << ld; yy++)
 337           line[yy] = scale_pixel_maps[0][(yy * _height_at_scale[0]) / (_height_at_scale[s] << ld)];
 338
 339         for(int y = 0; y < _height_at_scale[s]; y++) {
 340           unsigned char *dest_line = scale_pixel_maps[s][y];
 341           for(int x = 0; x < _width_at_scale[s]; x++) {
 342             int u = 0;
 343             for(int yy = (y << ld); yy < ((y+1) << ld); yy++)
 344               for(int xx = (x << ld); xx < ((x+1) << ld); xx++)
 345                 u += line[yy][delta_column[xx]];
 346             dest_line[x] = (u >> (2 * ld));
 347           }
 348         }
 349       }
 350     } else {
 351
 352       // Other scales are computed by halfing one of the already
 353       // computed scales
 354
 355       int r = s - global.nb_scales_per_power_of_two;
 356       for(int y = 0; y < _height_at_scale[s]; y++)
 357         for(int x = 0; x < _width_at_scale[s]; x++)
 358           scale_pixel_maps[s][y][x] =
 359             (int(scale_pixel_maps[r][y*2 + 0][x*2 + 0]) +
 360              int(scale_pixel_maps[r][y*2 + 0][x*2 + 1]) +
 361              int(scale_pixel_maps[r][y*2 + 1][x*2 + 0]) +
 362              int(scale_pixel_maps[r][y*2 + 1][x*2 + 1])) >> 2;
 363     }
 364
 365   }
 366
 367   // Allocate the memory for the edge maps at various scales and set
 368   // up the pointer arrays to speed-up accesses
 369
 370   unsigned int *sum_pixel_map = new unsigned int[_width * _height];
 371   unsigned int *sum_sq_pixel_map = new unsigned int[_width * _height];
 372
 373   for(int s = 0; s < _nb_scales; s++)
 374     compute_one_scale_edge_maps(_width_at_scale[s],
 375                                 _height_at_scale[s],
 376                                 _scale_edge_map[s],
 377                                 scale_pixel_maps[s][0],
 378                                 sum_pixel_map,
 379                                 sum_sq_pixel_map);
 380
 381   delete[] sum_pixel_map;
 382   delete[] sum_sq_pixel_map;
 383
 384   delete[] pixel_maps;
 385   delete[] line_pixel_maps;
 386   delete[] scale_pixel_maps;
 387 }
 388
 389 RichImage::RichImage() : Image() {
 390   _width_at_scale = 0;
 391   _height_at_scale = 0;
 392   _edge_map = 0;
 393   _line_edge_map = 0;
 394   _tag_edge_map = 0;
 395   _scale_edge_map = 0;
 396 }
 397
 398 RichImage::RichImage(int width, int height) : Image(width, height) {
 399   _width_at_scale = 0;
 400   _height_at_scale = 0;
 401   _edge_map = 0;
 402   _line_edge_map = 0;
 403   _tag_edge_map = 0;
 404   _scale_edge_map = 0;
 405 }
 406
 407 RichImage::~RichImage() {
 408   delete[] _edge_map;
 409   delete[] _line_edge_map;
 410   delete[] _tag_edge_map;
 411   delete[] _scale_edge_map;
 412   delete[] _width_at_scale;
 413   delete[] _height_at_scale;
 414 }
 415
 416 void RichImage::write(ostream *out) {
 417   Image::write(out);
 418 }
 419
 420 void RichImage::read(istream *in) {
 421   Image::read(in);
 422 }
 423
 424 void RichImage::write_tag_png(const char *filename) {
 425   const int nb_cols= 4;
 426   const int nb_rows = ((nb_tags() + nb_cols - 1) / nb_cols);
 427
 428   int height_total = 0;
 429   for(int s = 0; s < _nb_scales; s++) height_total += _height_at_scale[s];
 430
 431   RGBImage image(_width * nb_cols, height_total * nb_rows);
 432
 433   for(int y = 0; y < image.height(); y++) for(int x = 0; x < image.width(); x++)
 434                                             image.set_pixel(x, y, 0, 255, 0);
 435
 436   int sum_height = 0;
 437
 438   for(int s = 0; s < _nb_scales; s++) {
 439     for(int t = 0; t < nb_tags(); t++) {
 440       int x0 = (t%nb_cols) * _width, y0 = sum_height + (t/nb_cols) * _height_at_scale[s];
 441       for(int y = 0; y < _height_at_scale[s]; y++) for(int x = 0; x < _width_at_scale[s]; x++) {
 442           int c = 255 - min(255, max(0, int(255 * nb_tags_in_window(s, t, x, y, x+1, y+1))));
 443           image.set_pixel(x0 + x, y0 + y, c, c, c);
 444         }
 445     }
 446     sum_height += nb_rows * _height_at_scale[s];
 447   }
 448
 449   image.write_png(filename);
 450 }