From 243b48abf7922957fc3e0876bf9f12f893a8327f Mon Sep 17 00:00:00 2001 From: Franz Zotter Date: Wed, 17 Jun 2009 08:39:06 +0000 Subject: modified the chebyshev12 routine to put out readily normalized circular harmonics and added [mtx_circular_harmonics] to the file mtx_spherical_harmonics, as well as a corresponding helpfile. why: useful for circular (1D) ambisonics, e.g. ... svn path=/trunk/externals/iem/iemmatrix/; revision=11792 --- src/mtx_spherical_harmonics.c | 121 ++++++++++++++++++++- src/mtx_spherical_harmonics/chebyshev12.c | 10 +- src/mtx_spherical_harmonics/sharmonics.c | 2 +- .../sharmonics_normalization.c | 9 ++ 4 files changed, 135 insertions(+), 7 deletions(-) (limited to 'src') diff --git a/src/mtx_spherical_harmonics.c b/src/mtx_spherical_harmonics.c index b0c270a..8ecc9a4 100644 --- a/src/mtx_spherical_harmonics.c +++ b/src/mtx_spherical_harmonics.c @@ -35,7 +35,22 @@ struct _MTXSh_ size_t l; }; -static void allocMTXdata (MTXSh *x) +static t_class *mtx_circular_harmonics_class; + +typedef struct _MTXCh_ MTXCh; +struct _MTXCh_ +{ + t_object x_obj; + t_outlet *list_ch_out; + t_atom *list_ch; + double *phi; + Cheby12WorkSpace *wc; + size_t nmax; + size_t l; +}; + + +static void allocMTXShdata (MTXSh *x) { x->phi=(double*)calloc(x->l,sizeof(double)); x->theta=(double*)calloc(x->l,sizeof(double)); @@ -104,7 +119,7 @@ static void mTXShMatrix (MTXSh *x, t_symbol *s, if (x->l!=columns) { deleteMTXShdata(x); x->l=columns; - allocMTXdata(x); + allocMTXShdata(x); } for (n=0;nl;n++) { x->phi[n]=(double) atom_getfloat(argv+n); @@ -125,6 +140,92 @@ static void mTXShMatrix (MTXSh *x, t_symbol *s, } +} + +static void allocMTXChdata (MTXCh *x) +{ + x->phi=(double*)calloc(x->l,sizeof(double)); + x->wc=chebyshev12_alloc(x->nmax,x->l); + x->list_ch=(t_atom*)calloc(x->l*(2*x->nmax+1)+2,sizeof(t_atom)); +} + +static void deleteMTXChdata (MTXCh *x) +{ + if (x->phi!=0) + free(x->phi); + if (x->list_ch!=0) + free(x->list_ch); + chebyshev12_free(x->wc); + x->wc=0; + x->list_ch=0; + x->phi=0; +} + +static void *newMTXCh (t_symbol *s, int argc, t_atom *argv) +{ + int nmax; + MTXCh *x = (MTXCh *) pd_new (mtx_circular_harmonics_class); + x->list_ch_out = outlet_new (&x->x_obj, gensym("matrix")); + x->list_ch = 0; + x->phi = 0; + x->wc = 0; + x->l=0; + nmax=(int) atom_getfloat(argv); + if (nmax<0) + nmax=0; + x->nmax=nmax; + + return ((void *) x); +} + +static void mTXChBang (MTXCh *x) +{ + if (x->list_ch!=0) { + outlet_anything(x->list_ch_out, gensym("matrix"), x->l*(2*x->nmax+1)+2, x->list_ch); + } +} + +static void mTXChMatrix (MTXCh *x, t_symbol *s, + int argc, t_atom *argv) +{ + int rows = atom_getint (argv++); + int columns = atom_getint (argv++); + int size = rows * columns; + int in_size = argc-2; + int n; + + + /* size check */ + if (!size) + post("mtx_circular_harmonics: invalid dimensions"); + else if (in_sizel!=columns) { + deleteMTXChdata(x); + x->l=columns; + allocMTXChdata(x); + } + for (n=0;nl;n++) { + x->phi[n]=(double) atom_getfloat(argv+n); + } + + if (x->wc!=0) { + chebyshev12(x->phi, x->wc); + in_size=x->l*(2*x->nmax+1); + SETFLOAT(x->list_ch,(float)x->l); + SETFLOAT(x->list_ch+1,(float)(2*x->nmax+1)); + for (n=0;nlist_ch+n+2,(float)x->wc->t[n]); + mTXChBang(x); + } + else + post("mtx_circular_harmonics: memory error, no operation"); + } + + } void mtx_spherical_harmonics_setup (void) @@ -139,6 +240,22 @@ void mtx_spherical_harmonics_setup (void) class_addmethod (mtx_spherical_harmonics_class, (t_method) mTXShMatrix, gensym("matrix"), A_GIMME,0); } + +void mtx_circular_harmonics_setup (void) +{ + mtx_circular_harmonics_class = class_new + (gensym("mtx_circular_harmonics"), + (t_newmethod) newMTXCh, + (t_method) deleteMTXChdata, + sizeof (MTXCh), + CLASS_DEFAULT, A_GIMME, 0); + class_addbang (mtx_circular_harmonics_class, (t_method) mTXChBang); + class_addmethod (mtx_circular_harmonics_class, (t_method) mTXChMatrix, gensym("matrix"), A_GIMME,0); +} + void iemtx_spherical_harmonics_setup(void){ + mtx_circular_harmonics_setup(); mtx_spherical_harmonics_setup(); } + + diff --git a/src/mtx_spherical_harmonics/chebyshev12.c b/src/mtx_spherical_harmonics/chebyshev12.c index 462341b..64af232 100644 --- a/src/mtx_spherical_harmonics/chebyshev12.c +++ b/src/mtx_spherical_harmonics/chebyshev12.c @@ -1,5 +1,5 @@ /* - Evaluates all circular harmonics + Evaluates all fully normalized circular harmonics at the angles phi up to the order nmax. using the recurrence for the Chebyshev polynomials of the first and second kind @@ -42,6 +42,8 @@ void chebyshev12(double *phi, Cheby12WorkSpace *wc) { const int incr=2*wc->nmax+1; double *cosphi; double *sinphi; + const double oneoversqrt2pi=1.0/sqrt(2.0*M_PI); + const double oneoversqrtpi=1.0/sqrt(M_PI); // memory allocation if ((wc!=0)&&(phi!=0)) { if ((cosphi=(double*)calloc(wc->l,sizeof(double)))==0) { @@ -56,9 +58,9 @@ void chebyshev12(double *phi, Cheby12WorkSpace *wc) { 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]; + wc->t[l0]=oneoversqrt2pi; + wc->t[l0+1]=cosphi[l]*oneoversqrtpi; + wc->t[l0-1]=sinphi[l]*oneoversqrtpi; } // recurrence for n>1 for (n=2; n<=wc->nmax; n++) { diff --git a/src/mtx_spherical_harmonics/sharmonics.c b/src/mtx_spherical_harmonics/sharmonics.c index 5d6f853..dbfcb18 100644 --- a/src/mtx_spherical_harmonics/sharmonics.c +++ b/src/mtx_spherical_harmonics/sharmonics.c @@ -43,7 +43,7 @@ static void sharmonics_initlegendrenormlzd(SHWorkSpace *ws) { } } -// multiplying Chebyshev sin/cos to the preliminary result +// multiplying normalized Chebyshev sin/cos to the preliminary result // Y_n^m(phi,theta) = Y_n^m(theta) * T_m(phi) // ny0 and nt0 denote where the position (n,m)=(n,0) or m=0 is in the arrays // ly0 and lt0 denote the starting position for one vertex in the arrays diff --git a/src/mtx_spherical_harmonics/sharmonics_normalization.c b/src/mtx_spherical_harmonics/sharmonics_normalization.c index ca97292..d7fe00a 100644 --- a/src/mtx_spherical_harmonics/sharmonics_normalization.c +++ b/src/mtx_spherical_harmonics/sharmonics_normalization.c @@ -27,8 +27,14 @@ SHNorml *sharmonics_normalization_new (const size_t nmax) { wn=0; } else { + /* + deprecated: // computing N_n^m for m=0, wrongly normalized wn->n[0]=sqrt(1/(2*M_PI)); + */ + + // computing N_n^m for m=0, + wn->n[0]=oneoversqrt2; for (n=1,n0=1; n<=nmax; n++) { wn->n[n0]=wn->n[0] * sqrt(2*n+1); n0+=n+1; @@ -40,11 +46,14 @@ SHNorml *sharmonics_normalization_new (const size_t nmax) { } n0+=n+1; } + /* + deprecated: // correcting normalization of N_n^0 for (n=0,n0=0; n<=nmax; n++) { wn->n[n0]*=oneoversqrt2; n0+=n+1; } + */ } } return wn; -- cgit v1.2.1