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