spherical harmonics arithmetic

David Mayerich
1 parent 18f207ce
Showing 3 changed files with 51 additions and 3 deletions Show diff stats
stim/math/constants.h
stim/math/spharmonics.h
stim/visualization/gl_spharmonics.h
@@ -4,7 +4,6 @@
 #define STIM_PI		3.1415926535897932384626433832795028841971693993751058209749445923078164062862
 #define STIM_TAU	2 * STIM_PI
  
-#include "stim/cuda/cudatools/callable.h"
 namespace stim{
 	const double PI		=	3.1415926535897932384626433832795028841971693993751058209749445923078164062862;
 	const double TAU	=	2 * stim::PI;
@@ -5,7 +5,6 @@
 #include <boost/math/special_functions/spherical_harmonic.hpp>
 #include <vector>
  
-#define PI 3.14159
 #define WIRE_SCALE 1.001
 namespace stim{
  
@@ -32,6 +31,12 @@ protected:
  
  
 public:
+	spharmonics() {
+		mcN = 0;
+	}
+	spharmonics(size_t c) : spharmonics() {
+		resize(c);
+	}
  
 	void push(double c){
 		C.push_back(c);
@@ -163,6 +168,38 @@ public:
 		return fx;
 	}
  
+	//overload arithmetic operations
+	spharmonics<T> operator*(T rhs) const {
+		spharmonics<T> result(C.size());				//create a new spherical harmonics object
+		for (size_t c = 0; c < C.size(); c++)			//for each coefficient
+			result.C[c] = C[c] * rhs;					//calculate the factor and store the result in the new spharmonics object
+		return result;
+	}
+
+	spharmonics<T> operator+(spharmonics<T> rhs) {
+		size_t low = std::min(C.size(), rhs.C.size());				//store the number of coefficients in the lowest object
+		size_t high = std::max(C.size(), rhs.C.size());				//store the number of coefficients in the result
+		bool rhs_lowest = false;								//true if rhs has the lowest number of coefficients
+		if (rhs.C.size() < C.size()) rhs_lowest = true;			//if rhs has a lower number of coefficients, set the flag
+
+		spharmonics<T> result(high);								//create a new object
+		size_t c;
+		for (c = 0; c < low; c++)						//perform the first batch of additions
+			result.C[c] = C[c] + rhs.C[c];						//perform the addition
+
+		for (c = low; c < high; c++) {
+			if (rhs_lowest)
+				result.C[c] = C[c];
+			else
+				result.C[c] = rhs.C[c];
+		}
+		return result;
+	}
+
+	spharmonics<T> operator-(spharmonics<T> rhs) {
+		return (*this) + (rhs * (T)(-1));
+	}
+
 };		//end class sph_harmonics
  
  
@@ -163,6 +163,19 @@ protected:
 	}
  
 public:
+	gl_spharmonics<T>() {
+
+	}
+
+	gl_spharmonics<T>(const spharmonics<T> copy) : gl_spharmonics<T>() {
+		C = copy.C;
+	}
+
+	gl_spharmonics<T>& operator=(const spharmonics<T> rhs) {
+		gl_spharmonics<T> result(rhs.C.size());
+		result.C = spharmonics<T>::rhs.C;
+		return result;
+	}
  
 	void glRender(){
 		//set all OpenGL parameters required for drawing
@@ -177,7 +190,6 @@ public:
 		gen_function();
 	}
  
-
 };		//end gl_spharmonics