Merge branch 'develop' into master

README.md

@@ -3,14 +3,20 @@
 * git: [https://github.com/leomccormack/Spatial_Audio_Framework](https://github.com/leomccormack/Spatial_Audio_Framework)
 * doxygen: [https://leomccormack.github.io/Spatial_Audio_Framework/](https://leomccormack.github.io/Spatial_Audio_Framework/)
 
+# About
+
+The Spatial_Audio_Framework (SAF) is an open-source, cross-platform framework for developing spatial audio related applications in C/C++. While originally intended as a resource for researchers in the field, the framework has gradually grown into a comparatively large and comprehensively documented code-base comprising a number of distinct **modules**; with each module targeting a specific sub-field of spatial audio (e.g. Ambisonics encoding/decoding, spherical array processing, amplitude-panning, HRIR processing, room simulation, etc.). Several **examples** have also been included in the repository, which serve to demonstrate the functionality of the framework.
+
+Owing to its modular design, expanding the framework is also relatively straightforward; with contributions from researchers and developers of spatial audio technologies being actively encouraged! :-)
+
 # Prerequisites
 
-The Spatial_Audio_Framework (SAF) utilises the following external libraries:
+The framework utilises the following external libraries:
 * Any library (or libraries) conforming to the [CBLAS](https://en.wikipedia.org/wiki/Basic_Linear_Algebra_Subprograms#Implementations) and [LAPACK](https://en.wikipedia.org/wiki/LAPACK) standards
 * (**Optional**) [netCDF](https://www.unidata.ucar.edu/software/netcdf/) for reading [SOFA](https://www.sofaconventions.org/mediawiki/index.php/SOFA_(Spatially_Oriented_Format_for_Acoustics)) files
 * (**Optional**) Intel's [Integrated Performance Primitives (IPP)](https://software.intel.com/content/www/us/en/develop/tools/integrated-performance-primitives.html) for the FFT.
 
-In order to inform SAF which CBLAS/LAPACK supporting library/libraries you have linked to your project, simply add **one** of the following pre-processor definitions:
+In order to inform SAF which CBLAS/LAPACK supporting library/libraries you have linked to your project, simply add **one** of the following global pre-processor definitions:
 ```
 SAF_USE_INTEL_MKL             # great option, but only for x86 architectures    
 SAF_USE_OPEN_BLAS_AND_LAPACKE # good option, works on everything

examples/src/ambi_bin/ambi_bin.c

@@ -55,6 +55,8 @@ void ambi_bin_create
     *phAmbi = (void*)pData;
     int band;
 
+    printf(SAF_VERSION_LICENSE_STRING);
+
     /* default user parameters */
     for (band = 0; band<HYBRID_BANDS; band++)
         pData->EQ[band] = 1.0f;

examples/src/ambi_dec/ambi_dec.c

@@ -53,6 +53,8 @@ void ambi_dec_create
     *phAmbi = (void*)pData;
     int i, j, ch, band;
 
+    printf(SAF_VERSION_LICENSE_STRING);
+
     /* default user parameters */
     pData->masterOrder = pData->new_masterOrder = 1;
     for (band = 0; band<HYBRID_BANDS; band++)

examples/src/ambi_drc/ambi_drc.c

@@ -46,6 +46,8 @@ void ambi_drc_create
 {
     ambi_drc_data* pData = (ambi_drc_data*)malloc1d(sizeof(ambi_drc_data));
     *phAmbi = (void*)pData;
+
+    printf(SAF_VERSION_LICENSE_STRING);
 
     /* afSTFT stuff and audio buffers*/
     pData->hSTFT = NULL;

examples/src/ambi_enc/ambi_enc.c

@@ -33,6 +33,9 @@ void ambi_enc_create
     ambi_enc_data* pData = (ambi_enc_data*)malloc1d(sizeof(ambi_enc_data));
     *phAmbi = (void*)pData;
     int i;
+
+    printf(SAF_VERSION_LICENSE_STRING);
+
     pData->order = 1;
 
     /* default user parameters */

examples/src/ambi_roomsim/ambi_roomsim.c

@@ -33,6 +33,9 @@ void ambi_roomsim_create
     ambi_roomsim_data* pData = (ambi_roomsim_data*)malloc1d(sizeof(ambi_roomsim_data));
     *phAmbi = (void*)pData;
     int i;
+
+    printf(SAF_VERSION_LICENSE_STRING);
+
     pData->order = 1;
 
     /* default user parameters */ 

examples/src/array2sh/array2sh.c

@@ -56,6 +56,8 @@ void array2sh_create
     array2sh_data* pData = (array2sh_data*)malloc1d(sizeof(array2sh_data));
     *phA2sh = (void*)pData;
 
+    printf(SAF_VERSION_LICENSE_STRING);
+
     /* defualt parameters */
     array2sh_createArray(&(pData->arraySpecs)); 
     pData->filterType = FILTER_TIKHONOV;

examples/src/beamformer/beamformer.c

@@ -34,6 +34,8 @@ void beamformer_create
     *phBeam = (void*)pData;
     int i, ch;
 
+    printf(SAF_VERSION_LICENSE_STRING);
+
     /* default user parameters */
     pData->beamOrder = 1;
     for(i=0; i<MAX_NUM_BEAMS; i++){

examples/src/binauraliser/binauraliser.c

@@ -40,6 +40,8 @@ void binauraliser_create
     *phBin = (void*)pData;
     int ch;
 
+    printf(SAF_VERSION_LICENSE_STRING);
+
     /* user parameters */
     binauraliser_loadPreset(SOURCE_CONFIG_PRESET_DEFAULT, pData->src_dirs_deg, &(pData->new_nSources), &(pData->input_nDims)); /*check setStateInformation if you change default preset*/
     pData->nSources = pData->new_nSources;

examples/src/decorrelator/decorrelator.c

@@ -32,6 +32,8 @@ void decorrelator_create
     decorrelator_data* pData = (decorrelator_data*)malloc1d(sizeof(decorrelator_data));
     *phDecor = (void*)pData;
 
+    printf(SAF_VERSION_LICENSE_STRING);
+
     /* default user parameters */
     pData->nChannels = 1;
     pData->enableTransientDucker = 0;

examples/src/dirass/dirass.c

@@ -45,6 +45,8 @@ void dirass_create
     *phDir = (void*)pData;
     int i;
 
+    printf(SAF_VERSION_LICENSE_STRING);
+
     /* Default user parameters */
     pData->inputOrder = pData->new_inputOrder = SH_ORDER_FIRST;
     pData->beamType = STATIC_BEAM_TYPE_HYPERCARDIOID;

examples/src/matrixconv/matrixconv.c

@@ -32,6 +32,8 @@ void matrixconv_create
     matrixconv_data* pData = (matrixconv_data*)malloc1d(sizeof(matrixconv_data));
     *phMCnv = (void*)pData;
 
+    printf(SAF_VERSION_LICENSE_STRING);
+
     /* Default user parameters */
     pData->nInputChannels = 1;
     pData->enablePartitionedConv = 0;

examples/src/multiconv/multiconv.c

@@ -32,6 +32,8 @@ void multiconv_create
     multiconv_data* pData = (multiconv_data*)malloc1d(sizeof(multiconv_data));
     *phMCnv = (void*)pData;
 
+    printf(SAF_VERSION_LICENSE_STRING);
+
     /* Default user parameters */
     pData->nChannels = 1;
     pData->enablePartitionedConv = 0;

examples/src/panner/panner.c

@@ -51,6 +51,8 @@ void panner_create
     *phPan = (void*)pData;
     int ch, dummy;
 
+    printf(SAF_VERSION_LICENSE_STRING);
+
     /* default user parameters */
     panner_loadSourcePreset(SOURCE_CONFIG_PRESET_DEFAULT, pData->src_dirs_deg, &(pData->new_nSources), &(dummy)); /*check setStateInformation if you change default preset*/
     pData->nSources = pData->new_nSources;

examples/src/pitch_shifter/pitch_shifter.c

@@ -32,7 +32,9 @@ void pitch_shifter_create
 {
     pitch_shifter_data* pData = (pitch_shifter_data*)malloc1d(sizeof(pitch_shifter_data));
     *phPS = (void*)pData;
-
+
+    printf(SAF_VERSION_LICENSE_STRING);
+
     /* Default user parameters */
     pData->new_nChannels = pData->nChannels = 1;
     pData->pitchShift_factor = 0.5f; 

examples/src/powermap/powermap.c

@@ -40,6 +40,8 @@ void powermap_create
     *phPm = (void*)pData;
     int n, i, band;
 
+    printf(SAF_VERSION_LICENSE_STRING);
+
     /* Default user parameters */
     pData->masterOrder = pData->new_masterOrder = SH_ORDER_FIRST;
     for(band=0; band<HYBRID_BANDS; band++){

examples/src/rotator/rotator.c

@@ -40,6 +40,8 @@ void rotator_create
 {
     rotator_data* pData = (rotator_data*)malloc1d(sizeof(rotator_data));
     *phRot = (void*)pData;
+
+    printf(SAF_VERSION_LICENSE_STRING);
 
     pData->recalc_M_rotFLAG = 1;
 

examples/src/sldoa/sldoa.c

@@ -55,6 +55,8 @@ void sldoa_create
     *phSld = (void*)pData;
     int i, j, band;
 
+    printf(SAF_VERSION_LICENSE_STRING);
+
     /* Default user parameters */
     pData->new_masterOrder = pData->masterOrder = 1;
     for(band=0; band<HYBRID_BANDS; band++){

framework/include/saf.h

@@ -90,7 +90,7 @@
 
 #define SAF_VERSION_MAJOR 1        /**< Major version */
 #define SAF_VERSION_MINOR 1        /**< Minor version */
-#define SAF_VERSION_PATCH 4        /**< Patch version */
+#define SAF_VERSION_PATCH 5        /**< Patch version */
 #define SAF_VERSION_SPECIAL ""     /**< Append text ("alpha", "beta", "") */
 #define MKSTRING_(s) #s            /**< Stringify */
 #define MKSTRING(s) MKSTRING_(s)   /**< Stringify */
@@ -104,6 +104,11 @@
                            SAF_VERSION_MINOR) "." MKSTRING(SAF_VERSION_PATCH) \
                            "" SAF_VERSION_SPECIAL
 
+
+/** The Spatial_Audio_Framework Version and License as a string */
+#define SAF_VERSION_LICENSE_STRING "SAF Version: " SAF_VERSION_STRING \
+                                   ", License: "   SAF_LICENSE_STRING "\n"
+
 /** The Spatial_Audio_Framework Version and License as a banner */
 #define SAF_VERSION_BANNER \
   "     _____     _____     ______                                        \n" \

framework/modules/saf_hoa/saf_hoa.c

@@ -158,6 +158,9 @@ void getRSH
     int i, nSH;
     float scale;
     float* dirs_rad;
+
+    if(nDirs<1)
+        return;
 
     nSH = (N+1)*(N+1);
     scale = sqrtf(4.0f*SAF_PI);
@@ -191,6 +194,9 @@ void getRSH_recur
     float sleg_n[8], sleg_n_1[8], sleg_n_2[8], ssin_el, sfactorials_n[15];
     float* leg_n, *leg_n_1, *leg_n_2, *sin_el, *factorials_n;
 
+    if(nDirs<1)
+        return;
+
     if(N<=7 && nDirs == 1){
         /* Single direction optimisation for up to 7th order */
         leg_n = sleg_n;

framework/modules/saf_sh/saf_sh.c

@@ -198,6 +198,9 @@ void getSHreal
     double* Lnm, *CosSin;
     double *p_nm, *cos_incl;
     double *norm_real;
+
+    if(nDirs<1)
+        return;
 
     Lnm = malloc1d((2*order+1)*nDirs*sizeof(double));
     norm_real = malloc1d((2*order+1)*sizeof(double));
@@ -261,6 +264,9 @@ void getSHreal_recur
     float sleg_n[8], sleg_n_1[8], sleg_n_2[8], scos_incl, sfactorials_n[15];
     float* leg_n, *leg_n_1, *leg_n_2, *cos_incl, *factorials_n;
 
+    if(nDirs<1)
+        return;
+
     if(N<=7 && nDirs == 1){
         /* Single direction optimisation for up to 7th order */
         leg_n = sleg_n;

framework/modules/saf_tracker/saf_tracker.c

@@ -63,23 +63,23 @@ void tracker3d_create
     pData->tpars.init_birth = CLAMP(pData->tpars.init_birth, 0.0f, 0.99f);
     pData->tpars.alpha_death = CLAMP(pData->tpars.alpha_death, 1.0f, 20.0f);
     pData->tpars.beta_death = CLAMP(pData->tpars.beta_death, 1.0f, 20.0f);
-    pData->tpars.dt = MAX(pData->tpars.dt, 0.00001f);
-    pData->tpars.cd = MAX(pData->tpars.cd, 0.00001f);
+    pData->tpars.dt = MAX(pData->tpars.dt, 0.0001f);
+    pData->tpars.cd = MAX(pData->tpars.cd, 0.0001f);
     pData->tpars.W_avg_coeff = CLAMP(pData->tpars.W_avg_coeff, 0.0f, 0.99f);
-    pData->tpars.noiseSpecDen = MAX(pData->tpars.noiseSpecDen, 0.00001f);
+    pData->tpars.noiseSpecDen = MAX(pData->tpars.noiseSpecDen, 0.0001f);
     pData->tpars.noiseLikelihood = CLAMP(pData->tpars.noiseLikelihood, 0.0f, 0.99f);
-    pData->tpars.measNoiseSD = MAX(pData->tpars.measNoiseSD, 0.00001f);
+    pData->tpars.measNoiseSD = MAX(pData->tpars.measNoiseSD, 0.001f);
 
     /* Measurement noise PRIORs along the x,y,z axis, respectively  */
-    sd_xyz = tpars.measNoiseSD;
+    sd_xyz = pData->tpars.measNoiseSD;
     memset(pData->R, 0, 3*3*sizeof(float));
     pData->R[0][0] = powf(sd_xyz,2.0f);
     pData->R[1][1] = powf(sd_xyz,2.0f);
     pData->R[2][2] = powf(sd_xyz,2.0f);
 
     /* Noise spectral density along x, y, z axis qx,y,z which, (in combination
      * with sd_xyz), dictates how smooth the target tracks are. */
-    q_xyz = tpars.noiseSpecDen;
+    q_xyz = pData->tpars.noiseSpecDen;
 
     /* Dynamic and measurement models */
     const float F[6][6] =
@@ -93,25 +93,25 @@ void tracker3d_create
     Qc[3][3] = q_xyz;
     Qc[4][4] = q_xyz;
     Qc[5][5] = q_xyz;
-    lti_disc((float*)F, 6, 6, NULL, (float*)Qc, tpars.dt, (float*)pData->A, (float*)pData->Q);
+    lti_disc((float*)F, 6, 6, NULL, (float*)Qc, pData->tpars.dt, (float*)pData->A, (float*)pData->Q);
     memset(pData->H, 0, 3*6*sizeof(float));
     pData->H[0][0] = 1.0f;
     pData->H[1][1] = 1.0f;
     pData->H[2][2] = 1.0f;
 
     /* Create particles */
-    pData->SS = malloc1d(tpars.Np * sizeof(voidPtr));
-    pData->SS_resamp = malloc1d(tpars.Np * sizeof(voidPtr));
-    pData->W0 = 1.0f/(float)tpars.Np;
-    for(i=0; i<tpars.Np; i++){
-        tracker3d_particleCreate(&(pData->SS[i]), pData->W0, tpars.dt);
-        tracker3d_particleCreate(&(pData->SS_resamp[i]), pData->W0, tpars.dt);
+    pData->SS = malloc1d(pData->tpars.Np * sizeof(voidPtr));
+    pData->SS_resamp = malloc1d(pData->tpars.Np * sizeof(voidPtr));
+    pData->W0 = 1.0f/(float)pData->tpars.Np;
+    for(i=0; i<pData->tpars.Np; i++){
+        tracker3d_particleCreate(&(pData->SS[i]), pData->W0, pData->tpars.dt);
+        tracker3d_particleCreate(&(pData->SS_resamp[i]), pData->W0, pData->tpars.dt);
     }
 
     /* Event starting values */
     for(i=0; i<TRACKER3D_MAX_NUM_EVENTS; i++){
         pData->evta[i] = NULL;
-        tracker3d_particleCreate(&(pData->str[i]), pData->W0, tpars.dt);
+        tracker3d_particleCreate(&(pData->str[i]), pData->W0, pData->tpars.dt);
     }
     pData->incrementTime = 0;
 }

framework/modules/saf_utilities/saf_utility_qmf.c

@@ -542,11 +542,12 @@ void qmf_clearBuffers
     int i;
 
     /* flush analysis buffers */
-    for(i=0; i<h->nCHin; i++)
+    for(i=0; i<h->nCHin; i++){
         memset(h->buffer_ana[i], 0, h->hopsize * 10 * sizeof(float));
-    if(h->hybridmode){
-        memset(FLATTEN3D(h->qmfDelayBuffer), 0, h->nCHin*(h->hopsize-QMF_NBANDS_2_SUBDIVIDE)* ((QMF_HYBRID_FILTER_LENGTH-1)/2 + 1)*sizeof(float_complex));
-        memset(FLATTEN3D(h->hybBuffer), 0, h->nCHin*QMF_NBANDS_2_SUBDIVIDE*QMF_HYBRID_FILTER_LENGTH*sizeof(float_complex));
+        if(h->hybridmode){
+            memset(FLATTEN3D(h->qmfDelayBuffer), 0, h->nCHin*(h->hopsize-QMF_NBANDS_2_SUBDIVIDE)* ((QMF_HYBRID_FILTER_LENGTH-1)/2 + 1)*sizeof(float_complex));
+            memset(FLATTEN3D(h->hybBuffer), 0, h->nCHin*QMF_NBANDS_2_SUBDIVIDE*QMF_HYBRID_FILTER_LENGTH*sizeof(float_complex));
+        }
     }
 
     /* flush synthesis buffers */

framework/modules/saf_utilities/saf_utility_veclib.c

@@ -1749,6 +1749,61 @@ void utility_sglslvt
     free(b);
 }
 
+void utility_cglslvt
+(
+    const float_complex* A,
+    const int dim,
+    float_complex* B,
+    int nCol,
+    float_complex* X
+)
+{
+    veclib_int i;
+    veclib_int n = dim, nrhs = nCol, lda = dim, ldb = dim, info;
+    veclib_int* IPIV;
+    float_complex* a, *b;
+
+    assert(0); /* This function needs checking */
+
+    IPIV = malloc1d(dim*sizeof(veclib_int));
+    a = malloc1d(dim*dim*sizeof(float_complex));
+    b = malloc1d(dim*nrhs*sizeof(float_complex));
+
+    /* store locally - Hermitian */
+    cblas_ccopy(dim*dim, A, 1, a, 1);
+    for(i=0; i<dim*dim; i++)
+        a[i] = conj(a[i]);
+    cblas_ccopy(dim*nCol, B, 1, b, 1);
+    for(i=0; i<dim*nCol; i++)
+        b[i] = conj(b[i]);
+
+    /* solve Ax = b for each column in b (b is replaced by the solution: x) */
+#ifdef VECLIB_USE_CLAPACK_INTERFACE
+    info = clapack_cgesv(CblasColMajor, n, nrhs, a, ldb, IPIV, b, lda);
+#elif defined(VECLIB_USE_LAPACKE_INTERFACE)
+    info = LAPACKE_cgesv(CblasColMajor, n, nrhs, a, ldb, IPIV, b, lda );
+#elif defined(VECLIB_USE_LAPACK_FORTRAN_INTERFACE)
+    cgesv_( &n, &nrhs, (veclib_float_complex*)a, &ldb, IPIV, (veclib_float_complex*)b, &lda, &info );
+#endif
+
+    if(info!=0){
+        /* A is singular, solution not possible */
+        memset(X, 0, dim*nCol*sizeof(float_complex));
+#ifndef NDEBUG
+        saf_error_print(SAF_WARNING__FAILED_TO_SOLVE_LINEAR_EQUATION);
+#endif
+    }
+    else{
+        cblas_ccopy(dim*nCol, b, 1, X, 1);
+        for(i=0; i<dim*nCol; i++)
+            X[i] = conj(b[i]);
+    }
+
+    free(IPIV);
+    free(a);
+    free(b);
+}
+
 
 /* ========================================================================== */
 /*                      Symmetric Linear Solver (?slslv)                      */