git.cworth.org Git - vogl/blob - src/voglcore/vogl_math.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_math.cpp
  28 #include "vogl_core.h"
  29
  30 namespace vogl
  31 {
  32     namespace math
  33     {
  34         uint g_bitmasks[32] =
  35             {
  36                 1U << 0U, 1U << 1U, 1U << 2U, 1U << 3U,
  37                 1U << 4U, 1U << 5U, 1U << 6U, 1U << 7U,
  38                 1U << 8U, 1U << 9U, 1U << 10U, 1U << 11U,
  39                 1U << 12U, 1U << 13U, 1U << 14U, 1U << 15U,
  40                 1U << 16U, 1U << 17U, 1U << 18U, 1U << 19U,
  41                 1U << 20U, 1U << 21U, 1U << 22U, 1U << 23U,
  42                 1U << 24U, 1U << 25U, 1U << 26U, 1U << 27U,
  43                 1U << 28U, 1U << 29U, 1U << 30U, 1U << 31U
  44             };
  45
  46         bool is_prime(uint n)
  47         {
  48             if (n <= 2)
  49                 return n == 2;
  50
  51             if ((n & 1) == 0)
  52                 return false;
  53
  54             const uint end = (uint)sqrt(static_cast<double>(n));
  55
  56             for (uint i = 3; i <= end; i += 2)
  57                 if ((n % i) == 0)
  58                     return false;
  59
  60             return true;
  61         }
  62
  63         // Find the smallest prime >= n.
  64         uint get_prime(uint n)
  65         {
  66             while (!is_prime(n))
  67                 n++;
  68             return n;
  69         }
  70
  71         double compute_entropy(const uint8 *p, uint n)
  72         {
  73             uint hist[256];
  74             utils::zero_object(hist);
  75
  76             for (uint i = 0; i < n; i++)
  77                 hist[*p++]++;
  78
  79             double entropy = 0.0f;
  80
  81             const double invln2 = 1.0f / log(2.0f);
  82             for (uint i = 0; i < 256; i++)
  83             {
  84                 if (!hist[i])
  85                     continue;
  86
  87                 double prob = static_cast<double>(hist[i]) / n;
  88                 entropy += (-log(prob) * invln2) * hist[i];
  89             }
  90
  91             return entropy;
  92         }
  93
  94         void compute_lower_pow2_dim(int &width, int &height)
  95         {
  96             const int tex_width = width;
  97             const int tex_height = height;
  98
  99             width = 1;
 100             for (;;)
 101             {
 102                 if ((width * 2) > tex_width)
 103                     break;
 104                 width *= 2;
 105             }
 106
 107             height = 1;
 108             for (;;)
 109             {
 110                 if ((height * 2) > tex_height)
 111                     break;
 112                 height *= 2;
 113             }
 114         }
 115
 116         void compute_upper_pow2_dim(int &width, int &height)
 117         {
 118             if (!math::is_power_of_2((uint32)width))
 119                 width = math::next_pow2((uint32)width);
 120
 121             if (!math::is_power_of_2((uint32)height))
 122                 height = math::next_pow2((uint32)height);
 123         }
 124
 125         float gauss(int x, int y, float sigma_sqr)
 126         {
 127             float pow = expf(-((x * x + y * y) / (2.0f * sigma_sqr)));
 128             float g = (1.0f / (sqrtf(2.0f * cPi * sigma_sqr))) * pow;
 129             return g;
 130         }
 131
 132         // size_x/y should be odd
 133         void compute_gaussian_kernel(float *pDst, int size_x, int size_y, float sigma_sqr, uint flags)
 134         {
 135             VOGL_ASSERT(size_x & size_y & 1);
 136
 137             if (!(size_x | size_y))
 138                 return;
 139
 140             int mid_x = size_x / 2;
 141             int mid_y = size_y / 2;
 142
 143             double sum = 0;
 144             for (int x = 0; x < size_x; x++)
 145             {
 146                 for (int y = 0; y < size_y; y++)
 147                 {
 148                     float g;
 149                     if ((x > mid_x) && (y < mid_y))
 150                         g = pDst[(size_x - x - 1) + y * size_x];
 151                     else if ((x < mid_x) && (y > mid_y))
 152                         g = pDst[x + (size_y - y - 1) * size_x];
 153                     else if ((x > mid_x) && (y > mid_y))
 154                         g = pDst[(size_x - x - 1) + (size_y - y - 1) * size_x];
 155                     else
 156                         g = gauss(x - mid_x, y - mid_y, sigma_sqr);
 157
 158                     pDst[x + y * size_x] = g;
 159                     sum += g;
 160                 }
 161             }
 162
 163             if (flags & cComputeGaussianFlagNormalizeCenterToOne)
 164             {
 165                 sum = pDst[mid_x + mid_y * size_x];
 166             }
 167
 168             if (flags & (cComputeGaussianFlagNormalizeCenterToOne | cComputeGaussianFlagNormalize))
 169             {
 170                 double one_over_sum = 1.0f / sum;
 171                 for (int i = 0; i < size_x * size_y; i++)
 172                     pDst[i] = static_cast<float>(pDst[i] * one_over_sum);
 173
 174                 if (flags & cComputeGaussianFlagNormalizeCenterToOne)
 175                     pDst[mid_x + mid_y * size_x] = 1.0f;
 176             }
 177
 178             if (flags & cComputeGaussianFlagPrint)
 179             {
 180                 printf("{\n");
 181                 for (int y = 0; y < size_y; y++)
 182                 {
 183                     printf("  ");
 184                     for (int x = 0; x < size_x; x++)
 185                     {
 186                         printf("%f, ", pDst[x + y * size_x]);
 187                     }
 188                     printf("\n");
 189                 }
 190                 printf("}");
 191             }
 192         }
 193
 194     } // namespace math
 195 } // namespace vogl