rich_image.cc

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