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