misc/tools/fft-normalmap-to-heightmap.c

   1 /*
   2  *  FFT based normalmap to heightmap converter
   3  *  Copyright (C) 2010  Rudolf Polzer
   4  *
   5  *  This program is free software; you can redistribute it and/or modify
   6  *  it under the terms of the GNU General Public License as published by
   7  *  the Free Software Foundation; either version 2 of the License, or
   8  *  (at your option) any later version.
   9  *
  10  *  This program is distributed in the hope that it will be useful,
  11  *  but WITHOUT ANY WARRANTY; without even the implied warranty of
  12  *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  13  *  GNU General Public License for more details.
  14  *
  15  *  You should have received a copy of the GNU General Public License
  16  *  along with this program; if not, write to the Free Software
  17  *  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
  18  */
  19
  20 #undef C99
  21 #if __STDC_VERSION__ >= 199901L
  22 #define C99
  23 #endif
  24
  25 #ifdef C99
  26 #include <complex.h>
  27 #endif
  28
  29 #include <stdio.h>
  30 #include <stdlib.h>
  31 #include <string.h>
  32 #include <math.h>
  33
  34 #include <fftw3.h>
  35
  36 #define TWO_PI (4*atan2(1,1) * 2)
  37
  38 void nmap_to_hmap(unsigned char *map, const unsigned char *refmap, int w, int h, double scale, double offset, const double *filter, int filterw, int filterh, int renormalize, double highpass)
  39 {
  40         int x, y;
  41         int i, j;
  42         double fx, fy;
  43         double ffx, ffy;
  44         double nx, ny, nz;
  45         double v, vmin, vmax;
  46 #ifndef C99
  47         double save;
  48 #endif
  49
  50         fftw_complex *imgspace1 = fftw_malloc(w*h * sizeof(fftw_complex));
  51         fftw_complex *imgspace2 = fftw_malloc(w*h * sizeof(fftw_complex));
  52         fftw_complex *freqspace1 = fftw_malloc(w*h * sizeof(fftw_complex));
  53         fftw_complex *freqspace2 = fftw_malloc(w*h * sizeof(fftw_complex));
  54         fftw_plan i12f1 = fftw_plan_dft_2d(h, w, imgspace1, freqspace1, FFTW_FORWARD, FFTW_ESTIMATE);
  55         fftw_plan i22f2 = fftw_plan_dft_2d(h, w, imgspace2, freqspace2, FFTW_FORWARD, FFTW_ESTIMATE);
  56         fftw_plan f12i1 = fftw_plan_dft_2d(h, w, freqspace1, imgspace1, FFTW_BACKWARD, FFTW_ESTIMATE);
  57
  58         for(y = 0; y < h; ++y)
  59         for(x = 0; x < w; ++x)
  60         {
  61                 /*
  62                  * unnormalized normals:
  63                  * n_x = -dh/dx
  64                  * n_y = -dh/dy
  65                  * n_z = -dh/dh = -1
  66                  * BUT: darkplaces uses inverted normals, n_y actually is dh/dy by image pixel coordinates
  67                  */
  68                 nx = ((int)map[(w*y+x)*4+2] - 127.5) / 128;
  69                 ny = ((int)map[(w*y+x)*4+1] - 127.5) / 128;
  70                 nz = ((int)map[(w*y+x)*4+0] - 127.5) / 128;
  71
  72                 /* reconstruct the derivatives from here */
  73 #ifdef C99
  74                 imgspace1[(w*y+x)] =  nx / nz * w; /* = dz/dx */
  75                 imgspace2[(w*y+x)] = -ny / nz * h; /* = dz/dy */
  76 #else
  77                 imgspace1[(w*y+x)][0] =  nx / nz * w; /* = dz/dx */
  78                 imgspace1[(w*y+x)][1] = 0;
  79                 imgspace2[(w*y+x)][0] = -ny / nz * h; /* = dz/dy */
  80                 imgspace2[(w*y+x)][1] = 0;
  81 #endif
  82
  83                 if(renormalize)
  84                 {
  85                         double v = nx * nx + ny * ny + nz * nz;
  86                         if(v > 0)
  87                         {
  88                                 v = sqrt(v);
  89                                 nx /= v;
  90                                 ny /= v;
  91                                 nz /= v;
  92                                 map[(w*y+x)*4+2] = floor(nx * 127.5 + 128);
  93                                 map[(w*y+x)*4+1] = floor(ny * 127.5 + 128);
  94                                 map[(w*y+x)*4+0] = floor(nz * 127.5 + 128);
  95                         }
  96                 }
  97         }
  98
  99         /* see http://www.gamedev.net/community/forums/topic.asp?topic_id=561430 */
 100
 101         fftw_execute(i12f1);
 102         fftw_execute(i22f2);
 103
 104         for(y = 0; y < h; ++y)
 105         for(x = 0; x < w; ++x)
 106         {
 107                 fx = x * 1.0 / w;
 108                 fy = y * 1.0 / h;
 109                 if(fx > 0.5)
 110                         fx -= 1;
 111                 if(fy > 0.5)
 112                         fy -= 1;
 113                 if(filter)
 114                 {
 115                         // discontinous case
 116                         // we must invert whatever "filter" would do on (x, y)!
 117 #ifdef C99
 118                         fftw_complex response_x = 0;
 119                         fftw_complex response_y = 0;
 120                         double sum;
 121                         for(i = -filterh / 2; i <= filterh / 2; ++i)
 122                                 for(j = -filterw / 2; j <= filterw / 2; ++j)
 123                                 {
 124                                         response_x += filter[(i + filterh / 2) * filterw + j + filterw / 2] * cexp(-_Complex_I * TWO_PI * (j * fx + i * fy));
 125                                         response_y += filter[(i + filterh / 2) * filterw + j + filterw / 2] * cexp(-_Complex_I * TWO_PI * (i * fx + j * fy));
 126                                 }
 127
 128                         // we know:
 129                         //   fourier(df/dx)_xy = fourier(f)_xy * response_x
 130                         //   fourier(df/dy)_xy = fourier(f)_xy * response_y
 131                         // mult by conjugate of response_x, response_y:
 132                         //   conj(response_x) * fourier(df/dx)_xy = fourier(f)_xy * |response_x^2|
 133                         //   conj(response_y) * fourier(df/dy)_xy = fourier(f)_xy * |response_y^2|
 134                         // and
 135                         //   fourier(f)_xy = (conj(response_x) * fourier(df/dx)_xy + conj(response_y) * fourier(df/dy)_xy) / (|response_x|^2 + |response_y|^2)
 136
 137                         sum = cabs(response_x) * cabs(response_x) + cabs(response_y) * cabs(response_y);
 138
 139                         if(sum > 0)
 140                                 freqspace1[(w*y+x)] = (conj(response_x) * freqspace1[(w*y+x)] + conj(response_y) * freqspace2[(w*y+x)]) / sum;
 141                         else
 142                                 freqspace1[(w*y+x)] = 0;
 143 #else
 144                         fftw_complex response_x = {0, 0};
 145                         fftw_complex response_y = {0, 0};
 146                         double sum;
 147                         for(i = -filterh / 2; i <= filterh / 2; ++i)
 148                                 for(j = -filterw / 2; j <= filterw / 2; ++j)
 149                                 {
 150                                         response_x[0] += filter[(i + filterh / 2) * filterw + j + filterw / 2] * cos(-TWO_PI * (j * fx + i * fy));
 151                                         response_x[1] += filter[(i + filterh / 2) * filterw + j + filterw / 2] * sin(-TWO_PI * (j * fx + i * fy));
 152                                         response_y[0] += filter[(i + filterh / 2) * filterw + j + filterw / 2] * cos(-TWO_PI * (i * fx + j * fy));
 153                                         response_y[1] += filter[(i + filterh / 2) * filterw + j + filterw / 2] * sin(-TWO_PI * (i * fx + j * fy));
 154                                 }
 155
 156                         sum = response_x[0] * response_x[0] + response_x[1] * response_x[1]
 157                             + response_y[0] * response_y[0] + response_y[1] * response_y[1];
 158
 159                         if(sum > 0)
 160                         {
 161                                 double s = freqspace1[(w*y+x)][0];
 162                                 freqspace1[(w*y+x)][0] = (response_x[0] * s                      + response_x[1] * freqspace1[(w*y+x)][1] + response_y[0] * freqspace2[(w*y+x)][0] + response_y[1] * freqspace2[(w*y+x)][1]) / sum;
 163                                 freqspace1[(w*y+x)][1] = (response_x[0] * freqspace1[(w*y+x)][1] - response_x[1] * s                      + response_y[0] * freqspace2[(w*y+x)][1] - response_y[1] * freqspace2[(w*y+x)][0]) / sum;
 164                         }
 165                         else
 166                         {
 167                                 freqspace1[(w*y+x)][0] = 0;
 168                                 freqspace1[(w*y+x)][1] = 0;
 169                         }
 170 #endif
 171                 }
 172                 else
 173                 {
 174                         // continuous integration case
 175                         /* these must have the same sign as fx and fy (so ffx*fx + ffy*fy is nonzero), otherwise do not matter */
 176                         /* it basically decides how artifacts are distributed */
 177                         ffx = fx;
 178                         ffy = fy;
 179 #ifdef C99
 180                         if(fx||fy)
 181                                 freqspace1[(w*y+x)] = _Complex_I * (ffx * freqspace1[(w*y+x)] + ffy * freqspace2[(w*y+x)]) / (ffx*fx + ffy*fy) / TWO_PI;
 182                         else
 183                                 freqspace1[(w*y+x)] = 0;
 184 #else
 185                         if(fx||fy)
 186                         {
 187                                 save = freqspace1[(w*y+x)][0];
 188                                 freqspace1[(w*y+x)][0] = -(ffx * freqspace1[(w*y+x)][1] + ffy * freqspace2[(w*y+x)][1]) / (ffx*fx + ffy*fy) / TWO_PI;
 189                                 freqspace1[(w*y+x)][1] =  (ffx * save + ffy * freqspace2[(w*y+x)][0]) / (ffx*fx + ffy*fy) / TWO_PI;
 190                         }
 191                         else
 192                         {
 193                                 freqspace1[(w*y+x)][0] = 0;
 194                                 freqspace1[(w*y+x)][1] = 0;
 195                         }
 196 #endif
 197                 }
 198                 if(highpass > 0)
 199                 {
 200                         double f1 = (fabs(fx)*highpass);
 201                         double f2 = (fabs(fy)*highpass);
 202                         // if either of them is < 1, phase out (min at 0.5)
 203                         double f =
 204                                 (f1 <= 0.5 ? 0 : (f1 >= 1 ? 1 : ((f1 - 0.5) * 2.0)))
 205                                 *
 206                                 (f2 <= 0.5 ? 0 : (f2 >= 1 ? 1 : ((f2 - 0.5) * 2.0)));
 207 #ifdef C99
 208                         freqspace1[(w*y+x)] *= f;
 209 #else
 210                         freqspace1[(w*y+x)][0] *= f;
 211                         freqspace1[(w*y+x)][1] *= f;
 212 #endif
 213                 }
 214         }
 215
 216         fftw_execute(f12i1);
 217
 218         /* renormalize, find min/max */
 219         vmin = vmax = 0;
 220         for(y = 0; y < h; ++y)
 221         for(x = 0; x < w; ++x)
 222         {
 223 #ifdef C99
 224                 v = creal(imgspace1[(w*y+x)] /= pow(w*h, 1.5));
 225 #else
 226                 v = (imgspace1[(w*y+x)][0] /= pow(w*h, 1.5));
 227                 // imgspace1[(w*y+x)][1] /= pow(w*h, 1.5);
 228                 // this value is never used
 229 #endif
 230                 if(v < vmin || (x == 0 && y == 0))
 231                         vmin = v;
 232                 if(v > vmax || (x == 0 && y == 0))
 233                         vmax = v;
 234         }
 235
 236         if(refmap)
 237         {
 238                 double f, a;
 239                 double o, s;
 240                 double sa, sfa, sffa, sfva, sva;
 241                 double mi, ma;
 242                 sa = sfa = sffa = sfva = sva = 0;
 243                 mi = 1;
 244                 ma = -1;
 245                 for(y = 0; y < h; ++y)
 246                 for(x = 0; x < w; ++x)
 247                 {
 248                         a = (int)refmap[(w*y+x)*4+3];
 249                         v = (refmap[(w*y+x)*4+0]*0.114 + refmap[(w*y+x)*4+1]*0.587 + refmap[(w*y+x)*4+2]*0.299);
 250                         v = (v - 128.0) / 127.0;
 251 #ifdef C99
 252                         f = creal(imgspace1[(w*y+x)]);
 253 #else
 254                         f = imgspace1[(w*y+x)][0];
 255 #endif
 256                         if(a <= 0)
 257                                 continue;
 258                         if(v < mi)
 259                                 mi = v;
 260                         if(v > ma)
 261                                 ma = v;
 262                         sa += a;
 263                         sfa += f*a;
 264                         sffa += f*f*a;
 265                         sfva += f*v*a;
 266                         sva += v*a;
 267                 }
 268                 if(mi < ma)
 269                 {
 270                         /* linear regression ftw */
 271                         o = (sfa*sfva - sffa*sva) / (sfa*sfa-sa*sffa);
 272                         s = (sfa*sva - sa*sfva) / (sfa*sfa-sa*sffa);
 273                 }
 274                 else /* all values of v are equal, so we cannot get scale; we can still get offset */
 275                 {
 276                         o = ((sva - sfa) / sa);
 277                         s = 1;
 278                 }
 279
 280                 /*
 281                  * now apply user-given offset and scale to these values
 282                  * (x * s + o) * scale + offset
 283                  * x * s * scale + o * scale + offset
 284                  */
 285                 offset += o * scale;
 286                 scale *= s;
 287         }
 288         else if(scale == 0)
 289         {
 290                 /*
 291                  * map vmin to -1
 292                  * map vmax to +1
 293                  */
 294                 scale = 2 / (vmax - vmin);
 295                 offset = -(vmax + vmin) / (vmax - vmin);
 296         }
 297
 298         printf("Min: %f\nAvg: %f\nMax: %f\nScale: %f\nOffset: %f\nScaled-Min: %f\nScaled-Avg: %f\nScaled-Max: %f\n",
 299                 vmin, 0.0, vmax, scale, offset, vmin * scale + offset, offset, vmax * scale + offset);
 300
 301         for(y = 0; y < h; ++y)
 302         for(x = 0; x < w; ++x)
 303         {
 304 #ifdef C99
 305                 v = creal(imgspace1[(w*y+x)]);
 306 #else
 307                 v = imgspace1[(w*y+x)][0];
 308 #endif
 309                 v = v * scale + offset;
 310                 if(v < -1)
 311                         v = -1;
 312                 if(v > 1)
 313                         v = 1;
 314                 map[(w*y+x)*4+3] = floor(128.5 + 127 * v);
 315         }
 316
 317         fftw_destroy_plan(i12f1);
 318         fftw_destroy_plan(i22f2);
 319         fftw_destroy_plan(f12i1);
 320
 321         fftw_free(freqspace2);
 322         fftw_free(freqspace1);
 323         fftw_free(imgspace2);
 324         fftw_free(imgspace1);
 325 }
 326
 327 void hmap_to_nmap(unsigned char *map, int w, int h, int src_chan, double scale)
 328 {
 329         int x, y;
 330         double fx, fy;
 331         double nx, ny, nz;
 332         double v;
 333 #ifndef C99
 334         double save;
 335 #endif
 336
 337         fftw_complex *imgspace1 = fftw_malloc(w*h * sizeof(fftw_complex));
 338         fftw_complex *imgspace2 = fftw_malloc(w*h * sizeof(fftw_complex));
 339         fftw_complex *freqspace1 = fftw_malloc(w*h * sizeof(fftw_complex));
 340         fftw_complex *freqspace2 = fftw_malloc(w*h * sizeof(fftw_complex));
 341         fftw_plan i12f1 = fftw_plan_dft_2d(h, w, imgspace1, freqspace1, FFTW_FORWARD, FFTW_ESTIMATE);
 342         fftw_plan f12i1 = fftw_plan_dft_2d(h, w, freqspace1, imgspace1, FFTW_BACKWARD, FFTW_ESTIMATE);
 343         fftw_plan f22i2 = fftw_plan_dft_2d(h, w, freqspace2, imgspace2, FFTW_BACKWARD, FFTW_ESTIMATE);
 344
 345         for(y = 0; y < h; ++y)
 346         for(x = 0; x < w; ++x)
 347         {
 348                 switch(src_chan)
 349                 {
 350                         case 0:
 351                         case 1:
 352                         case 2:
 353                         case 3:
 354                                 v = map[(w*y+x)*4+src_chan];
 355                                 break;
 356                         case 4:
 357                                 v = (map[(w*y+x)*4+0] + map[(w*y+x)*4+1] + map[(w*y+x)*4+2]) / 3;
 358                                 break;
 359                         default:
 360                         case 5:
 361                                 v = (map[(w*y+x)*4+0]*0.114 + map[(w*y+x)*4+1]*0.587 + map[(w*y+x)*4+2]*0.299);
 362                                 break;
 363                 }
 364 #ifdef C99
 365                 imgspace1[(w*y+x)] = (v - 128.0) / 127.0;
 366 #else
 367                 imgspace1[(w*y+x)][0] = (v - 128.0) / 127.0;
 368                 imgspace1[(w*y+x)][1] = 0;
 369 #endif
 370                 if(v < 1)
 371                         v = 1; /* do not write alpha zero */
 372                 map[(w*y+x)*4+3] = floor(v + 0.5);
 373         }
 374
 375         /* see http://www.gamedev.net/community/forums/topic.asp?topic_id=561430 */
 376
 377         fftw_execute(i12f1);
 378
 379         for(y = 0; y < h; ++y)
 380         for(x = 0; x < w; ++x)
 381         {
 382                 fx = x;
 383                 fy = y;
 384                 if(fx > w/2)
 385                         fx -= w;
 386                 if(fy > h/2)
 387                         fy -= h;
 388 #ifdef DISCONTINUOUS
 389                 fx = sin(fx * TWO_PI / w);
 390                 fy = sin(fy * TWO_PI / h);
 391 #else
 392 #ifdef C99
 393                 /* a lowpass to prevent the worst */
 394                 freqspace1[(w*y+x)] *= 1 - pow(abs(fx) / (double)(w/2), 1);
 395                 freqspace1[(w*y+x)] *= 1 - pow(abs(fy) / (double)(h/2), 1);
 396 #else
 397                 /* a lowpass to prevent the worst */
 398                 freqspace1[(w*y+x)][0] *= 1 - pow(abs(fx) / (double)(w/2), 1);
 399                 freqspace1[(w*y+x)][1] *= 1 - pow(abs(fx) / (double)(w/2), 1);
 400                 freqspace1[(w*y+x)][0] *= 1 - pow(abs(fy) / (double)(h/2), 1);
 401                 freqspace1[(w*y+x)][1] *= 1 - pow(abs(fy) / (double)(h/2), 1);
 402 #endif
 403 #endif
 404 #ifdef C99
 405                 freqspace2[(w*y+x)] = TWO_PI*_Complex_I * fy * freqspace1[(w*y+x)]; /* y derivative */
 406                 freqspace1[(w*y+x)] = TWO_PI*_Complex_I * fx * freqspace1[(w*y+x)]; /* x derivative */
 407 #else
 408                 freqspace2[(w*y+x)][0] = -TWO_PI * fy * freqspace1[(w*y+x)][1]; /* y derivative */
 409                 freqspace2[(w*y+x)][1] =  TWO_PI * fy * freqspace1[(w*y+x)][0];
 410                 save = freqspace1[(w*y+x)][0];
 411                 freqspace1[(w*y+x)][0] = -TWO_PI * fx * freqspace1[(w*y+x)][1]; /* x derivative */
 412                 freqspace1[(w*y+x)][1] =  TWO_PI * fx * save;
 413 #endif
 414         }
 415
 416         fftw_execute(f12i1);
 417         fftw_execute(f22i2);
 418
 419         scale /= (w*h);
 420
 421         for(y = 0; y < h; ++y)
 422         for(x = 0; x < w; ++x)
 423         {
 424 #ifdef C99
 425                 nx = creal(imgspace1[(w*y+x)]);
 426                 ny = creal(imgspace2[(w*y+x)]);
 427 #else
 428                 nx = imgspace1[(w*y+x)][0];
 429                 ny = imgspace2[(w*y+x)][0];
 430 #endif
 431                 nx /= w;
 432                 ny /= h;
 433                 nz = -1 / scale;
 434                 v = -sqrt(nx*nx + ny*ny + nz*nz);
 435                 nx /= v;
 436                 ny /= v;
 437                 nz /= v;
 438                 ny = -ny; /* DP inverted normals */
 439                 map[(w*y+x)*4+2] = floor(128 + 127.5 * nx);
 440                 map[(w*y+x)*4+1] = floor(128 + 127.5 * ny);
 441                 map[(w*y+x)*4+0] = floor(128 + 127.5 * nz);
 442         }
 443
 444         fftw_destroy_plan(i12f1);
 445         fftw_destroy_plan(f12i1);
 446         fftw_destroy_plan(f22i2);
 447
 448         fftw_free(freqspace2);
 449         fftw_free(freqspace1);
 450         fftw_free(imgspace2);
 451         fftw_free(imgspace1);
 452 }
 453
 454 void hmap_to_nmap_local(unsigned char *map, int w, int h, int src_chan, double scale, const double *filter, int filterw, int filterh)
 455 {
 456         int x, y;
 457         double nx, ny, nz;
 458         double v;
 459         int i, j;
 460         double *img_reduced = malloc(w*h * sizeof(double));
 461
 462         for(y = 0; y < h; ++y)
 463         for(x = 0; x < w; ++x)
 464         {
 465                 switch(src_chan)
 466                 {
 467                         case 0:
 468                         case 1:
 469                         case 2:
 470                         case 3:
 471                                 v = map[(w*y+x)*4+src_chan];
 472                                 break;
 473                         case 4:
 474                                 v = (map[(w*y+x)*4+0] + map[(w*y+x)*4+1] + map[(w*y+x)*4+2]) / 3;
 475                                 break;
 476                         default:
 477                         case 5:
 478                                 v = (map[(w*y+x)*4+0]*0.114 + map[(w*y+x)*4+1]*0.587 + map[(w*y+x)*4+2]*0.299);
 479                                 break;
 480                 }
 481                 img_reduced[(w*y+x)] = (v - 128.0) / 127.0;
 482                 if(v < 1)
 483                         v = 1; /* do not write alpha zero */
 484                 map[(w*y+x)*4+3] = floor(v + 0.5);
 485         }
 486
 487         for(y = 0; y < h; ++y)
 488         for(x = 0; x < w; ++x)
 489         {
 490                 nz = -1 / scale;
 491                 nx = ny = 0;
 492
 493                 for(i = -filterh / 2; i <= filterh / 2; ++i)
 494                         for(j = -filterw / 2; j <= filterw / 2; ++j)
 495                         {
 496                                 nx += img_reduced[w*((y+i+h)%h)+(x+j+w)%w] * filter[(i + filterh / 2) * filterw + j + filterw / 2];
 497                                 ny += img_reduced[w*((y+j+h)%h)+(x+i+w)%w] * filter[(i + filterh / 2) * filterw + j + filterw / 2];
 498                         }
 499
 500                 v = -sqrt(nx*nx + ny*ny + nz*nz);
 501                 nx /= v;
 502                 ny /= v;
 503                 nz /= v;
 504                 ny = -ny; /* DP inverted normals */
 505                 map[(w*y+x)*4+2] = floor(128 + 127.5 * nx);
 506                 map[(w*y+x)*4+1] = floor(128 + 127.5 * ny);
 507                 map[(w*y+x)*4+0] = floor(128 + 127.5 * nz);
 508         }
 509
 510         free(img_reduced);
 511 }
 512
 513 unsigned char *FS_LoadFile(const char *fn, int *len)
 514 {
 515         unsigned char *buf = NULL;
 516         int n;
 517         FILE *f = fopen(fn, "rb");
 518         *len = 0;
 519         if(!f)
 520                 return NULL;
 521         for(;;)
 522         {
 523                 buf = realloc(buf, *len + 65536);
 524                 if(!buf)
 525                 {
 526                         fclose(f);
 527                         free(buf);
 528                         *len = 0;
 529                         return NULL;
 530                 }
 531                 n = fread(buf + *len, 1, 65536, f);
 532                 if(n < 0)
 533                 {
 534                         fclose(f);
 535                         free(buf);
 536                         *len = 0;
 537                         return NULL;
 538                 }
 539                 *len += n;
 540                 if(n < 65536)
 541                         break;
 542         }
 543         return buf;
 544 }
 545
 546 int FS_WriteFile(const char *fn, unsigned char *data, int len)
 547 {
 548         FILE *f = fopen(fn, "wb");
 549         if(!f)
 550                 return 0;
 551         if(fwrite(data, len, 1, f) != 1)
 552         {
 553                 fclose(f);
 554                 return 0;
 555         }
 556         if(fclose(f))
 557                 return 0;
 558         return 1;
 559 }
 560
 561 /* START stuff that originates from image.c in DarkPlaces */
 562 int image_width, image_height;
 563
 564 typedef struct _TargaHeader
 565 {
 566         unsigned char   id_length, colormap_type, image_type;
 567         unsigned short  colormap_index, colormap_length;
 568         unsigned char   colormap_size;
 569         unsigned short  x_origin, y_origin, width, height;
 570         unsigned char   pixel_size, attributes;
 571 }
 572 TargaHeader;
 573
 574 void PrintTargaHeader(TargaHeader *t)
 575 {
 576         printf("TargaHeader:\nuint8 id_length = %i;\nuint8 colormap_type = %i;\nuint8 image_type = %i;\nuint16 colormap_index = %i;\nuint16 colormap_length = %i;\nuint8 colormap_size = %i;\nuint16 x_origin = %i;\nuint16 y_origin = %i;\nuint16 width = %i;\nuint16 height = %i;\nuint8 pixel_size = %i;\nuint8 attributes = %i;\n", t->id_length, t->colormap_type, t->image_type, t->colormap_index, t->colormap_length, t->colormap_size, t->x_origin, t->y_origin, t->width, t->height, t->pixel_size, t->attributes);
 577 }
 578
 579 unsigned char *LoadTGA_BGRA (const unsigned char *f, int filesize)
 580 {
 581         int x, y, pix_inc, row_inci, runlen, alphabits;
 582         unsigned char *image_buffer;
 583         unsigned int *pixbufi;
 584         const unsigned char *fin, *enddata;
 585         TargaHeader targa_header;
 586         unsigned int palettei[256];
 587         union
 588         {
 589                 unsigned int i;
 590                 unsigned char b[4];
 591         }
 592         bgra;
 593
 594         if (filesize < 19)
 595                 return NULL;
 596
 597         enddata = f + filesize;
 598
 599         targa_header.id_length = f[0];
 600         targa_header.colormap_type = f[1];
 601         targa_header.image_type = f[2];
 602
 603         targa_header.colormap_index = f[3] + f[4] * 256;
 604         targa_header.colormap_length = f[5] + f[6] * 256;
 605         targa_header.colormap_size = f[7];
 606         targa_header.x_origin = f[8] + f[9] * 256;
 607         targa_header.y_origin = f[10] + f[11] * 256;
 608         targa_header.width = image_width = f[12] + f[13] * 256;
 609         targa_header.height = image_height = f[14] + f[15] * 256;
 610         targa_header.pixel_size = f[16];
 611         targa_header.attributes = f[17];
 612
 613         if (image_width > 32768 || image_height > 32768 || image_width <= 0 || image_height <= 0)
 614         {
 615                 printf("LoadTGA: invalid size\n");
 616                 PrintTargaHeader(&targa_header);
 617                 return NULL;
 618         }
 619
 620         /* advance to end of header */
 621         fin = f + 18;
 622
 623         /* skip TARGA image comment (usually 0 bytes) */
 624         fin += targa_header.id_length;
 625
 626         /* read/skip the colormap if present (note: according to the TARGA spec it */
 627         /* can be present even on 1color or greyscale images, just not used by */
 628         /* the image data) */
 629         if (targa_header.colormap_type)
 630         {
 631                 if (targa_header.colormap_length > 256)
 632                 {
 633                         printf("LoadTGA: only up to 256 colormap_length supported\n");
 634                         PrintTargaHeader(&targa_header);
 635                         return NULL;
 636                 }
 637                 if (targa_header.colormap_index)
 638                 {
 639                         printf("LoadTGA: colormap_index not supported\n");
 640                         PrintTargaHeader(&targa_header);
 641                         return NULL;
 642                 }
 643                 if (targa_header.colormap_size == 24)
 644                 {
 645                         for (x = 0;x < targa_header.colormap_length;x++)
 646                         {
 647                                 bgra.b[0] = *fin++;
 648                                 bgra.b[1] = *fin++;
 649                                 bgra.b[2] = *fin++;
 650                                 bgra.b[3] = 255;
 651                                 palettei[x] = bgra.i;
 652                         }
 653                 }
 654                 else if (targa_header.colormap_size == 32)
 655                 {
 656                         memcpy(palettei, fin, targa_header.colormap_length*4);
 657                         fin += targa_header.colormap_length * 4;
 658                 }
 659                 else
 660                 {
 661                         printf("LoadTGA: Only 32 and 24 bit colormap_size supported\n");
 662                         PrintTargaHeader(&targa_header);
 663                         return NULL;
 664                 }
 665         }
 666
 667         /* check our pixel_size restrictions according to image_type */
 668         switch (targa_header.image_type & ~8)
 669         {
 670         case 2:
 671                 if (targa_header.pixel_size != 24 && targa_header.pixel_size != 32)
 672                 {
 673                         printf("LoadTGA: only 24bit and 32bit pixel sizes supported for type 2 and type 10 images\n");
 674                         PrintTargaHeader(&targa_header);
 675                         return NULL;
 676                 }
 677                 break;
 678         case 3:
 679                 /* set up a palette to make the loader easier */
 680                 for (x = 0;x < 256;x++)
 681                 {
 682                         bgra.b[0] = bgra.b[1] = bgra.b[2] = x;
 683                         bgra.b[3] = 255;
 684                         palettei[x] = bgra.i;
 685                 }
 686                 /* fall through to colormap case */
 687         case 1:
 688                 if (targa_header.pixel_size != 8)
 689                 {
 690                         printf("LoadTGA: only 8bit pixel size for type 1, 3, 9, and 11 images supported\n");
 691                         PrintTargaHeader(&targa_header);
 692                         return NULL;
 693                 }
 694                 break;
 695         default:
 696                 printf("LoadTGA: Only type 1, 2, 3, 9, 10, and 11 targa RGB images supported, image_type = %i\n", targa_header.image_type);
 697                 PrintTargaHeader(&targa_header);
 698                 return NULL;
 699         }
 700
 701         if (targa_header.attributes & 0x10)
 702         {
 703                 printf("LoadTGA: origin must be in top left or bottom left, top right and bottom right are not supported\n");
 704                 return NULL;
 705         }
 706
 707         /* number of attribute bits per pixel, we only support 0 or 8 */
 708         alphabits = targa_header.attributes & 0x0F;
 709         if (alphabits != 8 && alphabits != 0)
 710         {
 711                 printf("LoadTGA: only 0 or 8 attribute (alpha) bits supported\n");
 712                 return NULL;
 713         }
 714
 715         image_buffer = (unsigned char *)malloc(image_width * image_height * 4);
 716         if (!image_buffer)
 717         {
 718                 printf("LoadTGA: not enough memory for %i by %i image\n", image_width, image_height);
 719                 return NULL;
 720         }
 721
 722         /* If bit 5 of attributes isn't set, the image has been stored from bottom to top */
 723         if ((targa_header.attributes & 0x20) == 0)
 724         {
 725                 pixbufi = (unsigned int*)image_buffer + (image_height - 1)*image_width;
 726                 row_inci = -image_width*2;
 727         }
 728         else
 729         {
 730                 pixbufi = (unsigned int*)image_buffer;
 731                 row_inci = 0;
 732         }
 733
 734         x = 0;
 735         y = 0;
 736         pix_inc = 1;
 737         if ((targa_header.image_type & ~8) == 2)
 738                 pix_inc = (targa_header.pixel_size + 7) / 8;
 739         switch (targa_header.image_type)
 740         {
 741         case 1: /* colormapped, uncompressed */
 742         case 3: /* greyscale, uncompressed */
 743                 if (fin + image_width * image_height * pix_inc > enddata)
 744                         break;
 745                 for (y = 0;y < image_height;y++, pixbufi += row_inci)
 746                         for (x = 0;x < image_width;x++)
 747                                 *pixbufi++ = palettei[*fin++];
 748                 break;
 749         case 2:
 750                 /* BGR or BGRA, uncompressed */
 751                 if (fin + image_width * image_height * pix_inc > enddata)
 752                         break;
 753                 if (targa_header.pixel_size == 32 && alphabits)
 754                 {
 755                         for (y = 0;y < image_height;y++)
 756                                 memcpy(pixbufi + y * (image_width + row_inci), fin + y * image_width * pix_inc, image_width*4);
 757                 }
 758                 else
 759                 {
 760                         for (y = 0;y < image_height;y++, pixbufi += row_inci)
 761                         {
 762                                 for (x = 0;x < image_width;x++, fin += pix_inc)
 763                                 {
 764                                         bgra.b[0] = fin[0];
 765                                         bgra.b[1] = fin[1];
 766                                         bgra.b[2] = fin[2];
 767                                         bgra.b[3] = 255;
 768                                         *pixbufi++ = bgra.i;
 769                                 }
 770                         }
 771                 }
 772                 break;
 773         case 9: /* colormapped, RLE */
 774         case 11: /* greyscale, RLE */
 775                 for (y = 0;y < image_height;y++, pixbufi += row_inci)
 776                 {
 777                         for (x = 0;x < image_width;)
 778                         {
 779                                 if (fin >= enddata)
 780                                         break; /* error - truncated file */
 781                                 runlen = *fin++;
 782                                 if (runlen & 0x80)
 783                                 {
 784                                         /* RLE - all pixels the same color */
 785                                         runlen += 1 - 0x80;
 786                                         if (fin + pix_inc > enddata)
 787                                                 break; /* error - truncated file */
 788                                         if (x + runlen > image_width)
 789                                                 break; /* error - line exceeds width */
 790                                         bgra.i = palettei[*fin++];
 791                                         for (;runlen--;x++)
 792                                                 *pixbufi++ = bgra.i;
 793                                 }
 794                                 else
 795                                 {
 796                                         /* uncompressed - all pixels different color */
 797                                         runlen++;
 798                                         if (fin + pix_inc * runlen > enddata)
 799                                                 break; /* error - truncated file */
 800                                         if (x + runlen > image_width)
 801                                                 break; /* error - line exceeds width */
 802                                         for (;runlen--;x++)
 803                                                 *pixbufi++ = palettei[*fin++];
 804                                 }
 805                         }
 806
 807                         if (x != image_width)
 808                         {
 809                                 /* pixbufi is useless now */
 810                                 printf("LoadTGA: corrupt file\n");
 811                                 break;
 812                         }
 813                 }
 814                 break;
 815         case 10:
 816                 /* BGR or BGRA, RLE */
 817                 if (targa_header.pixel_size == 32 && alphabits)
 818                 {
 819                         for (y = 0;y < image_height;y++, pixbufi += row_inci)
 820                         {
 821                                 for (x = 0;x < image_width;)
 822                                 {
 823                                         if (fin >= enddata)
 824                                                 break; /* error - truncated file */
 825                                         runlen = *fin++;
 826                                         if (runlen & 0x80)
 827                                         {
 828                                                 /* RLE - all pixels the same color */
 829                                                 runlen += 1 - 0x80;
 830                                                 if (fin + pix_inc > enddata)
 831                                                         break; /* error - truncated file */
 832                                                 if (x + runlen > image_width)
 833                                                         break; /* error - line exceeds width */
 834                                                 bgra.b[0] = fin[0];
 835                                                 bgra.b[1] = fin[1];
 836                                                 bgra.b[2] = fin[2];
 837                                                 bgra.b[3] = fin[3];
 838                                                 fin += pix_inc;
 839                                                 for (;runlen--;x++)
 840                                                         *pixbufi++ = bgra.i;
 841                                         }
 842                                         else
 843                                         {
 844                                                 /* uncompressed - all pixels different color */
 845                                                 runlen++;
 846                                                 if (fin + pix_inc * runlen > enddata)
 847                                                         break; /* error - truncated file */
 848                                                 if (x + runlen > image_width)
 849                                                         break; /* error - line exceeds width */
 850                                                 for (;runlen--;x++)
 851                                                 {
 852                                                         bgra.b[0] = fin[0];
 853                                                         bgra.b[1] = fin[1];
 854                                                         bgra.b[2] = fin[2];
 855                                                         bgra.b[3] = fin[3];
 856                                                         fin += pix_inc;
 857                                                         *pixbufi++ = bgra.i;
 858                                                 }
 859                                         }
 860                                 }
 861
 862                                 if (x != image_width)
 863                                 {
 864                                         /* pixbufi is useless now */
 865                                         printf("LoadTGA: corrupt file\n");
 866                                         break;
 867                                 }
 868                         }
 869                 }
 870                 else
 871                 {
 872                         for (y = 0;y < image_height;y++, pixbufi += row_inci)
 873                         {
 874                                 for (x = 0;x < image_width;)
 875                                 {
 876                                         if (fin >= enddata)
 877                                                 break; /* error - truncated file */
 878                                         runlen = *fin++;
 879                                         if (runlen & 0x80)
 880                                         {
 881                                                 /* RLE - all pixels the same color */
 882                                                 runlen += 1 - 0x80;
 883                                                 if (fin + pix_inc > enddata)
 884                                                         break; /* error - truncated file */
 885                                                 if (x + runlen > image_width)
 886                                                         break; /* error - line exceeds width */
 887                                                 bgra.b[0] = fin[0];
 888                                                 bgra.b[1] = fin[1];
 889                                                 bgra.b[2] = fin[2];
 890                                                 bgra.b[3] = 255;
 891                                                 fin += pix_inc;
 892                                                 for (;runlen--;x++)
 893                                                         *pixbufi++ = bgra.i;
 894                                         }
 895                                         else
 896                                         {
 897                                                 /* uncompressed - all pixels different color */
 898                                                 runlen++;
 899                                                 if (fin + pix_inc * runlen > enddata)
 900                                                         break; /* error - truncated file */
 901                                                 if (x + runlen > image_width)
 902                                                         break; /* error - line exceeds width */
 903                                                 for (;runlen--;x++)
 904                                                 {
 905                                                         bgra.b[0] = fin[0];
 906                                                         bgra.b[1] = fin[1];
 907                                                         bgra.b[2] = fin[2];
 908                                                         bgra.b[3] = 255;
 909                                                         fin += pix_inc;
 910                                                         *pixbufi++ = bgra.i;
 911                                                 }
 912                                         }
 913                                 }
 914
 915                                 if (x != image_width)
 916                                 {
 917                                         /* pixbufi is useless now */
 918                                         printf("LoadTGA: corrupt file\n");
 919                                         break;
 920                                 }
 921                         }
 922                 }
 923                 break;
 924         default:
 925                 /* unknown image_type */
 926                 break;
 927         }
 928
 929         return image_buffer;
 930 }
 931
 932 int Image_WriteTGABGRA (const char *filename, int width, int height, const unsigned char *data)
 933 {
 934         int y;
 935         unsigned char *buffer, *out;
 936         const unsigned char *in, *end;
 937         int ret;
 938
 939         buffer = (unsigned char *)malloc(width*height*4 + 18);
 940
 941         memset (buffer, 0, 18);
 942         buffer[2] = 2;          /* uncompressed type */
 943         buffer[12] = (width >> 0) & 0xFF;
 944         buffer[13] = (width >> 8) & 0xFF;
 945         buffer[14] = (height >> 0) & 0xFF;
 946         buffer[15] = (height >> 8) & 0xFF;
 947
 948         for (y = 3;y < width*height*4;y += 4)
 949                 if (data[y] < 255)
 950                         break;
 951
 952         if (y < width*height*4)
 953         {
 954                 /* save the alpha channel */
 955                 buffer[16] = 32;        /* pixel size */
 956                 buffer[17] = 8; /* 8 bits of alpha */
 957
 958                 /* flip upside down */
 959                 out = buffer + 18;
 960                 for (y = height - 1;y >= 0;y--)
 961                 {
 962                         memcpy(out, data + y * width * 4, width * 4);
 963                         out += width*4;
 964                 }
 965         }
 966         else
 967         {
 968                 /* save only the color channels */
 969                 buffer[16] = 24;        /* pixel size */
 970                 buffer[17] = 0; /* 8 bits of alpha */
 971
 972                 /* truncate bgra to bgr and flip upside down */
 973                 out = buffer + 18;
 974                 for (y = height - 1;y >= 0;y--)
 975                 {
 976                         in = data + y * width * 4;
 977                         end = in + width * 4;
 978                         for (;in < end;in += 4)
 979                         {
 980                                 *out++ = in[0];
 981                                 *out++ = in[1];
 982                                 *out++ = in[2];
 983                         }
 984                 }
 985         }
 986         ret = FS_WriteFile (filename, buffer, out - buffer);
 987
 988         free(buffer);
 989
 990         return ret;
 991 }
 992 /* START stuff that originates from image.c in DarkPlaces */
 993
 994 int usage(const char *me)
 995 {
 996         printf("Usage: %s <infile_norm.tga> <outfile_normandheight.tga> filtertype [<scale> [<offset> [<infile_ref.tga>]]] (get heightmap from normalmap)\n", me);
 997         printf("or:    %s <infile_height.tga> <outfile_normandheight.tga> filtertype -1 [<scale>] (read from B)\n", me);
 998         printf("or:    %s <infile_height.tga> <outfile_normandheight.tga> filtertype -2 [<scale>] (read from G)\n", me);
 999         printf("or:    %s <infile_height.tga> <outfile_normandheight.tga> filtertype -3 [<scale>] (read from R)\n", me);
1000         printf("or:    %s <infile_height.tga> <outfile_normandheight.tga> filtertype -4 [<scale>] (read from A)\n", me);
1001         printf("or:    %s <infile_height.tga> <outfile_normandheight.tga> filtertype -5 [<scale>] (read from (R+G+B)/3)\n", me);
1002         printf("or:    %s <infile_height.tga> <outfile_normandheight.tga> filtertype -6 [<scale>] (read from Y)\n", me);
1003         return 1;
1004 }
1005
1006 static const double filter_scharr3[3][3] = {
1007         {  -3/32.0, 0,  3/32.0 },
1008         { -10/32.0, 0, 10/32.0 },
1009         {  -3/32.0, 0,  3/32.0 }
1010 };
1011
1012 static const double filter_prewitt3[3][3] = {
1013         { -1/6.0, 0, 1/6.0 },
1014         { -1/6.0, 0, 1/6.0 },
1015         { -1/6.0, 0, 1/6.0 }
1016 };
1017
1018 // pathologic for inverting
1019 static const double filter_sobel3[3][3] = {
1020         { -1/8.0, 0, 1/8.0 },
1021         { -2/8.0, 0, 2/8.0 },
1022         { -1/8.0, 0, 1/8.0 }
1023 };
1024
1025 // pathologic for inverting
1026 static const double filter_sobel5[5][5] = {
1027         { -1/128.0,  -2/128.0, 0,  2/128.0, 1/128.0 },
1028         { -4/128.0,  -8/128.0, 0,  8/128.0, 4/128.0 },
1029         { -6/128.0, -12/128.0, 0, 12/128.0, 6/128.0 },
1030         { -4/128.0,  -8/128.0, 0,  8/128.0, 4/128.0 },
1031         { -1/128.0,  -2/128.0, 0,  2/128.0, 1/128.0 }
1032 };
1033
1034 // pathologic for inverting
1035 static const double filter_prewitt5[5][5] = {
1036         { -1/40.0, -2/40.0, 0, 2/40.0, 1/40.0 },
1037         { -1/40.0, -2/40.0, 0, 2/40.0, 1/40.0 },
1038         { -1/40.0, -2/40.0, 0, 2/40.0, 1/40.0 },
1039         { -1/40.0, -2/40.0, 0, 2/40.0, 1/40.0 },
1040         { -1/40.0, -2/40.0, 0, 2/40.0, 1/40.0 }
1041 };
1042
1043 static const double filter_trivial[1][3] = {
1044         { -0.5, 0, 0.5 }
1045 };
1046
1047 int main(int argc, char **argv)
1048 {
1049         const char *infile, *outfile, *reffile;
1050         double scale, offset;
1051         int nmaplen, w, h;
1052         int renormalize = 0;
1053         double highpass = 0;
1054         unsigned char *nmapdata, *nmap, *refmap;
1055         const char *filtertype;
1056         const double *filter = NULL;
1057         int filterw = 0, filterh = 0;
1058 #define USE_FILTER(f) \
1059         do \
1060         { \
1061                 filterw = sizeof(*(f)) / sizeof(**(f)); \
1062                 filterh = sizeof((f)) / sizeof(*(f)); \
1063                 filter = &(f)[0][0]; \
1064         } \
1065         while(0)
1066
1067         if(argc > 1)
1068                 infile = argv[1];
1069         else
1070                 return usage(*argv);
1071
1072         if(argc > 2)
1073                 outfile = argv[2];
1074         else
1075                 return usage(*argv);
1076
1077         if(argc > 3)
1078                 filtertype = argv[3];
1079         else
1080                 return usage(*argv);
1081
1082         if(argc > 4)
1083                 scale = atof(argv[4]);
1084         else
1085                 scale = 0;
1086
1087         if(argc > 5)
1088                 offset = atof(argv[5]);
1089         else
1090                 offset = (scale<0) ? 1 : 0;
1091
1092         if(argc > 6)
1093                 reffile = argv[6];
1094         else
1095                 reffile = NULL;
1096
1097         if(getenv("FFT_NORMALMAP_TO_HEIGHTMAP_RENORMALIZE"))
1098                 renormalize = atoi(getenv("FFT_NORMALMAP_TO_HEIGHTMAP_RENORMALIZE"));
1099         if(getenv("FFT_NORMALMAP_TO_HEIGHTMAP_HIGHPASS"))
1100                 highpass = atof(getenv("FFT_NORMALMAP_TO_HEIGHTMAP_HIGHPASS"));
1101
1102         nmapdata = FS_LoadFile(infile, &nmaplen);
1103         if(!nmapdata)
1104         {
1105                 printf("FS_LoadFile failed\n");
1106                 return 2;
1107         }
1108         nmap = LoadTGA_BGRA(nmapdata, nmaplen);
1109         free(nmapdata);
1110         if(!nmap)
1111         {
1112                 printf("LoadTGA_BGRA failed\n");
1113                 return 2;
1114         }
1115         w = image_width;
1116         h = image_height;
1117
1118         if(reffile)
1119         {
1120                 nmapdata = FS_LoadFile(reffile, &nmaplen);
1121                 if(!nmapdata)
1122                 {
1123                         printf("FS_LoadFile failed\n");
1124                         return 2;
1125                 }
1126                 refmap = LoadTGA_BGRA(nmapdata, nmaplen);
1127                 free(nmapdata);
1128                 if(!refmap)
1129                 {
1130                         printf("LoadTGA_BGRA failed\n");
1131                         return 2;
1132                 }
1133                 if(image_width != w || image_height != h)
1134                 {
1135                         printf("reference map must have same size as input normalmap\n");
1136                         return 2;
1137                 }
1138         }
1139         else
1140                 refmap = NULL;
1141
1142         if(!strcmp(filtertype, "trivial"))
1143                 USE_FILTER(filter_trivial);
1144         if(!strcmp(filtertype, "prewitt3"))
1145                 USE_FILTER(filter_prewitt3);
1146         if(!strcmp(filtertype, "scharr3"))
1147                 USE_FILTER(filter_scharr3);
1148         if(!strcmp(filtertype, "sobel3"))
1149                 USE_FILTER(filter_sobel3);
1150         if(!strcmp(filtertype, "prewitt5"))
1151                 USE_FILTER(filter_prewitt5);
1152         if(!strcmp(filtertype, "sobel5"))
1153                 USE_FILTER(filter_sobel5);
1154
1155         if(scale < 0)
1156         {
1157                 if(filter)
1158                         hmap_to_nmap_local(nmap, image_width, image_height, -scale-1, offset, filter, filterw, filterh);
1159                 else
1160                         hmap_to_nmap(nmap, image_width, image_height, -scale-1, offset);
1161         }
1162         else
1163                 nmap_to_hmap(nmap, refmap, image_width, image_height, scale, offset, filter, filterw, filterh, renormalize, highpass);
1164
1165         if(!Image_WriteTGABGRA(outfile, image_width, image_height, nmap))
1166         {
1167                 printf("Image_WriteTGABGRA failed\n");
1168                 free(nmap);
1169                 return 2;
1170         }
1171         free(nmap);
1172         return 0;
1173 }