Line data Source code
1 : /******************************************************************************************************
2 :
3 : (C) 2022-2025 IVAS codec Public Collaboration with portions copyright Dolby International AB, Ericsson AB,
4 : Fraunhofer-Gesellschaft zur Foerderung der angewandten Forschung e.V., Huawei Technologies Co. LTD.,
5 : Koninklijke Philips N.V., Nippon Telegraph and Telephone Corporation, Nokia Technologies Oy, Orange,
6 : Panasonic Holdings Corporation, Qualcomm Technologies, Inc., VoiceAge Corporation, and other
7 : contributors to this repository. All Rights Reserved.
8 :
9 : This software is protected by copyright law and by international treaties.
10 : The IVAS codec Public Collaboration consisting of Dolby International AB, Ericsson AB,
11 : Fraunhofer-Gesellschaft zur Foerderung der angewandten Forschung e.V., Huawei Technologies Co. LTD.,
12 : Koninklijke Philips N.V., Nippon Telegraph and Telephone Corporation, Nokia Technologies Oy, Orange,
13 : Panasonic Holdings Corporation, Qualcomm Technologies, Inc., VoiceAge Corporation, and other
14 : contributors to this repository retain full ownership rights in their respective contributions in
15 : the software. This notice grants no license of any kind, including but not limited to patent
16 : license, nor is any license granted by implication, estoppel or otherwise.
17 :
18 : Contributors are required to enter into the IVAS codec Public Collaboration agreement before making
19 : contributions.
20 :
21 : This software is provided "AS IS", without any express or implied warranties. The software is in the
22 : development stage. It is intended exclusively for experts who have experience with such software and
23 : solely for the purpose of inspection. All implied warranties of non-infringement, merchantability
24 : and fitness for a particular purpose are hereby disclaimed and excluded.
25 :
26 : Any dispute, controversy or claim arising under or in relation to providing this software shall be
27 : submitted to and settled by the final, binding jurisdiction of the courts of Munich, Germany in
28 : accordance with the laws of the Federal Republic of Germany excluding its conflict of law rules and
29 : the United Nations Convention on Contracts on the International Sales of Goods.
30 :
31 : *******************************************************************************************************/
32 :
33 : #include <stdint.h>
34 : #include <assert.h>
35 : #include "options.h"
36 : #include "prot.h"
37 : #include <math.h>
38 : #include "ivas_prot_rend.h"
39 : #include "ivas_rom_rend.h"
40 : #include "ivas_cnst.h"
41 : #include "wmc_auto.h"
42 :
43 :
44 : /*---------------------------------------------------------------------*
45 : * Local function prototypes
46 : *---------------------------------------------------------------------*/
47 :
48 : static void getPeriodicBSplineSampVec( float *BfVec, int16_t *AzIdx, const int16_t NumBFs, const float t, int16_t *num_az_idx, const float knot_interval, const float azimKSeq_0, const int16_t azimSegSamples, const float *azimBsShape, const int16_t subSampFactor );
49 : static void getStandardBSplineSampVec( float *BfVec, int16_t *NzIdx, int16_t *num_idx, const int16_t NumBFs, const float t, const float *KSeq, const int16_t SegSamples, const int16_t *BsLen, const int16_t *BsStart, const float *BsShape );
50 : static void GenerateFilter( const float elev, float azim, ModelParams_t *model, ModelEval_t *modelEval );
51 : static void GenerateITD( const float elev, float azim, ModelParamsITD_t *model, ModelEval_t *modelEval );
52 : static void SkipSmallest_ValueIndex( int16_t *use_inds, const ValueIndex_t *VI, const int16_t N, const int16_t n_smallest );
53 :
54 :
55 : /*-------------------------------------------------------------------*
56 : * TDREND_REND_RenderSourceHRFilt()
57 : *
58 : * Renders each object using the HR filters
59 : --------------------------------------------------------------------*/
60 :
61 3461327 : ivas_error TDREND_REND_RenderSourceHRFilt(
62 : TDREND_SRC_t *Src_p, /* i/o: The source to be rendered */
63 : const float *hrf_left_delta, /* i : Left filter interpolation delta */
64 : const float *hrf_right_delta, /* i : Right filter interpolation delta */
65 : const int16_t intp_count, /* i : Interpolation count */
66 : float output_buf[][L_SPATIAL_SUBFR_48k], /* o : Output buffer */
67 : const int16_t subframe_length /* i : Subframe length in use */
68 : )
69 : {
70 : float LeftOutputFrame[L_SPATIAL_SUBFR_48k];
71 : float RightOutputFrame[L_SPATIAL_SUBFR_48k];
72 :
73 3461327 : TDREND_Apply_ITD( Src_p->InputFrame_p, LeftOutputFrame, RightOutputFrame, &Src_p->previtd, Src_p->itd, Src_p->mem_itd, subframe_length );
74 3461327 : TDREND_firfilt( LeftOutputFrame, Src_p->hrf_left_prev, hrf_left_delta, intp_count, Src_p->mem_hrf_left, subframe_length, Src_p->filterlength, Src_p->Gain, Src_p->prevGain );
75 3461327 : TDREND_firfilt( RightOutputFrame, Src_p->hrf_right_prev, hrf_right_delta, intp_count, Src_p->mem_hrf_right, subframe_length, Src_p->filterlength, Src_p->Gain, Src_p->prevGain );
76 3461327 : Src_p->prevGain = Src_p->Gain;
77 :
78 : /* Copy to accumulative output frame */
79 3461327 : v_add( LeftOutputFrame, output_buf[0], output_buf[0], subframe_length );
80 3461327 : v_add( RightOutputFrame, output_buf[1], output_buf[1], subframe_length );
81 :
82 :
83 3461327 : return IVAS_ERR_OK;
84 : }
85 :
86 :
87 : /*-------------------------------------------------------------------*
88 : * GetFilterFromAngle()
89 : *
90 : * Obtain an HR filter corresponding to the input angle pair.
91 : * This version uses the HR filter model.
92 : --------------------------------------------------------------------*/
93 :
94 2136738 : void GetFilterFromAngle(
95 : TDREND_HRFILT_FiltSet_t *HrFiltSet_p, /* i/o: HR filter set structure */
96 : const float Elev, /* i : Elevation, degrees */
97 : float Azim, /* i : Azimuth, degrees */
98 : const int16_t filterlength, /* i : Filter length */
99 : float *hrf_left, /* o : Left HR filter */
100 : float *hrf_right, /* o : Right HR filter */
101 : int16_t *itd /* o : ITD value */
102 : )
103 : {
104 :
105 2136738 : GenerateFilter( Elev, Azim, &HrFiltSet_p->ModelParams, &HrFiltSet_p->ModelEval );
106 :
107 2136738 : mvr2r( HrFiltSet_p->ModelEval.hrfModL, hrf_left, filterlength );
108 2136738 : mvr2r( HrFiltSet_p->ModelEval.hrfModR, hrf_right, filterlength );
109 :
110 : /* 4. Evaluate the ITD */
111 2136738 : if ( HrFiltSet_p->ModelParams.UseItdModel )
112 : {
113 2136738 : GenerateITD( Elev, Azim, &HrFiltSet_p->ModelParamsITD, &HrFiltSet_p->ModelEval );
114 2136738 : *itd = (int16_t) HrFiltSet_p->ModelEval.itdMod;
115 : }
116 : else
117 : {
118 0 : *itd = 0;
119 : }
120 :
121 2136738 : return;
122 : }
123 :
124 :
125 : /*-------------------------------------------------------------------*
126 : * GenerateFilter()
127 : *
128 : * Generate an HR filter using the B Spline model.
129 : --------------------------------------------------------------------*/
130 :
131 2136738 : static void GenerateFilter(
132 : const float elev, /* i : Elevation angle, degrees */
133 : float azim, /* i : Azimuth angle, degrees */
134 : ModelParams_t *model, /* i : Model parameters structure */
135 : ModelEval_t *modelEval /* i/o: Model evaluation structure */
136 : )
137 : {
138 : int16_t qp, p, k, i;
139 : int32_t index;
140 : int16_t AzIdx[HRTF_MODEL_BSPLINE_NUM_COEFFS][HRTF_MODEL_BSPLINE_NUM_COEFFS], EvIdx[HRTF_MODEL_BSPLINE_NUM_COEFFS]; /* non-zero basis functions */
141 : int16_t num_az_idx[HRTF_MODEL_BSPLINE_NUM_COEFFS];
142 : int16_t num_ev_idx;
143 : int16_t BM_idx[HRTF_MODEL_BSPLINE_NUM_COEFFS_SQ];
144 : float knot_interval;
145 : float ETotL, ETotR, ESynL, ESynR;
146 : float ScaleL, ScaleR;
147 : int16_t iSec;
148 :
149 2136738 : getStandardBSplineSampVec( modelEval->elevBfVec, EvIdx, &num_ev_idx, model->elevDim3, elev, model->elevKSeq,
150 2136738 : model->elevSegSamples, model->elevBsLen, model->elevBsStart, model->elevBsShape );
151 :
152 9625344 : for ( p = 0; p < num_ev_idx; p++ )
153 : {
154 : /* Wrap the requested azimuth to the range of the BSplines */
155 7488606 : azim = fmodf( azim + model->azimKSeq[p][0], 360 );
156 7488606 : if ( azim < 0 )
157 : {
158 1044807 : azim += 360.0f;
159 : }
160 7488606 : azim = azim - model->azimKSeq[p][0];
161 :
162 7488606 : if ( model->azimDim3[EvIdx[p]] == 1 ) /* Constant basis function */
163 : {
164 28648 : num_az_idx[p] = 1;
165 28648 : AzIdx[p][0] = 0;
166 28648 : modelEval->azimBfVec[p][0] = 1.0f;
167 : }
168 : else
169 : {
170 7459958 : k = EvIdx[p];
171 7459958 : knot_interval = ( model->azimKSeq[k][model->azimDim3[k]] - model->azimKSeq[k][0] ) / model->azimDim3[k];
172 :
173 7459958 : getPeriodicBSplineSampVec( modelEval->azimBfVec[p], AzIdx[p], model->azimDim3[k], azim, &num_az_idx[p],
174 7459958 : knot_interval, model->azimKSeq[k][0], model->azimSegSamples[model->azimShapeIdx[k]],
175 7459958 : model->azimBsShape[model->azimShapeIdx[k]], model->azimShapeSampFactor[k] );
176 : }
177 : }
178 :
179 : /* Compute BM */
180 2136738 : qp = 0;
181 9625344 : for ( p = 0; p < num_ev_idx; p++ )
182 : {
183 36279990 : for ( i = 0; i < num_az_idx[p]; i++ )
184 : {
185 28791384 : modelEval->BM[qp + i] = modelEval->elevBfVec[p] * modelEval->azimBfVec[p][i];
186 28791384 : BM_idx[qp + i] = model->azim_start_idx[EvIdx[p]] + AzIdx[p][i];
187 : }
188 7488606 : qp = qp + num_az_idx[p];
189 : }
190 :
191 : /* Compute HR filters, approximate optimized model evaluation */
192 8546952 : for ( iSec = 0; iSec < HRTF_MODEL_N_SECTIONS; iSec++ )
193 : {
194 6410214 : ETotL = 0.0f;
195 6410214 : ETotR = 0.0f;
196 6410214 : ESynL = 0.0f;
197 6410214 : ESynR = 0.0f;
198 :
199 : /* Energy is precalculated part updated with square of BM value. Store index for sorting */
200 92784366 : for ( i = 0; i < qp; i++ )
201 : {
202 86374152 : modelEval->BMEnergiesL[i].val = modelEval->BM[i] * modelEval->BM[i] * model->EL[iSec * model->AlphaN + BM_idx[i]];
203 86374152 : modelEval->BMEnergiesR[i].val = modelEval->BM[i] * modelEval->BM[i] * model->ER[iSec * model->AlphaN + BM_idx[i]];
204 86374152 : modelEval->BMEnergiesL[i].i = i;
205 86374152 : modelEval->BMEnergiesR[i].i = i;
206 86374152 : ETotL += modelEval->BMEnergiesL[i].val;
207 86374152 : ETotR += modelEval->BMEnergiesR[i].val;
208 : }
209 :
210 : /* Number of basis components actually used. */
211 6410214 : p = min( HRTF_MODEL_N_CPTS_VAR[iSec], qp );
212 6410214 : SkipSmallest_ValueIndex( modelEval->UseIndsL, modelEval->BMEnergiesL, qp, qp - p );
213 6410214 : SkipSmallest_ValueIndex( modelEval->UseIndsR, modelEval->BMEnergiesR, qp, qp - p );
214 :
215 : /* Account for lost energy */
216 80950320 : for ( i = 0; i < p; i++ )
217 : {
218 74540106 : ESynL += modelEval->BMEnergiesL[modelEval->UseIndsL[i]].val;
219 74540106 : ESynR += modelEval->BMEnergiesR[modelEval->UseIndsR[i]].val;
220 : }
221 6410214 : ScaleL = sqrtf( ETotL / ( ESynL + EPSILON ) );
222 6410214 : ScaleR = sqrtf( ETotR / ( ESynR + EPSILON ) );
223 :
224 : /* Build using only the most energetic components. */
225 273197394 : for ( k = model->iSecFirst[iSec]; k <= model->iSecLast[iSec]; k++ )
226 : {
227 266787180 : modelEval->hrfModL[k] = 0.0;
228 266787180 : modelEval->hrfModR[k] = 0.0;
229 3365487156 : for ( i = 0; i < p; i++ )
230 : {
231 3098699976 : index = BM_idx[modelEval->BMEnergiesL[modelEval->UseIndsL[i]].i] + k * model->AlphaN;
232 3098699976 : modelEval->hrfModL[k] += modelEval->BM[modelEval->BMEnergiesL[modelEval->UseIndsL[i]].i] * model->AlphaL[index];
233 3098699976 : index = BM_idx[modelEval->BMEnergiesR[modelEval->UseIndsR[i]].i] + k * model->AlphaN;
234 3098699976 : modelEval->hrfModR[k] += modelEval->BM[modelEval->BMEnergiesR[modelEval->UseIndsR[i]].i] * model->AlphaR[index];
235 : }
236 : /* Account for lost energy */
237 266787180 : modelEval->hrfModL[k] *= ScaleL;
238 266787180 : modelEval->hrfModR[k] *= ScaleR;
239 : }
240 : }
241 :
242 2136738 : return;
243 : }
244 :
245 :
246 : /*-------------------------------------------------------------------*
247 : * GenerateITD()
248 : *
249 : * Generates an ITD value from the B Spline model.
250 : --------------------------------------------------------------------*/
251 :
252 2136738 : static void GenerateITD(
253 : const float elev, /* i : Elevation angle, degrees */
254 : float azim, /* i : Azimuth angle, degrees */
255 : ModelParamsITD_t *model, /* i : ITD Model parameters structure */
256 : ModelEval_t *modelEval /* i/o: Model evaluation structure */
257 : )
258 : {
259 : int16_t qp, p, i;
260 : int32_t index;
261 : int16_t AlphaN;
262 : int16_t elev_offset;
263 : float azim_itd;
264 : int16_t AzIdx[HRTF_MODEL_BSPLINE_NUM_COEFFS], EvIdx[HRTF_MODEL_BSPLINE_NUM_COEFFS]; /* non-zero basis functions */
265 : int16_t num_az_idx, num_ev_idx;
266 : int16_t BM_idx[HRTF_MODEL_BSPLINE_NUM_COEFFS_SQ];
267 :
268 : /* Wrap the requested azimuth to the range of the BSplines */
269 2136738 : azim = fmodf( azim + model->azimKSeq[0], 360 );
270 2136738 : if ( azim < 0 )
271 : {
272 1044807 : azim += 360.0f;
273 : }
274 2136738 : azim = azim - model->azimKSeq[0];
275 :
276 2136738 : if ( (int16_t) abs( (int16_t) elev ) != 90 )
277 : {
278 2136738 : getStandardBSplineSampVec( modelEval->elevBfVecITD, EvIdx, &num_ev_idx, model->elevDim3, elev, model->elevKSeq,
279 2136738 : model->elevSegSamples, model->elevBsLen, model->elevBsStart, model->elevBsShape );
280 2136738 : azim_itd = azim;
281 2136738 : if ( azim > 180 )
282 : {
283 : /* Flip spline functions around 180 deg */
284 1044074 : azim_itd = 360 - azim;
285 : }
286 2136738 : getStandardBSplineSampVec( modelEval->azimBfVecITD, AzIdx, &num_az_idx, ( model->azimDim3 + 1 ) / 2, azim_itd, model->azimKSeq,
287 2136738 : model->azimSegSamples, model->azimBsLen, model->azimBsStart, model->azimBsShape );
288 2136738 : if ( azim > 180 )
289 : {
290 : /* Flip spline functions around 180 deg */
291 5220370 : for ( i = 0; i < HRTF_MODEL_BSPLINE_NUM_COEFFS; i++ )
292 : {
293 4176296 : AzIdx[i] = model->azimDim3 - 1 - AzIdx[i];
294 : }
295 : }
296 : }
297 : else
298 : {
299 0 : num_az_idx = 1;
300 0 : num_ev_idx = 1;
301 0 : modelEval->elevBfVecITD[0] = 1.0f;
302 0 : if ( elev < 0 )
303 : {
304 0 : EvIdx[0] = 0;
305 : }
306 : else
307 : {
308 0 : EvIdx[0] = model->elevDim3 - 1;
309 : }
310 : }
311 :
312 : /* Compute BM_ITD */
313 2136738 : elev_offset = 0;
314 2136738 : if ( model->elevKSeq[0] == -90 )
315 : {
316 2136738 : elev_offset = -( model->azimDim3 - 1 );
317 : }
318 2136738 : qp = 0;
319 9957061 : for ( p = 0; p < num_ev_idx; p++ )
320 : {
321 7820323 : if ( EvIdx[p] == 0 && (int16_t) ( model->elevKSeq[EvIdx[p]] ) == -90 )
322 : {
323 10332 : modelEval->BM_ITD[qp] = modelEval->elevBfVecITD[p];
324 10332 : BM_idx[qp] = EvIdx[p] * model->azimDim3;
325 10332 : qp = qp + 1;
326 : }
327 7809991 : else if ( EvIdx[p] == ( model->elevDim3 - 1 ) && (int16_t) ( model->elevKSeq[EvIdx[p] - 2] ) == 90 )
328 : {
329 : /* NB: -2 in if() condition above as number of knot points is numBF-2 */
330 12034 : modelEval->BM_ITD[qp] = modelEval->elevBfVecITD[p];
331 12034 : BM_idx[qp] = EvIdx[p] * model->azimDim3 + elev_offset;
332 12034 : qp = qp + 1;
333 : }
334 : else
335 : {
336 37873860 : for ( i = 0; i < num_az_idx; i++ )
337 : {
338 30075903 : modelEval->BM_ITD[qp + i] = modelEval->elevBfVecITD[p] * modelEval->azimBfVecITD[i];
339 30075903 : BM_idx[qp + i] = EvIdx[p] * model->azimDim3 + elev_offset + AzIdx[i];
340 : }
341 7797957 : qp = qp + num_az_idx;
342 : }
343 : }
344 :
345 : /* Compute ITD */
346 2136738 : AlphaN = model->elevDim3 * model->azimDim3 + elev_offset;
347 2136738 : if ( model->elevKSeq[model->elevDim3 - 3] == 90 ) /* Constant azimuth basis function */
348 : {
349 2136738 : AlphaN = AlphaN - ( model->azimDim3 - 1 );
350 : }
351 :
352 : /* Matrix multiplcation (row x column) */
353 2136738 : modelEval->itdMod = 0.0f;
354 32235007 : for ( i = 0; i < qp; i++ )
355 : {
356 30098269 : index = BM_idx[i];
357 30098269 : modelEval->itdMod += modelEval->BM_ITD[i] * model->W[index];
358 : }
359 :
360 2136738 : modelEval->itdMod = -(float) round_f( modelEval->itdMod * model->resamp_factor ); /* Resample and reverse sign of ITD as convention in model and renderer are opposite */
361 :
362 2136738 : return;
363 : }
364 :
365 :
366 : /*-------------------------------------------------------------------*
367 : * getStandardBSplineSampVec()
368 : *
369 : * Obtain a periodic sampled B Spline basis vector.
370 : --------------------------------------------------------------------*/
371 :
372 7459958 : static void getPeriodicBSplineSampVec(
373 : float *BfVec, /* i/o: values for non-zero basis functions */
374 : int16_t *AzIdx, /* i/o: indices of non-zero basis functions */
375 : const int16_t NumBFs, /* i : the number of basis functions = third index of Bf. */
376 : const float t, /* i : azimuth */
377 : int16_t *num_az_idx, /* i : Number of azimuth indices */
378 : const float knot_interval, /* i : The knot interval */
379 : const float azimKSeq_0, /* i : Knot sequence */
380 : const int16_t azimSegSamples, /* i : Samples per segment */
381 : const float *azimBsShape, /* i : Basis shape */
382 : const int16_t subSampFactor /* i : Subsampling factor */
383 : )
384 : {
385 : int16_t i, nI, d0, d;
386 :
387 7459958 : int16_t SegSamples = azimSegSamples / subSampFactor;
388 :
389 : /* index of closest sample point */
390 7459958 : d0 = (int16_t) round_f( ( t - azimKSeq_0 ) / ( knot_interval / SegSamples ) );
391 :
392 : /* find segment */
393 7459958 : nI = d0 / SegSamples;
394 :
395 7459958 : *num_az_idx = HRTF_MODEL_BSPLINE_NUM_COEFFS;
396 :
397 7459958 : if ( d0 % SegSamples == 0 )
398 : {
399 1077096 : ( *num_az_idx )--; /* on the knot points, the last basis function is zero */
400 : }
401 :
402 36222694 : for ( i = 0; i < *num_az_idx; i++ )
403 : {
404 28762736 : d = d0 - ( i + nI - 1 ) * SegSamples; /* offset of knot_interval */
405 28762736 : BfVec[i] = azimBsShape[(int16_t) abs( d ) * subSampFactor];
406 28762736 : AzIdx[i] = ( nI + i ) % NumBFs;
407 : }
408 :
409 7459958 : return;
410 : }
411 :
412 :
413 : /*-------------------------------------------------------------------*
414 : * getStandardBSplineSampVec()
415 : *
416 : * Obtain a sampled B Spline basis vector.
417 : --------------------------------------------------------------------*/
418 :
419 6410214 : static void getStandardBSplineSampVec(
420 : float *BfVec, /* i/o: values for non-zero basis functions */
421 : int16_t *NzIdx, /* i/o: indices of non-zero basis functions */
422 : int16_t *num_idx, /* i/o: number of non-zero indices */
423 : const int16_t NumBFs, /* i : the number of basis functions = third index of Bf */
424 : const float t, /* i : estimation point */
425 : const float *KSeq, /* i : knot sequence including multiplicities */
426 : const int16_t SegSamples, /* i : samples per segment */
427 : const int16_t *BsLen, /* i : lengths of basis shapes */
428 : const int16_t *BsStart, /* i : start of basis shapes */
429 : const float *BsShape /* i : basis shapes */
430 : )
431 : {
432 : int16_t i, nI;
433 : float knot_interval;
434 : int16_t d0, d, shape_idx, start_idx;
435 :
436 6410214 : knot_interval = ( KSeq[NumBFs - 3] - KSeq[0] ) / ( NumBFs - 3 ); /* assuming triple knot at the first knot */
437 : /* index of closest sample point */
438 6410214 : d0 = (int16_t) round_f( ( t - KSeq[0] ) / ( knot_interval / SegSamples ) );
439 : /* find segment */
440 6410214 : nI = d0 / SegSamples;
441 :
442 6410214 : *num_idx = HRTF_MODEL_BSPLINE_NUM_COEFFS;
443 :
444 6410214 : if ( d0 % SegSamples == 0 )
445 : {
446 2090466 : ( *num_idx )--; /* on the knot points, the last basis function is zero */
447 : }
448 :
449 :
450 29960604 : for ( i = 0; i < *num_idx; i++ )
451 : {
452 23550390 : start_idx = max( 0, i + nI - 3 );
453 23550390 : shape_idx = min( i + nI, min( 3, NumBFs - 1 - ( i + nI ) ) );
454 23550390 : d = d0 - start_idx * SegSamples;
455 :
456 23550390 : if ( i + nI > NumBFs - 4 ) /* reverse full shape */
457 : {
458 1171017 : d = BsLen[shape_idx] - 1 - d;
459 : }
460 22379373 : else if ( d > BsLen[shape_idx] - 1 ) /* reverse half shape */
461 : {
462 9632079 : d = 2 * ( BsLen[shape_idx] - 1 ) - d;
463 : }
464 23550390 : BfVec[i] = BsShape[BsStart[shape_idx] + (int16_t) abs( d )]; /*TT, verify if abs is needed */
465 23550390 : NzIdx[i] = nI + i;
466 : }
467 :
468 6410214 : return;
469 : }
470 :
471 :
472 : /*-------------------------------------------------------------------*
473 : * HRTF_model_precalc()
474 : *
475 : * Set precalculated parameters for HRTF model
476 : --------------------------------------------------------------------*/
477 :
478 1403 : void HRTF_model_precalc(
479 : ModelParams_t *model /* i/o: HRTF model parameters */
480 : )
481 : {
482 : int16_t sec_length;
483 : int16_t i;
484 :
485 1403 : sec_length = model->K / HRTF_MODEL_N_SECTIONS;
486 :
487 5612 : for ( i = 0; i < HRTF_MODEL_N_SECTIONS; i++ )
488 : {
489 4209 : model->iSecFirst[i] = i * sec_length;
490 : }
491 4209 : for ( i = 0; i < HRTF_MODEL_N_SECTIONS - 1; i++ )
492 : {
493 2806 : model->iSecLast[i] = ( i + 1 ) * sec_length - 1;
494 : }
495 1403 : model->iSecLast[HRTF_MODEL_N_SECTIONS - 1] = model->K - 1; /* Final section is longer if (K % nSec) > 0 */
496 :
497 1403 : maximum_s( model->azimDim3, model->elevDim3, &model->azimDim3Max );
498 :
499 1403 : return;
500 : }
501 :
502 :
503 : /*-------------------------------------------------------------------*
504 : * BSplineModelEvalDealloc()
505 : *
506 : * Deallocate BSpline HR Filter model
507 : --------------------------------------------------------------------*/
508 :
509 1232 : void BSplineModelEvalDealloc(
510 : ModelParams_t *model, /* i : Model parameters */
511 : ModelEval_t *modelEval /* i : Model evaluation structure */
512 : )
513 : {
514 : /* Allocated in LoadBSplineBinary() */
515 : int16_t i;
516 :
517 1232 : if ( model->modelROM )
518 : {
519 :
520 1232 : free( (void *) model->azimBsShape ); /* void* cast needed to please both gcc and Visual studio compilers. Deallocating const float** should be fine and gcc agrees, but Visual studio complains. */
521 :
522 19712 : for ( i = 0; i < model->elevDim3; i++ )
523 : {
524 18480 : free( model->azimKSeq[i] );
525 : }
526 1232 : free( model->azimKSeq );
527 :
528 1232 : if ( modelEval != NULL )
529 : {
530 1232 : free( modelEval->hrfModL );
531 1232 : free( modelEval->hrfModR );
532 : }
533 : }
534 :
535 1232 : return;
536 : }
537 :
538 :
539 : /*-------------------------------------------------------------------*
540 : * SkipSmallest_ValueIndex()
541 : *
542 : * Returns indices to the largest values in a ValueIndex list,
543 : * unordered (i.e. skip the n_smallest values, return the remainder).
544 : --------------------------------------------------------------------*/
545 :
546 12820428 : static void SkipSmallest_ValueIndex(
547 : int16_t *use_inds, /* i/o: List of indices to use */
548 : const ValueIndex_t *VI, /* i : List of value-index items */
549 : const int16_t N, /* i : Length of list */
550 : const int16_t n_smallest /* i : Number of items to skip */
551 : )
552 : {
553 : int16_t i, j, k, flag;
554 : int16_t skip_inds[HRTF_MODEL_BSPLINE_NUM_COEFFS_SQ];
555 : float candidate_max; /* Stores the max value in the shortlist (next to be kicked off) and its index */
556 : int16_t candidate_max_i;
557 :
558 : /* Initialise with first entries of VI */
559 12820428 : candidate_max = 0.0;
560 12820428 : candidate_max_i = 0;
561 36488520 : for ( j = 0; j < n_smallest; j++ )
562 : {
563 23668092 : skip_inds[j] = j;
564 23668092 : if ( candidate_max < VI[j].val )
565 : {
566 16433884 : candidate_max = VI[j].val;
567 16433884 : candidate_max_i = j;
568 : }
569 : }
570 :
571 : /* Look in the remainder of the list for smaller values */
572 161900640 : for ( i = n_smallest; i < N; i++ )
573 : {
574 353875186 : for ( j = 0; j < n_smallest; j++ )
575 : {
576 229285142 : if ( VI[i].val < VI[skip_inds[j]].val )
577 : {
578 : /* Found a smaller value, so it goes into the list, replacing candidate_max. */
579 24490168 : skip_inds[candidate_max_i] = i;
580 24490168 : candidate_max = VI[i].val;
581 : /* Update candidate_max */
582 114371689 : for ( k = 0; k < n_smallest; k++ )
583 : {
584 89881521 : if ( VI[skip_inds[k]].val > candidate_max )
585 : {
586 24148164 : candidate_max = VI[skip_inds[k]].val;
587 24148164 : candidate_max_i = k;
588 : }
589 : }
590 24490168 : break;
591 : }
592 : }
593 : }
594 :
595 : /* Build the list of indices that will not be skipped */
596 12820428 : k = 0;
597 185568732 : for ( j = 0; j < N; j++ )
598 : {
599 172748304 : flag = 1;
600 484848696 : for ( i = 0; i < n_smallest; i++ )
601 : {
602 335768484 : if ( skip_inds[i] == j )
603 : {
604 23668092 : flag = 0;
605 23668092 : break;
606 : }
607 : }
608 172748304 : if ( flag )
609 : {
610 149080212 : use_inds[k] = j;
611 149080212 : k++;
612 : }
613 : }
614 :
615 12820428 : return;
616 : }
|
