Commit 39d843e3d8634588fe028cdb5948b0421d5aeec8
1 parent
761ebaa9
fixed spharmonics and extened gl_spharmonics to allow for interpolation between …
…different harmonic functions
Showing
2 changed files
with
54 additions
and
4 deletions
Show diff stats
stim/math/spharmonics.h
1 | #ifndef STIM_SPH_HARMONICS | 1 | #ifndef STIM_SPH_HARMONICS |
2 | #define STIM_SPH_HARMONICS | 2 | #define STIM_SPH_HARMONICS |
3 | 3 | ||
4 | +#include <complex> | ||
4 | #include <stim/math/vector.h> | 5 | #include <stim/math/vector.h> |
5 | #include <boost/math/special_functions/spherical_harmonic.hpp> | 6 | #include <boost/math/special_functions/spherical_harmonic.hpp> |
7 | +#include <stim/math/constants.h> | ||
6 | #include <vector> | 8 | #include <vector> |
7 | 9 | ||
8 | #define WIRE_SCALE 1.001 | 10 | #define WIRE_SCALE 1.001 |
@@ -11,16 +13,20 @@ namespace stim{ | @@ -11,16 +13,20 @@ namespace stim{ | ||
11 | template<class T> | 13 | template<class T> |
12 | class spharmonics{ | 14 | class spharmonics{ |
13 | 15 | ||
16 | +public: | ||
17 | + std::vector<T> C; //list of SH coefficients | ||
18 | + | ||
14 | protected: | 19 | protected: |
15 | 20 | ||
16 | - std::vector<T> C; //list of SH coefficients | ||
17 | 21 | ||
18 | unsigned int mcN; //number of Monte-Carlo samples | 22 | unsigned int mcN; //number of Monte-Carlo samples |
19 | 23 | ||
20 | //calculate the value of the SH basis function (l, m) at (theta, phi) | 24 | //calculate the value of the SH basis function (l, m) at (theta, phi) |
21 | //here, theta = [0, PI], phi = [0, 2*PI] | 25 | //here, theta = [0, PI], phi = [0, 2*PI] |
22 | double SH(int l, int m, double theta, double phi){ | 26 | double SH(int l, int m, double theta, double phi){ |
23 | - return boost::math::spherical_harmonic_r(l, m, phi, theta); | 27 | + std::complex<T> result = boost::math::spherical_harmonic(l, m, phi, theta); |
28 | + return result.imag() + result.real(); | ||
29 | +// return boost::math::spherical_harmonic_i(l, m, phi, theta); | ||
24 | } | 30 | } |
25 | 31 | ||
26 | unsigned int coeff_1d(unsigned int l, int m){ | 32 | unsigned int coeff_1d(unsigned int l, int m){ |
@@ -55,6 +61,16 @@ public: | @@ -55,6 +61,16 @@ public: | ||
55 | C[c] = value; | 61 | C[c] = value; |
56 | } | 62 | } |
57 | 63 | ||
64 | + unsigned int getSize() const{ | ||
65 | + return C.size(); | ||
66 | + } | ||
67 | + | ||
68 | + std::vector<T> getC() const{ | ||
69 | + return C; | ||
70 | + } | ||
71 | + | ||
72 | + | ||
73 | + | ||
58 | /// Initialize Monte-Carlo sampling of a function using N spherical harmonics coefficients | 74 | /// Initialize Monte-Carlo sampling of a function using N spherical harmonics coefficients |
59 | 75 | ||
60 | /// @param N is the number of spherical harmonics coefficients used to represent the user function | 76 | /// @param N is the number of spherical harmonics coefficients used to represent the user function |
@@ -167,7 +183,7 @@ public: | @@ -167,7 +183,7 @@ public: | ||
167 | 183 | ||
168 | return fx; | 184 | return fx; |
169 | } | 185 | } |
170 | - | 186 | +/* |
171 | //overload arithmetic operations | 187 | //overload arithmetic operations |
172 | spharmonics<T> operator*(T rhs) const { | 188 | spharmonics<T> operator*(T rhs) const { |
173 | spharmonics<T> result(C.size()); //create a new spherical harmonics object | 189 | spharmonics<T> result(C.size()); //create a new spherical harmonics object |
@@ -199,7 +215,7 @@ public: | @@ -199,7 +215,7 @@ public: | ||
199 | spharmonics<T> operator-(spharmonics<T> rhs) { | 215 | spharmonics<T> operator-(spharmonics<T> rhs) { |
200 | return (*this) + (rhs * (T)(-1)); | 216 | return (*this) + (rhs * (T)(-1)); |
201 | } | 217 | } |
202 | - | 218 | +*/ |
203 | }; //end class sph_harmonics | 219 | }; //end class sph_harmonics |
204 | 220 | ||
205 | 221 |
stim/visualization/gl_spharmonics.h
@@ -190,6 +190,40 @@ public: | @@ -190,6 +190,40 @@ public: | ||
190 | gen_function(); | 190 | gen_function(); |
191 | } | 191 | } |
192 | 192 | ||
193 | + //overload arithmetic operations | ||
194 | + gl_spharmonics<T> operator*(T rhs) const { | ||
195 | + gl_spharmonics<T> result(C.size()); //create a new sph erical harmonics object | ||
196 | + for (size_t c = 0; c < C.size(); c++) //for each coefficient | ||
197 | + result.C[c] = C[c] * rhs; //calculat e the factor and store the result in the new spharmonics object | ||
198 | + return result; | ||
199 | + } | ||
200 | + | ||
201 | + gl_spharmonics<T> operator+(gl_spharmonics<T> rhs) { | ||
202 | + size_t low = std::min(C.size(), rhs.C.size()); //store th e number of coefficients in the lowest object | ||
203 | + size_t high = std::max(C.size(), rhs.C.size()); //store th e number of coefficients in the result | ||
204 | + bool rhs_lowest = false; //true if rhs has the lowest number of coefficients | ||
205 | + if (rhs.C.size() < C.size()) rhs_lowest = true; //if rhs has a low er number of coefficients, set the flag | ||
206 | + | ||
207 | + gl_spharmonics<T> result(high); //create a new object | ||
208 | + size_t c; | ||
209 | + for (c = 0; c < low; c++) //perform the first batch of additions | ||
210 | + result.C[c] = C[c] + rhs.C[c]; // perform the addition | ||
211 | + | ||
212 | + for (c = low; c < high; c++) { | ||
213 | + if (rhs_lowest) | ||
214 | + result.C[c] = C[c]; | ||
215 | + else | ||
216 | + result.C[c] = rhs.C[c]; | ||
217 | + } | ||
218 | + return result; | ||
219 | + } | ||
220 | + | ||
221 | + gl_spharmonics<T> operator-(gl_spharmonics<T> rhs) { | ||
222 | + return (*this) + (rhs * (T)(-1)); | ||
223 | + } | ||
224 | + | ||
225 | + | ||
226 | + | ||
193 | }; //end gl_spharmonics | 227 | }; //end gl_spharmonics |
194 | 228 | ||
195 | 229 |