1 /*
2  * loudness.c
3  *
4  *  Created on: 17.01.2017
5  *      Author: hoene
6  */
7 
8 #include "../resampler/speex_resampler.h"
9 #include "mysofa.h"
10 #include "mysofa_export.h"
11 #include "tools.h"
12 #include <assert.h>
13 #include <float.h>
14 #include <math.h>
15 #include <stdio.h>
16 #include <stdlib.h>
17 #include <string.h>
18 
mysofa_loudness(struct MYSOFA_HRTF * hrtf)19 MYSOFA_EXPORT float mysofa_loudness(struct MYSOFA_HRTF *hrtf) {
20   float c[3], factor;
21   float min = FLT_MAX;
22   int radius = 0;
23   unsigned int i, index = 0;
24   int cartesian =
25       verifyAttribute(hrtf->SourcePosition.attributes, "Type", "cartesian");
26 
27   /*
28    * find frontal source position
29    */
30   for (i = 0; i + 2 < hrtf->SourcePosition.elements; i += hrtf->C) {
31     c[0] = hrtf->SourcePosition.values[i];
32     c[1] = hrtf->SourcePosition.values[i + 1];
33     c[2] = hrtf->SourcePosition.values[i + 2];
34 
35     if (cartesian)
36       mysofa_c2s(c);
37 
38     if (min > c[0] + c[1]) {
39       min = c[0] + c[1];
40       radius = c[2];
41       index = i;
42     } else if (min == c[0] + c[1] && radius < c[2]) {
43       radius = c[2];
44       index = i;
45     }
46   }
47 
48   /* get loudness of frontal fir filter, for both channels */
49   factor = loudness(hrtf->DataIR.values + (index / hrtf->C) * hrtf->N * hrtf->R,
50                     hrtf->N * hrtf->R);
51   factor = sqrtf(2 / factor);
52   if (fequals(factor, 1.f))
53     return 1.f;
54 
55   scaleArray(hrtf->DataIR.values, hrtf->DataIR.elements, factor);
56 
57   return factor;
58 }
59