d2/trans_abstract.h

   1 // Copyright 2002, 2004, 2007 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  * trans_abstract.h: Abstract transformation superclass.
  23  */
  24
  25 #ifndef __trans_abstract_h__
  26 #define __trans_abstract_h__
  27
  28 #include "image.h"
  29 #include "point.h"
  30 #include "pixel.h"
  31
  32 #ifndef M_PI
  33 #define M_PI 3.14159265358979323846
  34 #endif
  35
  36 /*
  37  * Number of coefficients used in correcting barrel distortion.
  38  */
  39
  40 #define BARREL_DEGREE 5
  41
  42 /*
  43  * Acceptable error for inverse barrel distortion, measured in scaled output
  44  * pixels.
  45  */
  46
  47 #define BARREL_INV_ERROR 0.01
  48
  49 struct trans_abstract {
  50 private:
  51         ale_pos bdc[BARREL_DEGREE];                             // barrel-dist. coeffs.
  52         unsigned int bdcnum;                                    // number of bdcs
  53
  54 protected:
  55         ale_pos scale_factor;
  56         unsigned int input_height, input_width;
  57
  58         virtual void specific_rescale(ale_pos factor) = 0;
  59         virtual void reset_memos() = 0;
  60         virtual void specific_set_dimensions(const image *im) = 0;
  61
  62 public:
  63
  64         struct elem_bounds_int_t {
  65                 unsigned int imin, imax, jmin, jmax;
  66
  67                 int satisfies_min_dim(unsigned int min_dimension) {
  68                         if (imax - imin < min_dimension
  69                          || jmax - jmin < min_dimension)
  70                                 return 0;
  71
  72                         return 1;
  73                 }
  74         };
  75
  76         trans_abstract() {
  77                 bdcnum = 0;
  78         }
  79
  80         trans_abstract &operator=(const trans_abstract &ta) {
  81                 scale_factor = ta.scale_factor;
  82                 input_width = ta.input_width;
  83                 input_height = ta.input_height;
  84
  85                 bdcnum = ta.bdcnum;
  86
  87                 assert (bdcnum < BARREL_DEGREE);
  88
  89                 for (unsigned int d = 0; d < bdcnum; d++)
  90                         bdc[d] = ta.bdc[d];
  91
  92                 return *this;
  93         }
  94
  95         trans_abstract (const trans_abstract &ta) {
  96                 operator=(ta);
  97         }
  98
  99         /*
 100          * Returns non-zero if the transformation might be non-Euclidean.
 101          */
 102         virtual int is_projective() const = 0;
 103
 104         /*
 105          * Get scale factor.
 106          */
 107
 108         ale_pos scale() const {
 109                 return scale_factor;
 110         }
 111
 112         /*
 113          * Get width of input image.
 114          */
 115         ale_pos scaled_width() const {
 116                 return (input_width * scale_factor);
 117         }
 118
 119         /*
 120          * Get unscaled width of input image.
 121          */
 122         unsigned int unscaled_width() const {
 123                 return (unsigned int) input_width;
 124         }
 125
 126         /*
 127          * Get height of input image;
 128          */
 129         ale_pos scaled_height() const {
 130                return (input_height * scale_factor);
 131         }
 132
 133         /*
 134          * Get unscaled height of input image.
 135          */
 136         unsigned int unscaled_height() const {
 137                 return (unsigned int) input_height;
 138         }
 139
 140         /*
 141          * Barrel distortion radial component.
 142          */
 143         ale_pos bdr(ale_pos r) const {
 144                 assert (bdcnum < BARREL_DEGREE);
 145                 ale_pos s = r;
 146                 for (unsigned int d = 0; d < bdcnum; d++)
 147                         s += bdc[d] * (pow(r, d + 2) - r);
 148                 return s;
 149         }
 150
 151         /*
 152          * Derivative of the barrel distortion radial component.
 153          */
 154         ale_pos bdrd(ale_pos r) const {
 155                 assert (bdcnum < BARREL_DEGREE);
 156                 ale_pos s = 1;
 157                 for (unsigned int d = 0; d < bdcnum; d++)
 158                         s += bdc[d] * (pow(r, d + 1) - 1);
 159                 return s;
 160         }
 161
 162         /*
 163          * Barrel distortion.
 164          */
 165         struct point bd(struct point p) const {
 166                 if (bdcnum > 0) {
 167                         point half_diag = point(unscaled_height(), unscaled_width()) / 2;
 168
 169                         p -= half_diag;
 170
 171                         ale_pos r = p.norm() / half_diag.norm();
 172
 173                         if (r > 0.00001)
 174                                 p *= bdr(r)/r;
 175
 176                         p += half_diag;
 177                 }
 178
 179                 return p;
 180         }
 181
 182         /*
 183          * Barrel distortion inverse.
 184          */
 185         struct point bdi(struct point p) const {
 186                 if (bdcnum > 0) {
 187                         point half_diag = point(unscaled_height(), unscaled_width()) / 2;
 188
 189                         p -= half_diag;
 190
 191                         ale_pos r = p.norm() / half_diag.norm();
 192                         ale_pos s = r;
 193
 194                         while (fabs(r - bdr(s)) * half_diag.norm() > BARREL_INV_ERROR)
 195                                 s += (r - bdr(s)) / bdrd(s);
 196
 197                         if (r > 0.0001)
 198                                 p *= s / r;
 199
 200                         p += half_diag;
 201                 }
 202
 203                 assert (!isnan(p[0]) && !isnan(p[1]));
 204
 205                 return p;
 206         }
 207
 208         /*
 209          * Transformation sans barrel distortion
 210          */
 211         virtual struct point pe(struct point p) const = 0;
 212
 213         /*
 214          * Transformation inverse sans barrel distortion
 215          */
 216         virtual struct point pei(struct point p) const = 0;
 217
 218         /*
 219          * Map unscaled point p.
 220          */
 221         struct point transform_unscaled(struct point p) const {
 222                 return pe(bdi(p));
 223         }
 224
 225         /*
 226          * Transform point p.
 227          *
 228          * Barrel distortion correction followed by a projective/euclidean
 229          * transformation.
 230          */
 231         struct point transform_scaled(struct point p) const {
 232                 return transform_unscaled(p / scale_factor);
 233         }
 234
 235 #if 0
 236         /*
 237          * operator() is the transformation operator.
 238          */
 239         struct point operator()(struct point p) {
 240                 return transform(p);
 241         }
 242 #endif
 243
 244         /*
 245          * Map point p using the inverse of the transform into
 246          * the unscaled image space.
 247          */
 248         struct point unscaled_inverse_transform(struct point p) const {
 249                 return bd(pei(p));
 250         }
 251
 252         /*
 253          * Map point p using the inverse of the transform.
 254          *
 255          * Projective/euclidean inverse followed by barrel distortion.
 256          */
 257         struct point scaled_inverse_transform(struct point p) const {
 258                 assert (p.defined());
 259                 point q = unscaled_inverse_transform(p);
 260
 261                 q[0] *= scale_factor;
 262                 q[1] *= scale_factor;
 263
 264                 return q;
 265         }
 266
 267         /*
 268          * Calculate projective transformation parameters from a euclidean
 269          * transformation.
 270          */
 271         virtual void eu_to_gpt() = 0;
 272
 273         /*
 274          * Set the tonal multiplier.
 275          */
 276         virtual void set_tonal_multiplier(pixel p) = 0;
 277
 278         /*
 279          * Get the tonal multiplier.
 280          */
 281         virtual pixel get_tonal_multiplier(struct point p) const = 0;
 282         virtual pixel get_inverse_tonal_multiplier(struct point p) const = 0;
 283
 284         /*
 285          * Modify a euclidean transform in the indicated manner.
 286          */
 287         virtual void eu_modify(int i1, ale_pos diff) = 0;
 288
 289         /*
 290          * Rotate about a given point in the original reference frame.
 291          */
 292         virtual void eu_rotate_about_scaled(point center, ale_pos diff) = 0;
 293
 294         /*
 295          * Modify all euclidean parameters at once.
 296          */
 297         virtual void eu_set(ale_pos eu[3]) = 0;
 298
 299         /*
 300          * Get the specified euclidean parameter
 301          */
 302         virtual ale_pos eu_get(int param) const = 0;
 303
 304         /*
 305          * Modify a projective transform in the indicated manner.
 306          */
 307         virtual void gpt_modify(int i1, int i2, ale_pos diff) = 0;
 308
 309         /*
 310          * Modify a projective transform according to the group operation.
 311          */
 312         virtual void gr_modify(int i1, int i2, ale_pos diff) = 0;
 313
 314         /*
 315          * Modify all projective parameters at once.
 316          */
 317         virtual void gpt_set(point x[4]) = 0;
 318
 319         virtual void gpt_set(point x1, point x2, point x3, point x4) = 0;
 320
 321         /*
 322          * Snap positional parameters to the specified resolution.
 323          */
 324
 325         virtual void snap(ale_pos interval) = 0;
 326
 327         /*
 328          * Get the specified projective parameter
 329          */
 330         virtual point gpt_get(int point) const = 0;
 331
 332         /*
 333          * Get the specified projective parameter
 334          */
 335         virtual ale_pos gpt_get(int point, int dim) = 0;
 336
 337         /*
 338          * Check equality of transformation parameters.
 339          */
 340         virtual int operator==(const trans_abstract &t) const {
 341                 /*
 342                  * Small tolerances (< 10^-6?) can cause odd errors,
 343                  * possibly due to float<->double conversion issues.
 344                  */
 345                 double zero_tolerance = 0.01;
 346
 347                 if (scale() != t.scale())
 348                         return 0;
 349
 350                 if (is_projective() != t.is_projective())
 351                         return 0;
 352
 353                 if (is_projective()) {
 354                         assert (t.is_projective());
 355                         for (int i = 0; i < 4; i++)
 356                         for (int d = 0; d < 2; d++) {
 357                                 double abs_difference = fabs(gpt_get(i)[d] - t.gpt_get(i)[d]);
 358
 359                                 if (abs_difference > zero_tolerance)
 360                                         return 0;
 361                         }
 362                 } else {
 363                         assert (!t.is_projective());
 364                         for (int i = 0; i < 3; i++) {
 365                                 double abs_difference = fabs(eu_get(i) - t.eu_get(i));
 366
 367                                 if (abs_difference > zero_tolerance)
 368                                         return 0;
 369                         }
 370                 }
 371
 372                 return 1;
 373         }
 374
 375         virtual int operator!=(const trans_abstract &t) const {
 376                 return !(operator==(t));
 377         }
 378
 379
 380         /*
 381          * Translate by a given amount
 382          */
 383         virtual void translate(point p) = 0;
 384
 385         /*
 386          * Rotate by a given amount about a given point.
 387          */
 388         virtual void rotate(point p, ale_pos degrees) = 0;
 389
 390         /*
 391          * Set the specified barrel distortion parameter.
 392          */
 393         void bd_set(unsigned int degree, ale_pos value) {
 394                 assert (degree < bdcnum);
 395                 bdc[degree] = value;
 396         }
 397
 398         /*
 399          * Set all barrel distortion parameters.
 400          */
 401         void bd_set(unsigned int degree, ale_pos values[BARREL_DEGREE]) {
 402                 assert (degree <= BARREL_DEGREE);
 403                 bdcnum = degree;
 404                 for (unsigned int d = 0; d < degree; d++)
 405                         bdc[d] = values[d];
 406         }
 407
 408         /*
 409          * Get all barrel distortion parameters.
 410          */
 411         void bd_get(ale_pos result[BARREL_DEGREE]) {
 412                 for (unsigned int d = 0; d < bdcnum; d++)
 413                         result[d] = bdc[d];
 414         }
 415
 416         /*
 417          * Get the specified barrel distortion parameter.
 418          */
 419         ale_pos bd_get(unsigned int degree) {
 420                 assert (degree < bdcnum);
 421                 return bdc[degree];
 422         }
 423
 424         /*
 425          * Get the number of barrel distortion parameters.
 426          */
 427         unsigned int bd_count() {
 428                 return bdcnum;
 429         }
 430
 431         /*
 432          * Get the maximum allowable number of barrel distortion parameters.
 433          */
 434         unsigned int bd_max() {
 435                 return BARREL_DEGREE;
 436         }
 437
 438         /*
 439          * Modify the specified barrel distortion parameter.
 440          */
 441         void bd_modify(unsigned int degree, ale_pos diff) {
 442                 assert (degree < bdcnum);
 443                 bd_set(degree, bd_get(degree) + diff);
 444         }
 445
 446         /*
 447          * Rescale a transform with a given factor.
 448          */
 449         void rescale(ale_pos factor) {
 450                 specific_rescale(factor);
 451                 scale_factor *= factor;
 452         }
 453
 454         /*
 455          * Set a new domain.
 456          */
 457
 458         void set_domain(unsigned int new_height, unsigned int new_width) {
 459                 reset_memos();
 460                 input_width = new_width;
 461                 input_height = new_height;
 462         }
 463
 464         /*
 465          * Set the dimensions of the image.
 466          */
 467         void set_dimensions(const image *im) {
 468
 469                 int new_height = (int) im->height();
 470                 int new_width  = (int) im->width();
 471
 472                 reset_memos();
 473                 specific_set_dimensions(im);
 474                 input_height = new_height;
 475                 input_width  = new_width;
 476         }
 477
 478         /*
 479          * Get the position and dimensions of a pixel P mapped from one
 480          * coordinate system to another, using the forward transformation.
 481          * This function uses scaled input coordinates.
 482          */
 483         virtual void map_area(point p, point *q, ale_pos d[2]) {
 484
 485                 /*
 486                  * Determine the coordinates in the target frame for the source
 487                  * image pixel P and two adjacent source pixels.
 488                  */
 489
 490                     (*q) = transform_scaled(p);
 491                 point q0 = transform_scaled(point(p[0] + 1, p[1]));
 492                 point q1 = transform_scaled(point(p[0], p[1] + 1));
 493
 494                 /*
 495                  * Calculate the distance between source image pixel and
 496                  * adjacent source pixels, measured in the coordinate system of
 497                  * the target frame.
 498                  */
 499
 500                 ale_pos ui = fabs(q0[0] - (*q)[0]);
 501                 ale_pos uj = fabs(q0[1] - (*q)[1]);
 502                 ale_pos vi = fabs(q1[0] - (*q)[0]);
 503                 ale_pos vj = fabs(q1[1] - (*q)[1]);
 504
 505                 /*
 506                  * We map the area of the source image pixel P onto the target
 507                  * frame as a rectangular area oriented on the target frame's
 508                  * axes.  Note that this results in an area that may be the
 509                  * wrong shape or orientation.
 510                  *
 511                  * We define two estimates of the rectangle's dimensions below.
 512                  * For rotations of 0, 90, 180, or 270 degrees, max and sum are
 513                  * identical.  For other orientations, sum is too large and max
 514                  * is too small.  We use the mean of max and sum, which we then
 515                  * divide by two to obtain the distance between the center and
 516                  * the edge.
 517                  */
 518
 519                 ale_pos maxi = (ui > vi) ? ui : vi;
 520                 ale_pos maxj = (uj > vj) ? uj : vj;
 521                 ale_pos sumi = ui + vi;
 522                 ale_pos sumj = uj + vj;
 523
 524                 d[0] = (maxi + sumi) / 4;
 525                 d[1] = (maxj + sumj) / 4;
 526         }
 527
 528         /*
 529          * Get the position and dimensions of a pixel P mapped from one
 530          * coordinate system to another, using the forward transformation.
 531          * This function uses unscaled input coordinates.
 532          */
 533         virtual void map_area_unscaled(point p, point *q, ale_pos d[2]) {
 534
 535                 /*
 536                  * Determine the coordinates in the target frame for the source
 537                  * image pixel P and two adjacent source pixels.
 538                  */
 539
 540                     (*q) = transform_unscaled(p);
 541                 point q0 = transform_unscaled(point(p[0] + 1, p[1]));
 542                 point q1 = transform_unscaled(point(p[0], p[1] + 1));
 543
 544                 /*
 545                  * Calculate the distance between source image pixel and
 546                  * adjacent source pixels, measured in the coordinate system of
 547                  * the target frame.
 548                  */
 549
 550                 ale_pos ui = fabs(q0[0] - (*q)[0]);
 551                 ale_pos uj = fabs(q0[1] - (*q)[1]);
 552                 ale_pos vi = fabs(q1[0] - (*q)[0]);
 553                 ale_pos vj = fabs(q1[1] - (*q)[1]);
 554
 555                 /*
 556                  * We map the area of the source image pixel P onto the target
 557                  * frame as a rectangular area oriented on the target frame's
 558                  * axes.  Note that this results in an area that may be the
 559                  * wrong shape or orientation.
 560                  *
 561                  * We define two estimates of the rectangle's dimensions below.
 562                  * For rotations of 0, 90, 180, or 270 degrees, max and sum are
 563                  * identical.  For other orientations, sum is too large and max
 564                  * is too small.  We use the mean of max and sum, which we then
 565                  * divide by two to obtain the distance between the center and
 566                  * the edge.
 567                  */
 568
 569                 ale_pos maxi = (ui > vi) ? ui : vi;
 570                 ale_pos maxj = (uj > vj) ? uj : vj;
 571                 ale_pos sumi = ui + vi;
 572                 ale_pos sumj = uj + vj;
 573
 574                 d[0] = (maxi + sumi) / 4;
 575                 d[1] = (maxj + sumj) / 4;
 576         }
 577
 578         /*
 579          * Get the position and dimensions of a pixel P mapped from one
 580          * coordinate system to another, using the inverse transformation.  If
 581          * SCALE_FACTOR is not equal to one, divide out the scale factor to
 582          * obtain unscaled coordinates.  This method is very similar to the
 583          * map_area method above.
 584          */
 585         virtual void unscaled_map_area_inverse(point p, point *q, ale_pos d[2]) {
 586
 587                 /*
 588                  * Determine the coordinates in the target frame for the source
 589                  * image pixel P and two adjacent source pixels.
 590                  */
 591
 592                     (*q) = scaled_inverse_transform(p);
 593                 point q0 = scaled_inverse_transform(point(p[0] + 1, p[1]));
 594                 point q1 = scaled_inverse_transform(point(p[0], p[1] + 1));
 595
 596
 597                 /*
 598                  * Calculate the distance between source image pixel and
 599                  * adjacent source pixels, measured in the coordinate system of
 600                  * the target frame.
 601                  */
 602
 603                 ale_pos ui = fabs(q0[0] - (*q)[0]);
 604                 ale_pos uj = fabs(q0[1] - (*q)[1]);
 605                 ale_pos vi = fabs(q1[0] - (*q)[0]);
 606                 ale_pos vj = fabs(q1[1] - (*q)[1]);
 607
 608                 /*
 609                  * We map the area of the source image pixel P onto the target
 610                  * frame as a rectangular area oriented on the target frame's
 611                  * axes.  Note that this results in an area that may be the
 612                  * wrong shape or orientation.
 613                  *
 614                  * We define two estimates of the rectangle's dimensions below.
 615                  * For rotations of 0, 90, 180, or 270 degrees, max and sum are
 616                  * identical.  For other orientations, sum is too large and max
 617                  * is too small.  We use the mean of max and sum, which we then
 618                  * divide by two to obtain the distance between the center and
 619                  * the edge.
 620                  */
 621
 622                 ale_pos maxi = (ui > vi) ? ui : vi;
 623                 ale_pos maxj = (uj > vj) ? uj : vj;
 624                 ale_pos sumi = ui + vi;
 625                 ale_pos sumj = uj + vj;
 626
 627                 d[0] = (maxi + sumi) / 4;
 628                 d[1] = (maxj + sumj) / 4;
 629
 630                 if (scale_factor != 1) {
 631                         d[0] /= scale_factor;
 632                         d[1] /= scale_factor;
 633                         (*q)[0] /= scale_factor;
 634                         (*q)[1] /= scale_factor;
 635                 }
 636         }
 637
 638         /*
 639          * Modify all projective parameters at once.  Accommodate bugs in the
 640          * version 0 transformation file handler (ALE versions 0.4.0p1 and
 641          * earlier).  This code is only called when using a transformation data
 642          * file created with an old version of ALE.
 643          */
 644         virtual void gpt_v0_set(point x[4]) = 0;
 645
 646         /*
 647          * Modify all euclidean parameters at once.  Accommodate bugs in the
 648          * version 0 transformation file handler (ALE versions 0.4.0p1 and
 649          * earlier).  This code is only called when using a transformation data
 650          * file created with an old version of ALE.
 651          */
 652         virtual void eu_v0_set(ale_pos eu[3]) = 0;
 653
 654         virtual void debug_output() = 0;
 655
 656         virtual ~trans_abstract() {
 657         }
 658 };
 659
 660 #endif