d2/render/psf/rasterizer.h

   1 // Copyright 2003 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 #ifndef __psf_rasterize_h__
  22 #define __psf_rasterize_h__
  23
  24 #include "../../point.h"
  25 #include "raster.h"
  26 #include "psf.h"
  27
  28 /*
  29  * Point-spread function rasterizer.
  30  *
  31  * These operations rasterize a PSF to a multiple of the resolution of the
  32  * rendering grid for a given frame.
  33  */
  34
  35 class rasterizer : public raster {
  36         psf *input;
  37
  38 public:
  39         unsigned int varieties() const {
  40                 return input->varieties();
  41         }
  42
  43         unsigned int select(unsigned int i, unsigned int j) const {
  44                 return input->select(i, j);
  45         }
  46
  47         rasterizer (psf *input, transformation t) {
  48                 this->input = input;
  49
  50                 _height = -input->min_i();
  51
  52                 if (input->max_i() > _height)
  53                         _height = input->max_i();
  54
  55                 _width = -input->min_j();
  56
  57                 if (input->max_j() > _width)
  58                         _width = input->max_j();
  59
  60                 /*
  61                  * Approximate the desired resolution.
  62                  *
  63                  * Assume that maximum resolution is reached at (at least) one
  64                  * of the corners of the image.  (This should be true for
  65                  * projective, Euclidean, and translational transformations,
  66                  * but it would be worthwhile to check/prove this for the
  67                  * projective case at some point, since it's a bit less
  68                  * obvious.)
  69                  */
  70
  71                 point min_diff;
  72
  73                 /*
  74                  * XXX: this loop breaks when height <= 1 or width <= 1.
  75                  */
  76
  77                 for (unsigned int i = 0; i < t.unscaled_height(); i += (t.unscaled_height() - 1))
  78                 for (unsigned int j = 0; j < t.unscaled_width();  j += (t.unscaled_width()  - 1)) {
  79                         point corner = point(i, j);
  80                         point delta1 = corner - t.scaled_inverse_transform(t.transform_scaled(corner) + point(1, 0));
  81                         point delta2 = corner - t.scaled_inverse_transform(t.transform_scaled(corner) + point(0, 1));
  82
  83                         for (int index = 0; index < 2; index++) {
  84                                 ale_real d1 = fabs(delta1[index]);
  85                                 ale_real d2 = fabs(delta2[index]);
  86
  87                                 /*
  88                                  * Take the largest change in each direction.
  89                                  */
  90
  91                                 ale_real delta = (d1 > d2) ? d1 : d2;
  92
  93                                 if ((i == 0 && j == 0) || delta < min_diff[index])
  94                                         min_diff[index] = delta;
  95                         }
  96                 }
  97
  98                 ale_real resolution_multiplier = 3;  /* Arbitrary */
  99
 100                 _filter_dim_i = (int) ceil(2 * _height * resolution_multiplier / min_diff[0]);
 101                 _filter_dim_j = (int) ceil(2 * _width * resolution_multiplier / min_diff[1]);
 102
 103                 /*
 104                  * Ensure that the array has an odd number of elements in each
 105                  * direction.  This allows us to move the center to the right
 106                  * place when using FFTW.
 107                  */
 108
 109                 if (_filter_dim_i % 2 == 0)
 110                         _filter_dim_i++;
 111                 if (_filter_dim_j % 2 == 0)
 112                         _filter_dim_j++;
 113
 114                 /*
 115                  * Determine the number of arrays to create.
 116                  */
 117
 118                 num_arrays = input->varieties();
 119
 120                 /*
 121                  * Create arrays
 122                  */
 123
 124                 response_arrays = (ale_real **)malloc(num_arrays * sizeof(ale_real *));
 125
 126                 if (!response_arrays) {
 127                         fprintf(stderr, "Could not allocate in rasterizer.\n");
 128                         exit(1);
 129                 }
 130
 131                 ale_real stepsize_i = (2 * _height) / _filter_dim_i;
 132                 ale_real stepsize_j = (2 * _width) / _filter_dim_j;
 133                 ale_real divisor = stepsize_i * stepsize_j;
 134
 135                 for (unsigned int n = 0; n < num_arrays; n++) {
 136                         response_arrays[n] = (ale_real *)malloc(_filter_dim_i * _filter_dim_j * 3
 137                                                             * sizeof(ale_real));
 138
 139                         if (!response_arrays[n]) {
 140                                 fprintf(stderr, "Could not allocate in rasterizer.\n");
 141                                 exit(1);
 142                         }
 143
 144                         for (unsigned int i = 0; i < _filter_dim_i; i++)
 145                         for (unsigned int j = 0; j < _filter_dim_j; j++) {
 146                                 psf_result r
 147                                         = (*input)(-_height + stepsize_i * i,
 148                                                    -_height + stepsize_i * (i + 1),
 149                                                    -_width  + stepsize_j * j,
 150                                                    -_width  + stepsize_j * (j + 1), n);
 151
 152                                 for (unsigned int k = 0; k < 3; k++) {
 153                                         response_arrays[n][i * _filter_dim_j * 3 + j * 3 + k]
 154                                                 = r.matrix(k, k) / divisor;
 155                                 }
 156                         }
 157                 }
 158
 159 #if 0
 160                 avg_response = (ale_real *)malloc(_filter_dim_i * _filter_dim_j * 3
 161                                               * sizeof(ale_real));
 162
 163                 if (!avg_response) {
 164                         fprintf(stderr, "Could not allocate in rasterizer.\n");
 165                         exit(1);
 166                 }
 167
 168                 for (unsigned int i = 0; i < _filter_dim_i; i++)
 169                 for (unsigned int j = 0; j < _filter_dim_j; j++) {
 170                         psf::psf_result r
 171                                 = (*input)(-_height + stepsize_i * i,
 172                                            -_height + stepsize_i * (i + 1),
 173                                            -_width  + stepsize_j * j,
 174                                            -_width  + stepsize_j * (j + 1));
 175
 176                         for (unsigned int k = 0; k < 3; k++)
 177                                 avg_response[i * _filter_dim_j * 3 + j * 3 + k]
 178                                         = r.matrix(k, k) / divisor;
 179                 }
 180 #endif
 181
 182                 compute_integrals();
 183
 184 //              fprintf(stderr, "(w=%f h=%f we=%d he=%d [", _width, _height, _filter_dim_j, _filter_dim_i);
 185 //              for (unsigned int i = 0; i < _filter_dim_i; i++)
 186 //              for (unsigned int j = 0; j < _filter_dim_j; j++)
 187 //                      fprintf(stderr, "%f ", response_arrays[0][i * _filter_dim_j * 3 + j * 3 + 0]);
 188 //              fprintf(stderr, "])");
 189         }
 190 };
 191
 192 #endif