git.cworth.org Git - vogl/blob - src/voglcore/vogl_image_utils.cpp

   1 /**************************************************************************
   2  *
   3  * Copyright 2013-2014 RAD Game Tools and Valve Software
   4  * Copyright 2010-2014 Rich Geldreich and Tenacious Software LLC
   5  * All Rights Reserved.
   6  *
   7  * Permission is hereby granted, free of charge, to any person obtaining a copy
   8  * of this software and associated documentation files (the "Software"), to deal
   9  * in the Software without restriction, including without limitation the rights
  10  * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
  11  * copies of the Software, and to permit persons to whom the Software is
  12  * furnished to do so, subject to the following conditions:
  13  *
  14  * The above copyright notice and this permission notice shall be included in
  15  * all copies or substantial portions of the Software.
  16  *
  17  * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
  18  * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
  19  * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
  20  * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
  21  * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
  22  * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
  23  * THE SOFTWARE.
  24  *
  25  **************************************************************************/
  26
  27 // File: vogl_image_utils.cpp
  28 #include "vogl_core.h"
  29 #include "vogl_image_utils.h"
  30 #include "vogl_console.h"
  31 #include "vogl_resampler.h"
  32 #include "vogl_threaded_resampler.h"
  33 #include "vogl_strutils.h"
  34 #include "vogl_file_utils.h"
  35 #include "vogl_threading.h"
  36 #include "vogl_miniz.h"
  37 #include "vogl_jpge.h"
  38 #include "vogl_cfile_stream.h"
  39 #include "vogl_mipmapped_texture.h"
  40 #include "vogl_buffer_stream.h"
  41 #include "vogl_vector2d.h"
  42
  43 #include "vogl_stb_image.h"
  44
  45 #include "vogl_jpgd.h"
  46
  47 #include "vogl_pixel_format.h"
  48
  49 namespace vogl
  50 {
  51     const float cInfinitePSNR = 999999.0f;
  52     const uint VOGL_LARGEST_SUPPORTED_IMAGE_DIMENSION = 16384;
  53
  54     namespace image_utils
  55     {
  56         bool read_from_stream_stb(data_stream_serializer &serializer, image_u8 &img)
  57         {
  58             uint8_vec buf;
  59             if (!serializer.read_entire_file(buf))
  60                 return false;
  61
  62             int x = 0, y = 0, n = 0;
  63             unsigned char *pData = stbi_load_from_memory(buf.get_ptr(), buf.size_in_bytes(), &x, &y, &n, 4);
  64
  65             if (!pData)
  66                 return false;
  67
  68             if ((x > (int)VOGL_LARGEST_SUPPORTED_IMAGE_DIMENSION) || (y > (int)VOGL_LARGEST_SUPPORTED_IMAGE_DIMENSION))
  69             {
  70                 stbi_image_free(pData);
  71                 return false;
  72             }
  73
  74             const bool has_alpha = ((n == 2) || (n == 4));
  75
  76             img.crop(x, y);
  77
  78             bool grayscale = true;
  79
  80             for (int py = 0; py < y; py++)
  81             {
  82                 const color_quad_u8 *pSrc = reinterpret_cast<const color_quad_u8 *>(pData) + (py * x);
  83                 color_quad_u8 *pDst = img.get_scanline(py);
  84                 color_quad_u8 *pDst_end = pDst + x;
  85
  86                 while (pDst != pDst_end)
  87                 {
  88                     color_quad_u8 c(*pSrc++);
  89                     if (!has_alpha)
  90                         c.a = 255;
  91
  92                     if (!c.is_grayscale())
  93                         grayscale = false;
  94
  95                     *pDst++ = c;
  96                 }
  97             }
  98
  99             stbi_image_free(pData);
 100
 101             img.reset_comp_flags();
 102             img.set_grayscale(grayscale);
 103             img.set_component_valid(3, has_alpha);
 104
 105             return true;
 106         }
 107
 108         bool read_from_stream_jpgd(data_stream_serializer &serializer, image_u8 &img)
 109         {
 110             uint8_vec buf;
 111             if (!serializer.read_entire_file(buf))
 112                 return false;
 113
 114             int width = 0, height = 0, actual_comps = 0;
 115             unsigned char *pSrc_img = jpgd::decompress_jpeg_image_from_memory(buf.get_ptr(), buf.size_in_bytes(), &width, &height, &actual_comps, 4);
 116             if (!pSrc_img)
 117                 return false;
 118
 119             if (math::maximum(width, height) > (int)VOGL_LARGEST_SUPPORTED_IMAGE_DIMENSION)
 120             {
 121                 vogl_free(pSrc_img);
 122                 return false;
 123             }
 124
 125             if (!img.grant_ownership(reinterpret_cast<color_quad_u8 *>(pSrc_img), width, height))
 126             {
 127                 vogl_free(pSrc_img);
 128                 return false;
 129             }
 130
 131             img.reset_comp_flags();
 132             img.set_grayscale(actual_comps == 1);
 133             img.set_component_valid(3, false);
 134
 135             return true;
 136         }
 137
 138         bool read_from_stream(image_u8 &dest, data_stream_serializer &serializer, uint read_flags)
 139         {
 140             if (read_flags > cReadFlagsAllFlags)
 141             {
 142                 VOGL_ASSERT_ALWAYS;
 143                 return false;
 144             }
 145
 146             if (!serializer.get_stream())
 147             {
 148                 VOGL_ASSERT_ALWAYS;
 149                 return false;
 150             }
 151
 152             dynamic_string ext(serializer.get_name());
 153             file_utils::get_extension(ext);
 154
 155             if ((ext == "jpg") || (ext == "jpeg"))
 156             {
 157                 // Use my jpeg decoder by default because it supports progressive jpeg's.
 158                 if ((read_flags & cReadFlagForceSTB) == 0)
 159                 {
 160                     return image_utils::read_from_stream_jpgd(serializer, dest);
 161                 }
 162             }
 163
 164             return image_utils::read_from_stream_stb(serializer, dest);
 165         }
 166
 167         bool read_from_file(image_u8 &dest, const char *pFilename, uint read_flags)
 168         {
 169             if (read_flags > cReadFlagsAllFlags)
 170             {
 171                 VOGL_ASSERT_ALWAYS;
 172                 return false;
 173             }
 174
 175             cfile_stream file_stream;
 176             if (!file_stream.open(pFilename))
 177                 return false;
 178
 179             data_stream_serializer serializer(file_stream);
 180             return read_from_stream(dest, serializer, read_flags);
 181         }
 182
 183         bool write_to_file(const char *pFilename, const image_u8 &img, uint write_flags, int grayscale_comp_index)
 184         {
 185             if ((grayscale_comp_index < -1) || (grayscale_comp_index > 3))
 186             {
 187                 VOGL_ASSERT_ALWAYS;
 188                 return false;
 189             }
 190
 191             if (!img.get_width())
 192             {
 193                 VOGL_ASSERT_ALWAYS;
 194                 return false;
 195             }
 196
 197             dynamic_string ext(pFilename);
 198             bool is_jpeg = false;
 199             if (file_utils::get_extension(ext))
 200             {
 201                 is_jpeg = ((ext == "jpg") || (ext == "jpeg"));
 202
 203                 if ((ext != "png") && (ext != "bmp") && (ext != "tga") && (!is_jpeg))
 204                 {
 205                     console::error("vogl::image_utils::write_to_file: Can only write .BMP, .TGA, .PNG, or .JPG files!\n");
 206                     return false;
 207                 }
 208             }
 209
 210             vogl::vector<uint8> temp;
 211             uint num_src_chans = 0;
 212             const void *pSrc_img = NULL;
 213
 214             if (is_jpeg)
 215             {
 216                 write_flags |= cWriteFlagIgnoreAlpha;
 217             }
 218
 219             if ((img.get_comp_flags() & pixel_format_helpers::cCompFlagGrayscale) || (write_flags & image_utils::cWriteFlagGrayscale))
 220             {
 221                 VOGL_ASSERT(grayscale_comp_index < 4);
 222                 if (grayscale_comp_index > 3)
 223                     grayscale_comp_index = 3;
 224
 225                 temp.resize(img.get_total_pixels());
 226
 227                 for (uint y = 0; y < img.get_height(); y++)
 228                 {
 229                     const color_quad_u8 *pSrc = img.get_scanline(y);
 230                     const color_quad_u8 *pSrc_end = pSrc + img.get_width();
 231                     uint8 *pDst = &temp[y * img.get_width()];
 232
 233                     if (img.get_comp_flags() & pixel_format_helpers::cCompFlagGrayscale)
 234                     {
 235                         while (pSrc != pSrc_end)
 236                             *pDst++ = (*pSrc++)[1];
 237                     }
 238                     else if (grayscale_comp_index < 0)
 239                     {
 240                         while (pSrc != pSrc_end)
 241                             *pDst++ = static_cast<uint8>((*pSrc++).get_luma());
 242                     }
 243                     else
 244                     {
 245                         while (pSrc != pSrc_end)
 246                             *pDst++ = (*pSrc++)[grayscale_comp_index];
 247                     }
 248                 }
 249
 250                 pSrc_img = &temp[0];
 251                 num_src_chans = 1;
 252             }
 253             else if ((!img.is_component_valid(3)) || (write_flags & cWriteFlagIgnoreAlpha))
 254             {
 255                 temp.resize(img.get_total_pixels() * 3);
 256
 257                 for (uint y = 0; y < img.get_height(); y++)
 258                 {
 259                     const color_quad_u8 *pSrc = img.get_scanline(y);
 260                     const color_quad_u8 *pSrc_end = pSrc + img.get_width();
 261                     uint8 *pDst = &temp[y * img.get_width() * 3];
 262
 263                     while (pSrc != pSrc_end)
 264                     {
 265                         const color_quad_u8 c(*pSrc++);
 266
 267                         pDst[0] = c.r;
 268                         pDst[1] = c.g;
 269                         pDst[2] = c.b;
 270
 271                         pDst += 3;
 272                     }
 273                 }
 274
 275                 num_src_chans = 3;
 276                 pSrc_img = &temp[0];
 277             }
 278             else
 279             {
 280                 num_src_chans = 4;
 281                 pSrc_img = img.get_ptr();
 282             }
 283
 284             bool success = false;
 285             if (ext == "png")
 286             {
 287                 size_t png_image_size = 0;
 288                 void *pPNG_image_data = tdefl_write_image_to_png_file_in_memory_ex(pSrc_img, img.get_width(), img.get_height(), num_src_chans, &png_image_size, MZ_DEFAULT_LEVEL, MZ_FALSE);
 289                 if (!pPNG_image_data)
 290                     return false;
 291                 success = file_utils::write_buf_to_file(pFilename, pPNG_image_data, png_image_size);
 292                 mz_free(pPNG_image_data);
 293             }
 294             else if (is_jpeg)
 295             {
 296                 jpge::params params;
 297                 if (write_flags & cWriteFlagJPEGQualityLevelMask)
 298                     params.m_quality = math::clamp<uint>((write_flags & cWriteFlagJPEGQualityLevelMask) >> cWriteFlagJPEGQualityLevelShift, 1U, 100U);
 299                 params.m_two_pass_flag = (write_flags & cWriteFlagJPEGTwoPass) != 0;
 300                 params.m_no_chroma_discrim_flag = (write_flags & cWriteFlagJPEGNoChromaDiscrim) != 0;
 301
 302                 if (write_flags & cWriteFlagJPEGH1V1)
 303                     params.m_subsampling = jpge::H1V1;
 304                 else if (write_flags & cWriteFlagJPEGH2V1)
 305                     params.m_subsampling = jpge::H2V1;
 306                 else if (write_flags & cWriteFlagJPEGH2V2)
 307                     params.m_subsampling = jpge::H2V2;
 308
 309                 success = jpge::compress_image_to_jpeg_file(pFilename, img.get_width(), img.get_height(), num_src_chans, (const jpge::uint8 *)pSrc_img, params);
 310             }
 311             else
 312             {
 313                 success = ((ext == "bmp" ? stbi_write_bmp : stbi_write_tga)(pFilename, img.get_width(), img.get_height(), num_src_chans, pSrc_img) == VOGL_TRUE);
 314             }
 315             return success;
 316         }
 317
 318         bool write_to_file(const char *pFilename, const image_f &img, const vec4F &l, const vec4F &h, uint write_flags, int grayscale_comp_index)
 319         {
 320             image_u8 temp(img.get_width(), img.get_height());
 321
 322             vec4F scale(h - l);
 323             for (uint i = 0; i < 4; i++)
 324                 scale[i] = (scale[i] != 0.0f) ? (255.0f / scale[i]) : 0;
 325
 326             for (uint y = 0; y < img.get_height(); y++)
 327             {
 328                 for (uint x = 0; x < img.get_width(); x++)
 329                 {
 330                     const color_quad_f &o = img(x, y);
 331
 332                     color_quad_u8 c;
 333                     for (uint i = 0; i < 4; i++)
 334                         c[i] = static_cast<uint8>(math::clamp(math::float_to_int_nearest((o[i] - l[i]) * scale[i]), 0, 255));
 335
 336                     temp(x, y) = c;
 337                 }
 338             }
 339
 340             return write_to_file(pFilename, temp, write_flags, grayscale_comp_index);
 341         }
 342
 343         bool has_alpha(const image_u8 &img)
 344         {
 345             for (uint y = 0; y < img.get_height(); y++)
 346                 for (uint x = 0; x < img.get_width(); x++)
 347                     if (img(x, y).a < 255)
 348                         return true;
 349
 350             return false;
 351         }
 352
 353         void renorm_normal_map(image_u8 &img)
 354         {
 355             for (uint y = 0; y < img.get_height(); y++)
 356             {
 357                 for (uint x = 0; x < img.get_width(); x++)
 358                 {
 359                     color_quad_u8 &c = img(x, y);
 360                     if ((c.r == 128) && (c.g == 128) && (c.b == 128))
 361                         continue;
 362
 363                     vec3F v(c.r, c.g, c.b);
 364                     v *= 1.0f / 255.0f;
 365                     v *= 2.0f;
 366                     v -= vec3F(1.0f);
 367                     v.clamp(-1.0f, 1.0f);
 368
 369                     float length = v.length();
 370                     if (length < .077f)
 371                         c.set(128, 128, 128, c.a);
 372                     else if (fabs(length - 1.0f) > .077f)
 373                     {
 374                         if (length)
 375                             v /= length;
 376
 377                         for (uint i = 0; i < 3; i++)
 378                             c[i] = static_cast<uint8>(math::clamp<float>(floor((v[i] + 1.0f) * .5f * 255.0f + .5f), 0.0f, 255.0f));
 379
 380                         if ((c.r == 128) && (c.g == 128))
 381                         {
 382                             if (c.b < 128)
 383                                 c.b = 0;
 384                             else
 385                                 c.b = 255;
 386                         }
 387                     }
 388                 }
 389             }
 390         }
 391
 392         bool is_normal_map(const image_u8 &img, const char *pFilename)
 393         {
 394             float score = 0.0f;
 395
 396             uint num_invalid_pixels = 0;
 397
 398             // TODO: Derive better score from pixel mean, eigenvecs/vals
 399             //vogl::vector<vec3F> pixels;
 400
 401             for (uint y = 0; y < img.get_height(); y++)
 402             {
 403                 for (uint x = 0; x < img.get_width(); x++)
 404                 {
 405                     const color_quad_u8 &c = img(x, y);
 406
 407                     if (c.b < 123)
 408                     {
 409                         num_invalid_pixels++;
 410                         continue;
 411                     }
 412                     else if ((c.r != 128) || (c.g != 128) || (c.b != 128))
 413                     {
 414                         vec3F v(c.r, c.g, c.b);
 415                         v -= vec3F(128.0f);
 416                         v /= vec3F(127.0f);
 417                         //pixels.push_back(v);
 418                         v.clamp(-1.0f, 1.0f);
 419
 420                         float norm = v.norm();
 421                         if ((norm < 0.83f) || (norm > 1.29f))
 422                             num_invalid_pixels++;
 423                     }
 424                 }
 425             }
 426
 427             score -= math::clamp(float(num_invalid_pixels) / (img.get_width() * img.get_height()) - .026f, 0.0f, 1.0f) * 5.0f;
 428
 429             if (pFilename)
 430             {
 431                 dynamic_string str(pFilename);
 432                 str.tolower();
 433
 434                 if (str.contains("normal") || str.contains("local") || str.contains("nmap"))
 435                     score += 1.0f;
 436
 437                 if (str.contains("diffuse") || str.contains("spec") || str.contains("gloss"))
 438                     score -= 1.0f;
 439             }
 440
 441             return score >= 0.0f;
 442         }
 443
 444         bool resample_single_thread(const image_u8 &src, image_u8 &dst, const resample_params &params)
 445         {
 446             const uint src_width = src.get_width();
 447             const uint src_height = src.get_height();
 448
 449             if (math::maximum(src_width, src_height) > VOGL_RESAMPLER_MAX_DIMENSION)
 450             {
 451                 printf("Image is too large!\n");
 452                 return EXIT_FAILURE;
 453             }
 454
 455             const int cMaxComponents = 4;
 456             if (((int)params.m_num_comps < 1) || ((int)params.m_num_comps > (int)cMaxComponents))
 457                 return false;
 458
 459             const uint dst_width = params.m_dst_width;
 460             const uint dst_height = params.m_dst_height;
 461
 462             if ((math::minimum(dst_width, dst_height) < 1) || (math::maximum(dst_width, dst_height) > VOGL_RESAMPLER_MAX_DIMENSION))
 463             {
 464                 printf("Image is too large!\n");
 465                 return EXIT_FAILURE;
 466             }
 467
 468             if ((src_width == dst_width) && (src_height == dst_height))
 469             {
 470                 dst = src;
 471                 return true;
 472             }
 473
 474             dst.clear();
 475             dst.crop(params.m_dst_width, params.m_dst_height);
 476
 477             // Partial gamma correction looks better on mips. Set to 1.0 to disable gamma correction.
 478             const float source_gamma = params.m_source_gamma; //1.75f;
 479
 480             float srgb_to_linear[256];
 481             if (params.m_srgb)
 482             {
 483                 for (int i = 0; i < 256; ++i)
 484                     srgb_to_linear[i] = (float)pow(i * 1.0f / 255.0f, source_gamma);
 485             }
 486
 487             const int linear_to_srgb_table_size = 8192;
 488             unsigned char linear_to_srgb[linear_to_srgb_table_size];
 489
 490             const float inv_linear_to_srgb_table_size = 1.0f / linear_to_srgb_table_size;
 491             const float inv_source_gamma = 1.0f / source_gamma;
 492
 493             if (params.m_srgb)
 494             {
 495                 for (int i = 0; i < linear_to_srgb_table_size; ++i)
 496                 {
 497                     int k = (int)(255.0f * pow(i * inv_linear_to_srgb_table_size, inv_source_gamma) + .5f);
 498                     if (k < 0)
 499                         k = 0;
 500                     else if (k > 255)
 501                         k = 255;
 502                     linear_to_srgb[i] = (unsigned char)k;
 503                 }
 504             }
 505
 506             Resampler *resamplers[cMaxComponents];
 507             vogl::vector<float> samples[cMaxComponents];
 508
 509             resamplers[0] = vogl_new(Resampler, src_width, src_height, dst_width, dst_height,
 510                                        params.m_wrapping ? Resampler::BOUNDARY_WRAP : Resampler::BOUNDARY_CLAMP, 0.0f, 1.0f,
 511                                        params.m_pFilter, (Resampler::Contrib_List *)NULL, (Resampler::Contrib_List *)NULL, params.m_filter_scale, params.m_filter_scale, params.m_x_ofs, params.m_y_ofs);
 512             samples[0].resize(src_width);
 513
 514             for (uint i = 1; i < params.m_num_comps; i++)
 515             {
 516                 resamplers[i] = vogl_new(Resampler, src_width, src_height, dst_width, dst_height,
 517                                            params.m_wrapping ? Resampler::BOUNDARY_WRAP : Resampler::BOUNDARY_CLAMP, 0.0f, 1.0f,
 518                                            params.m_pFilter, resamplers[0]->get_clist_x(), resamplers[0]->get_clist_y(), params.m_filter_scale, params.m_filter_scale, params.m_x_ofs, params.m_y_ofs);
 519                 samples[i].resize(src_width);
 520             }
 521
 522             uint dst_y = 0;
 523
 524             for (uint src_y = 0; src_y < src_height; src_y++)
 525             {
 526                 const color_quad_u8 *pSrc = src.get_scanline(src_y);
 527
 528                 for (uint x = 0; x < src_width; x++)
 529                 {
 530                     for (uint c = 0; c < params.m_num_comps; c++)
 531                     {
 532                         const uint comp_index = params.m_first_comp + c;
 533                         const uint8 v = (*pSrc)[comp_index];
 534
 535                         if (!params.m_srgb || (comp_index == 3))
 536                             samples[c][x] = v * (1.0f / 255.0f);
 537                         else
 538                             samples[c][x] = srgb_to_linear[v];
 539                     }
 540
 541                     pSrc++;
 542                 }
 543
 544                 for (uint c = 0; c < params.m_num_comps; c++)
 545                 {
 546                     if (!resamplers[c]->put_line(&samples[c][0]))
 547                     {
 548                         for (uint i = 0; i < params.m_num_comps; i++)
 549                             vogl_delete(resamplers[i]);
 550                         return false;
 551                     }
 552                 }
 553
 554                 for (;;)
 555                 {
 556                     uint c;
 557                     for (c = 0; c < params.m_num_comps; c++)
 558                     {
 559                         const uint comp_index = params.m_first_comp + c;
 560
 561                         const float *pOutput_samples = resamplers[c]->get_line();
 562                         if (!pOutput_samples)
 563                             break;
 564
 565                         const bool linear = !params.m_srgb || (comp_index == 3);
 566                         VOGL_ASSERT(dst_y < dst_height);
 567                         color_quad_u8 *pDst = dst.get_scanline(dst_y);
 568
 569                         for (uint x = 0; x < dst_width; x++)
 570                         {
 571                             if (linear)
 572                             {
 573                                 int c = (int)(255.0f * pOutput_samples[x] + .5f);
 574                                 if (c < 0)
 575                                     c = 0;
 576                                 else if (c > 255)
 577                                     c = 255;
 578                                 (*pDst)[comp_index] = (unsigned char)c;
 579                             }
 580                             else
 581                             {
 582                                 int j = (int)(linear_to_srgb_table_size * pOutput_samples[x] + .5f);
 583                                 if (j < 0)
 584                                     j = 0;
 585                                 else if (j >= linear_to_srgb_table_size)
 586                                     j = linear_to_srgb_table_size - 1;
 587                                 (*pDst)[comp_index] = linear_to_srgb[j];
 588                             }
 589
 590                             pDst++;
 591                         }
 592                     }
 593                     if (c < params.m_num_comps)
 594                         break;
 595
 596                     dst_y++;
 597                 }
 598             }
 599
 600             for (uint i = 0; i < params.m_num_comps; i++)
 601                 vogl_delete(resamplers[i]);
 602
 603             return true;
 604         }
 605
 606         bool resample(const image_f &src, image_f &dst, const resample_params &params)
 607         {
 608             const uint src_width = src.get_width();
 609             const uint src_height = src.get_height();
 610
 611             if (math::maximum(src_width, src_height) > VOGL_RESAMPLER_MAX_DIMENSION)
 612             {
 613                 printf("Image is too large!\n");
 614                 return EXIT_FAILURE;
 615             }
 616
 617             const int cMaxComponents = 4;
 618             if (((int)params.m_num_comps < 1) || ((int)params.m_num_comps > (int)cMaxComponents))
 619                 return false;
 620
 621             const uint dst_width = params.m_dst_width;
 622             const uint dst_height = params.m_dst_height;
 623
 624             if ((math::minimum(dst_width, dst_height) < 1) || (math::maximum(dst_width, dst_height) > VOGL_RESAMPLER_MAX_DIMENSION))
 625             {
 626                 printf("Image is too large!\n");
 627                 return EXIT_FAILURE;
 628             }
 629
 630             if ((src_width == dst_width) && (src_height == dst_height))
 631             {
 632                 dst = src;
 633                 return true;
 634             }
 635
 636             dst.clear();
 637             dst.crop(params.m_dst_width, params.m_dst_height);
 638
 639             Resampler *resamplers[cMaxComponents];
 640             vogl::vector<float> samples[cMaxComponents];
 641
 642             resamplers[0] = vogl_new(Resampler, src_width, src_height, dst_width, dst_height,
 643                                        params.m_wrapping ? Resampler::BOUNDARY_WRAP : Resampler::BOUNDARY_CLAMP, 0.0f, 1.0f,
 644                                        params.m_pFilter, (Resampler::Contrib_List *)NULL, (Resampler::Contrib_List *)NULL, params.m_filter_scale, params.m_filter_scale, params.m_x_ofs, params.m_y_ofs);
 645             samples[0].resize(src_width);
 646
 647             for (uint i = 1; i < params.m_num_comps; i++)
 648             {
 649                 resamplers[i] = vogl_new(Resampler, src_width, src_height, dst_width, dst_height,
 650                                            params.m_wrapping ? Resampler::BOUNDARY_WRAP : Resampler::BOUNDARY_CLAMP, 0.0f, 1.0f,
 651                                            params.m_pFilter, resamplers[0]->get_clist_x(), resamplers[0]->get_clist_y(), params.m_filter_scale, params.m_filter_scale, params.m_x_ofs, params.m_y_ofs);
 652                 samples[i].resize(src_width);
 653             }
 654
 655             uint dst_y = 0;
 656
 657             for (uint src_y = 0; src_y < src_height; src_y++)
 658             {
 659                 const color_quad_f *pSrc = src.get_scanline(src_y);
 660
 661                 for (uint x = 0; x < src_width; x++)
 662                 {
 663                     for (uint c = 0; c < params.m_num_comps; c++)
 664                     {
 665                         const uint comp_index = params.m_first_comp + c;
 666                         samples[c][x] = (*pSrc)[comp_index];
 667                     }
 668
 669                     pSrc++;
 670                 }
 671
 672                 for (uint c = 0; c < params.m_num_comps; c++)
 673                 {
 674                     if (!resamplers[c]->put_line(&samples[c][0]))
 675                     {
 676                         for (uint i = 0; i < params.m_num_comps; i++)
 677                             vogl_delete(resamplers[i]);
 678                         return false;
 679                     }
 680                 }
 681
 682                 for (;;)
 683                 {
 684                     uint c;
 685                     for (c = 0; c < params.m_num_comps; c++)
 686                     {
 687                         const uint comp_index = params.m_first_comp + c;
 688
 689                         const float *pOutput_samples = resamplers[c]->get_line();
 690                         if (!pOutput_samples)
 691                             break;
 692
 693                         VOGL_ASSERT(dst_y < dst_height);
 694                         color_quad_f *pDst = dst.get_scanline(dst_y);
 695
 696                         for (uint x = 0; x < dst_width; x++)
 697                         {
 698                             (*pDst)[comp_index] = pOutput_samples[x];
 699                             pDst++;
 700                         }
 701                     }
 702
 703                     if (c < params.m_num_comps)
 704                         break;
 705
 706                     dst_y++;
 707                 }
 708             }
 709
 710             for (uint i = 0; i < params.m_num_comps; i++)
 711                 vogl_delete(resamplers[i]);
 712
 713             return true;
 714         }
 715
 716         bool resample_multithreaded(const image_u8 &src, image_u8 &dst, const resample_params &params)
 717         {
 718             const uint src_width = src.get_width();
 719             const uint src_height = src.get_height();
 720
 721             if (math::maximum(src_width, src_height) > VOGL_RESAMPLER_MAX_DIMENSION)
 722             {
 723                 printf("Image is too large!\n");
 724                 return EXIT_FAILURE;
 725             }
 726
 727             const int cMaxComponents = 4;
 728             if (((int)params.m_num_comps < 1) || ((int)params.m_num_comps > (int)cMaxComponents))
 729                 return false;
 730
 731             const uint dst_width = params.m_dst_width;
 732             const uint dst_height = params.m_dst_height;
 733
 734             if ((math::minimum(dst_width, dst_height) < 1) || (math::maximum(dst_width, dst_height) > VOGL_RESAMPLER_MAX_DIMENSION))
 735             {
 736                 printf("Image is too large!\n");
 737                 return EXIT_FAILURE;
 738             }
 739
 740             if ((src_width == dst_width) && (src_height == dst_height))
 741             {
 742                 dst = src;
 743                 return true;
 744             }
 745
 746             dst.clear();
 747
 748             // Partial gamma correction looks better on mips. Set to 1.0 to disable gamma correction.
 749             const float source_gamma = params.m_source_gamma; //1.75f;
 750
 751             float srgb_to_linear[256];
 752             if (params.m_srgb)
 753             {
 754                 for (int i = 0; i < 256; ++i)
 755                     srgb_to_linear[i] = (float)pow(i * 1.0f / 255.0f, source_gamma);
 756             }
 757
 758             const int linear_to_srgb_table_size = 8192;
 759             unsigned char linear_to_srgb[linear_to_srgb_table_size];
 760
 761             const float inv_linear_to_srgb_table_size = 1.0f / linear_to_srgb_table_size;
 762             const float inv_source_gamma = 1.0f / source_gamma;
 763
 764             if (params.m_srgb)
 765             {
 766                 for (int i = 0; i < linear_to_srgb_table_size; ++i)
 767                 {
 768                     int k = (int)(255.0f * pow(i * inv_linear_to_srgb_table_size, inv_source_gamma) + .5f);
 769                     if (k < 0)
 770                         k = 0;
 771                     else if (k > 255)
 772                         k = 255;
 773                     linear_to_srgb[i] = (unsigned char)k;
 774                 }
 775             }
 776
 777             task_pool tp;
 778             tp.init(g_number_of_processors - 1);
 779
 780             threaded_resampler resampler(tp);
 781             threaded_resampler::params p;
 782             p.m_src_width = src_width;
 783             p.m_src_height = src_height;
 784             p.m_dst_width = dst_width;
 785             p.m_dst_height = dst_height;
 786             p.m_sample_low = 0.0f;
 787             p.m_sample_high = 1.0f;
 788             p.m_boundary_op = params.m_wrapping ? Resampler::BOUNDARY_WRAP : Resampler::BOUNDARY_CLAMP;
 789             p.m_Pfilter_name = params.m_pFilter;
 790             p.m_filter_x_scale = params.m_filter_scale;
 791             p.m_filter_y_scale = params.m_filter_scale;
 792             p.m_x_ofs = params.m_x_ofs;
 793             p.m_y_ofs = params.m_y_ofs;
 794
 795             uint resampler_comps = 4;
 796             if (params.m_num_comps == 1)
 797             {
 798                 p.m_fmt = threaded_resampler::cPF_Y_F32;
 799                 resampler_comps = 1;
 800             }
 801             else if (params.m_num_comps <= 3)
 802                 p.m_fmt = threaded_resampler::cPF_RGBX_F32;
 803             else
 804                 p.m_fmt = threaded_resampler::cPF_RGBA_F32;
 805
 806             vogl::vector<float> src_samples;
 807             vogl::vector<float> dst_samples;
 808
 809             if (!src_samples.try_resize(src_width * src_height * resampler_comps))
 810                 return false;
 811
 812             if (!dst_samples.try_resize(dst_width * dst_height * resampler_comps))
 813                 return false;
 814
 815             p.m_pSrc_pixels = src_samples.get_ptr();
 816             p.m_src_pitch = src_width * resampler_comps * sizeof(float);
 817             p.m_pDst_pixels = dst_samples.get_ptr();
 818             p.m_dst_pitch = dst_width * resampler_comps * sizeof(float);
 819
 820             for (uint src_y = 0; src_y < src_height; src_y++)
 821             {
 822                 const color_quad_u8 *pSrc = src.get_scanline(src_y);
 823                 float *pDst = src_samples.get_ptr() + src_width * resampler_comps * src_y;
 824
 825                 for (uint x = 0; x < src_width; x++)
 826                 {
 827                     for (uint c = 0; c < params.m_num_comps; c++)
 828                     {
 829                         const uint comp_index = params.m_first_comp + c;
 830                         const uint8 v = (*pSrc)[comp_index];
 831
 832                         if (!params.m_srgb || (comp_index == 3))
 833                             pDst[c] = v * (1.0f / 255.0f);
 834                         else
 835                             pDst[c] = srgb_to_linear[v];
 836                     }
 837
 838                     pSrc++;
 839                     pDst += resampler_comps;
 840                 }
 841             }
 842
 843             if (!resampler.resample(p))
 844                 return false;
 845
 846             src_samples.clear();
 847
 848             if (!dst.crop(params.m_dst_width, params.m_dst_height))
 849                 return false;
 850
 851             for (uint dst_y = 0; dst_y < dst_height; dst_y++)
 852             {
 853                 const float *pSrc = dst_samples.get_ptr() + dst_width * resampler_comps * dst_y;
 854                 color_quad_u8 *pDst = dst.get_scanline(dst_y);
 855
 856                 for (uint x = 0; x < dst_width; x++)
 857                 {
 858                     color_quad_u8 dst(0, 0, 0, 255);
 859
 860                     for (uint c = 0; c < params.m_num_comps; c++)
 861                     {
 862                         const uint comp_index = params.m_first_comp + c;
 863                         const float v = pSrc[c];
 864
 865                         if ((!params.m_srgb) || (comp_index == 3))
 866                         {
 867                             int c = static_cast<int>(255.0f * v + .5f);
 868                             if (c < 0)
 869                                 c = 0;
 870                             else if (c > 255)
 871                                 c = 255;
 872                             dst[comp_index] = (unsigned char)c;
 873                         }
 874                         else
 875                         {
 876                             int j = static_cast<int>(linear_to_srgb_table_size * v + .5f);
 877                             if (j < 0)
 878                                 j = 0;
 879                             else if (j >= linear_to_srgb_table_size)
 880                                 j = linear_to_srgb_table_size - 1;
 881                             dst[comp_index] = linear_to_srgb[j];
 882                         }
 883                     }
 884
 885                     *pDst++ = dst;
 886
 887                     pSrc += resampler_comps;
 888                 }
 889             }
 890
 891             return true;
 892         }
 893
 894         bool resample(const image_u8 &src, image_u8 &dst, const resample_params &params)
 895         {
 896             if ((params.m_multithreaded) && (g_number_of_processors > 1))
 897                 return resample_multithreaded(src, dst, params);
 898             else
 899                 return resample_single_thread(src, dst, params);
 900         }
 901
 902         bool compute_delta(image_u8 &dest, const image_u8 &a, const image_u8 &b, uint scale)
 903         {
 904             if ((a.get_width() != b.get_width()) || (a.get_height() != b.get_height()))
 905                 return false;
 906
 907             dest.crop(a.get_width(), b.get_height());
 908
 909             for (uint y = 0; y < a.get_height(); y++)
 910             {
 911                 for (uint x = 0; x < a.get_width(); x++)
 912                 {
 913                     const color_quad_u8 &ca = a(x, y);
 914                     const color_quad_u8 &cb = b(x, y);
 915
 916                     color_quad_u8 cd;
 917                     for (uint c = 0; c < 4; c++)
 918                     {
 919                         int d = (ca[c] - cb[c]) * scale + 128;
 920                         d = math::clamp(d, 0, 255);
 921                         cd[c] = static_cast<uint8>(d);
 922                     }
 923
 924                     dest(x, y) = cd;
 925                 }
 926             }
 927
 928             return true;
 929         }
 930
 931         // See page 605 of http://www.cns.nyu.edu/pub/eero/wang03-reprint.pdf
 932         // Image Quality Assessment: From Error Visibility to Structural Similarity
 933         double compute_block_ssim(uint t, const uint8 *pX, const uint8 *pY)
 934         {
 935             double ave_x = 0.0f;
 936             double ave_y = 0.0f;
 937             for (uint i = 0; i < t; i++)
 938             {
 939                 ave_x += pX[i];
 940                 ave_y += pY[i];
 941             }
 942
 943             ave_x /= t;
 944             ave_y /= t;
 945
 946             double var_x = 0.0f;
 947             double var_y = 0.0f;
 948             for (uint i = 0; i < t; i++)
 949             {
 950                 var_x += math::square(pX[i] - ave_x);
 951                 var_y += math::square(pY[i] - ave_y);
 952             }
 953
 954             var_x = sqrt(var_x / (t - 1));
 955             var_y = sqrt(var_y / (t - 1));
 956
 957             double covar_xy = 0.0f;
 958             for (uint i = 0; i < t; i++)
 959                 covar_xy += (pX[i] - ave_x) * (pY[i] - ave_y);
 960
 961             covar_xy /= (t - 1);
 962
 963             const double c1 = 6.5025;  //(255*.01)^2
 964             const double c2 = 58.5225; //(255*.03)^2
 965
 966             double n = (2.0f * ave_x * ave_y + c1) * (2.0f * covar_xy + c2);
 967             double d = (ave_x * ave_x + ave_y * ave_y + c1) * (var_x * var_x + var_y * var_y + c2);
 968
 969             return n / d;
 970         }
 971
 972         double compute_ssim(const image_u8 &a, const image_u8 &b, int channel_index)
 973         {
 974             const uint N = 8;
 975             uint8 sx[N * N], sy[N * N];
 976
 977             double total_ssim = 0.0f;
 978             uint total_blocks = 0;
 979
 980             for (uint y = 0; y < a.get_height(); y += N)
 981             {
 982                 for (uint x = 0; x < a.get_width(); x += N)
 983                 {
 984                     for (uint iy = 0; iy < N; iy++)
 985                     {
 986                         for (uint ix = 0; ix < N; ix++)
 987                         {
 988                             if (channel_index < 0)
 989                                 sx[ix + iy * N] = (uint8)a.get_clamped(x + ix, y + iy).get_luma();
 990                             else
 991                                 sx[ix + iy * N] = (uint8)a.get_clamped(x + ix, y + iy)[channel_index];
 992
 993                             if (channel_index < 0)
 994                                 sy[ix + iy * N] = (uint8)b.get_clamped(x + ix, y + iy).get_luma();
 995                             else
 996                                 sy[ix + iy * N] = (uint8)b.get_clamped(x + ix, y + iy)[channel_index];
 997                         }
 998                     }
 999
1000                     double ssim = compute_block_ssim(N * N, sx, sy);
1001                     total_ssim += ssim;
1002                     total_blocks++;
1003
1004                     //uint ssim_c = (uint)math::clamp<double>(ssim * 127.0f + 128.0f, 0, 255);
1005                     //yimg(x / N, y / N).set(ssim_c, ssim_c, ssim_c, 255);
1006                 }
1007             }
1008
1009             if (!total_blocks)
1010                 return 0.0f;
1011
1012             //save_to_file_stb_or_miniz("ssim.tga", yimg, cWriteFlagGrayscale);
1013
1014             return total_ssim / total_blocks;
1015         }
1016
1017         void print_ssim(const image_u8 &src_img, const image_u8 &dst_img)
1018         {
1019             double y_ssim = compute_ssim(src_img, dst_img, -1);
1020             console::printf("Luma MSSIM: %f, Scaled [.8,1]: %f", y_ssim, (y_ssim - .8f) / .2f);
1021
1022             double r_ssim = compute_ssim(src_img, dst_img, 0);
1023             console::printf("   R MSSIM: %f", r_ssim);
1024
1025             double g_ssim = compute_ssim(src_img, dst_img, 1);
1026             console::printf("   G MSSIM: %f", g_ssim);
1027
1028             double b_ssim = compute_ssim(src_img, dst_img, 2);
1029             console::printf("   B MSSIM: %f", b_ssim);
1030         }
1031
1032         void error_metrics::print(const char *pName) const
1033         {
1034             if (mPeakSNR >= cInfinitePSNR)
1035                 console::printf("%s Error: Max: %3u, Mean: %3.3f, MSE: %3.3f, RMSE: %3.3f, PSNR: Infinite", pName, mMax, mMean, mMeanSquared, mRootMeanSquared);
1036             else
1037                 console::printf("%s Error: Max: %3u, Mean: %3.3f, MSE: %3.3f, RMSE: %3.3f, PSNR: %3.3f, SSIM: %1.6f", pName, mMax, mMean, mMeanSquared, mRootMeanSquared, mPeakSNR, mSSIM);
1038         }
1039
1040         bool error_metrics::compute(const image_u8 &a, const image_u8 &b, uint first_channel, uint num_channels, bool average_component_error)
1041         {
1042             //if ( (!a.get_width()) || (!b.get_height()) || (a.get_width() != b.get_width()) || (a.get_height() != b.get_height()) )
1043             //   return false;
1044
1045             const uint width = math::minimum(a.get_width(), b.get_width());
1046             const uint height = math::minimum(a.get_height(), b.get_height());
1047
1048             VOGL_ASSERT((first_channel < 4U) && (first_channel + num_channels <= 4U));
1049
1050             // Histogram approach due to Charles Bloom.
1051             double hist[256];
1052             utils::zero_object(hist);
1053
1054             for (uint y = 0; y < height; y++)
1055             {
1056                 for (uint x = 0; x < width; x++)
1057                 {
1058                     const color_quad_u8 &ca = a(x, y);
1059                     const color_quad_u8 &cb = b(x, y);
1060
1061                     if (!num_channels)
1062                         hist[labs(ca.get_luma() - cb.get_luma())]++;
1063                     else
1064                     {
1065                         for (uint c = 0; c < num_channels; c++)
1066                             hist[labs(ca[first_channel + c] - cb[first_channel + c])]++;
1067                     }
1068                 }
1069             }
1070
1071             mMax = 0;
1072             double sum = 0.0f, sum2 = 0.0f;
1073             for (uint i = 0; i < 256; i++)
1074             {
1075                 if (!hist[i])
1076                     continue;
1077
1078                 mMax = math::maximum(mMax, i);
1079
1080                 double x = i * hist[i];
1081
1082                 sum += x;
1083                 sum2 += i * x;
1084             }
1085
1086             // See http://bmrc.berkeley.edu/courseware/cs294/fall97/assignment/psnr.html
1087             double total_values = width * height;
1088
1089             if (average_component_error)
1090                 total_values *= math::clamp<uint>(num_channels, 1, 4);
1091
1092             mMean = math::clamp<double>(sum / total_values, 0.0f, 255.0f);
1093             mMeanSquared = math::clamp<double>(sum2 / total_values, 0.0f, 255.0f * 255.0f);
1094
1095             mRootMeanSquared = sqrt(mMeanSquared);
1096
1097             if (!mRootMeanSquared)
1098                 mPeakSNR = cInfinitePSNR;
1099             else
1100                 mPeakSNR = math::clamp<double>(log10(255.0f / mRootMeanSquared) * 20.0f, 0.0f, 500.0f);
1101
1102             if (!a.has_rgb() || !b.has_rgb())
1103                 mSSIM = 0;
1104             else
1105                 mSSIM = compute_ssim(a, b, (num_channels != 1) ? -1 : first_channel);
1106
1107             return true;
1108         }
1109
1110         void print_image_metrics(const image_u8 &src_img, const image_u8 &dst_img)
1111         {
1112             if ((!src_img.get_width()) || (!dst_img.get_height()) || (src_img.get_width() != dst_img.get_width()) || (src_img.get_height() != dst_img.get_height()))
1113                 console::printf("print_image_metrics: Image resolutions don't match exactly (%ux%u) vs. (%ux%u)", src_img.get_width(), src_img.get_height(), dst_img.get_width(), dst_img.get_height());
1114
1115             image_utils::error_metrics error_metrics;
1116
1117             if (src_img.has_rgb() || dst_img.has_rgb())
1118             {
1119                 error_metrics.compute(src_img, dst_img, 0, 3, false);
1120                 error_metrics.print("RGB Total  ");
1121
1122                 error_metrics.compute(src_img, dst_img, 0, 3, true);
1123                 error_metrics.print("RGB Average");
1124
1125                 error_metrics.compute(src_img, dst_img, 0, 0);
1126                 error_metrics.print("Luma       ");
1127
1128                 error_metrics.compute(src_img, dst_img, 0, 1);
1129                 error_metrics.print("Red        ");
1130
1131                 error_metrics.compute(src_img, dst_img, 1, 1);
1132                 error_metrics.print("Green      ");
1133
1134                 error_metrics.compute(src_img, dst_img, 2, 1);
1135                 error_metrics.print("Blue       ");
1136             }
1137
1138             if (src_img.has_alpha() || dst_img.has_alpha())
1139             {
1140                 error_metrics.compute(src_img, dst_img, 3, 1);
1141                 error_metrics.print("Alpha      ");
1142             }
1143         }
1144
1145         static uint8 regen_z(uint x, uint y)
1146         {
1147             float vx = math::clamp((x - 128.0f) * 1.0f / 127.0f, -1.0f, 1.0f);
1148             float vy = math::clamp((y - 128.0f) * 1.0f / 127.0f, -1.0f, 1.0f);
1149             float vz = sqrt(math::clamp(1.0f - vx * vx - vy * vy, 0.0f, 1.0f));
1150
1151             vz = vz * 127.0f + 128.0f;
1152
1153             if (vz < 128.0f)
1154                 vz -= .5f;
1155             else
1156                 vz += .5f;
1157
1158             int ib = math::float_to_int(vz);
1159
1160             return static_cast<uint8>(math::clamp(ib, 0, 255));
1161         }
1162
1163         void convert_image(image_u8 &img, image_utils::conversion_type conv_type)
1164         {
1165             switch (conv_type)
1166             {
1167                 case image_utils::cConversion_To_CCxY:
1168                 {
1169                     img.set_comp_flags(static_cast<pixel_format_helpers::component_flags>(pixel_format_helpers::cCompFlagRValid | pixel_format_helpers::cCompFlagGValid | pixel_format_helpers::cCompFlagAValid | pixel_format_helpers::cCompFlagLumaChroma));
1170                     break;
1171                 }
1172                 case image_utils::cConversion_From_CCxY:
1173                 {
1174                     img.set_comp_flags(static_cast<pixel_format_helpers::component_flags>(pixel_format_helpers::cCompFlagRValid | pixel_format_helpers::cCompFlagGValid | pixel_format_helpers::cCompFlagBValid));
1175                     break;
1176                 }
1177                 case image_utils::cConversion_To_xGxR:
1178                 {
1179                     img.set_comp_flags(static_cast<pixel_format_helpers::component_flags>(pixel_format_helpers::cCompFlagGValid | pixel_format_helpers::cCompFlagAValid | pixel_format_helpers::cCompFlagNormalMap));
1180                     break;
1181                 }
1182                 case image_utils::cConversion_From_xGxR:
1183                 {
1184                     img.set_comp_flags(static_cast<pixel_format_helpers::component_flags>(pixel_format_helpers::cCompFlagRValid | pixel_format_helpers::cCompFlagGValid | pixel_format_helpers::cCompFlagBValid | pixel_format_helpers::cCompFlagNormalMap));
1185                     break;
1186                 }
1187                 case image_utils::cConversion_To_xGBR:
1188                 {
1189                     img.set_comp_flags(static_cast<pixel_format_helpers::component_flags>(pixel_format_helpers::cCompFlagGValid | pixel_format_helpers::cCompFlagBValid | pixel_format_helpers::cCompFlagAValid | pixel_format_helpers::cCompFlagNormalMap));
1190                     break;
1191                 }
1192                 case image_utils::cConversion_To_AGBR:
1193                 {
1194                     img.set_comp_flags(static_cast<pixel_format_helpers::component_flags>(pixel_format_helpers::cCompFlagRValid | pixel_format_helpers::cCompFlagGValid | pixel_format_helpers::cCompFlagBValid | pixel_format_helpers::cCompFlagAValid | pixel_format_helpers::cCompFlagNormalMap));
1195                     break;
1196                 }
1197                 case image_utils::cConversion_From_xGBR:
1198                 {
1199                     img.set_comp_flags(static_cast<pixel_format_helpers::component_flags>(pixel_format_helpers::cCompFlagRValid | pixel_format_helpers::cCompFlagGValid | pixel_format_helpers::cCompFlagBValid | pixel_format_helpers::cCompFlagNormalMap));
1200                     break;
1201                 }
1202                 case image_utils::cConversion_From_AGBR:
1203                 {
1204                     img.set_comp_flags(static_cast<pixel_format_helpers::component_flags>(pixel_format_helpers::cCompFlagRValid | pixel_format_helpers::cCompFlagGValid | pixel_format_helpers::cCompFlagBValid | pixel_format_helpers::cCompFlagAValid | pixel_format_helpers::cCompFlagNormalMap));
1205                     break;
1206                 }
1207                 case image_utils::cConversion_XY_to_XYZ:
1208                 {
1209                     img.set_comp_flags(static_cast<pixel_format_helpers::component_flags>(pixel_format_helpers::cCompFlagRValid | pixel_format_helpers::cCompFlagGValid | pixel_format_helpers::cCompFlagBValid | pixel_format_helpers::cCompFlagNormalMap));
1210                     break;
1211                 }
1212                 case cConversion_Y_To_A:
1213                 {
1214                     img.set_comp_flags(static_cast<pixel_format_helpers::component_flags>(img.get_comp_flags() | pixel_format_helpers::cCompFlagAValid));
1215                     break;
1216                 }
1217                 case cConversion_A_To_RGBA:
1218                 {
1219                     img.set_comp_flags(static_cast<pixel_format_helpers::component_flags>(pixel_format_helpers::cCompFlagRValid | pixel_format_helpers::cCompFlagGValid | pixel_format_helpers::cCompFlagBValid | pixel_format_helpers::cCompFlagAValid));
1220                     break;
1221                 }
1222                 case cConversion_Y_To_RGB:
1223                 {
1224                     img.set_comp_flags(static_cast<pixel_format_helpers::component_flags>(pixel_format_helpers::cCompFlagRValid | pixel_format_helpers::cCompFlagGValid | pixel_format_helpers::cCompFlagBValid | pixel_format_helpers::cCompFlagGrayscale | (img.has_alpha() ? pixel_format_helpers::cCompFlagAValid : 0)));
1225                     break;
1226                 }
1227                 case cConversion_To_Y:
1228                 {
1229                     img.set_comp_flags(static_cast<pixel_format_helpers::component_flags>(img.get_comp_flags() | pixel_format_helpers::cCompFlagGrayscale));
1230                     break;
1231                 }
1232                 default:
1233                 {
1234                     VOGL_ASSERT(false);
1235                     return;
1236                 }
1237             }
1238
1239             for (uint y = 0; y < img.get_height(); y++)
1240             {
1241                 for (uint x = 0; x < img.get_width(); x++)
1242                 {
1243                     color_quad_u8 src(img(x, y));
1244                     color_quad_u8 dst;
1245
1246                     switch (conv_type)
1247                     {
1248                         case image_utils::cConversion_To_CCxY:
1249                         {
1250                             color::RGB_to_YCC(dst, src);
1251                             break;
1252                         }
1253                         case image_utils::cConversion_From_CCxY:
1254                         {
1255                             color::YCC_to_RGB(dst, src);
1256                             break;
1257                         }
1258                         case image_utils::cConversion_To_xGxR:
1259                         {
1260                             dst.r = 0;
1261                             dst.g = src.g;
1262                             dst.b = 0;
1263                             dst.a = src.r;
1264                             break;
1265                         }
1266                         case image_utils::cConversion_From_xGxR:
1267                         {
1268                             dst.r = src.a;
1269                             dst.g = src.g;
1270                             // This is kinda iffy, we're assuming the image is a normal map here.
1271                             dst.b = regen_z(src.a, src.g);
1272                             dst.a = 255;
1273                             break;
1274                         }
1275                         case image_utils::cConversion_To_xGBR:
1276                         {
1277                             dst.r = 0;
1278                             dst.g = src.g;
1279                             dst.b = src.b;
1280                             dst.a = src.r;
1281                             break;
1282                         }
1283                         case image_utils::cConversion_To_AGBR:
1284                         {
1285                             dst.r = src.a;
1286                             dst.g = src.g;
1287                             dst.b = src.b;
1288                             dst.a = src.r;
1289                             break;
1290                         }
1291                         case image_utils::cConversion_From_xGBR:
1292                         {
1293                             dst.r = src.a;
1294                             dst.g = src.g;
1295                             dst.b = src.b;
1296                             dst.a = 255;
1297                             break;
1298                         }
1299                         case image_utils::cConversion_From_AGBR:
1300                         {
1301                             dst.r = src.a;
1302                             dst.g = src.g;
1303                             dst.b = src.b;
1304                             dst.a = src.r;
1305                             break;
1306                         }
1307                         case image_utils::cConversion_XY_to_XYZ:
1308                         {
1309                             dst.r = src.r;
1310                             dst.g = src.g;
1311                             // This is kinda iffy, we're assuming the image is a normal map here.
1312                             dst.b = regen_z(src.r, src.g);
1313                             dst.a = 255;
1314                             break;
1315                         }
1316                         case image_utils::cConversion_Y_To_A:
1317                         {
1318                             dst.r = src.r;
1319                             dst.g = src.g;
1320                             dst.b = src.b;
1321                             dst.a = static_cast<uint8>(src.get_luma());
1322                             break;
1323                         }
1324                         case image_utils::cConversion_Y_To_RGB:
1325                         {
1326                             uint8 y = static_cast<uint8>(src.get_luma());
1327                             dst.r = y;
1328                             dst.g = y;
1329                             dst.b = y;
1330                             dst.a = src.a;
1331                             break;
1332                         }
1333                         case image_utils::cConversion_A_To_RGBA:
1334                         {
1335                             dst.r = src.a;
1336                             dst.g = src.a;
1337                             dst.b = src.a;
1338                             dst.a = src.a;
1339                             break;
1340                         }
1341                         case image_utils::cConversion_To_Y:
1342                         {
1343                             uint8 y = static_cast<uint8>(src.get_luma());
1344                             dst.r = y;
1345                             dst.g = y;
1346                             dst.b = y;
1347                             dst.a = src.a;
1348                             break;
1349                         }
1350                         default:
1351                         {
1352                             VOGL_ASSERT(false);
1353                             dst = src;
1354                             break;
1355                         }
1356                     }
1357
1358                     img(x, y) = dst;
1359                 }
1360             }
1361         }
1362
1363         image_utils::conversion_type get_conversion_type(bool cooking, pixel_format fmt)
1364         {
1365             image_utils::conversion_type conv_type = image_utils::cConversion_Invalid;
1366
1367             if (cooking)
1368             {
1369                 switch (fmt)
1370                 {
1371                     case PIXEL_FMT_DXT5_CCxY:
1372                     {
1373                         conv_type = image_utils::cConversion_To_CCxY;
1374                         break;
1375                     }
1376                     case PIXEL_FMT_DXT5_xGxR:
1377                     {
1378                         conv_type = image_utils::cConversion_To_xGxR;
1379                         break;
1380                     }
1381                     case PIXEL_FMT_DXT5_xGBR:
1382                     {
1383                         conv_type = image_utils::cConversion_To_xGBR;
1384                         break;
1385                     }
1386                     case PIXEL_FMT_DXT5_AGBR:
1387                     {
1388                         conv_type = image_utils::cConversion_To_AGBR;
1389                         break;
1390                     }
1391                     default:
1392                         break;
1393                 }
1394             }
1395             else
1396             {
1397                 switch (fmt)
1398                 {
1399                     case PIXEL_FMT_3DC:
1400                     case PIXEL_FMT_DXN:
1401                     {
1402                         conv_type = image_utils::cConversion_XY_to_XYZ;
1403                         break;
1404                     }
1405                     case PIXEL_FMT_DXT5_CCxY:
1406                     {
1407                         conv_type = image_utils::cConversion_From_CCxY;
1408                         break;
1409                     }
1410                     case PIXEL_FMT_DXT5_xGxR:
1411                     {
1412                         conv_type = image_utils::cConversion_From_xGxR;
1413                         break;
1414                     }
1415                     case PIXEL_FMT_DXT5_xGBR:
1416                     {
1417                         conv_type = image_utils::cConversion_From_xGBR;
1418                         break;
1419                     }
1420                     case PIXEL_FMT_DXT5_AGBR:
1421                     {
1422                         conv_type = image_utils::cConversion_From_AGBR;
1423                         break;
1424                     }
1425                     default:
1426                         break;
1427                 }
1428             }
1429
1430             return conv_type;
1431         }
1432
1433         image_utils::conversion_type get_image_conversion_type_from_vogl_format(vogl_format fmt)
1434         {
1435             switch (fmt)
1436             {
1437                 case cCRNFmtDXT5_CCxY:
1438                     return image_utils::cConversion_To_CCxY;
1439                 case cCRNFmtDXT5_xGxR:
1440                     return image_utils::cConversion_To_xGxR;
1441                 case cCRNFmtDXT5_xGBR:
1442                     return image_utils::cConversion_To_xGBR;
1443                 case cCRNFmtDXT5_AGBR:
1444                     return image_utils::cConversion_To_AGBR;
1445                 default:
1446                     break;
1447             }
1448             return image_utils::cConversion_Invalid;
1449         }
1450
1451         double compute_std_dev(uint n, const color_quad_u8 *pPixels, uint first_channel, uint num_channels)
1452         {
1453             if (!n)
1454                 return 0.0f;
1455
1456             double sum = 0.0f;
1457             double sum2 = 0.0f;
1458
1459             for (uint i = 0; i < n; i++)
1460             {
1461                 const color_quad_u8 &cp = pPixels[i];
1462
1463                 if (!num_channels)
1464                 {
1465                     uint l = cp.get_luma();
1466                     sum += l;
1467                     sum2 += l * l;
1468                 }
1469                 else
1470                 {
1471                     for (uint c = 0; c < num_channels; c++)
1472                     {
1473                         uint l = cp[first_channel + c];
1474                         sum += l;
1475                         sum2 += l * l;
1476                     }
1477                 }
1478             }
1479
1480             double w = math::maximum(1U, num_channels) * n;
1481             sum /= w;
1482             sum2 /= w;
1483
1484             double var = sum2 - sum * sum;
1485             var = math::maximum<double>(var, 0.0f);
1486
1487             return sqrt(var);
1488         }
1489
1490         uint8 *read_image_from_memory(const uint8 *pImage, int nSize, int *pWidth, int *pHeight, int *pActualComps, int req_comps, const char *pFilename)
1491         {
1492             *pWidth = 0;
1493             *pHeight = 0;
1494             *pActualComps = 0;
1495
1496             if ((req_comps < 1) || (req_comps > 4))
1497                 return NULL;
1498
1499             mipmapped_texture tex;
1500
1501             buffer_stream buf_stream(pImage, nSize);
1502             buf_stream.set_name(pFilename);
1503             data_stream_serializer serializer(buf_stream);
1504
1505             if (!tex.read_from_stream(serializer))
1506                 return NULL;
1507
1508             if (tex.is_packed())
1509             {
1510                 if (!tex.unpack_from_dxt(true))
1511                     return NULL;
1512             }
1513
1514             image_u8 img;
1515             image_u8 *pImg = tex.get_level_image(0, 0, img);
1516             if (!pImg)
1517                 return NULL;
1518
1519             *pWidth = tex.get_width();
1520             *pHeight = tex.get_height();
1521
1522             if (pImg->has_alpha())
1523                 *pActualComps = 4;
1524             else if (pImg->is_grayscale())
1525                 *pActualComps = 1;
1526             else
1527                 *pActualComps = 3;
1528
1529             uint8 *pDst = NULL;
1530             if (req_comps == 4)
1531             {
1532                 pDst = (uint8 *)vogl_malloc(tex.get_total_pixels() * sizeof(uint32));
1533                 uint8 *pSrc = (uint8 *)pImg->get_ptr();
1534                 memcpy(pDst, pSrc, tex.get_total_pixels() * sizeof(uint32));
1535             }
1536             else
1537             {
1538                 image_u8 luma_img;
1539                 if (req_comps == 1)
1540                 {
1541                     luma_img = *pImg;
1542                     luma_img.convert_to_grayscale();
1543                     pImg = &luma_img;
1544                 }
1545
1546                 pixel_packer packer(req_comps, 8);
1547                 uint32 n;
1548                 pDst = image_utils::pack_image(*pImg, packer, n);
1549             }
1550
1551             return pDst;
1552         }
1553
1554         // http://lodev.org/cgtutor/fourier.html
1555         // TODO: Rewrite this - the license is too restrictive
1556         // N=width, M=height
1557         // image in memory is: [x][y][c], C bytes per pixel
1558         template <uint C>
1559         void FFT2D(int n, int m, bool inverse, const float *gRe, const float *gIm, float *GRe, float *GIm)
1560         {
1561             int l2n = 0, p = 1; //l2n will become log_2(n)
1562             while (p < n)
1563             {
1564                 p *= 2;
1565                 l2n++;
1566             }
1567             int l2m = 0;
1568             p = 1; //l2m will become log_2(m)
1569             while (p < m)
1570             {
1571                 p *= 2;
1572                 l2m++;
1573             }
1574
1575             m = 1 << l2m;
1576             n = 1 << l2n; //Make sure m and n will be powers of 2, otherwise you'll get in an infinite loop
1577
1578             //Erase all history of this array
1579             for (int x = 0; x < m; x++) //for each column
1580             {
1581                 for (int y = 0; y < m; y++) //for each row
1582                 {
1583                     for (int c = 0; c < static_cast<int>(C); c++) //for each color component
1584                     {
1585                         GRe[C * m * x + C * y + c] = gRe[C * m * x + C * y + c];
1586                         GIm[C * m * x + C * y + c] = gIm[C * m * x + C * y + c];
1587                     }
1588                 }
1589             }
1590
1591             //Bit reversal of each row
1592             int j;
1593             for (int y = 0; y < m; y++) //for each row
1594             {
1595                 for (int c = 0; c < static_cast<int>(C); c++) //for each color component
1596                 {
1597                     j = 0;
1598                     for (int i = 0; i < n - 1; i++)
1599                     {
1600                         GRe[C * m * i + C * y + c] = gRe[C * m * j + C * y + c];
1601                         GIm[C * m * i + C * y + c] = gIm[C * m * j + C * y + c];
1602                         int k = n / 2;
1603                         while (k <= j)
1604                         {
1605                             j -= k;
1606                             k /= 2;
1607                         }
1608                         j += k;
1609                     }
1610                 }
1611             }
1612
1613             //Bit reversal of each column
1614             float tx = 0, ty = 0;
1615             for (int x = 0; x < n; x++) //for each column
1616             {
1617                 for (int c = 0; c < static_cast<int>(C); c++) //for each color component
1618                 {
1619                     j = 0;
1620                     for (int i = 0; i < m - 1; i++)
1621                     {
1622                         if (i < j)
1623                         {
1624                             tx = GRe[C * m * x + C * i + c];
1625                             ty = GIm[C * m * x + C * i + c];
1626                             GRe[C * m * x + C * i + c] = GRe[C * m * x + C * j + c];
1627                             GIm[C * m * x + C * i + c] = GIm[C * m * x + C * j + c];
1628                             GRe[C * m * x + C * j + c] = tx;
1629                             GIm[C * m * x + C * j + c] = ty;
1630                         }
1631                         int k = m / 2;
1632                         while (k <= j)
1633                         {
1634                             j -= k;
1635                             k /= 2;
1636                         }
1637                         j += k;
1638                     }
1639                 }
1640             }
1641
1642             //Calculate the FFT of the columns
1643             for (int x = 0; x < n; x++) //for each column
1644             {
1645                 for (int c = 0; c < static_cast<int>(C); c++) //for each color component
1646                 {
1647                     //This is the 1D FFT:
1648                     float ca = -1.0;
1649                     float sa = 0.0;
1650                     int l1 = 1, l2 = 1;
1651                     for (int l = 0; l < l2n; l++)
1652                     {
1653                         l1 = l2;
1654                         l2 *= 2;
1655                         float u1 = 1.0;
1656                         float u2 = 0.0;
1657                         for (int j = 0; j < l1; j++)
1658                         {
1659                             for (int i = j; i < n; i += l2)
1660                             {
1661                                 int i1 = i + l1;
1662                                 float t1 = u1 * GRe[C * m * x + C * i1 + c] - u2 * GIm[C * m * x + C * i1 + c];
1663                                 float t2 = u1 * GIm[C * m * x + C * i1 + c] + u2 * GRe[C * m * x + C * i1 + c];
1664                                 GRe[C * m * x + C * i1 + c] = GRe[C * m * x + C * i + c] - t1;
1665                                 GIm[C * m * x + C * i1 + c] = GIm[C * m * x + C * i + c] - t2;
1666                                 GRe[C * m * x + C * i + c] += t1;
1667                                 GIm[C * m * x + C * i + c] += t2;
1668                             }
1669                             float z = u1 * ca - u2 * sa;
1670                             u2 = u1 * sa + u2 * ca;
1671                             u1 = z;
1672                         }
1673                         sa = sqrtf((1.0f - ca) / 2.0f);
1674                         if (!inverse)
1675                             sa = -sa;
1676                         ca = sqrtf((1.0f + ca) / 2.0f);
1677                     }
1678                 }
1679             }
1680
1681             //Calculate the FFT of the rows
1682             for (int y = 0; y < m; y++) //for each row
1683             {
1684                 for (int c = 0; c < static_cast<int>(C); c++) //for each color component
1685                 {
1686                     //This is the 1D FFT:
1687                     float ca = -1.0;
1688                     float sa = 0.0;
1689                     int l1 = 1, l2 = 1;
1690                     for (int l = 0; l < l2m; l++)
1691                     {
1692                         l1 = l2;
1693                         l2 *= 2;
1694                         float u1 = 1.0;
1695                         float u2 = 0.0;
1696                         for (int j = 0; j < l1; j++)
1697                         {
1698                             for (int i = j; i < n; i += l2)
1699                             {
1700                                 int i1 = i + l1;
1701                                 float t1 = u1 * GRe[C * m * i1 + C * y + c] - u2 * GIm[C * m * i1 + C * y + c];
1702                                 float t2 = u1 * GIm[C * m * i1 + C * y + c] + u2 * GRe[C * m * i1 + C * y + c];
1703                                 GRe[C * m * i1 + C * y + c] = GRe[C * m * i + C * y + c] - t1;
1704                                 GIm[C * m * i1 + C * y + c] = GIm[C * m * i + C * y + c] - t2;
1705                                 GRe[C * m * i + C * y + c] += t1;
1706                                 GIm[C * m * i + C * y + c] += t2;
1707                             }
1708                             float z = u1 * ca - u2 * sa;
1709                             u2 = u1 * sa + u2 * ca;
1710                             u1 = z;
1711                         }
1712                         sa = sqrtf((1.0f - ca) / 2.0f);
1713                         if (!inverse)
1714                             sa = -sa;
1715                         ca = sqrtf((1.0f + ca) / 2.0f);
1716                     }
1717                 }
1718             }
1719
1720             int d;
1721             if (inverse)
1722                 d = n;
1723             else
1724                 d = m;
1725             float inv_d = 1.0f / d;
1726
1727             for (int x = 0; x < n; x++)
1728             {
1729                 for (int y = 0; y < m; y++)
1730                 {
1731                     for (int c = 0; c < static_cast<int>(C); c++) //for every value of the buffers
1732                     {
1733                         GRe[C * m * x + C * y + c] *= inv_d;
1734                         GIm[C * m * x + C * y + c] *= inv_d;
1735                     }
1736                 }
1737             }
1738         }
1739
1740         // N=width, M=height
1741         // image in memory is: [x][y][c], C bytes per pixel
1742         //void FFT2D(int n, int m, bool inverse, const float *gRe, const float *gIm, float *GRe, float *GIm)
1743         bool forward_fourier_transform(const image_u8 &src, image_u8 &dst, uint comp_mask)
1744         {
1745             const uint width = src.get_width();
1746             const uint height = src.get_height();
1747             const uint half_width = src.get_width() / 2;
1748             const uint half_height = src.get_height() / 2;
1749
1750             if (!math::is_power_of_2(width) || !math::is_power_of_2(height))
1751             {
1752                 VOGL_ASSERT_ALWAYS;
1753                 return false;
1754             }
1755
1756             dst.resize(width * 2, height);
1757
1758             vogl::vector2D<float> re0(height, width), im0(height, width);
1759             vogl::vector2D<float> re1(height, width), im1(height, width);
1760
1761             for (uint c = 0; c < 4; c++)
1762             {
1763                 if ((comp_mask & (1 << c)) == 0)
1764                     continue;
1765
1766                 for (uint y = 0; y < height; y++)
1767                     for (uint x = 0; x < width; x++)
1768                         re0(y, x) = src(x, y)[c];
1769
1770                 FFT2D<1>(width, height, false, re0.get_ptr(), im0.get_ptr(), re1.get_ptr(), im1.get_ptr());
1771
1772                 for (uint y = 0; y < height; y++)
1773                 {
1774                     for (uint x = 0; x < width; x++)
1775                     {
1776                         dst(x, y)[c] = static_cast<uint8>(math::clamp<float>(128.5f + re1.at_wrapped(y + half_height, x + half_width), 0.0f, 255.0f));
1777                         dst(x + width, y)[c] = static_cast<uint8>(math::clamp<float>(128.5f + im1.at_wrapped(y + half_height, x + half_width), 0.0f, 255.0f));
1778                     }
1779                 }
1780             }
1781
1782             return true;
1783         }
1784
1785         void gaussian_filter(image_u8 &dst, const image_u8 &orig_img, uint width_divisor, uint height_divisor, uint odd_filter_width, float sigma_sqr, bool wrapping)
1786         {
1787             VOGL_ASSERT(odd_filter_width && (odd_filter_width & 1));
1788             odd_filter_width |= 1;
1789
1790             vector2D<float> kernel(odd_filter_width, odd_filter_width);
1791             math::compute_gaussian_kernel(kernel.get_ptr(), odd_filter_width, odd_filter_width, sigma_sqr, math::cComputeGaussianFlagNormalize);
1792
1793             const int dst_width = orig_img.get_width() / width_divisor;
1794             const int dst_height = orig_img.get_height() / height_divisor;
1795
1796             const int H = odd_filter_width / 2;
1797             const int L = -H;
1798
1799             dst.crop(dst_width, dst_height);
1800
1801             for (int oy = 0; oy < dst_height; oy++)
1802             {
1803                 for (int ox = 0; ox < dst_width; ox++)
1804                 {
1805                     vec4F c(0.0f);
1806
1807                     for (int yd = L; yd <= H; yd++)
1808                     {
1809                         int y = oy * height_divisor + (height_divisor >> 1) + yd;
1810
1811                         for (int xd = L; xd <= H; xd++)
1812                         {
1813                             int x = ox * width_divisor + (width_divisor >> 1) + xd;
1814
1815                             const color_quad_u8 &p = orig_img.get_clamped_or_wrapped(x, y, wrapping, wrapping);
1816
1817                             float w = kernel(xd + H, yd + H);
1818                             c[0] += p[0] * w;
1819                             c[1] += p[1] * w;
1820                             c[2] += p[2] * w;
1821                             c[3] += p[3] * w;
1822                         }
1823                     }
1824
1825                     dst(ox, oy).set(math::float_to_int_nearest(c[0]), math::float_to_int_nearest(c[1]), math::float_to_int_nearest(c[2]), math::float_to_int_nearest(c[3]));
1826                 }
1827             }
1828         }
1829
1830         void convolution_filter(image_u8 &dst, const image_u8 &src, const float *pWeights, uint M, uint N, bool wrapping)
1831         {
1832             if (((M & 1) == 0) || ((N & 1) == 0))
1833             {
1834                 VOGL_ASSERT_ALWAYS;
1835                 return;
1836             }
1837
1838             const int width = src.get_width();
1839             const int height = src.get_height();
1840
1841             dst.crop(width, height);
1842
1843             const int HM = M / 2;
1844             const int HN = N / 2;
1845
1846             VOGL_ASSERT((HM * 2 + 1) == static_cast<int>(M));
1847             VOGL_ASSERT((HN * 2 + 1) == static_cast<int>(N));
1848
1849             for (int dst_y = 0; dst_y < height; dst_y++)
1850             {
1851                 for (int dst_x = 0; dst_x < width; dst_x++)
1852                 {
1853                     vec4F c(0.0f);
1854
1855                     for (int yd = -HM; yd <= HM; yd++)
1856                     {
1857                         int src_y = wrapping ? math::posmod(dst_y + yd, height) : math::clamp<int>(dst_y + yd, 0, height - 1);
1858
1859                         for (int xd = -HN; xd <= HN; xd++)
1860                         {
1861                             int src_x = wrapping ? math::posmod(dst_x + xd, width) : math::clamp<int>(dst_x + xd, 0, width - 1);
1862                             float w = pWeights[(yd + HM) * N + (xd + HN)];
1863
1864                             const color_quad_u8 &src_c = src(src_x, src_y);
1865                             c[0] += src_c[0] * w;
1866                             c[1] += src_c[1] * w;
1867                             c[2] += src_c[2] * w;
1868                             c[3] += src_c[3] * w;
1869                         }
1870                     }
1871
1872                     dst(dst_x, dst_y).set(static_cast<int>(c[0] + .5f), static_cast<int>(c[1] + .5f), static_cast<int>(c[2] + .5f), static_cast<int>(c[3] + .5f));
1873                 }
1874             }
1875         }
1876
1877         void convolution_filter(image_f &dst, const image_f &src, const float *pWeights, uint M, uint N, bool wrapping)
1878         {
1879             if (((M & 1) == 0) || ((N & 1) == 0))
1880             {
1881                 VOGL_ASSERT_ALWAYS;
1882                 return;
1883             }
1884
1885             const int width = src.get_width();
1886             const int height = src.get_height();
1887
1888             dst.crop(width, height);
1889
1890             const int HM = M / 2;
1891             const int HN = N / 2;
1892
1893             VOGL_ASSERT((HM * 2 + 1) == static_cast<int>(M));
1894             VOGL_ASSERT((HN * 2 + 1) == static_cast<int>(N));
1895
1896             for (int dst_y = 0; dst_y < height; dst_y++)
1897             {
1898                 for (int dst_x = 0; dst_x < width; dst_x++)
1899                 {
1900                     vec4F c(0.0f);
1901
1902                     for (int yd = -HM; yd <= HM; yd++)
1903                     {
1904                         int src_y = wrapping ? math::posmod(dst_y + yd, height) : math::clamp<int>(dst_y + yd, 0, height - 1);
1905
1906                         for (int xd = -HN; xd <= HN; xd++)
1907                         {
1908                             int src_x = wrapping ? math::posmod(dst_x + xd, width) : math::clamp<int>(dst_x + xd, 0, width - 1);
1909                             float w = pWeights[(yd + HM) * N + (xd + HN)];
1910
1911                             const color_quad_u8 &src_c = src(src_x, src_y);
1912                             c[0] += src_c[0] * w;
1913                             c[1] += src_c[1] * w;
1914                             c[2] += src_c[2] * w;
1915                             c[3] += src_c[3] * w;
1916                         }
1917                     }
1918
1919                     dst(dst_x, dst_y).set(c[0], c[1], c[2], c[3]);
1920                 }
1921             }
1922         }
1923
1924     } // namespace image_utils
1925
1926 } // namespace vogl