aboutsummaryrefslogtreecommitdiffstats
path: root/utils/makehrtf.c
diff options
context:
space:
mode:
Diffstat (limited to 'utils/makehrtf.c')
-rw-r--r--utils/makehrtf.c113
1 files changed, 98 insertions, 15 deletions
diff --git a/utils/makehrtf.c b/utils/makehrtf.c
index ea56ffcd..f9da9429 100644
--- a/utils/makehrtf.c
+++ b/utils/makehrtf.c
@@ -2,7 +2,7 @@
* HRTF utility for producing and demonstrating the process of creating an
* OpenAL Soft compatible HRIR data set.
*
- * Copyright (C) 2011-2012 Christopher Fitzgerald
+ * Copyright (C) 2011-2014 Christopher Fitzgerald
*
* This program is free software; you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
@@ -58,7 +58,7 @@
* 1999
*/
-/* Needed for 64-bit unsigned integer. */
+// Needed for 64-bit unsigned integer.
#include "config.h"
#include <stdio.h>
@@ -153,6 +153,10 @@
#define MIN_TRUNCSIZE (8)
#define MAX_TRUNCSIZE (128)
+// The limits to the custom head radius on the command line.
+#define MIN_CUSTOM_RADIUS (0.05)
+#define MAX_CUSTOM_RADIUS (0.15)
+
// The truncation window size must be a multiple of the below value to allow
// for vectorized convolution.
#define MOD_TRUNCSIZE (8)
@@ -162,6 +166,8 @@
#define DEFAULT_SURFACE (1)
#define DEFAULT_LIMIT (24.0)
#define DEFAULT_TRUNCSIZE (32)
+#define DEFAULT_HEAD_MODEL (HM_DATASET)
+#define DEFAULT_CUSTOM_RADIUS (0.0)
// The four-character-codes for RIFF/RIFX WAVE file chunks.
#define FOURCC_RIFF (0x46464952) // 'RIFF'
@@ -208,6 +214,13 @@ enum ElementTypeT {
ET_FP // Floating-point elements.
};
+// Head model used for calculating the impulse delays.
+enum HeadModelT {
+ HM_NONE = 0,
+ HM_DATASET , // Measure the onset from the dataset.
+ HM_SPHERE // Calculate the onset using a spherical head model.
+};
+
// Desired output format from the command line.
enum OutputFormatT {
OF_NONE = 0,
@@ -239,6 +252,7 @@ typedef unsigned long long uint8;
typedef enum ByteOrderT ByteOrderT;
typedef enum SourceFormatT SourceFormatT;
typedef enum ElementTypeT ElementTypeT;
+typedef enum HeadModelT HeadModelT;
typedef enum OutputFormatT OutputFormatT;
typedef struct TokenReaderT TokenReaderT;
@@ -724,8 +738,8 @@ static int StrSubst (const char * in, const char * pat, const char * rep, const
return (! truncated);
}
-// Provide missing math routines for MSVC.
-#ifdef _MSC_VER
+// Provide missing math routines for MSVC versions < 1800 (Visual Studio 2013).
+#if defined(_MSC_VER) && _MSC_VER < 1800
static double round (double val) {
if (val < 0.0)
return (ceil (val - 0.5));
@@ -1683,6 +1697,25 @@ static int LoadSource (SourceRefT * src, const uint hrirRate, const uint n, doub
return (result);
}
+// Calculate the onset time of an HRIR and average it with any existing
+// timing for its elevation and azimuth.
+static void AverageHrirOnset (const double * hrir, const double f, const uint ei, const uint ai, const HrirDataT * hData) {
+ double mag;
+ uint n, i, j;
+
+ mag = 0.0;
+ n = hData -> mIrPoints;
+ for (i = 0; i < n; i ++)
+ mag = fmax (fabs (hrir [i]), mag);
+ mag *= 0.15;
+ for (i = 0; i < n; i ++) {
+ if (fabs (hrir [i]) >= mag)
+ break;
+ }
+ j = hData -> mEvOffset [ei] + ai;
+ hData -> mHrtds [j] = Lerp (hData -> mHrtds [j], ((double) i) / hData -> mIrRate, f);
+}
+
// Calculate the magnitude response of an HRIR and average it with any
// existing responses for its elevation and azimuth.
static void AverageHrirMagnitude (const double * hrir, const double f, const uint ei, const uint ai, const HrirDataT * hData) {
@@ -1867,6 +1900,26 @@ static void CalcAzIndices (const HrirDataT * hData, const uint ei, const double
(* jf) = af;
}
+// Synthesize any missing onset timings at the bottom elevations. This just
+// blends between slightly exaggerated known onsets. Not an accurate model.
+static void SynthesizeOnsets (HrirDataT * hData) {
+ uint oi, e, a, j0, j1;
+ double t, of, jf;
+
+ oi = hData -> mEvStart;
+ t = 0.0;
+ for (a = 0; a < hData -> mAzCount [oi]; a ++)
+ t += hData -> mHrtds [hData -> mEvOffset [oi] + a];
+ hData -> mHrtds [0] = 1.32e-4 + (t / hData -> mAzCount [oi]);
+ for (e = 1; e < hData -> mEvStart; e ++) {
+ of = ((double) e) / hData -> mEvStart;
+ for (a = 0; a < hData -> mAzCount [e]; a ++) {
+ CalcAzIndices (hData, oi, a * 2.0 * M_PI / hData -> mAzCount [e], & j0, & j1, & jf);
+ hData -> mHrtds [hData -> mEvOffset [e] + a] = Lerp (hData -> mHrtds [0], Lerp (hData -> mHrtds [j0], hData -> mHrtds [j1], jf), of);
+ }
+ }
+}
+
/* Attempt to synthesize any missing HRIRs at the bottom elevations. Right
* now this just blends the lowest elevation HRIRs together and applies some
* attenuation and high frequency damping. It is a simple, if inaccurate
@@ -1966,9 +2019,9 @@ static double CalcLTD (const double ev, const double az, const double rad, const
return (dlp / 343.3);
}
-// Calculate the effective head-related time delays for the each HRIR, now
-// that they are minimum-phase.
-static void CalculateHrtds (HrirDataT * hData) {
+// Calculate the effective head-related time delays for each minimum-phase
+// HRIR.
+static void CalculateHrtds (const HeadModelT model, const double radius, HrirDataT * hData) {
double minHrtd, maxHrtd;
uint e, a, j;
double t;
@@ -1978,9 +2031,13 @@ static void CalculateHrtds (HrirDataT * hData) {
for (e = 0; e < hData -> mEvCount; e ++) {
for (a = 0; a < hData -> mAzCount [e]; a ++) {
j = hData -> mEvOffset [e] + a;
- t = CalcLTD ((-90.0 + (e * 180.0 / (hData -> mEvCount - 1))) * M_PI / 180.0,
- (a * 360.0 / hData -> mAzCount [e]) * M_PI / 180.0,
- hData -> mRadius, hData -> mDistance);
+ if (model == HM_DATASET) {
+ t = hData -> mHrtds [j] * radius / hData -> mRadius;
+ } else {
+ t = CalcLTD ((-90.0 + (e * 180.0 / (hData -> mEvCount - 1))) * M_PI / 180.0,
+ (a * 360.0 / hData -> mAzCount [e]) * M_PI / 180.0,
+ radius, hData -> mDistance);
+ }
hData -> mHrtds [j] = t;
maxHrtd = fmax (t, maxHrtd);
minHrtd = fmin (t, minHrtd);
@@ -2372,7 +2429,7 @@ static int ReadSourceRef (TokenReaderT * tr, SourceRefT * src) {
}
// Process the list of sources in the data set definition.
-static int ProcessSources (TokenReaderT * tr, HrirDataT * hData) {
+static int ProcessSources (const HeadModelT model, TokenReaderT * tr, HrirDataT * hData) {
uint * setCount = NULL, * setFlag = NULL;
double * hrir = NULL;
uint line, col, ei, ai;
@@ -2393,6 +2450,8 @@ static int ProcessSources (TokenReaderT * tr, HrirDataT * hData) {
for (;;) {
if (ReadSourceRef (tr, & src)) {
if (LoadSource (& src, hData -> mIrRate, hData -> mIrPoints, hrir)) {
+ if (model == HM_DATASET)
+ AverageHrirOnset (hrir, 1.0 / factor, ei, ai, hData);
AverageHrirMagnitude (hrir, 1.0 / factor, ei, ai, hData);
factor += 1.0;
if (! TrIsOperator (tr, "+"))
@@ -2452,7 +2511,7 @@ static int ProcessSources (TokenReaderT * tr, HrirDataT * hData) {
* resulting data set as desired. If the input name is NULL it will read
* from standard input.
*/
-static int ProcessDefinition (const char * inName, const uint outRate, const uint fftSize, const int equalize, const int surface, const double limit, const uint truncSize, const OutputFormatT outFormat, const char * outName) {
+static int ProcessDefinition (const char * inName, const uint outRate, const uint fftSize, const int equalize, const int surface, const double limit, const uint truncSize, const HeadModelT model, const double radius, const OutputFormatT outFormat, const char * outName) {
FILE * fp = NULL;
TokenReaderT tr;
HrirDataT hData;
@@ -2486,7 +2545,7 @@ static int ProcessDefinition (const char * inName, const uint outRate, const uin
}
hData . mHrirs = CreateArray (hData . mIrCount * hData . mIrSize);
hData . mHrtds = CreateArray (hData . mIrCount);
- if (! ProcessSources (& tr, & hData)) {
+ if (! ProcessSources (model, & tr, & hData)) {
DestroyArray (hData . mHrtds);
DestroyArray (hData . mHrirs);
if (inName != NULL)
@@ -2512,11 +2571,13 @@ static int ProcessDefinition (const char * inName, const uint outRate, const uin
fprintf (stdout, "Truncating minimum-phase HRIRs...\n");
hData . mIrPoints = truncSize;
fprintf (stdout, "Synthesizing missing elevations...\n");
+ if (model == HM_DATASET)
+ SynthesizeOnsets (& hData);
SynthesizeHrirs (& hData);
fprintf (stdout, "Normalizing final HRIRs...\n");
NormalizeHrirs (& hData);
fprintf (stdout, "Calculating impulse delays...\n");
- CalculateHrtds (& hData);
+ CalculateHrtds (model, (radius > DEFAULT_CUSTOM_RADIUS) ? radius : hData . mRadius, & hData);
snprintf (rateStr, 8, "%u", hData . mIrRate);
StrSubst (outName, "%r", rateStr, MAX_PATH_LEN, expName);
switch (outFormat) {
@@ -2547,6 +2608,8 @@ int main (const int argc, const char * argv []) {
int equalize, surface;
double limit;
uint truncSize;
+ HeadModelT model;
+ double radius;
char * end = NULL;
if (argc < 2) {
@@ -2573,6 +2636,9 @@ int main (const int argc, const char * argv []) {
fprintf (stdout, " average (default: %.2f).\n", DEFAULT_LIMIT);
fprintf (stdout, " -w=<points> Specify the size of the truncation window that's applied\n");
fprintf (stdout, " after minimum-phase reconstruction (default: %u).\n", DEFAULT_TRUNCSIZE);
+ fprintf (stdout, " -d={dataset| Specify the model used for calculating the head-delay timing\n");
+ fprintf (stdout, " sphere} values (default: %s).\n", ((DEFAULT_HEAD_MODEL == HM_DATASET) ? "dataset" : "sphere"));
+ fprintf (stdout, " -c=<size> Use a customized head radius measured ear-to-ear in meters.\n");
fprintf (stdout, " -i=<filename> Specify an HRIR definition file to use (defaults to stdin).\n");
fprintf (stdout, " -o=<filename> Specify an output file. Overrides command-selected default.\n");
fprintf (stdout, " Use of '%%r' will be substituted with the data set sample rate.\n");
@@ -2601,6 +2667,8 @@ int main (const int argc, const char * argv []) {
surface = DEFAULT_SURFACE;
limit = DEFAULT_LIMIT;
truncSize = DEFAULT_TRUNCSIZE;
+ model = DEFAULT_HEAD_MODEL;
+ radius = DEFAULT_CUSTOM_RADIUS;
while (argi < argc) {
if (strncmp (argv [argi], "-r=", 3) == 0) {
outRate = strtoul (& argv [argi] [3], & end, 10);
@@ -2648,6 +2716,21 @@ int main (const int argc, const char * argv []) {
fprintf (stderr, "Error: Expected a value from %u to %u in multiples of %u for '-w'.\n", MIN_TRUNCSIZE, MAX_TRUNCSIZE, MOD_TRUNCSIZE);
return (-1);
}
+ } else if (strncmp (argv [argi], "-d=", 3) == 0) {
+ if (strcmp (& argv [argi] [3], "dataset") == 0) {
+ model = HM_DATASET;
+ } else if (strcmp (& argv [argi] [3], "sphere") == 0) {
+ model = HM_SPHERE;
+ } else {
+ fprintf (stderr, "Error: Expected 'dataset' or 'sphere' for '-d'.\n");
+ return (-1);
+ }
+ } else if (strncmp (argv [argi], "-c=", 3) == 0) {
+ radius = strtod (& argv [argi] [3], & end);
+ if ((end [0] != '\0') || (radius < MIN_CUSTOM_RADIUS) || (radius > MAX_CUSTOM_RADIUS)) {
+ fprintf (stderr, "Error: Expected a value from %.2f to %.2f for '-c'.\n", MIN_CUSTOM_RADIUS, MAX_CUSTOM_RADIUS);
+ return (-1);
+ }
} else if (strncmp (argv [argi], "-i=", 3) == 0) {
inName = & argv [argi] [3];
} else if (strncmp (argv [argi], "-o=", 3) == 0) {
@@ -2658,7 +2741,7 @@ int main (const int argc, const char * argv []) {
}
argi ++;
}
- if (! ProcessDefinition (inName, outRate, fftSize, equalize, surface, limit, truncSize, outFormat, outName))
+ if (! ProcessDefinition (inName, outRate, fftSize, equalize, surface, limit, truncSize, model, radius, outFormat, outName))
return (-1);
fprintf (stdout, "Operation completed.\n");
return (0);
generated by cgit v1.2.3 (git 2.39.1) at 2025-04-23 12:30:14 +0000