aboutsummaryrefslogtreecommitdiff
path: root/src/mtx_sh/chebyshev12.c
diff options
context:
space:
mode:
authorFranz Zotter <fzotter@users.sourceforge.net>2009-01-12 17:09:13 +0000
committerFranz Zotter <fzotter@users.sourceforge.net>2009-01-12 17:09:13 +0000
commit55253e4e15b093a0826fc10b376f7724cea2ba7e (patch)
tree1e1c06129db9f86a12daddcf7c985d530d0872d9 /src/mtx_sh/chebyshev12.c
parentbacae787abd380690126e6adf82ab0cd19f5abec (diff)
added [mtx_sh] for spherical harmonics evaluation.
svn path=/trunk/externals/iem/iemmatrix/; revision=10515
Diffstat (limited to 'src/mtx_sh/chebyshev12.c')
-rw-r--r--src/mtx_sh/chebyshev12.c75
1 files changed, 75 insertions, 0 deletions
diff --git a/src/mtx_sh/chebyshev12.c b/src/mtx_sh/chebyshev12.c
new file mode 100644
index 0000000..2466ac0
--- /dev/null
+++ b/src/mtx_sh/chebyshev12.c
@@ -0,0 +1,75 @@
+/*
+ Evaluates all circular harmonics
+ at the angles phi up to the order nmax.
+ using the recurrence for the Chebyshev
+ polynomials of the first and second kind
+ T has the dimensions length(phi) x 2nmax+1
+
+ Implementation by Franz Zotter, Institute of Electronic Music and Acoustics
+ (IEM), University of Music and Dramatic Arts (KUG), Graz, Austria
+ http://iem.at/Members/zotter, 2008.
+
+ This code is published under the Gnu General Public License, see
+ "LICENSE.txt"
+*/
+
+#include "mtx_sh/chebyshev12.h"
+
+Cheby12WorkSpace *chebyshev12_alloc(const size_t nmax, const size_t l) {
+ Cheby12WorkSpace *wc;
+ // memory allocation
+ if ((wc=(Cheby12WorkSpace*)calloc(1,sizeof(Cheby12WorkSpace)))!=0) {
+ wc->l=l;
+ wc->nmax=nmax;
+ if ((wc->t=(double*)calloc(l*(2*nmax+1),sizeof(double)))==0) {
+ free(wc);
+ return 0;
+ }
+ return wc;
+ }
+ return 0;
+}
+
+void chebyshev12_free(Cheby12WorkSpace *wc) {
+ if (wc!=0) {
+ free(wc->t);
+ free(wc);
+ }
+}
+
+void chebyshev12(double *phi, Cheby12WorkSpace *wc) {
+ int l,l0,n;
+ const int incr=2*wc->nmax+1;
+ double *cosphi;
+ double *sinphi;
+ // memory allocation
+ if ((wc!=0)&&(phi!=0)) {
+ if ((cosphi=(double*)calloc(wc->l,sizeof(double)))==0) {
+ return;
+ }
+ if ((sinphi=(double*)calloc(wc->l,sizeof(double)))==0) {
+ free(cosphi);
+ return;
+ }
+ // constants and initialization
+ for (l=0, l0=wc->nmax; l<wc->l; l++, l0+=incr) {
+ cosphi[l]=cos(phi[l]);
+ sinphi[l]=sin(phi[l]);
+ // initial value T_0=1
+ wc->t[l0]=1;
+ wc->t[l0+1]=cosphi[l];
+ wc->t[l0-1]=sinphi[l];
+ }
+ // recurrence for n>1
+ for (n=2; n<=wc->nmax; n++) {
+ for (l=0, l0=wc->nmax; l<wc->l; l++, l0+=incr) {
+ wc->t[l0+n]=cosphi[l]* wc->t[l0+n-1] - sinphi[l]* wc->t[l0-n+1];
+ wc->t[l0-n]=sinphi[l]* wc->t[l0+n-1] + cosphi[l]* wc->t[l0-n+1];
+ }
+ }
+ free(cosphi);
+ free(sinphi);
+ }
+}
+
+