SAF
Loading...
Searching...
No Matches
spreader.c
Go to the documentation of this file.
1/*
2 * Copyright 2021 Leo McCormack
3 *
4 * Permission to use, copy, modify, and/or distribute this software for any
5 * purpose with or without fee is hereby granted, provided that the above
6 * copyright notice and this permission notice appear in all copies.
7 *
8 * THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL WARRANTIES WITH
9 * REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY
10 * AND FITNESS. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT,
11 * INDIRECT, OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM
12 * LOSS OF USE, DATA OR PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR
13 * OTHER TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR
14 * PERFORMANCE OF THIS SOFTWARE.
15 */
16
31
32#include "spreader_internal.h"
33
35(
36 void ** const phSpr
37)
38{
40 *phSpr = (void*)pData;
41 int band, t, src;
42
43 /* user parameters */
44 pData->sofa_filepath = NULL;
45 pData->nSources = 1;
47 pData->useDefaultHRIRsFLAG = 1;
48 pData->covAvgCoeff = 0.85f;
49 memset(pData->src_spread, 0, SPREADER_MAX_NUM_SOURCES*sizeof(float));
50 memset(pData->src_dirs_deg, 0, SPREADER_MAX_NUM_SOURCES*2*sizeof(float));
51
52 /* time-frequency transform + buffers */
53 pData->fs = 48000.0f;
54 pData->hSTFT = NULL;
55 pData->inputFrameTD = (float**)malloc2d(MAX_NUM_INPUTS, SPREADER_FRAME_SIZE, sizeof(float));
56 pData->outframeTD = (float**)malloc2d(MAX_NUM_OUTPUTS, SPREADER_FRAME_SIZE, sizeof(float));
57 pData->inputframeTF = (float_complex***)malloc3d(HYBRID_BANDS, MAX_NUM_INPUTS, TIME_SLOTS, sizeof(float_complex));
58 pData->protoframeTF = (float_complex***)malloc3d(HYBRID_BANDS, MAX_NUM_OUTPUTS, TIME_SLOTS, sizeof(float_complex));
59 pData->decorframeTF = (float_complex***)malloc3d(HYBRID_BANDS, MAX_NUM_OUTPUTS, TIME_SLOTS, sizeof(float_complex));
60 pData->spreadframeTF = (float_complex***)malloc3d(HYBRID_BANDS, MAX_NUM_OUTPUTS, TIME_SLOTS, sizeof(float_complex));
61 pData->outputframeTF = (float_complex***)malloc3d(HYBRID_BANDS, MAX_NUM_OUTPUTS, TIME_SLOTS, sizeof(float_complex));
62
63 /* Internal */
64 pData->Q = pData->nGrid = pData->h_len = 0;
65 pData->h_fs = 0.0f;
66 pData->h_grid = NULL;
67 pData->H_grid = NULL;
68 for(band=0; band<HYBRID_BANDS; band++)
69 pData->HHH[band] = NULL;
70 pData->grid_dirs_deg = NULL;
71 pData->grid_dirs_xyz = NULL;
72 pData->weights = NULL;
73 pData->angles = NULL;
74 for(src=0; src<SPREADER_MAX_NUM_SOURCES; src++){
75 pData->hDecor[src] = NULL;
76 pData->Cy[src] = NULL;
77 pData->Cproto[src] = NULL;
78 pData->prev_M[src] = NULL;
79 pData->prev_Mr[src] = NULL;
80 pData->dirActive[src] = NULL;
81 }
82 pData->new_M = NULL;
83 pData->new_Mr = NULL;
84 pData->interp_M = NULL;
85 pData->interp_Mr = NULL;
86 pData->interp_Mr_cmplx = NULL;
87 for(t=0; t<TIME_SLOTS; t++){
88 pData->interpolatorFadeIn[t] = ((float)t+1.0f)/(float)TIME_SLOTS;
89 pData->interpolatorFadeOut[t] = 1.0f - ((float)t+1.0f)/(float)TIME_SLOTS;
90 }
91
92 // Hotfix:
93 pData->_tmpFrame = malloc1d(MAX_NUM_CHANNELS*TIME_SLOTS *sizeof(float_complex));
94 pData->_H_tmp = malloc1d(MAX_NUM_CHANNELS *sizeof(float_complex));
95 pData->_Cy = malloc1d(MAX_NUM_CHANNELS*MAX_NUM_CHANNELS *sizeof(float_complex));
96 pData->_E_dir = malloc1d(MAX_NUM_CHANNELS*MAX_NUM_CHANNELS *sizeof(float_complex));
97 pData->_V = malloc1d(MAX_NUM_OUTPUTS*MAX_NUM_OUTPUTS *sizeof(float_complex));
98 pData->_D = malloc1d(MAX_NUM_OUTPUTS*MAX_NUM_OUTPUTS *sizeof(float_complex));
99 pData->_Cproto = malloc1d(MAX_NUM_OUTPUTS*MAX_NUM_OUTPUTS *sizeof(float_complex));
100
101 /* Optimal mixing */
102 pData->hCdf = NULL;
103 pData->hCdf_res = NULL;
104 pData->Qmix = NULL;
105 pData->Qmix_cmplx = NULL;
106 pData->Cr = NULL;
107 pData->Cr_cmplx = NULL;
108
109 /* flags/status */
110 pData->new_procMode = pData->procMode;
111 pData->new_nSources = pData->nSources;
112 pData->progressBar0_1 = 0.0f;
114 strcpy(pData->progressBarText,"");
117}
118
120(
121 void ** const phSpr
122)
123{
124 spreader_data *pData = (spreader_data*)(*phSpr);
125 int band, src;
126
127 if (pData != NULL) {
128 /* not safe to free memory during intialisation/processing loop */
129 while (pData->codecStatus == CODEC_STATUS_INITIALISING ||
131 SAF_SLEEP(10);
132 }
133
134 free(pData->sofa_filepath);
135
136 /* free afSTFT and buffers */
137 if(pData->hSTFT !=NULL)
138 afSTFT_destroy(&(pData->hSTFT));
139 free(pData->inputFrameTD);
140 free(pData->outframeTD);
141 free(pData->inputframeTF);
142 free(pData->decorframeTF);
143 free(pData->spreadframeTF);
144 free(pData->outputframeTF);
145
146 /* internal */
147 free(pData->h_grid);
148 free(pData->H_grid);
149 for(band=0; band<HYBRID_BANDS; band++)
150 free(pData->HHH[band]);
151 free(pData->grid_dirs_deg);
152 free(pData->grid_dirs_xyz);
153 free(pData->weights);
154 free(pData->angles);
155 for(src=0; src<SPREADER_MAX_NUM_SOURCES; src++){
156 latticeDecorrelator_destroy(&(pData->hDecor[src]));
157 free(pData->Cy[src]);
158 free(pData->Cproto[src]);
159 free(pData->prev_M[src]);
160 free(pData->prev_Mr[src]);
161 free(pData->dirActive[src]);
162 }
163 free(pData->new_M);
164 free(pData->new_Mr);
165 free(pData->interp_M);
166 free(pData->interp_Mr);
167 free(pData->interp_Mr_cmplx);
168
169 // Hotfix:
170 free(pData->_tmpFrame);
171 free(pData->_H_tmp);
172 free(pData->_Cy);
173 free(pData->_E_dir);
174 free(pData->_V);
175 free(pData->_D);
176 free(pData->_Cproto);
177
178 /* Optimal mixing */
179 cdf4sap_cmplx_destroy(&(pData->hCdf));
180 cdf4sap_destroy(&(pData->hCdf_res));
181 free(pData->Qmix);
182 free(pData->Qmix_cmplx);
183 free(pData->Cr);
184 free(pData->Cr_cmplx);
185
186 free(pData->progressBarText);
187
188 free(pData);
189 pData = NULL;
190 *phSpr = NULL;
191 }
192}
193
195(
196 void * const hSpr,
197 int sampleRate
198)
199{
200 spreader_data *pData = (spreader_data*)(hSpr);
201
202 /* define frequency vector */
203 pData->fs = sampleRate;
204 afSTFT_getCentreFreqs(pData->hSTFT, (float)sampleRate, HYBRID_BANDS, pData->freqVector);
205}
206
208(
209 void* const hSpr
210)
211{
212 spreader_data *pData = (spreader_data*)(hSpr);
213 int q, band, ng, nSources, src;
214 float_complex scaleC;
215#ifdef SAF_ENABLE_SOFA_READER_MODULE
218#endif
219 SPREADER_PROC_MODES procMode;
220 float_complex H_tmp[MAX_NUM_CHANNELS];
221 const float_complex calpha = cmplxf(1.0f, 0.0f), cbeta = cmplxf(0.0f, 0.0f);
222
224 return; /* re-init not required, or already happening */
225 while (pData->procStatus == PROC_STATUS_ONGOING){
226 /* re-init required, but we need to wait for the current processing loop to end */
227 pData->codecStatus = CODEC_STATUS_INITIALISING; /* indicate that we want to init */
228 SAF_SLEEP(10);
229 }
230
231 nSources = pData->new_nSources;
232 procMode = pData->new_procMode;
233
234 /* for progress bar */
236 strcpy(pData->progressBarText,"Initialising");
237 pData->progressBar0_1 = 0.0f;
238
239 /* Load measurements (e.g. HRIRs, Microphone array IRs etc.) */
240#ifndef SAF_ENABLE_SOFA_READER_MODULE
241 pData->useDefaultHRIRsFLAG = 1;
242#endif
243 if(pData->useDefaultHRIRsFLAG){
244 /* Load default HRIR data */
245 pData->Q = NUM_EARS;
247 pData->h_len = __default_hrir_len;
248 pData->h_fs = (float)__default_hrir_fs;
249 pData->h_grid = realloc1d(pData->h_grid, pData->nGrid * (pData->Q) * (pData->h_len) * sizeof(float));
250 memcpy(pData->h_grid, (float*)__default_hrirs, pData->nGrid * (pData->Q) * (pData->h_len) * sizeof(float));
251 pData->grid_dirs_deg = realloc1d(pData->grid_dirs_deg, pData->nGrid * 2 * sizeof(float));
252 memcpy(pData->grid_dirs_deg, (float*)__default_hrir_dirs_deg, pData->nGrid * 2 * sizeof(float));
253 }
254#ifdef SAF_ENABLE_SOFA_READER_MODULE
255 else{
256 /* Use sofa loader */
258 if(error!=SAF_SOFA_OK){
259 /* if failed, then load default data instead */
260 pData->useDefaultHRIRsFLAG = 1;
261 saf_print_warning("Unable to load the specified SOFA file. Using default HRIR data instead");
262 spreader_initCodec(hSpr);
263 }
264 pData->h_fs = sofa.DataSamplingRate;
265 pData->h_len = sofa.DataLengthIR;
266 pData->nGrid = sofa.nSources;
267 pData->h_grid = realloc1d(pData->h_grid, pData->nGrid*(pData->Q)*(pData->h_len)*sizeof(float));
268 memcpy(pData->h_grid, sofa.DataIR, pData->nGrid*(pData->Q)*(pData->h_len)*sizeof(float));
269 pData->grid_dirs_deg = realloc1d(pData->grid_dirs_deg, pData->nGrid*2*sizeof(float));
270 cblas_scopy(pData->nGrid, sofa.SourcePosition, 3, pData->grid_dirs_deg, 2); /* azi */
271 cblas_scopy(pData->nGrid, &sofa.SourcePosition[1], 3, &pData->grid_dirs_deg[1], 2); /* elev */
272 saf_sofa_close(&sofa);
273 }
274#endif
275
276 /* Convert from the 0..360 convention, to -180..180, and pre-compute unit Cartesian vectors */
278 pData->grid_dirs_xyz = realloc1d(pData->grid_dirs_xyz, pData->nGrid*3*sizeof(float));
279 unitSph2cart(pData->grid_dirs_deg, pData->nGrid, 1, pData->grid_dirs_xyz);
280
281 /* Initialise time-frequency transform and decorrelators */
282 afSTFT_destroy(&(pData->hSTFT));
283 afSTFT_create(&(pData->hSTFT), nSources, pData->Q, HOP_SIZE, 0, 1, AFSTFT_BANDS_CH_TIME);
284 int orders[4] = {20, 15, 6, 6}; /* 20th order up to 700Hz, 15th->2.4kHz, 6th->4kHz, 3rd->12kHz, NONE(only delays)->Nyquist */
285 //float freqCutoffs[4] = {600.0f, 2.6e3f, 4.5e3f, 12e3f};
286 float freqCutoffs[4] = {900.0f, 6.8e3f, 12e3f, 24e3f};
287 const int maxDelay = 12;
288 for(src=0; src<SPREADER_MAX_NUM_SOURCES; src++){
289 latticeDecorrelator_destroy(&(pData->hDecor[src]));
290 latticeDecorrelator_create(&(pData->hDecor[src]), (float)pData->fs, HOP_SIZE, pData->freqVector, HYBRID_BANDS, pData->Q, orders, freqCutoffs, 4, maxDelay, 0, 0.75f);
291 }
292
293 /* Convert to filterbank coefficients and pre-compute outer products */
294 pData->H_grid = realloc1d(pData->H_grid, HYBRID_BANDS*(pData->Q)*pData->nGrid*sizeof(float_complex));
295 afSTFT_FIRtoFilterbankCoeffs(pData->h_grid, pData->nGrid, pData->Q, pData->h_len, HOP_SIZE, 0, 1, pData->H_grid);
296 pData->weights = realloc1d(pData->weights, pData->nGrid*sizeof(float));
297 getVoronoiWeights(pData->grid_dirs_deg, pData->nGrid, 0, pData->weights);
298 cblas_sscal(pData->nGrid, 1.0f/FOURPI, pData->weights, 1);
299 for(band=0; band<HYBRID_BANDS; band++){
300 pData->HHH[band] = (float_complex**)realloc2d((void**)pData->HHH[band], pData->nGrid, pData->Q * (pData->Q), sizeof(float_complex));
301 for(ng=0; ng<pData->nGrid; ng++){
302 for(q=0; q<pData->Q; q++)
303 H_tmp[q] = pData->H_grid[band*(pData->Q)*pData->nGrid + q*pData->nGrid + ng];
304 cblas_cgemm(CblasRowMajor, CblasNoTrans, CblasConjTrans, pData->Q, pData->Q, 1, &calpha,
305 H_tmp, 1,
306 H_tmp, 1, &cbeta,
307 pData->HHH[band][ng], pData->Q);
308 scaleC = cmplxf(pData->weights[ng], 0.0f);
309 cblas_cscal(pData->Q * (pData->Q), &scaleC, pData->HHH[band][ng], 1);
310 }
311 }
312 pData->angles = realloc1d(pData->angles, pData->nGrid*sizeof(float));
313
314 /* OM structures */
315 cdf4sap_cmplx_destroy(&(pData->hCdf));
316 cdf4sap_cmplx_create(&(pData->hCdf), pData->Q, pData->Q);
317 cdf4sap_destroy(&(pData->hCdf_res));
318 cdf4sap_create(&(pData->hCdf_res), pData->Q, pData->Q);
319 pData->Qmix = realloc1d(pData->Qmix, pData->Q*(pData->Q)*sizeof(float));
320 memset(pData->Qmix, 0, pData->Q*(pData->Q)*sizeof(float));
321 pData->Qmix_cmplx = realloc1d(pData->Qmix_cmplx, pData->Q*(pData->Q)*sizeof(float_complex));
322 memset(pData->Qmix_cmplx, 0, pData->Q*(pData->Q)*sizeof(float_complex));
323 for(q=0; q<pData->Q; q++){
324 pData->Qmix[q*(pData->Q)+q] = 1.0f;
325 pData->Qmix_cmplx[q*(pData->Q)+q] = cmplxf(1.0f, 0.0f);
326 }
327 pData->Cr = realloc1d(pData->Cr, pData->Q*(pData->Q)*sizeof(float));
328 pData->Cr_cmplx = realloc1d(pData->Cr_cmplx, pData->Q*(pData->Q)*sizeof(float_complex));
329
330 /* mixing matrices and buffers */
331 for(src=0; src<SPREADER_MAX_NUM_SOURCES; src++){
332 pData->Cy[src] = (float_complex**)realloc2d((void**)pData->Cy[src], HYBRID_BANDS, (pData->Q)*(pData->Q), sizeof(float_complex));
333 memset(FLATTEN2D(pData->Cy[src]), 0, HYBRID_BANDS * (pData->Q)*(pData->Q) * sizeof(float_complex));
334 pData->Cproto[src] = (float_complex**)realloc2d((void**)pData->Cproto[src], HYBRID_BANDS, (pData->Q)*(pData->Q), sizeof(float_complex));
335 memset(FLATTEN2D(pData->Cproto[src]), 0, HYBRID_BANDS * (pData->Q)*(pData->Q) * sizeof(float_complex));
336 pData->prev_M[src] = (float_complex**)realloc2d((void**)pData->prev_M[src], HYBRID_BANDS, (pData->Q)*(pData->Q), sizeof(float_complex));
337 memset(FLATTEN2D(pData->prev_M[src]), 0, HYBRID_BANDS * (pData->Q)*(pData->Q) * sizeof(float_complex));
338 pData->prev_Mr[src] = (float**)realloc2d((void**)pData->prev_Mr[src], HYBRID_BANDS, (pData->Q)*(pData->Q), sizeof(float));
339 memset(FLATTEN2D(pData->prev_Mr[src]), 0, HYBRID_BANDS * (pData->Q)*(pData->Q) * sizeof(float));
340 pData->dirActive[src] = realloc1d(pData->dirActive[src], pData->nGrid * sizeof(int));
341 memset(pData->dirActive[src], 0, pData->nGrid*sizeof(int));
342 }
343 pData->new_M = (float_complex**)realloc2d((void**)pData->new_M, HYBRID_BANDS, (pData->Q)*(pData->Q), sizeof(float_complex));
344 pData->new_Mr = (float**)realloc2d((void**)pData->new_Mr, HYBRID_BANDS, (pData->Q)*(pData->Q), sizeof(float));
345 pData->interp_M = realloc1d(pData->interp_M, (pData->Q)*(pData->Q) * sizeof(float_complex));
346 pData->interp_Mr = realloc1d(pData->interp_Mr, (pData->Q)*(pData->Q) * sizeof(float));
347 pData->interp_Mr_cmplx = realloc1d(pData->interp_Mr_cmplx, (pData->Q)*(pData->Q) * sizeof(float_complex));
348 memset(pData->interp_Mr_cmplx, 0, (pData->Q)*(pData->Q) * sizeof(float_complex));
349
350 /* New config */
351 pData->nSources = nSources;
352 pData->procMode = procMode;
353
354 /* done! */
355 strcpy(pData->progressBarText,"Done!");
356 pData->progressBar0_1 = 1.0f;
358}
359
361(
362 void * const hSpr,
363 const float *const * inputs,
364 float* const* const outputs,
365 int nInputs,
366 int nOutputs,
367 int nSamples
368)
369{
370 spreader_data *pData = (spreader_data*)(hSpr);
371 int q, src, ng, ch, i, j, band, t, nSources, Q, centre_ind, nSpread;
372 float trace, Ey, Eproto, Gcomp;
373 float src_dirs_deg[SPREADER_MAX_NUM_SOURCES][2], src_dir_xyz[3], CprotoDiag[MAX_NUM_OUTPUTS*MAX_NUM_OUTPUTS], src_spread[MAX_NUM_OUTPUTS];
374 float_complex scaleC, tmp;
375#if 0
376 float_complex Cx[MAX_NUM_OUTPUTS*MAX_NUM_OUTPUTS];
377 float CxDiag[MAX_NUM_OUTPUTS*MAX_NUM_OUTPUTS];
378#endif
379 SPREADER_PROC_MODES procMode;
380 const float_complex calpha = cmplxf(1.0f, 0.0f), cbeta = cmplxf(0.0f, 0.0f);
381
382 /* copy user parameters to local variables */
383 procMode = pData->procMode;
384 nSources = pData->nSources;
385 Q = pData->Q;
386 for(i=0; i<nSources; i++){
387 src_dirs_deg[i][0] = pData->src_dirs_deg[i][0];
388 src_dirs_deg[i][1] = pData->src_dirs_deg[i][1];
389 src_spread[i] = pData->src_spread[i];
390 }
391
392 /* apply binaural panner */
393 if ((nSamples == SPREADER_FRAME_SIZE) && (pData->codecStatus==CODEC_STATUS_INITIALISED) ){
395
396 /* Load time-domain data */
397 for(i=0; i < SAF_MIN(nSources,nInputs); i++)
398 utility_svvcopy(inputs[i], SPREADER_FRAME_SIZE, pData->inputFrameTD[i]);
399 for(; i<nSources; i++)
400 memset(pData->inputFrameTD[i], 0, SPREADER_FRAME_SIZE * sizeof(float));
401
402 /* Apply time-frequency transform (TFT) */
404
405 /* Zero output buffer */
406 for(band=0; band<HYBRID_BANDS; band++)
407 memset(FLATTEN2D(pData->outputframeTF[band]), 0, Q*TIME_SLOTS*sizeof(float_complex));
408
409 /* Loop over sources */
410 for(src=0; src<nSources; src++){
411 /* Find the "spread" indices */
412 unitSph2cart(src_dirs_deg[src], 1, 1, src_dir_xyz);
413 cblas_sgemm(CblasRowMajor, CblasNoTrans, CblasNoTrans, pData->nGrid, 1, 3, 1.0f,
414 pData->grid_dirs_xyz, 3,
415 src_dir_xyz, 1, 0.0f,
416 pData->angles, 1);
417 for(i=0; i<pData->nGrid; i++)
418 pData->angles[i] = acosf(SAF_MIN(pData->angles[i], 0.9999999f))*180.0f/SAF_PI;
419 utility_siminv(pData->angles, pData->nGrid, &centre_ind);
420
421 /* Define Prototype signals */
422 switch(procMode){
423 case SPREADER_MODE_NAIVE: /* fall through */
424 case SPREADER_MODE_OM:
425 for(band=0; band<HYBRID_BANDS; band++){
426 if(pData->freqVector[band]<MAX_SPREAD_FREQ){
427 /* Loop over all angles, and sum the H_grid's within the spreading area */
428 memset(pData->_H_tmp, 0, Q*sizeof(float_complex));
429 for(ng=0,nSpread=0; ng<pData->nGrid; ng++){
430 if(pData->angles[ng] <= (src_spread[src]/2.0f)){
431 for(q=0; q<Q; q++)
432 pData->_H_tmp[q] = ccaddf(pData->_H_tmp[q], pData->H_grid[band*Q*pData->nGrid + q*pData->nGrid + ng]);
433 nSpread++;
434 pData->dirActive[src][ng] = 1;
435 }
436 else
437 pData->dirActive[src][ng] = 0;
438 }
439 }
440 else
441 nSpread = 0;
442
443 /* If no directions found in the spread area, then just include the nearest one */
444 if(nSpread==0){
445 for(q=0; q<Q; q++)
446 pData->_H_tmp[q] = pData->H_grid[band*Q*pData->nGrid + q*pData->nGrid + centre_ind];
447 nSpread=1;
448 }
449
450 /* Apply */
451 cblas_cgemm(CblasRowMajor, CblasNoTrans, CblasNoTrans, Q, TIME_SLOTS, 1, &calpha,
452 pData->_H_tmp, 1,
453 pData->inputframeTF[band][src], TIME_SLOTS, &cbeta,
454 FLATTEN2D(pData->protoframeTF[band]), TIME_SLOTS);
455
456 /* Scale by number of spreading directions */
457 cblas_sscal(/*re+im*/2*Q*TIME_SLOTS, 1.0f/(float)nSpread, (float*)FLATTEN2D(pData->protoframeTF[band]), 1);
458 }
459 break;
460#if 0
461 case SPREADER_MODE_OM:
462 /* Use the centre direction as the prototype */
463 for(band=0; band<HYBRID_BANDS; band++){
464 for(q=0; q<Q; q++)
465 H_tmp[q] = pData->H_grid[band*Q*pData->nGrid + q*pData->nGrid + centre_ind];
466 cblas_cgemm(CblasRowMajor, CblasNoTrans, CblasNoTrans, Q, TIME_SLOTS, 1, &calpha,
467 H_tmp, 1,
468 pData->inputframeTF[band][src], TIME_SLOTS, &cbeta,
469 FLATTEN2D(pData->protoframeTF[band]), TIME_SLOTS);
470 }
471 break;
472#endif
473
475 /* Replicate the mono signal for all Q channels */
476 for(band=0; band<HYBRID_BANDS; band++)
477 for(q=0; q<Q; q++)
478 memcpy(pData->protoframeTF[band][q], pData->inputframeTF[band][src], TIME_SLOTS*sizeof(float_complex));
479 break;
480 }
481
482 /* Main processing */
483 if(procMode==SPREADER_MODE_NAIVE) {
484 /* If naive mode, then we're already done... */
485 for(band=0; band<HYBRID_BANDS; band++)
486 memcpy(FLATTEN2D(pData->spreadframeTF[band]), FLATTEN2D(pData->protoframeTF[band]), Q*TIME_SLOTS*sizeof(float_complex));
487 }
488 else{
489 /* Apply decorrelation of prototype signals */
491
492 /* Compute prototype covariance matrix and average over time */
493 for(band=0; band<HYBRID_BANDS; band++){
494 cblas_cgemm(CblasRowMajor, CblasNoTrans, CblasConjTrans, Q, Q, TIME_SLOTS, &calpha,
495 FLATTEN2D(pData->protoframeTF[band]), TIME_SLOTS,
496 FLATTEN2D(pData->protoframeTF[band]), TIME_SLOTS, &cbeta,
497 pData->_Cproto, Q);
498 cblas_sscal(/*re+im*/2*Q*Q, pData->covAvgCoeff, (float*)pData->Cproto[src][band], 1);
499 cblas_saxpy(/*re+im*/2*Q*Q, 1.0f-pData->covAvgCoeff, (float*)pData->_Cproto, 1, (float*)pData->Cproto[src][band], 1);
500 }
501
502 /* Define target covariance matrices */
503 for(band=0; band<HYBRID_BANDS; band++){
504 /* Sum the H_array outer product matrices for the whole spreading area */
505 if(pData->freqVector[band]<MAX_SPREAD_FREQ){
506 memset(pData->_Cy, 0, Q*Q*sizeof(float_complex));
507 memset(pData->_H_tmp, 0, Q*sizeof(float_complex));
508 for(ng=0, nSpread=0; ng<pData->nGrid; ng++){
509 if(pData->angles[ng] <= (src_spread[src]/2.0f)){
510 cblas_caxpy(Q*Q, &calpha, pData->HHH[band][ng], 1, pData->_Cy, 1);
511 for(q=0; q<Q; q++)
512 pData->_H_tmp[q] = ccaddf(pData->_H_tmp[q], pData->H_grid[band*Q*pData->nGrid + q*pData->nGrid + ng]);
513 nSpread++;
514 pData->dirActive[src][ng] = 1;
515 }
516 else
517 pData->dirActive[src][ng] = 0;
518 }
519 }
520 else
521 nSpread = 0;
522
523 /* If no directions found in the spread area, then just include the nearest one */
524 if(nSpread==0) {
525 cblas_caxpy(Q*Q, &calpha, pData->HHH[band][centre_ind], 1, pData->_Cy, 1);
526 for(q=0; q<Q; q++)
527 pData->_H_tmp[q] = pData->H_grid[band*Q*pData->nGrid + q*pData->nGrid + centre_ind];
528 nSpread++;
529 }
530#if 1
531 /* Impose target energies too */
532 if (procMode == SPREADER_MODE_OM && pData->freqVector[band]<MAX_SPREAD_FREQ){
533 /* Normalise Cy */
534 trace = 0.0f;
535 for(q=0; q<Q; q++)
536 trace += crealf(pData->_Cy[q*Q+q]);
537 cblas_sscal(/*re+im*/2*Q*Q, 1.0f/(trace+2.23e-9f), (float*)pData->_Cy, 1);
538
539 /* Compute signals for the centre of the spread */
540 for(q=0; q<Q; q++)
541 pData->_H_tmp[q] = pData->H_grid[band*Q*pData->nGrid + q*pData->nGrid + centre_ind];
542 cblas_cgemm(CblasRowMajor, CblasNoTrans, CblasNoTrans, Q, TIME_SLOTS, 1, &calpha,
543 pData->_H_tmp, 1,
544 pData->inputframeTF[band][src], TIME_SLOTS, &cbeta,
545 pData->_tmpFrame, TIME_SLOTS);
546
547 /* Introduce their channel energies into the target */
548 cblas_cgemm(CblasRowMajor, CblasNoTrans, CblasConjTrans, Q, Q, TIME_SLOTS, &calpha,
549 pData->_tmpFrame, TIME_SLOTS,
550 pData->_tmpFrame, TIME_SLOTS, &cbeta,
551 pData->_E_dir, Q);
552 trace = 0.0f;
553 for(q=0; q<Q; q++)
554 trace += crealf(pData->_E_dir[q*Q+q]);
555 cblas_sscal(/*re+im*/2*Q*Q, trace, (float*)pData->_Cy, 1);
556 }
557#endif
558 /* Average over time */
559 cblas_sscal(/*re+im*/2*Q*Q, pData->covAvgCoeff, (float*)pData->Cy[src][band], 1);
560 cblas_saxpy(/*re+im*/2*Q*Q, 1.0f-pData->covAvgCoeff, (float*)pData->_Cy, 1, (float*)pData->Cy[src][band], 1);
561 }
562
563 /* Formulate mixing matrices */
564 switch(procMode){
565 case SPREADER_MODE_NAIVE: saf_print_error("Shouldn't have gotten this far?"); break;
567 /* For normalising the level of Cy */
568 Ey = Eproto = 0.0f;
569 for(band=0; band<HYBRID_BANDS; band++){
570 for(i=0; i<Q; i++){
571 Ey += crealf(pData->Cy[src][band][i*Q+i]);
572 Eproto += crealf(pData->Cproto[src][band][i*Q+i])+0.000001f;
573 }
574 }
575 Gcomp = sqrtf(Eproto/(Ey+2.23e-9f));
576
577 /* Compute mixing matrix per band */
578 for(band=0; band<HYBRID_BANDS; band++){
579 memcpy(pData->_Cy, pData->Cy[src][band], Q*Q*sizeof(float_complex));
580 cblas_sscal(/*re+im*/2*Q*Q, Gcomp, (float*)pData->_Cy, 1);
581 utility_cseig(NULL, pData->_Cy, Q, 1, pData->_V, pData->_D, NULL);
582 for(i=0; i<Q; i++)
583 for(j=0; j<Q; j++)
584 pData->_D[i*Q+j] = i==j ? csqrtf(pData->_D[i*Q+j]) : cmplxf(0.0f, 0.0f);
585 cblas_cgemm(CblasRowMajor, CblasNoTrans, CblasNoTrans, Q, Q, Q, &calpha,
586 pData->_V, Q,
587 pData->_D, Q, &cbeta,
588 pData->new_M[band], Q);
589 }
590 break;
591
592 case SPREADER_MODE_OM:
593 for(band=0; band<HYBRID_BANDS; band++){
594 if(pData->freqVector[band]<MAX_SPREAD_FREQ){
595#if 1
596 /* Diagonalise and diagonally load the Cproto matrices */
597 cblas_ccopy(Q*Q, pData->Cproto[src][band], 1, pData->_Cproto, 1);
598 for(i=0; i<Q; i++){
599 for(j=0; j<Q; j++){
600 if(i==j)
601 pData->_Cproto[i*Q+i] = craddf(pData->_Cproto[i*Q+i], 0.00001f);
602 CprotoDiag[i*Q+j] = i==j ? crealf(pData->_Cproto[i*Q+i]) : 0.0f;
603 }
604 }
605
606 /* Compute mixing matrices */
607 formulate_M_and_Cr_cmplx(pData->hCdf, pData->_Cproto, pData->Cy[src][band], pData->Qmix_cmplx, 0, 0.2f, pData->new_M[band], pData->Cr_cmplx);
608 for(i=0; i<Q*Q; i++)
609 pData->Cr[i] = crealf(pData->Cr_cmplx[i]);
610 formulate_M_and_Cr(pData->hCdf_res, CprotoDiag, pData->Cr, pData->Qmix, 0, 0.2f, pData->new_Mr[band], NULL);
611#else
612 for(q=0; q<Q; q++)
613 H_tmp[q] = pData->H_grid[band*Q*pData->nGrid + q*pData->nGrid + centre_ind];
614 cblas_cgemm(CblasRowMajor, CblasNoTrans, CblasConjTrans, Q, Q, 1, &calpha,
615 H_tmp, 1,
616 H_tmp, 1, &cbeta,
617 Cx, Q);
618 memset(CxDiag, 0, Q*Q*sizeof(float));
619 for(i=0; i<Q; i++)
620 CxDiag[i*Q+i] = crealf(Cx[i*Q+i]);
621
622 /* Compute mixing matrices */
623 formulate_M_and_Cr_cmplx(pData->hCdf, Cx, pData->Cy[src][band], pData->Qmix_cmplx, 0, 0.2f, pData->new_M[band], pData->Cr_cmplx);
624 for(i=0; i<Q*Q; i++)
625 pData->Cr[i] = crealf(pData->Cr_cmplx[i]);
626 formulate_M_and_Cr(pData->hCdf_res, CxDiag, pData->Cr, pData->Qmix, 0, 0.2f, pData->new_Mr[band], NULL);
627#endif
628 }
629 else{
630 memcpy(pData->new_M[band], pData->Qmix_cmplx, Q*Q*sizeof(float_complex));
631 memset(pData->new_Mr[band], 0, Q*Q*sizeof(float));
632 }
633 }
634 break;
635 }
636
637 /* Apply mixing matrices */
638 for(band=0; band<HYBRID_BANDS; band++){
639 for(t=0; t<TIME_SLOTS; t++){
640 scaleC = cmplxf(pData->interpolatorFadeIn[t], 0.0f);
641 utility_cvsmul(pData->new_M[band], &scaleC, Q*Q, pData->interp_M);
642 cblas_saxpy(/*re+im*/2*Q*Q, pData->interpolatorFadeOut[t], (float*)pData->prev_M[src][band], 1, (float*)pData->interp_M, 1);
643 for(i=0; i<Q; i++) {
644 cblas_cdotu_sub(Q, (float_complex*)(&(pData->interp_M[i*Q])), 1,
645 FLATTEN2D((procMode == SPREADER_MODE_EVD ? pData->decorframeTF[band] : pData->protoframeTF[band])) + t,
646 TIME_SLOTS, &(pData->spreadframeTF[band][i][t]));
647 }
648 }
649
650 /* Also mix in the residual part */
651 if(procMode == SPREADER_MODE_OM){
652 if(pData->freqVector[band]<MAX_SPREAD_FREQ){
653 for(t=0; t<TIME_SLOTS; t++){
654 utility_svsmul(pData->new_Mr[band], &(pData->interpolatorFadeIn[t]), Q*Q, pData->interp_Mr);
655 cblas_saxpy(Q*Q, pData->interpolatorFadeOut[t], pData->prev_Mr[src][band], 1, pData->interp_Mr, 1);
656 cblas_scopy(Q*Q, pData->interp_Mr, 1, (float*)pData->interp_Mr_cmplx, 2);
657 for(i=0; i<Q; i++){
658 cblas_cdotu_sub(Q, (float_complex*)(&(pData->interp_Mr_cmplx[i*Q])), 1, FLATTEN2D(pData->decorframeTF[band]) + t, TIME_SLOTS, &tmp);
659 pData->spreadframeTF[band][i][t] = ccaddf(pData->spreadframeTF[band][i][t], tmp);
660 }
661 }
662 }
663 }
664 }
665 }
666
667 /* Add the spread frame to the output frame, then move onto the next source... */
668 for(band=0; band<HYBRID_BANDS; band++)
669 cblas_saxpy(/*re+im*/2*Q*TIME_SLOTS, 1.0f, (float*)FLATTEN2D(pData->spreadframeTF[band]), 1, (float*)FLATTEN2D(pData->outputframeTF[band]), 1);
670
671 /* For next frame */
672 cblas_ccopy(HYBRID_BANDS*Q*Q, FLATTEN2D(pData->new_M), 1, FLATTEN2D(pData->prev_M[src]), 1);
673 cblas_scopy(HYBRID_BANDS*Q*Q, FLATTEN2D(pData->new_Mr), 1, FLATTEN2D(pData->prev_Mr[src]), 1);
674 }
675
676 /* inverse-TFT */
678
679 /* Copy to output buffer */
680 for (ch = 0; ch < SAF_MIN(Q, nOutputs); ch++)
681 utility_svvcopy(pData->outframeTD[ch], SPREADER_FRAME_SIZE, outputs[ch]);
682 for (; ch < nOutputs; ch++)
683 memset(outputs[ch], 0, SPREADER_FRAME_SIZE*sizeof(float));
684 }
685 else{
686 for (ch=0; ch < nOutputs; ch++)
687 memset(outputs[ch],0, SPREADER_FRAME_SIZE*sizeof(float));
688 }
689
691}
692
693/* Set Functions */
694
699
700void spreader_setSpreadingMode(void* const hSpr, int newMode)
701{
702 spreader_data *pData = (spreader_data*)(hSpr);
703 pData->new_procMode = newMode;
705}
706
707void spreader_setAveragingCoeff(void* const hSpr, float newValue)
708{
709 spreader_data *pData = (spreader_data*)(hSpr);
710 pData->covAvgCoeff = newValue;
711}
712
713void spreader_setSourceAzi_deg(void* const hSpr, int index, float newAzi_deg)
714{
715 spreader_data *pData = (spreader_data*)(hSpr);
716 saf_assert(index<SPREADER_MAX_NUM_SOURCES, "index exceeds the maximum number of sources permitted");
717 if(newAzi_deg>180.0f)
718 newAzi_deg = -360.0f + newAzi_deg;
719 newAzi_deg = SAF_MAX(newAzi_deg, -180.0f);
720 newAzi_deg = SAF_MIN(newAzi_deg, 180.0f);
721 if(pData->src_dirs_deg[index][0]!=newAzi_deg)
722 pData->src_dirs_deg[index][0] = newAzi_deg;
723}
724
725void spreader_setSourceElev_deg(void* const hSpr, int index, float newElev_deg)
726{
727 spreader_data *pData = (spreader_data*)(hSpr);
728 saf_assert(index<SPREADER_MAX_NUM_SOURCES, "index exceeds the maximum number of sources permitted");
729 newElev_deg = SAF_MAX(newElev_deg, -90.0f);
730 newElev_deg = SAF_MIN(newElev_deg, 90.0f);
731 if(pData->src_dirs_deg[index][1] != newElev_deg)
732 pData->src_dirs_deg[index][1] = newElev_deg;
733}
734
735void spreader_setSourceSpread_deg(void* const hSpr, int index, float newSpread_deg)
736{
737 spreader_data *pData = (spreader_data*)(hSpr);
738 saf_assert(index<SPREADER_MAX_NUM_SOURCES, "index exceeds the maximum number of sources permitted");
739 newSpread_deg = SAF_MAX(newSpread_deg, 0.0f);
740 newSpread_deg = SAF_MIN(newSpread_deg, 360.0f);
741 if(pData->src_spread[index] != newSpread_deg)
742 pData->src_spread[index] = newSpread_deg;
743}
744
745void spreader_setNumSources(void* const hSpr, int new_nSources)
746{
747 spreader_data *pData = (spreader_data*)(hSpr);
748 pData->new_nSources = SAF_CLAMP(new_nSources, 1, SPREADER_MAX_NUM_SOURCES);
750}
751
752void spreader_setUseDefaultHRIRsflag(void* const hSpr, int newState)
753{
754 spreader_data *pData = (spreader_data*)(hSpr);
755 if((!pData->useDefaultHRIRsFLAG) && (newState)){
756 pData->useDefaultHRIRsFLAG = newState;
757 spreader_refreshSettings(hSpr); // re-init and re-calc
758 }
759}
760
761void spreader_setSofaFilePath(void* const hSpr, const char* path)
762{
763 spreader_data *pData = (spreader_data*)(hSpr);
764
765 pData->sofa_filepath = realloc1d(pData->sofa_filepath, strlen(path) + 1);
766 strcpy(pData->sofa_filepath, path);
767 pData->useDefaultHRIRsFLAG = 0;
768 spreader_refreshSettings(hSpr); // re-init and re-calc
769}
770
771
772/* Get Functions */
773
775{
776 return SPREADER_FRAME_SIZE;
777}
778
780{
781 spreader_data *pData = (spreader_data*)(hSpr);
782 return pData->codecStatus;
783}
784
785float spreader_getProgressBar0_1(void* const hSpr)
786{
787 spreader_data *pData = (spreader_data*)(hSpr);
788 return pData->progressBar0_1;
789}
790
791void spreader_getProgressBarText(void* const hSpr, char* text)
792{
793 spreader_data *pData = (spreader_data*)(hSpr);
794 memcpy(text, pData->progressBarText, PROGRESSBARTEXT_CHAR_LENGTH*sizeof(char));
795}
796
797int* spreader_getDirectionActivePtr(void* const hSpr, int index)
798{
799 spreader_data *pData = (spreader_data*)(hSpr);
800 return pData->dirActive[index];
801}
802
803int spreader_getSpreadingMode(void* const hSpr)
804{
805 spreader_data *pData = (spreader_data*)(hSpr);
806 return pData->new_procMode;
807}
808
809float spreader_getAveragingCoeff(void* const hSpr)
810{
811 spreader_data *pData = (spreader_data*)(hSpr);
812 return pData->covAvgCoeff;
813}
814
815float spreader_getSourceAzi_deg(void* const hSpr, int index)
816{
817 spreader_data *pData = (spreader_data*)(hSpr);
818 saf_assert(index<SPREADER_MAX_NUM_SOURCES, "index exceeds the maximum number of sources permitted");
819 return pData->src_dirs_deg[index][0];
820}
821
822float spreader_getSourceElev_deg(void* const hSpr, int index)
823{
824 spreader_data *pData = (spreader_data*)(hSpr);
825 saf_assert(index<SPREADER_MAX_NUM_SOURCES, "index exceeds the maximum number of sources permitted");
826 return pData->src_dirs_deg[index][1];
827}
828
829float spreader_getSourceSpread_deg(void* const hSpr, int index)
830{
831 spreader_data *pData = (spreader_data*)(hSpr);
832 saf_assert(index<SPREADER_MAX_NUM_SOURCES, "index exceeds the maximum number of sources permitted");
833 return pData->src_spread[index];
834}
835
836int spreader_getNumSources(void* const hSpr)
837{
838 spreader_data *pData = (spreader_data*)(hSpr);
839 return pData->new_nSources;
840}
841
846
847int spreader_getNumOutputs(void* const hSpr)
848{
849 spreader_data *pData = (spreader_data*)(hSpr);
850 return pData->Q;
851}
852
853int spreader_getNDirs(void* const hSpr)
854{
855 spreader_data *pData = (spreader_data*)(hSpr);
856 return pData->nGrid;
857}
858
859float spreader_getIRAzi_deg(void* const hSpr, int index)
860{
861 spreader_data *pData = (spreader_data*)(hSpr);
862 if(pData->codecStatus == CODEC_STATUS_INITIALISED && pData->grid_dirs_deg!=NULL)
863 return pData->grid_dirs_deg[index*2+0];
864 else
865 return 0.0f;
866}
867
868float spreader_getIRElev_deg(void* const hSpr, int index)
869{
870 spreader_data *pData = (spreader_data*)(hSpr);
871 if(pData->codecStatus == CODEC_STATUS_INITIALISED && pData->grid_dirs_deg!=NULL)
872 return pData->grid_dirs_deg[index*2+1];
873 else
874 return 0.0f;
875}
876
877int spreader_getIRlength(void* const hSpr)
878{
879 spreader_data *pData = (spreader_data*)(hSpr);
880 return pData->h_len;
881}
882
883int spreader_getIRsamplerate(void* const hSpr)
884{
885 spreader_data *pData = (spreader_data*)(hSpr);
886 return (int)pData->h_fs;
887}
888
890{
891 spreader_data *pData = (spreader_data*)(hSpr);
892 return pData->useDefaultHRIRsFLAG;
893}
894
895char* spreader_getSofaFilePath(void* const hSpr)
896{
897 spreader_data *pData = (spreader_data*)(hSpr);
898 if(pData->sofa_filepath!=NULL)
899 return pData->sofa_filepath;
900 else
901 return "no_file";
902}
903
904
905int spreader_getDAWsamplerate(void* const hSpr)
906{
907 spreader_data *pData = (spreader_data*)(hSpr);
908 return pData->fs;
909}
910
912{
913 return 12*HOP_SIZE;
914}
#define MAX_NUM_CHANNELS
Maximum number of input/output channels supported.
Definition _common.h:234
#define MAX_NUM_INPUTS
Maximum number of input channels supported.
Definition _common.h:237
#define MAX_NUM_OUTPUTS
Maximum number of output channels supported.
Definition _common.h:240
#define PROGRESSBARTEXT_CHAR_LENGTH
Length of progress bar string.
Definition _common.h:231
@ PROC_STATUS_ONGOING
Codec is processing input audio, and should not be reinitialised at this time.
Definition _common.h:224
@ PROC_STATUS_NOT_ONGOING
Codec is not processing input audio, and may be reinitialised if needed.
Definition _common.h:226
CODEC_STATUS
Current status of the codec.
Definition _common.h:205
@ CODEC_STATUS_NOT_INITIALISED
Codec has not yet been initialised, or the codec configuration has changed.
Definition _common.h:208
@ CODEC_STATUS_INITIALISED
Codec is initialised and ready to process input audio.
Definition _common.h:206
@ CODEC_STATUS_INITIALISING
Codec is currently being initialised, input audio should not be processed.
Definition _common.h:211
void afSTFT_FIRtoFilterbankCoeffs(float *hIR, int N_dirs, int nCH, int ir_len, int hopSize, int LDmode, int hybridmode, float_complex *hFB)
Converts FIR filters into Filterbank Coefficients by passing them through the afSTFT filterbank.
Definition afSTFTlib.c:593
void afSTFT_create(void **const phSTFT, int nCHin, int nCHout, int hopsize, int lowDelayMode, int hybridmode, AFSTFT_FDDATA_FORMAT format)
Creates an instance of afSTFT.
Definition afSTFTlib.c:143
void afSTFT_backward_knownDimensions(void *const hSTFT, float_complex ***dataFD, int framesize, int dataFD_nCH, int dataFD_nHops, float **dataTD)
Performs backward afSTFT transform (dataFD dimensions are known)
Definition afSTFTlib.c:391
void afSTFT_forward_knownDimensions(void *const hSTFT, float **dataTD, int framesize, int dataFD_nCH, int dataFD_nHops, float_complex ***dataFD)
Performs forward afSTFT transform (dataFD dimensions are known)
Definition afSTFTlib.c:268
void afSTFT_getCentreFreqs(void *const hSTFT, float fs, int nBands, float *freqVector)
Returns current frequency vector.
Definition afSTFTlib.c:546
void afSTFT_destroy(void **const phSTFT)
Destroys an instance of afSTFT.
Definition afSTFTlib.c:199
@ AFSTFT_BANDS_CH_TIME
nBands x nChannels x nTimeHops
Definition afSTFTlib.h:80
#define TIME_SLOTS
Number of STFT timeslots.
#define HOP_SIZE
STFT hop size.
#define HYBRID_BANDS
Number of frequency bands.
void formulate_M_and_Cr_cmplx(void *const hCdf, float_complex *Cx, float_complex *Cy, float_complex *Q, int useEnergyFLAG, float reg, float_complex *M, float_complex *Cr)
Computes the optimal mixing matrices.
void cdf4sap_create(void **const phCdf, int nXcols, int nYcols)
Creates an instance of the Covariance Domain Framework.
Definition saf_cdf4sap.c:85
void cdf4sap_cmplx_destroy(void **const phCdf)
Destroys an instance of the Covariance Domain Framework.
void cdf4sap_destroy(void **const phCdf)
Destroys an instance of the Covariance Domain Framework.
void cdf4sap_cmplx_create(void **const phCdf, int nXcols, int nYcols)
Creates an instance of the Covariance Domain Framework.
void formulate_M_and_Cr(void *const hCdf, float *Cx, float *Cy, float *Q, int useEnergyFLAG, float reg, float *M, float *Cr)
Computes the optimal mixing matrices.
const int __default_N_hrir_dirs
The number of directions/measurements in the default HRIR dataset.
const float __default_hrirs[836][2][256]
The default HRIR data for SAF.
const float __default_hrir_dirs_deg[836][2]
The measurement directions used for the default HRIR dataset.
const int __default_hrir_len
The length of the filters, in samples, for the default HRIR dataset.
const int __default_hrir_fs
The samplerate used to measure the default HRIR filters.
void saf_sofa_close(saf_sofa_container *c)
Frees all SOFA data in a sofa_container.
SAF_SOFA_ERROR_CODES saf_sofa_open(saf_sofa_container *h, char *sofa_filepath, SAF_SOFA_READER_OPTIONS option)
Fills a 'sofa_container' with data found in a SOFA file (GeneralFIR or SimpleFreeFieldHRIR),...
SAF_SOFA_ERROR_CODES
SOFA loader error codes.
@ SAF_SOFA_READER_OPTION_DEFAULT
The default option is SAF_SOFA_READER_OPTION_LIBMYSOFA.
@ SAF_SOFA_OK
None of the error checks failed.
#define saf_print_error(message)
Macro to print a error message along with the filename and line number.
#define saf_assert(x, message)
Macro to make an assertion, along with a string explaining its purpose.
#define SAF_CLAMP(a, min, max)
Ensures value "a" is clamped between the "min" and "max" values.
#define SAF_PI
pi constant (single precision)
#define SAF_MAX(a, b)
Returns the maximum of the two values.
void convert_0_360To_m180_180(float *dirs_deg, int nDirs)
Wraps around any angles exeeding 180 degrees (e.g., 200-> -160)
#define NUM_EARS
2 (true for most humans)
void getVoronoiWeights(float *dirs_deg, int nDirs, int diagFLAG, float *weights)
Computes the integration weights, based on the areas of each face of the corresponding Voronoi diagra...
#define FOURPI
4pi (single precision)
void utility_svvcopy(const float *a, const int len, float *c)
Single-precision, vector-vector copy, i.e.
void utility_cseig(void *const hWork, const float_complex *A, const int dim, int sortDecFLAG, float_complex *V, float_complex *D, float *eig)
Eigenvalue decomposition of a SYMMETRIC/HERMITION matrix: single precision complex,...
void latticeDecorrelator_apply(void *hDecor, float_complex ***inFrame, int nTimeSlots, float_complex ***decorFrame)
Applies the lattice all-pass-filter-based multi-channel signal decorrelator.
#define SAF_MIN(a, b)
Returns the minimum of the two values.
void latticeDecorrelator_create(void **phDecor, float fs, int hopsize, float *freqVector, int nBands, int nCH, int *orders, float *freqCutoffs, int nCutoffs, int maxDelay, int lookupOffset, float enComp_coeff)
Creates an instance of the lattice all-pass-filter-based multi-channel signal decorrelator.
void unitSph2cart(float *dirs, int nDirs, int anglesInDegreesFLAG, float *dirs_xyz)
Converts spherical coordinates to Cartesian coordinates of unit length.
#define saf_print_warning(message)
Macro to print a warning message along with the filename and line number.
void utility_cvsmul(float_complex *a, const float_complex *s, const int len, float_complex *c)
Single-precision, complex, multiplies each element in vector 'a' with a scalar 's',...
void latticeDecorrelator_destroy(void **phDecor)
Destroys an instance of the lattice all-pass-filter-based multi-channel signal decorrelator.
void utility_svsmul(float *a, const float *s, const int len, float *c)
Single-precision, multiplies each element in vector 'a' with a scalar 's', i.e.
void utility_siminv(const float *a, const int len, int *index)
Single-precision, index of minimum absolute value in a vector, i.e.
void ** malloc2d(size_t dim1, size_t dim2, size_t data_size)
2-D malloc (contiguously allocated, so use free() as usual to deallocate)
Definition md_malloc.c:89
void ** realloc2d(void **ptr, size_t dim1, size_t dim2, size_t data_size)
2-D realloc which does NOT retain previous data order
Definition md_malloc.c:115
void * malloc1d(size_t dim1_data_size)
1-D malloc (same as malloc, but with error checking)
Definition md_malloc.c:59
void * realloc1d(void *ptr, size_t dim1_data_size)
1-D realloc (same as realloc, but with error checking)
Definition md_malloc.c:79
void *** malloc3d(size_t dim1, size_t dim2, size_t dim3, size_t data_size)
3-D malloc (contiguously allocated, so use free() as usual to deallocate)
Definition md_malloc.c:151
#define FLATTEN2D(A)
Use this macro when passing a 2-D dynamic multi-dimensional array to memset, memcpy or any other func...
Definition md_malloc.h:65
int spreader_getNumSources(void *const hSpr)
Returns the number of inputs/sources in the current config.
Definition spreader.c:836
int spreader_getProcessingDelay(void)
Returns the processing delay in samples (may be used for delay compensation purposes)
Definition spreader.c:911
int spreader_getIRsamplerate(void *const hSpr)
Returns the IR sample rate.
Definition spreader.c:883
int spreader_getMaxNumSources(void)
Returns the maximum number of input sources supported by spreader.
Definition spreader.c:842
int spreader_getNDirs(void *const hSpr)
Returns the number of directions in the currently used HRIR set.
Definition spreader.c:853
void spreader_process(void *const hSpr, const float *const *inputs, float *const *const outputs, int nInputs, int nOutputs, int nSamples)
Spatialises and spreads the input signals in the user specified directions.
Definition spreader.c:361
float spreader_getIRAzi_deg(void *const hSpr, int index)
Returns the IR/TF azimuth for a given index, in DEGREES.
Definition spreader.c:859
void spreader_setSourceSpread_deg(void *const hSpr, int index, float newSpread_deg)
Sets the source spread for a specific channel index, in DEGREES.
Definition spreader.c:735
int * spreader_getDirectionActivePtr(void *const hSpr, int index)
Returns the pointer to a vector describing which directions are currently being used for the spreadin...
Definition spreader.c:797
void spreader_setSpreadingMode(void *const hSpr, int newMode)
Sets the spreading mode (see SPREADER_PROC_MODES)
Definition spreader.c:700
float spreader_getSourceElev_deg(void *const hSpr, int index)
Returns the source elevation for a given source index, in DEGREES.
Definition spreader.c:822
void spreader_create(void **const phSpr)
Creates an instance of the spreader.
Definition spreader.c:35
int spreader_getFrameSize(void)
Returns the processing framesize (i.e., number of samples processed with every _process() call )
Definition spreader.c:774
int spreader_getIRlength(void *const hSpr)
Returns the length of IRs in time-domain samples.
Definition spreader.c:877
int spreader_getDAWsamplerate(void *const hSpr)
Returns the DAW/Host sample rate.
Definition spreader.c:905
void spreader_initCodec(void *const hSpr)
Intialises the codec variables, based on current global/user parameters.
Definition spreader.c:208
float spreader_getAveragingCoeff(void *const hSpr)
Returns the averaging coefficient [0..1].
Definition spreader.c:809
float spreader_getProgressBar0_1(void *const hSpr)
(Optional) Returns current intialisation/processing progress, between 0..1
Definition spreader.c:785
float spreader_getIRElev_deg(void *const hSpr, int index)
Returns the IR/TF elevation for a given index, in DEGREES.
Definition spreader.c:868
void spreader_setSofaFilePath(void *const hSpr, const char *path)
Sets the file path for a .sofa file, in order to employ a custom HRIR set for the decoding.
Definition spreader.c:761
float spreader_getSourceSpread_deg(void *const hSpr, int index)
Returns the source spread for a given source index, in DEGREES.
Definition spreader.c:829
void spreader_getProgressBarText(void *const hSpr, char *text)
(Optional) Returns current intialisation/processing progress text
Definition spreader.c:791
int spreader_getNumOutputs(void *const hSpr)
Returns the number of ears possessed by the average homo sapien.
Definition spreader.c:847
void spreader_setSourceElev_deg(void *const hSpr, int index, float newElev_deg)
Sets the panning elevation for a specific channel index, in DEGREES.
Definition spreader.c:725
void spreader_setNumSources(void *const hSpr, int new_nSources)
Sets the number of input channels/sources to binauralise.
Definition spreader.c:745
void spreader_refreshSettings(void *const hSpr)
Sets all intialisation flags to 1; re-initialising all settings/variables as spreader is currently co...
Definition spreader.c:695
int spreader_getUseDefaultHRIRsflag(void *const hSpr)
Returns the value of a flag used to dictate whether the default HRIRs in the Spatial_Audio_Framework ...
Definition spreader.c:889
void spreader_setSourceAzi_deg(void *const hSpr, int index, float newAzi_deg)
Sets the panning azimuth for a specific channel index, in DEGREES.
Definition spreader.c:713
float spreader_getSourceAzi_deg(void *const hSpr, int index)
Returns the source azimuth for a given source index, in DEGREES.
Definition spreader.c:815
void spreader_destroy(void **const phSpr)
Destroys an instance of the spreader.
Definition spreader.c:120
int spreader_getSpreadingMode(void *const hSpr)
Returns the spreading mode (see SPREADER_PROC_MODES)
Definition spreader.c:803
void spreader_setAveragingCoeff(void *const hSpr, float newValue)
Sets the averaging coefficient [0..1].
Definition spreader.c:707
char * spreader_getSofaFilePath(void *const hSpr)
Returns the file path for a .sofa file.
Definition spreader.c:895
void spreader_setUseDefaultHRIRsflag(void *const hSpr, int newState)
Sets flag to dictate whether the default HRIRs in the Spatial_Audio_Framework should be used (1),...
Definition spreader.c:752
void spreader_init(void *const hSpr, int sampleRate)
Initialises an instance of spreader with default settings.
Definition spreader.c:195
CODEC_STATUS spreader_getCodecStatus(void *const hSpr)
Returns current codec status codec status (see CODEC_STATUS enum)
Definition spreader.c:779
SPREADER_PROC_MODES
Available processing modes.
Definition spreader.h:55
@ SPREADER_MODE_NAIVE
Simple coherent copies of the input signal(s) areassigned to the spreading areas.
Definition spreader.h:56
@ SPREADER_MODE_EVD
Basic solution based on an Eigenvalue decomposition.
Definition spreader.h:59
@ SPREADER_MODE_OM
Optimal mixing solution.
Definition spreader.h:58
#define SPREADER_MAX_NUM_SOURCES
Maximum number of sources supported by the spreader example.
Definition spreader.h:52
void spreader_setCodecStatus(void *const hSpr, CODEC_STATUS newStatus)
Sets codec status (see CODEC_STATUS enum)
An arbitrary array panner (HRIRs, microphone array IRs, etc.) with coherent and incoherent spreading ...
#define SPREADER_FRAME_SIZE
Framesize, in time-domain samples.
#define MAX_SPREAD_FREQ
Maximum spread frequency, above which no spreading occurs.
SOFA container struct comprising all possible data that can be extracted from SOFA 1....
int DataLengthIR
Length of the IRs, in samples.
int nSources
Number of source/measurement positions.
float * SourcePosition
Source positions (refer to SourcePositionType & SourcePositionUnits for the convention and units); FL...
float DataSamplingRate
Sampling rate used to measure the IRs.
float * DataIR
The impulse response (IR) Data; FLAT:nSources x nReceivers x DataLengthIR.
Main structure for spreader.
_Atomic_FLOAT32 covAvgCoeff
Covariance matrix averaging coefficient, [0..1].
float interpolatorFadeOut[TIME_SLOTS]
Linear Interpolator - Fade out.
_Atomic_INT32 h_len
Length of time-domain filters, in samples.
void * hCdf_res
covariance domain framework handle for the residual
int fs
Host sampling rate, in Hz.
float_complex *** spreadframeTF
time-frequency domain spread frame; HYBRID_BANDS x MAX_NUM_OUTPUTS x TIME_SLOTS
_Atomic_CODEC_STATUS codecStatus
see CODEC_STATUS
float ** inputFrameTD
time-domain input frame; MAX_NUM_INPUTS x SPREADER_FRAME_SIZE
_Atomic_FLOAT32 h_fs
Sample rate used to measure the filters.
void * hSTFT
afSTFT handle
float_complex ** prev_M[SPREADER_MAX_NUM_SOURCES]
previous mixing matrices; HYBRID_BANDS x FLAT:(Q x Q)
float_complex *** inputframeTF
time-frequency domain input frame; HYBRID_BANDS x MAX_NUM_INPUTS x TIME_SLOTS
float_complex * interp_Mr_cmplx
Complex variant of interp_Mr.
float * h_grid
FLAT: nGrid x Q x h_len.
float * grid_dirs_xyz
Grid directions as unit-length Cartesian coordinates; FLAT: nGrid x 3.
float_complex ** Cy[SPREADER_MAX_NUM_SOURCES]
Target covariance matrices; HYBRID_BANDS x FLAT:(Q x Q)
float * Qmix
Identity; FLAT: Q x Q.
float_complex * Qmix_cmplx
Identity; FLAT: Q x Q.
float_complex ** HHH[HYBRID_BANDS]
Pre-computed array outer-products; HYBRID_BANDS x nGrid x FLAT: (Q x Q)
float ** prev_Mr[SPREADER_MAX_NUM_SOURCES]
previous residual mixing matrices; HYBRID_BANDS x FLAT:(Q x Q)
float_complex *** outputframeTF
time-frequency domain output frame; HYBRID_BANDS x MAX_NUM_OUTPUTS x TIME_SLOTS
void * hDecor[SPREADER_MAX_NUM_SOURCES]
handles for decorrelators
_Atomic_INT32 Q
Number of channels in the target playback setup; for example: 2 for binaural.
float_complex ** new_M
mixing matrices; HYBRID_BANDS x FLAT:(Q x Q)
float ** new_Mr
residual mixing matrices; HYBRID_BANDS x FLAT:(Q x Q)
int * dirActive[SPREADER_MAX_NUM_SOURCES]
1: IR direction currently used for spreading, 0: not
_Atomic_FLOAT32 src_spread[SPREADER_MAX_NUM_SOURCES]
Source spreading, in degrees.
float ** outframeTD
time-domain output frame; MAX_NUM_OUTPUTS x SPREADER_FRAME_SIZE
float_complex *** decorframeTF
time-frequency domain decorrelated frame; HYBRID_BANDS x MAX_NUM_OUTPUTS x TIME_SLOTS
char * progressBarText
Current (re)initialisation step, string.
_Atomic_FLOAT32 src_dirs_deg[SPREADER_MAX_NUM_SOURCES][2]
Source directions, in degrees.
char * sofa_filepath
SOFA file path.
float * Cr
Residual covariance; FLAT: Q x Q.
_Atomic_SPREADER_PROC_MODES new_procMode
See SPREADER_PROC_MODES (current value will be replaced by this after next re-init)
float_complex ** Cproto[SPREADER_MAX_NUM_SOURCES]
Current prototype covariance matrices; HYBRID_BANDS x FLAT:(Q x Q)
float freqVector[HYBRID_BANDS]
Frequency vector (filterbank centre frequencies)
float * angles
angles; nGrid x 1
_Atomic_INT32 nSources
Current number of input signals.
float * interp_Mr
Interpolated residual mixing matrix; FLAT:(Q x Q)
_Atomic_INT32 new_nSources
New number of input signals (current value will be replaced by this after next re-init)
_Atomic_INT32 useDefaultHRIRsFLAG
1: use default HRIRs in database, 0: use the measurements from SOFA file (can be anything,...
_Atomic_FLOAT32 progressBar0_1
Current (re)initialisation progress, between [0..1].
float_complex * H_grid
FLAT: HYBRID_BANDS x Q x nGrid.
void * hCdf
covariance domain framework handle
float * grid_dirs_deg
Grid directions, in degrees; FLAT: nGrid x 2.
_Atomic_PROC_STATUS procStatus
see PROC_STATUS
float interpolatorFadeIn[TIME_SLOTS]
Linear Interpolator - Fade in.
float * weights
Integration weights; nGrid x 1.
_Atomic_SPREADER_PROC_MODES procMode
See SPREADER_PROC_MODES.
_Atomic_INT32 nGrid
Number of directions/measurements/HRTFs etc.
float_complex * interp_M
Interpolated mixing matrix; FLAT:(Q x Q)
float_complex *** protoframeTF
time-frequency domain prototype frame; HYBRID_BANDS x MAX_NUM_OUTPUTS x TIME_SLOTS
float_complex * Cr_cmplx
Residual covariance; FLAT: Q x Q.