First release

This commit is contained in:
Paolo Cignoni 2008-05-28 08:55:04 +00:00
parent 0a0436a13c
commit a152d0e208
1 changed files with 155 additions and 0 deletions

View File

@ -0,0 +1,155 @@
/****************************************************************************
* VCGLib o o *
* Visual and Computer Graphics Library o o *
* _ O _ *
* Copyright(C) 2006 \/)\/ *
* Visual Computing Lab /\/| *
* ISTI - Italian National Research Council | *
* \ *
* All rights reserved. *
* *
* This program is free software; you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by *
* the Free Software Foundation; either version 2 of the License, or *
* (at your option) any later version. *
* *
* This program is distributed in the hope that it will be useful, *
* but WITHOUT ANY WARRANTY; without even the implied warranty of *
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
* GNU General Public License (http://www.gnu.org/licenses/gpl.txt) *
* for more details. *
* *
****************************************************************************/
#ifndef __VCGLIB_SPHERICAL_HARMONICS_H
#define __VCGLIB_SPHERICAL_HARMONICS_H
#include "vcg/math/base.h"
#include "vcg/math/legendre.h"
#include "vcg/math/factorial.h"
namespace vcg{
namespace math{
/**
* Although the Real Spherical Harmonic Function is correctly defined over any
* positive l and any -l <= m <= l, the two internal functions computing the
* imaginary and real parts of the Complex Spherical Harmonic Functions are defined
* for positive m only.
*/
template <typename ScalarType>
class SphericalHarmonics{
private :
inline static ScalarType scaling_factor(unsigned l, unsigned m)
{
return Sqrt( ( (2.0*l + 1.0) * Factorial<ScalarType>(l-m) ) / (4.0 * M_PI * Factorial<ScalarType>(l + m)) );;
}
inline static ScalarType complex_spherical_harmonic_re(unsigned l, unsigned m, ScalarType theta, ScalarType phi)
{
return scaling_factor(l, m) * Legendre<ScalarType>::AssociatedPolynomial(l, m, Cos(theta), Sin(theta)) * Cos(m * phi);
}
inline static ScalarType complex_spherical_harmonic_im(unsigned l, unsigned m, ScalarType theta, ScalarType phi)
{
return scaling_factor(l, m) * Legendre<ScalarType>::AssociatedPolynomial(l, m, Cos(theta), Sin(theta)) * Sin(m * phi);
}
ScalarType * coefficients;
int max_band;
public :
/**
* Returns the Real Spherical Harmonic Function
*
* l is any positive integer,
* m is such that -l <= m <= l
* theta is inside [0, PI]
* phi is inside [0, 2*PI]
*/
static ScalarType Real(unsigned l, int m, ScalarType theta, ScalarType phi)
{
assert((int)-l <= m && m <= (int)l && theta >= 0 && theta <= M_PI && phi >= 0 && phi <= 2 * M_PI);
if (m > 0) return SQRT_TWO * complex_spherical_harmonic_re(l, m, theta, phi);
else if (m == 0) return scaling_factor(l, 0) * Legendre<ScalarType>::Polynomial(l, Cos(theta));
else return SQRT_TWO * complex_spherical_harmonic_im(l, -m, theta, phi);
}
typedef ScalarType (*PolarFunction) (ScalarType theta, ScalarType phi);
static SphericalHarmonics Project(PolarFunction fun, unsigned max_band, unsigned n_samples)
{
const ScalarType weight = 4 * M_PI;
unsigned sqrt_n_samples = (unsigned int) Sqrt((int)n_samples);
unsigned actual_n_samples = sqrt_n_samples * sqrt_n_samples;
unsigned n_coeff = max_band * max_band;
ScalarType one_over_n = 1.0/(ScalarType)sqrt_n_samples;
RandomGenerator rand;
SphericalHarmonics sph;
sph.coefficients = new ScalarType[n_coeff];
sph.max_band = max_band;
int i = 0;
for (unsigned k = 0; k < n_coeff; k++ ) sph.coefficients[k] = 0;
for (unsigned a = 0; a < sqrt_n_samples; ++a )
{
for (unsigned b = 0; b < sqrt_n_samples; ++b)
{
ScalarType x = (a + rand(INT_MAX)/(ScalarType)INT_MAX) * one_over_n;
ScalarType y = (b + rand(INT_MAX)/(ScalarType)INT_MAX) * one_over_n;
ScalarType theta = 2.0 * Acos(Sqrt(1.0 - x));
ScalarType phi = 2.0 * M_PI * y;
for (int l = 0; l < (int)max_band; ++l)
{
for (int m = -l; m <= l; ++m)
{
int index = l * (l+1) + m;
sph.coefficients[index] += fun(theta, phi) * Real(l, m, theta, phi);
}
}
i++;
}
}
ScalarType factor = weight / actual_n_samples;
for(i = 0; i < (int)n_coeff; ++i)
{
sph.coefficients[i] *= factor;
}
return sph;
}
ScalarType operator()(ScalarType theta, ScalarType phi)
{
ScalarType f = 0;
for (int l = 0; l < max_band; ++l)
{
for (int m = -l; m <= l; ++m)
{
int index = l * (l+1) + m;
f += (coefficients[index] * Real(l, m, theta, phi));
}
}
return f;
}
};
}} //namespace vcg::math
#endif