diff options
Diffstat (limited to 'utils/makemhr')
| -rw-r--r-- | utils/makemhr/loaddef.cpp | 1976 | ||||
| -rw-r--r-- | utils/makemhr/loaddef.h | 30 | ||||
| -rw-r--r-- | utils/makemhr/loadsofa.cpp | 4 | ||||
| -rw-r--r-- | utils/makemhr/loadsofa.h | 6 | ||||
| -rw-r--r-- | utils/makemhr/makemhr.cpp | 2088 | ||||
| -rw-r--r-- | utils/makemhr/makemhr.h | 128 |
6 files changed, 2155 insertions, 2077 deletions
diff --git a/utils/makemhr/loaddef.cpp b/utils/makemhr/loaddef.cpp index e69de29b..d4e8db27 100644 --- a/utils/makemhr/loaddef.cpp +++ b/utils/makemhr/loaddef.cpp @@ -0,0 +1,1976 @@ +/* + * HRTF utility for producing and demonstrating the process of creating an + * OpenAL Soft compatible HRIR data set. + * + * Copyright (C) 2011-2019 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 + * the Free Software Foundation; either version 2 of the License, or + * (at your option) any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License along + * with this program; if not, write to the Free Software Foundation, Inc., + * 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. + * + * Or visit: http://www.gnu.org/licenses/old-licenses/gpl-2.0.html + */ + +#include <cstring> +#include <limits> + +#include "mysofa.h" + +#include "loaddef.h" + + +// The token reader's load interval in bytes. +#define TR_LOAD_SIZE (TR_RING_SIZE >> 2) + +// The maximum identifier length used when processing the data set +// definition. +#define MAX_IDENT_LEN (16) + +// The limits for the listener's head 'radius' in the data set definition. +#define MIN_RADIUS (0.05) +#define MAX_RADIUS (0.15) + +// The maximum number of channels that can be addressed for a WAVE file +// source listed in the data set definition. +#define MAX_WAVE_CHANNELS (65535) + +// The limits to the byte size for a binary source listed in the definition +// file. +#define MIN_BIN_SIZE (2) +#define MAX_BIN_SIZE (4) + +// The minimum number of significant bits for binary sources listed in the +// data set definition. The maximum is calculated from the byte size. +#define MIN_BIN_BITS (16) + +// The limits to the number of significant bits for an ASCII source listed in +// the data set definition. +#define MIN_ASCII_BITS (16) +#define MAX_ASCII_BITS (32) + +// The four-character-codes for RIFF/RIFX WAVE file chunks. +#define FOURCC_RIFF (0x46464952) // 'RIFF' +#define FOURCC_RIFX (0x58464952) // 'RIFX' +#define FOURCC_WAVE (0x45564157) // 'WAVE' +#define FOURCC_FMT (0x20746D66) // 'fmt ' +#define FOURCC_DATA (0x61746164) // 'data' +#define FOURCC_LIST (0x5453494C) // 'LIST' +#define FOURCC_WAVL (0x6C766177) // 'wavl' +#define FOURCC_SLNT (0x746E6C73) // 'slnt' + +// The supported wave formats. +#define WAVE_FORMAT_PCM (0x0001) +#define WAVE_FORMAT_IEEE_FLOAT (0x0003) +#define WAVE_FORMAT_EXTENSIBLE (0xFFFE) + + +enum ByteOrderT { + BO_NONE, + BO_LITTLE, + BO_BIG +}; + +// Source format for the references listed in the data set definition. +enum SourceFormatT { + SF_NONE, + SF_ASCII, // ASCII text file. + SF_BIN_LE, // Little-endian binary file. + SF_BIN_BE, // Big-endian binary file. + SF_WAVE, // RIFF/RIFX WAVE file. + SF_SOFA // Spatially Oriented Format for Accoustics (SOFA) file. +}; + +// Element types for the references listed in the data set definition. +enum ElementTypeT { + ET_NONE, + ET_INT, // Integer elements. + ET_FP // Floating-point elements. +}; + +// Source reference state used when loading sources. +struct SourceRefT { + SourceFormatT mFormat; + ElementTypeT mType; + uint mSize; + int mBits; + uint mChannel; + double mAzimuth; + double mElevation; + double mRadius; + uint mSkip; + uint mOffset; + char mPath[MAX_PATH_LEN+1]; +}; + + +/* Whitespace is not significant. It can process tokens as identifiers, numbers + * (integer and floating-point), strings, and operators. Strings must be + * encapsulated by double-quotes and cannot span multiple lines. + */ + +// Setup the reader on the given file. The filename can be NULL if no error +// output is desired. +void TrSetup(FILE *fp, const char *filename, TokenReaderT *tr) +{ + const char *name = nullptr; + + if(filename) + { + const char *slash = strrchr(filename, '/'); + if(slash) + { + const char *bslash = strrchr(slash+1, '\\'); + if(bslash) name = bslash+1; + else name = slash+1; + } + else + { + const char *bslash = strrchr(filename, '\\'); + if(bslash) name = bslash+1; + else name = filename; + } + } + + tr->mFile = fp; + tr->mName = name; + tr->mLine = 1; + tr->mColumn = 1; + tr->mIn = 0; + tr->mOut = 0; +} + +// Prime the reader's ring buffer, and return a result indicating that there +// is text to process. +static int TrLoad(TokenReaderT *tr) +{ + size_t toLoad, in, count; + + toLoad = TR_RING_SIZE - (tr->mIn - tr->mOut); + if(toLoad >= TR_LOAD_SIZE && !feof(tr->mFile)) + { + // Load TR_LOAD_SIZE (or less if at the end of the file) per read. + toLoad = TR_LOAD_SIZE; + in = tr->mIn&TR_RING_MASK; + count = TR_RING_SIZE - in; + if(count < toLoad) + { + tr->mIn += fread(&tr->mRing[in], 1, count, tr->mFile); + tr->mIn += fread(&tr->mRing[0], 1, toLoad-count, tr->mFile); + } + else + tr->mIn += fread(&tr->mRing[in], 1, toLoad, tr->mFile); + + if(tr->mOut >= TR_RING_SIZE) + { + tr->mOut -= TR_RING_SIZE; + tr->mIn -= TR_RING_SIZE; + } + } + if(tr->mIn > tr->mOut) + return 1; + return 0; +} + +// Error display routine. Only displays when the base name is not NULL. +static void TrErrorVA(const TokenReaderT *tr, uint line, uint column, const char *format, va_list argPtr) +{ + if(!tr->mName) + return; + fprintf(stderr, "\nError (%s:%u:%u): ", tr->mName, line, column); + vfprintf(stderr, format, argPtr); +} + +// Used to display an error at a saved line/column. +static void TrErrorAt(const TokenReaderT *tr, uint line, uint column, const char *format, ...) +{ + va_list argPtr; + + va_start(argPtr, format); + TrErrorVA(tr, line, column, format, argPtr); + va_end(argPtr); +} + +// Used to display an error at the current line/column. +static void TrError(const TokenReaderT *tr, const char *format, ...) +{ + va_list argPtr; + + va_start(argPtr, format); + TrErrorVA(tr, tr->mLine, tr->mColumn, format, argPtr); + va_end(argPtr); +} + +// Skips to the next line. +static void TrSkipLine(TokenReaderT *tr) +{ + char ch; + + while(TrLoad(tr)) + { + ch = tr->mRing[tr->mOut&TR_RING_MASK]; + tr->mOut++; + if(ch == '\n') + { + tr->mLine++; + tr->mColumn = 1; + break; + } + tr->mColumn ++; + } +} + +// Skips to the next token. +static int TrSkipWhitespace(TokenReaderT *tr) +{ + while(TrLoad(tr)) + { + char ch{tr->mRing[tr->mOut&TR_RING_MASK]}; + if(isspace(ch)) + { + tr->mOut++; + if(ch == '\n') + { + tr->mLine++; + tr->mColumn = 1; + } + else + tr->mColumn++; + } + else if(ch == '#') + TrSkipLine(tr); + else + return 1; + } + return 0; +} + +// Get the line and/or column of the next token (or the end of input). +static void TrIndication(TokenReaderT *tr, uint *line, uint *column) +{ + TrSkipWhitespace(tr); + if(line) *line = tr->mLine; + if(column) *column = tr->mColumn; +} + +// Checks to see if a token is (likely to be) an identifier. It does not +// display any errors and will not proceed to the next token. +static int TrIsIdent(TokenReaderT *tr) +{ + if(!TrSkipWhitespace(tr)) + return 0; + char ch{tr->mRing[tr->mOut&TR_RING_MASK]}; + return ch == '_' || isalpha(ch); +} + + +// Checks to see if a token is the given operator. It does not display any +// errors and will not proceed to the next token. +static int TrIsOperator(TokenReaderT *tr, const char *op) +{ + size_t out, len; + char ch; + + if(!TrSkipWhitespace(tr)) + return 0; + out = tr->mOut; + len = 0; + while(op[len] != '\0' && out < tr->mIn) + { + ch = tr->mRing[out&TR_RING_MASK]; + if(ch != op[len]) break; + len++; + out++; + } + if(op[len] == '\0') + return 1; + return 0; +} + +/* The TrRead*() routines obtain the value of a matching token type. They + * display type, form, and boundary errors and will proceed to the next + * token. + */ + +// Reads and validates an identifier token. +static int TrReadIdent(TokenReaderT *tr, const uint maxLen, char *ident) +{ + uint col, len; + char ch; + + col = tr->mColumn; + if(TrSkipWhitespace(tr)) + { + col = tr->mColumn; + ch = tr->mRing[tr->mOut&TR_RING_MASK]; + if(ch == '_' || isalpha(ch)) + { + len = 0; + do { + if(len < maxLen) + ident[len] = ch; + len++; + tr->mOut++; + if(!TrLoad(tr)) + break; + ch = tr->mRing[tr->mOut&TR_RING_MASK]; + } while(ch == '_' || isdigit(ch) || isalpha(ch)); + + tr->mColumn += len; + if(len < maxLen) + { + ident[len] = '\0'; + return 1; + } + TrErrorAt(tr, tr->mLine, col, "Identifier is too long.\n"); + return 0; + } + } + TrErrorAt(tr, tr->mLine, col, "Expected an identifier.\n"); + return 0; +} + +// Reads and validates (including bounds) an integer token. +static int TrReadInt(TokenReaderT *tr, const int loBound, const int hiBound, int *value) +{ + uint col, digis, len; + char ch, temp[64+1]; + + col = tr->mColumn; + if(TrSkipWhitespace(tr)) + { + col = tr->mColumn; + len = 0; + ch = tr->mRing[tr->mOut&TR_RING_MASK]; + if(ch == '+' || ch == '-') + { + temp[len] = ch; + len++; + tr->mOut++; + } + digis = 0; + while(TrLoad(tr)) + { + ch = tr->mRing[tr->mOut&TR_RING_MASK]; + if(!isdigit(ch)) break; + if(len < 64) + temp[len] = ch; + len++; + digis++; + tr->mOut++; + } + tr->mColumn += len; + if(digis > 0 && ch != '.' && !isalpha(ch)) + { + if(len > 64) + { + TrErrorAt(tr, tr->mLine, col, "Integer is too long."); + return 0; + } + temp[len] = '\0'; + *value = strtol(temp, nullptr, 10); + if(*value < loBound || *value > hiBound) + { + TrErrorAt(tr, tr->mLine, col, "Expected a value from %d to %d.\n", loBound, hiBound); + return 0; + } + return 1; + } + } + TrErrorAt(tr, tr->mLine, col, "Expected an integer.\n"); + return 0; +} + +// Reads and validates (including bounds) a float token. +static int TrReadFloat(TokenReaderT *tr, const double loBound, const double hiBound, double *value) +{ + uint col, digis, len; + char ch, temp[64+1]; + + col = tr->mColumn; + if(TrSkipWhitespace(tr)) + { + col = tr->mColumn; + len = 0; + ch = tr->mRing[tr->mOut&TR_RING_MASK]; + if(ch == '+' || ch == '-') + { + temp[len] = ch; + len++; + tr->mOut++; + } + + digis = 0; + while(TrLoad(tr)) + { + ch = tr->mRing[tr->mOut&TR_RING_MASK]; + if(!isdigit(ch)) break; + if(len < 64) + temp[len] = ch; + len++; + digis++; + tr->mOut++; + } + if(ch == '.') + { + if(len < 64) + temp[len] = ch; + len++; + tr->mOut++; + } + while(TrLoad(tr)) + { + ch = tr->mRing[tr->mOut&TR_RING_MASK]; + if(!isdigit(ch)) break; + if(len < 64) + temp[len] = ch; + len++; + digis++; + tr->mOut++; + } + if(digis > 0) + { + if(ch == 'E' || ch == 'e') + { + if(len < 64) + temp[len] = ch; + len++; + digis = 0; + tr->mOut++; + if(ch == '+' || ch == '-') + { + if(len < 64) + temp[len] = ch; + len++; + tr->mOut++; + } + while(TrLoad(tr)) + { + ch = tr->mRing[tr->mOut&TR_RING_MASK]; + if(!isdigit(ch)) break; + if(len < 64) + temp[len] = ch; + len++; + digis++; + tr->mOut++; + } + } + tr->mColumn += len; + if(digis > 0 && ch != '.' && !isalpha(ch)) + { + if(len > 64) + { + TrErrorAt(tr, tr->mLine, col, "Float is too long."); + return 0; + } + temp[len] = '\0'; + *value = strtod(temp, nullptr); + if(*value < loBound || *value > hiBound) + { + TrErrorAt(tr, tr->mLine, col, "Expected a value from %f to %f.\n", loBound, hiBound); + return 0; + } + return 1; + } + } + else + tr->mColumn += len; + } + TrErrorAt(tr, tr->mLine, col, "Expected a float.\n"); + return 0; +} + +// Reads and validates a string token. +static int TrReadString(TokenReaderT *tr, const uint maxLen, char *text) +{ + uint col, len; + char ch; + + col = tr->mColumn; + if(TrSkipWhitespace(tr)) + { + col = tr->mColumn; + ch = tr->mRing[tr->mOut&TR_RING_MASK]; + if(ch == '\"') + { + tr->mOut++; + len = 0; + while(TrLoad(tr)) + { + ch = tr->mRing[tr->mOut&TR_RING_MASK]; + tr->mOut++; + if(ch == '\"') + break; + if(ch == '\n') + { + TrErrorAt(tr, tr->mLine, col, "Unterminated string at end of line.\n"); + return 0; + } + if(len < maxLen) + text[len] = ch; + len++; + } + if(ch != '\"') + { + tr->mColumn += 1 + len; + TrErrorAt(tr, tr->mLine, col, "Unterminated string at end of input.\n"); + return 0; + } + tr->mColumn += 2 + len; + if(len > maxLen) + { + TrErrorAt(tr, tr->mLine, col, "String is too long.\n"); + return 0; + } + text[len] = '\0'; + return 1; + } + } + TrErrorAt(tr, tr->mLine, col, "Expected a string.\n"); + return 0; +} + +// Reads and validates the given operator. +static int TrReadOperator(TokenReaderT *tr, const char *op) +{ + uint col, len; + char ch; + + col = tr->mColumn; + if(TrSkipWhitespace(tr)) + { + col = tr->mColumn; + len = 0; + while(op[len] != '\0' && TrLoad(tr)) + { + ch = tr->mRing[tr->mOut&TR_RING_MASK]; + if(ch != op[len]) break; + len++; + tr->mOut++; + } + tr->mColumn += len; + if(op[len] == '\0') + return 1; + } + TrErrorAt(tr, tr->mLine, col, "Expected '%s' operator.\n", op); + return 0; +} + + +/************************* + *** File source input *** + *************************/ + +// Read a binary value of the specified byte order and byte size from a file, +// storing it as a 32-bit unsigned integer. +static int ReadBin4(FILE *fp, const char *filename, const ByteOrderT order, const uint bytes, uint32_t *out) +{ + uint8_t in[4]; + uint32_t accum; + uint i; + + if(fread(in, 1, bytes, fp) != bytes) + { + fprintf(stderr, "\nError: Bad read from file '%s'.\n", filename); + return 0; + } + accum = 0; + switch(order) + { + case BO_LITTLE: + for(i = 0;i < bytes;i++) + accum = (accum<<8) | in[bytes - i - 1]; + break; + case BO_BIG: + for(i = 0;i < bytes;i++) + accum = (accum<<8) | in[i]; + break; + default: + break; + } + *out = accum; + return 1; +} + +// Read a binary value of the specified byte order from a file, storing it as +// a 64-bit unsigned integer. +static int ReadBin8(FILE *fp, const char *filename, const ByteOrderT order, uint64_t *out) +{ + uint8_t in[8]; + uint64_t accum; + uint i; + + if(fread(in, 1, 8, fp) != 8) + { + fprintf(stderr, "\nError: Bad read from file '%s'.\n", filename); + return 0; + } + accum = 0ULL; + switch(order) + { + case BO_LITTLE: + for(i = 0;i < 8;i++) + accum = (accum<<8) | in[8 - i - 1]; + break; + case BO_BIG: + for(i = 0;i < 8;i++) + accum = (accum<<8) | in[i]; + break; + default: + break; + } + *out = accum; + return 1; +} + +/* Read a binary value of the specified type, byte order, and byte size from + * a file, converting it to a double. For integer types, the significant + * bits are used to normalize the result. The sign of bits determines + * whether they are padded toward the MSB (negative) or LSB (positive). + * Floating-point types are not normalized. + */ +static int ReadBinAsDouble(FILE *fp, const char *filename, const ByteOrderT order, const ElementTypeT type, const uint bytes, const int bits, double *out) +{ + union { + uint32_t ui; + int32_t i; + float f; + } v4; + union { + uint64_t ui; + double f; + } v8; + + *out = 0.0; + if(bytes > 4) + { + if(!ReadBin8(fp, filename, order, &v8.ui)) + return 0; + if(type == ET_FP) + *out = v8.f; + } + else + { + if(!ReadBin4(fp, filename, order, bytes, &v4.ui)) + return 0; + if(type == ET_FP) + *out = v4.f; + else + { + if(bits > 0) + v4.ui >>= (8*bytes) - (static_cast<uint>(bits)); + else + v4.ui &= (0xFFFFFFFF >> (32+bits)); + + if(v4.ui&static_cast<uint>(1<<(std::abs(bits)-1))) + v4.ui |= (0xFFFFFFFF << std::abs(bits)); + *out = v4.i / static_cast<double>(1<<(std::abs(bits)-1)); + } + } + return 1; +} + +/* Read an ascii value of the specified type from a file, converting it to a + * double. For integer types, the significant bits are used to normalize the + * result. The sign of the bits should always be positive. This also skips + * up to one separator character before the element itself. + */ +static int ReadAsciiAsDouble(TokenReaderT *tr, const char *filename, const ElementTypeT type, const uint bits, double *out) +{ + if(TrIsOperator(tr, ",")) + TrReadOperator(tr, ","); + else if(TrIsOperator(tr, ":")) + TrReadOperator(tr, ":"); + else if(TrIsOperator(tr, ";")) + TrReadOperator(tr, ";"); + else if(TrIsOperator(tr, "|")) + TrReadOperator(tr, "|"); + + if(type == ET_FP) + { + if(!TrReadFloat(tr, -std::numeric_limits<double>::infinity(), + std::numeric_limits<double>::infinity(), out)) + { + fprintf(stderr, "\nError: Bad read from file '%s'.\n", filename); + return 0; + } + } + else + { + int v; + if(!TrReadInt(tr, -(1<<(bits-1)), (1<<(bits-1))-1, &v)) + { + fprintf(stderr, "\nError: Bad read from file '%s'.\n", filename); + return 0; + } + *out = v / static_cast<double>((1<<(bits-1))-1); + } + return 1; +} + +// Read the RIFF/RIFX WAVE format chunk from a file, validating it against +// the source parameters and data set metrics. +static int ReadWaveFormat(FILE *fp, const ByteOrderT order, const uint hrirRate, SourceRefT *src) +{ + uint32_t fourCC, chunkSize; + uint32_t format, channels, rate, dummy, block, size, bits; + + chunkSize = 0; + do { + if(chunkSize > 0) + fseek(fp, static_cast<long>(chunkSize), SEEK_CUR); + if(!ReadBin4(fp, src->mPath, BO_LITTLE, 4, &fourCC) || + !ReadBin4(fp, src->mPath, order, 4, &chunkSize)) + return 0; + } while(fourCC != FOURCC_FMT); + if(!ReadBin4(fp, src->mPath, order, 2, &format) || + !ReadBin4(fp, src->mPath, order, 2, &channels) || + !ReadBin4(fp, src->mPath, order, 4, &rate) || + !ReadBin4(fp, src->mPath, order, 4, &dummy) || + !ReadBin4(fp, src->mPath, order, 2, &block)) + return 0; + block /= channels; + if(chunkSize > 14) + { + if(!ReadBin4(fp, src->mPath, order, 2, &size)) + return 0; + size /= 8; + if(block > size) + size = block; + } + else + size = block; + if(format == WAVE_FORMAT_EXTENSIBLE) + { + fseek(fp, 2, SEEK_CUR); + if(!ReadBin4(fp, src->mPath, order, 2, &bits)) + return 0; + if(bits == 0) + bits = 8 * size; + fseek(fp, 4, SEEK_CUR); + if(!ReadBin4(fp, src->mPath, order, 2, &format)) + return 0; + fseek(fp, static_cast<long>(chunkSize - 26), SEEK_CUR); + } + else + { + bits = 8 * size; + if(chunkSize > 14) + fseek(fp, static_cast<long>(chunkSize - 16), SEEK_CUR); + else + fseek(fp, static_cast<long>(chunkSize - 14), SEEK_CUR); + } + if(format != WAVE_FORMAT_PCM && format != WAVE_FORMAT_IEEE_FLOAT) + { + fprintf(stderr, "\nError: Unsupported WAVE format in file '%s'.\n", src->mPath); + return 0; + } + if(src->mChannel >= channels) + { + fprintf(stderr, "\nError: Missing source channel in WAVE file '%s'.\n", src->mPath); + return 0; + } + if(rate != hrirRate) + { + fprintf(stderr, "\nError: Mismatched source sample rate in WAVE file '%s'.\n", src->mPath); + return 0; + } + if(format == WAVE_FORMAT_PCM) + { + if(size < 2 || size > 4) + { + fprintf(stderr, "\nError: Unsupported sample size in WAVE file '%s'.\n", src->mPath); + return 0; + } + if(bits < 16 || bits > (8*size)) + { + fprintf(stderr, "\nError: Bad significant bits in WAVE file '%s'.\n", src->mPath); + return 0; + } + src->mType = ET_INT; + } + else + { + if(size != 4 && size != 8) + { + fprintf(stderr, "\nError: Unsupported sample size in WAVE file '%s'.\n", src->mPath); + return 0; + } + src->mType = ET_FP; + } + src->mSize = size; + src->mBits = static_cast<int>(bits); + src->mSkip = channels; + return 1; +} + +// Read a RIFF/RIFX WAVE data chunk, converting all elements to doubles. +static int ReadWaveData(FILE *fp, const SourceRefT *src, const ByteOrderT order, const uint n, double *hrir) +{ + int pre, post, skip; + uint i; + + pre = static_cast<int>(src->mSize * src->mChannel); + post = static_cast<int>(src->mSize * (src->mSkip - src->mChannel - 1)); + skip = 0; + for(i = 0;i < n;i++) + { + skip += pre; + if(skip > 0) + fseek(fp, skip, SEEK_CUR); + if(!ReadBinAsDouble(fp, src->mPath, order, src->mType, src->mSize, src->mBits, &hrir[i])) + return 0; + skip = post; + } + if(skip > 0) + fseek(fp, skip, SEEK_CUR); + return 1; +} + +// Read the RIFF/RIFX WAVE list or data chunk, converting all elements to +// doubles. +static int ReadWaveList(FILE *fp, const SourceRefT *src, const ByteOrderT order, const uint n, double *hrir) +{ + uint32_t fourCC, chunkSize, listSize, count; + uint block, skip, offset, i; + double lastSample; + + for(;;) + { + if(!ReadBin4(fp, src->mPath, BO_LITTLE, 4, &fourCC) || + !ReadBin4(fp, src->mPath, order, 4, &chunkSize)) + return 0; + + if(fourCC == FOURCC_DATA) + { + block = src->mSize * src->mSkip; + count = chunkSize / block; + if(count < (src->mOffset + n)) + { + fprintf(stderr, "\nError: Bad read from file '%s'.\n", src->mPath); + return 0; + } + fseek(fp, static_cast<long>(src->mOffset * block), SEEK_CUR); + if(!ReadWaveData(fp, src, order, n, &hrir[0])) + return 0; + return 1; + } + else if(fourCC == FOURCC_LIST) + { + if(!ReadBin4(fp, src->mPath, BO_LITTLE, 4, &fourCC)) + return 0; + chunkSize -= 4; + if(fourCC == FOURCC_WAVL) + break; + } + if(chunkSize > 0) + fseek(fp, static_cast<long>(chunkSize), SEEK_CUR); + } + listSize = chunkSize; + block = src->mSize * src->mSkip; + skip = src->mOffset; + offset = 0; + lastSample = 0.0; + while(offset < n && listSize > 8) + { + if(!ReadBin4(fp, src->mPath, BO_LITTLE, 4, &fourCC) || + !ReadBin4(fp, src->mPath, order, 4, &chunkSize)) + return 0; + listSize -= 8 + chunkSize; + if(fourCC == FOURCC_DATA) + { + count = chunkSize / block; + if(count > skip) + { + fseek(fp, static_cast<long>(skip * block), SEEK_CUR); + chunkSize -= skip * block; + count -= skip; + skip = 0; + if(count > (n - offset)) + count = n - offset; + if(!ReadWaveData(fp, src, order, count, &hrir[offset])) + return 0; + chunkSize -= count * block; + offset += count; + lastSample = hrir[offset - 1]; + } + else + { + skip -= count; + count = 0; + } + } + else if(fourCC == FOURCC_SLNT) + { + if(!ReadBin4(fp, src->mPath, order, 4, &count)) + return 0; + chunkSize -= 4; + if(count > skip) + { + count -= skip; + skip = 0; + if(count > (n - offset)) + count = n - offset; + for(i = 0; i < count; i ++) + hrir[offset + i] = lastSample; + offset += count; + } + else + { + skip -= count; + count = 0; + } + } + if(chunkSize > 0) + fseek(fp, static_cast<long>(chunkSize), SEEK_CUR); + } + if(offset < n) + { + fprintf(stderr, "\nError: Bad read from file '%s'.\n", src->mPath); + return 0; + } + return 1; +} + +// Load a source HRIR from an ASCII text file containing a list of elements +// separated by whitespace or common list operators (',', ';', ':', '|'). +static int LoadAsciiSource(FILE *fp, const SourceRefT *src, const uint n, double *hrir) +{ + TokenReaderT tr; + uint i, j; + double dummy; + + TrSetup(fp, nullptr, &tr); + for(i = 0;i < src->mOffset;i++) + { + if(!ReadAsciiAsDouble(&tr, src->mPath, src->mType, static_cast<uint>(src->mBits), &dummy)) + return 0; + } + for(i = 0;i < n;i++) + { + if(!ReadAsciiAsDouble(&tr, src->mPath, src->mType, static_cast<uint>(src->mBits), &hrir[i])) + return 0; + for(j = 0;j < src->mSkip;j++) + { + if(!ReadAsciiAsDouble(&tr, src->mPath, src->mType, static_cast<uint>(src->mBits), &dummy)) + return 0; + } + } + return 1; +} + +// Load a source HRIR from a binary file. +static int LoadBinarySource(FILE *fp, const SourceRefT *src, const ByteOrderT order, const uint n, double *hrir) +{ + uint i; + + fseek(fp, static_cast<long>(src->mOffset), SEEK_SET); + for(i = 0;i < n;i++) + { + if(!ReadBinAsDouble(fp, src->mPath, order, src->mType, src->mSize, src->mBits, &hrir[i])) + return 0; + if(src->mSkip > 0) + fseek(fp, static_cast<long>(src->mSkip), SEEK_CUR); + } + return 1; +} + +// Load a source HRIR from a RIFF/RIFX WAVE file. +static int LoadWaveSource(FILE *fp, SourceRefT *src, const uint hrirRate, const uint n, double *hrir) +{ + uint32_t fourCC, dummy; + ByteOrderT order; + + if(!ReadBin4(fp, src->mPath, BO_LITTLE, 4, &fourCC) || + !ReadBin4(fp, src->mPath, BO_LITTLE, 4, &dummy)) + return 0; + if(fourCC == FOURCC_RIFF) + order = BO_LITTLE; + else if(fourCC == FOURCC_RIFX) + order = BO_BIG; + else + { + fprintf(stderr, "\nError: No RIFF/RIFX chunk in file '%s'.\n", src->mPath); + return 0; + } + + if(!ReadBin4(fp, src->mPath, BO_LITTLE, 4, &fourCC)) + return 0; + if(fourCC != FOURCC_WAVE) + { + fprintf(stderr, "\nError: Not a RIFF/RIFX WAVE file '%s'.\n", src->mPath); + return 0; + } + if(!ReadWaveFormat(fp, order, hrirRate, src)) + return 0; + if(!ReadWaveList(fp, src, order, n, hrir)) + return 0; + return 1; +} + + + +// Load a Spatially Oriented Format for Accoustics (SOFA) file. +static MYSOFA_EASY* LoadSofaFile(SourceRefT *src, const uint hrirRate, const uint n) +{ + struct MYSOFA_EASY *sofa{mysofa_cache_lookup(src->mPath, (float)hrirRate)}; + if(sofa) return sofa; + + sofa = static_cast<MYSOFA_EASY*>(calloc(1, sizeof(*sofa))); + if(sofa == nullptr) + { + fprintf(stderr, "\nError: Out of memory.\n"); + return nullptr; + } + sofa->lookup = nullptr; + sofa->neighborhood = nullptr; + + int err; + sofa->hrtf = mysofa_load(src->mPath, &err); + if(!sofa->hrtf) + { + mysofa_close(sofa); + fprintf(stderr, "\nError: Could not load source file '%s'.\n", src->mPath); + return nullptr; + } + err = mysofa_check(sofa->hrtf); + if(err != MYSOFA_OK) +/* NOTE: Some valid SOFA files are failing this check. + { + mysofa_close(sofa); + fprintf(stderr, "\nError: Malformed source file '%s'.\n", src->mPath); + return nullptr; + }*/ + fprintf(stderr, "\nWarning: Supposedly malformed source file '%s'.\n", src->mPath); + if((src->mOffset + n) > sofa->hrtf->N) + { + mysofa_close(sofa); + fprintf(stderr, "\nError: Not enough samples in SOFA file '%s'.\n", src->mPath); + return nullptr; + } + if(src->mChannel >= sofa->hrtf->R) + { + mysofa_close(sofa); + fprintf(stderr, "\nError: Missing source receiver in SOFA file '%s'.\n", src->mPath); + return nullptr; + } + mysofa_tocartesian(sofa->hrtf); + sofa->lookup = mysofa_lookup_init(sofa->hrtf); + if(sofa->lookup == nullptr) + { + mysofa_close(sofa); + fprintf(stderr, "\nError: Out of memory.\n"); + return nullptr; + } + return mysofa_cache_store(sofa, src->mPath, (float)hrirRate); +} + +// Copies the HRIR data from a particular SOFA measurement. +static void ExtractSofaHrir(const MYSOFA_EASY *sofa, const uint index, const uint channel, const uint offset, const uint n, double *hrir) +{ + for(uint i{0u};i < n;i++) + hrir[i] = sofa->hrtf->DataIR.values[(index*sofa->hrtf->R + channel)*sofa->hrtf->N + offset + i]; +} + +// Load a source HRIR from a Spatially Oriented Format for Accoustics (SOFA) +// file. +static int LoadSofaSource(SourceRefT *src, const uint hrirRate, const uint n, double *hrir) +{ + struct MYSOFA_EASY *sofa; + float target[3]; + int nearest; + float *coords; + + sofa = LoadSofaFile(src, hrirRate, n); + if(sofa == nullptr) + return 0; + + /* NOTE: At some point it may be benficial or necessary to consider the + various coordinate systems, listener/source orientations, and + direciontal vectors defined in the SOFA file. + */ + target[0] = src->mAzimuth; + target[1] = src->mElevation; + target[2] = src->mRadius; + mysofa_s2c(target); + + nearest = mysofa_lookup(sofa->lookup, target); + if(nearest < 0) + { + fprintf(stderr, "\nError: Lookup failed in source file '%s'.\n", src->mPath); + return 0; + } + + coords = &sofa->hrtf->SourcePosition.values[3 * nearest]; + if(std::abs(coords[0] - target[0]) > 0.001 || std::abs(coords[1] - target[1]) > 0.001 || std::abs(coords[2] - target[2]) > 0.001) + { + fprintf(stderr, "\nError: No impulse response at coordinates (%.3fr, %.1fev, %.1faz) in file '%s'.\n", src->mRadius, src->mElevation, src->mAzimuth, src->mPath); + target[0] = coords[0]; + target[1] = coords[1]; + target[2] = coords[2]; + mysofa_c2s(target); + fprintf(stderr, " Nearest candidate at (%.3fr, %.1fev, %.1faz).\n", target[2], target[1], target[0]); + return 0; + } + + ExtractSofaHrir(sofa, nearest, src->mChannel, src->mOffset, n, hrir); + + return 1; +} + +// Load a source HRIR from a supported file type. +static int LoadSource(SourceRefT *src, const uint hrirRate, const uint n, double *hrir) +{ + FILE *fp{nullptr}; + if(src->mFormat != SF_SOFA) + { + if(src->mFormat == SF_ASCII) + fp = fopen(src->mPath, "r"); + else + fp = fopen(src->mPath, "rb"); + if(fp == nullptr) + { + fprintf(stderr, "\nError: Could not open source file '%s'.\n", src->mPath); + return 0; + } + } + int result; + switch(src->mFormat) + { + case SF_ASCII: + result = LoadAsciiSource(fp, src, n, hrir); + break; + case SF_BIN_LE: + result = LoadBinarySource(fp, src, BO_LITTLE, n, hrir); + break; + case SF_BIN_BE: + result = LoadBinarySource(fp, src, BO_BIG, n, hrir); + break; + case SF_WAVE: + result = LoadWaveSource(fp, src, hrirRate, n, hrir); + break; + case SF_SOFA: + result = LoadSofaSource(src, hrirRate, n, hrir); + break; + default: + result = 0; + } + if(fp) fclose(fp); + return result; +} + + +// Match the channel type from a given identifier. +static ChannelTypeT MatchChannelType(const char *ident) +{ + if(strcasecmp(ident, "mono") == 0) + return CT_MONO; + if(strcasecmp(ident, "stereo") == 0) + return CT_STEREO; + return CT_NONE; +} + + +// Process the data set definition to read and validate the data set metrics. +int ProcessMetrics(TokenReaderT *tr, const uint fftSize, const uint truncSize, HrirDataT *hData) +{ + int hasRate = 0, hasType = 0, hasPoints = 0, hasRadius = 0; + int hasDistance = 0, hasAzimuths = 0; + char ident[MAX_IDENT_LEN+1]; + uint line, col; + double fpVal; + uint points; + int intVal; + double distances[MAX_FD_COUNT]; + uint fdCount = 0; + uint evCounts[MAX_FD_COUNT]; + std::vector<uint> azCounts(MAX_FD_COUNT * MAX_EV_COUNT); + + TrIndication(tr, &line, &col); + while(TrIsIdent(tr)) + { + TrIndication(tr, &line, &col); + if(!TrReadIdent(tr, MAX_IDENT_LEN, ident)) + return 0; + if(strcasecmp(ident, "rate") == 0) + { + if(hasRate) + { + TrErrorAt(tr, line, col, "Redefinition of 'rate'.\n"); + return 0; + } + if(!TrReadOperator(tr, "=")) + return 0; + if(!TrReadInt(tr, MIN_RATE, MAX_RATE, &intVal)) + return 0; + hData->mIrRate = static_cast<uint>(intVal); + hasRate = 1; + } + else if(strcasecmp(ident, "type") == 0) + { + char type[MAX_IDENT_LEN+1]; + + if(hasType) + { + TrErrorAt(tr, line, col, "Redefinition of 'type'.\n"); + return 0; + } + if(!TrReadOperator(tr, "=")) + return 0; + + if(!TrReadIdent(tr, MAX_IDENT_LEN, type)) + return 0; + hData->mChannelType = MatchChannelType(type); + if(hData->mChannelType == CT_NONE) + { + TrErrorAt(tr, line, col, "Expected a channel type.\n"); + return 0; + } + hasType = 1; + } + else if(strcasecmp(ident, "points") == 0) + { + if(hasPoints) + { + TrErrorAt(tr, line, col, "Redefinition of 'points'.\n"); + return 0; + } + if(!TrReadOperator(tr, "=")) + return 0; + TrIndication(tr, &line, &col); + if(!TrReadInt(tr, MIN_POINTS, MAX_POINTS, &intVal)) + return 0; + points = static_cast<uint>(intVal); + if(fftSize > 0 && points > fftSize) + { + TrErrorAt(tr, line, col, "Value exceeds the overridden FFT size.\n"); + return 0; + } + if(points < truncSize) + { + TrErrorAt(tr, line, col, "Value is below the truncation size.\n"); + return 0; + } + hData->mIrPoints = points; + hData->mFftSize = fftSize; + hData->mIrSize = 1 + (fftSize / 2); + if(points > hData->mIrSize) + hData->mIrSize = points; + hasPoints = 1; + } + else if(strcasecmp(ident, "radius") == 0) + { + if(hasRadius) + { + TrErrorAt(tr, line, col, "Redefinition of 'radius'.\n"); + return 0; + } + if(!TrReadOperator(tr, "=")) + return 0; + if(!TrReadFloat(tr, MIN_RADIUS, MAX_RADIUS, &fpVal)) + return 0; + hData->mRadius = fpVal; + hasRadius = 1; + } + else if(strcasecmp(ident, "distance") == 0) + { + uint count = 0; + + if(hasDistance) + { + TrErrorAt(tr, line, col, "Redefinition of 'distance'.\n"); + return 0; + } + if(!TrReadOperator(tr, "=")) + return 0; + + for(;;) + { + if(!TrReadFloat(tr, MIN_DISTANCE, MAX_DISTANCE, &fpVal)) + return 0; + if(count > 0 && fpVal <= distances[count - 1]) + { + TrError(tr, "Distances are not ascending.\n"); + return 0; + } + distances[count++] = fpVal; + if(!TrIsOperator(tr, ",")) + break; + if(count >= MAX_FD_COUNT) + { + TrError(tr, "Exceeded the maximum of %d fields.\n", MAX_FD_COUNT); + return 0; + } + TrReadOperator(tr, ","); + } + if(fdCount != 0 && count != fdCount) + { + TrError(tr, "Did not match the specified number of %d fields.\n", fdCount); + return 0; + } + fdCount = count; + hasDistance = 1; + } + else if(strcasecmp(ident, "azimuths") == 0) + { + uint count = 0; + + if(hasAzimuths) + { + TrErrorAt(tr, line, col, "Redefinition of 'azimuths'.\n"); + return 0; + } + if(!TrReadOperator(tr, "=")) + return 0; + + evCounts[0] = 0; + for(;;) + { + if(!TrReadInt(tr, MIN_AZ_COUNT, MAX_AZ_COUNT, &intVal)) + return 0; + azCounts[(count * MAX_EV_COUNT) + evCounts[count]++] = static_cast<uint>(intVal); + if(TrIsOperator(tr, ",")) + { + if(evCounts[count] >= MAX_EV_COUNT) + { + TrError(tr, "Exceeded the maximum of %d elevations.\n", MAX_EV_COUNT); + return 0; + } + TrReadOperator(tr, ","); + } + else + { + if(evCounts[count] < MIN_EV_COUNT) + { + TrErrorAt(tr, line, col, "Did not reach the minimum of %d azimuth counts.\n", MIN_EV_COUNT); + return 0; + } + if(azCounts[count * MAX_EV_COUNT] != 1 || azCounts[(count * MAX_EV_COUNT) + evCounts[count] - 1] != 1) + { + TrError(tr, "Poles are not singular for field %d.\n", count - 1); + return 0; + } + count++; + if(!TrIsOperator(tr, ";")) + break; + + if(count >= MAX_FD_COUNT) + { + TrError(tr, "Exceeded the maximum number of %d fields.\n", MAX_FD_COUNT); + return 0; + } + evCounts[count] = 0; + TrReadOperator(tr, ";"); + } + } + if(fdCount != 0 && count != fdCount) + { + TrError(tr, "Did not match the specified number of %d fields.\n", fdCount); + return 0; + } + fdCount = count; + hasAzimuths = 1; + } + else + { + TrErrorAt(tr, line, col, "Expected a metric name.\n"); + return 0; + } + TrSkipWhitespace(tr); + } + if(!(hasRate && hasPoints && hasRadius && hasDistance && hasAzimuths)) + { + TrErrorAt(tr, line, col, "Expected a metric name.\n"); + return 0; + } + if(distances[0] < hData->mRadius) + { + TrError(tr, "Distance cannot start below head radius.\n"); + return 0; + } + if(hData->mChannelType == CT_NONE) + hData->mChannelType = CT_MONO; + if(!PrepareHrirData(fdCount, distances, evCounts, azCounts.data(), hData)) + { + fprintf(stderr, "Error: Out of memory.\n"); + exit(-1); + } + return 1; +} + +// Parse an index triplet from the data set definition. +static int ReadIndexTriplet(TokenReaderT *tr, const HrirDataT *hData, uint *fi, uint *ei, uint *ai) +{ + int intVal; + + if(hData->mFdCount > 1) + { + if(!TrReadInt(tr, 0, static_cast<int>(hData->mFdCount) - 1, &intVal)) + return 0; + *fi = static_cast<uint>(intVal); + if(!TrReadOperator(tr, ",")) + return 0; + } + else + { + *fi = 0; + } + if(!TrReadInt(tr, 0, static_cast<int>(hData->mFds[*fi].mEvCount) - 1, &intVal)) + return 0; + *ei = static_cast<uint>(intVal); + if(!TrReadOperator(tr, ",")) + return 0; + if(!TrReadInt(tr, 0, static_cast<int>(hData->mFds[*fi].mEvs[*ei].mAzCount) - 1, &intVal)) + return 0; + *ai = static_cast<uint>(intVal); + return 1; +} + +// Match the source format from a given identifier. +static SourceFormatT MatchSourceFormat(const char *ident) +{ + if(strcasecmp(ident, "ascii") == 0) + return SF_ASCII; + if(strcasecmp(ident, "bin_le") == 0) + return SF_BIN_LE; + if(strcasecmp(ident, "bin_be") == 0) + return SF_BIN_BE; + if(strcasecmp(ident, "wave") == 0) + return SF_WAVE; + if(strcasecmp(ident, "sofa") == 0) + return SF_SOFA; + return SF_NONE; +} + +// Match the source element type from a given identifier. +static ElementTypeT MatchElementType(const char *ident) +{ + if(strcasecmp(ident, "int") == 0) + return ET_INT; + if(strcasecmp(ident, "fp") == 0) + return ET_FP; + return ET_NONE; +} + +// Parse and validate a source reference from the data set definition. +static int ReadSourceRef(TokenReaderT *tr, SourceRefT *src) +{ + char ident[MAX_IDENT_LEN+1]; + uint line, col; + double fpVal; + int intVal; + + TrIndication(tr, &line, &col); + if(!TrReadIdent(tr, MAX_IDENT_LEN, ident)) + return 0; + src->mFormat = MatchSourceFormat(ident); + if(src->mFormat == SF_NONE) + { + TrErrorAt(tr, line, col, "Expected a source format.\n"); + return 0; + } + if(!TrReadOperator(tr, "(")) + return 0; + if(src->mFormat == SF_SOFA) + { + if(!TrReadFloat(tr, MIN_DISTANCE, MAX_DISTANCE, &fpVal)) + return 0; + src->mRadius = fpVal; + if(!TrReadOperator(tr, ",")) + return 0; + if(!TrReadFloat(tr, -90.0, 90.0, &fpVal)) + return 0; + src->mElevation = fpVal; + if(!TrReadOperator(tr, ",")) + return 0; + if(!TrReadFloat(tr, -360.0, 360.0, &fpVal)) + return 0; + src->mAzimuth = fpVal; + if(!TrReadOperator(tr, ":")) + return 0; + if(!TrReadInt(tr, 0, MAX_WAVE_CHANNELS, &intVal)) + return 0; + src->mType = ET_NONE; + src->mSize = 0; + src->mBits = 0; + src->mChannel = (uint)intVal; + src->mSkip = 0; + } + else if(src->mFormat == SF_WAVE) + { + if(!TrReadInt(tr, 0, MAX_WAVE_CHANNELS, &intVal)) + return 0; + src->mType = ET_NONE; + src->mSize = 0; + src->mBits = 0; + src->mChannel = static_cast<uint>(intVal); + src->mSkip = 0; + } + else + { + TrIndication(tr, &line, &col); + if(!TrReadIdent(tr, MAX_IDENT_LEN, ident)) + return 0; + src->mType = MatchElementType(ident); + if(src->mType == ET_NONE) + { + TrErrorAt(tr, line, col, "Expected a source element type.\n"); + return 0; + } + if(src->mFormat == SF_BIN_LE || src->mFormat == SF_BIN_BE) + { + if(!TrReadOperator(tr, ",")) + return 0; + if(src->mType == ET_INT) + { + if(!TrReadInt(tr, MIN_BIN_SIZE, MAX_BIN_SIZE, &intVal)) + return 0; + src->mSize = static_cast<uint>(intVal); + if(!TrIsOperator(tr, ",")) + src->mBits = static_cast<int>(8*src->mSize); + else + { + TrReadOperator(tr, ","); + TrIndication(tr, &line, &col); + if(!TrReadInt(tr, -2147483647-1, 2147483647, &intVal)) + return 0; + if(std::abs(intVal) < MIN_BIN_BITS || static_cast<uint>(std::abs(intVal)) > (8*src->mSize)) + { + TrErrorAt(tr, line, col, "Expected a value of (+/-) %d to %d.\n", MIN_BIN_BITS, 8*src->mSize); + return 0; + } + src->mBits = intVal; + } + } + else + { + TrIndication(tr, &line, &col); + if(!TrReadInt(tr, -2147483647-1, 2147483647, &intVal)) + return 0; + if(intVal != 4 && intVal != 8) + { + TrErrorAt(tr, line, col, "Expected a value of 4 or 8.\n"); + return 0; + } + src->mSize = static_cast<uint>(intVal); + src->mBits = 0; + } + } + else if(src->mFormat == SF_ASCII && src->mType == ET_INT) + { + if(!TrReadOperator(tr, ",")) + return 0; + if(!TrReadInt(tr, MIN_ASCII_BITS, MAX_ASCII_BITS, &intVal)) + return 0; + src->mSize = 0; + src->mBits = intVal; + } + else + { + src->mSize = 0; + src->mBits = 0; + } + + if(!TrIsOperator(tr, ";")) + src->mSkip = 0; + else + { + TrReadOperator(tr, ";"); + if(!TrReadInt(tr, 0, 0x7FFFFFFF, &intVal)) + return 0; + src->mSkip = static_cast<uint>(intVal); + } + } + if(!TrReadOperator(tr, ")")) + return 0; + if(TrIsOperator(tr, "@")) + { + TrReadOperator(tr, "@"); + if(!TrReadInt(tr, 0, 0x7FFFFFFF, &intVal)) + return 0; + src->mOffset = static_cast<uint>(intVal); + } + else + src->mOffset = 0; + if(!TrReadOperator(tr, ":")) + return 0; + if(!TrReadString(tr, MAX_PATH_LEN, src->mPath)) + return 0; + return 1; +} + +// Parse and validate a SOFA source reference from the data set definition. +static int ReadSofaRef(TokenReaderT *tr, SourceRefT *src) +{ + char ident[MAX_IDENT_LEN+1]; + uint line, col; + int intVal; + + TrIndication(tr, &line, &col); + if(!TrReadIdent(tr, MAX_IDENT_LEN, ident)) + return 0; + src->mFormat = MatchSourceFormat(ident); + if(src->mFormat != SF_SOFA) + { + TrErrorAt(tr, line, col, "Expected the SOFA source format.\n"); + return 0; + } + + src->mType = ET_NONE; + src->mSize = 0; + src->mBits = 0; + src->mChannel = 0; + src->mSkip = 0; + + if(TrIsOperator(tr, "@")) + { + TrReadOperator(tr, "@"); + if(!TrReadInt(tr, 0, 0x7FFFFFFF, &intVal)) + return 0; + src->mOffset = (uint)intVal; + } + else + src->mOffset = 0; + if(!TrReadOperator(tr, ":")) + return 0; + if(!TrReadString(tr, MAX_PATH_LEN, src->mPath)) + return 0; + return 1; +} + +// Match the target ear (index) from a given identifier. +static int MatchTargetEar(const char *ident) +{ + if(strcasecmp(ident, "left") == 0) + return 0; + if(strcasecmp(ident, "right") == 0) + return 1; + return -1; +} + +// Calculate the onset time of an HRIR and average it with any existing +// timing for its field, elevation, azimuth, and ear. +static double AverageHrirOnset(const uint rate, const uint n, const double *hrir, const double f, const double onset) +{ + std::vector<double> upsampled(10 * n); + { + ResamplerT rs; + ResamplerSetup(&rs, rate, 10 * rate); + ResamplerRun(&rs, n, hrir, 10 * n, upsampled.data()); + } + + double mag{0.0}; + for(uint i{0u};i < 10*n;i++) + mag = std::max(std::abs(upsampled[i]), mag); + + mag *= 0.15; + uint i{0u}; + for(;i < 10*n;i++) + { + if(std::abs(upsampled[i]) >= mag) + break; + } + return Lerp(onset, static_cast<double>(i) / (10*rate), f); +} + +// Calculate the magnitude response of an HRIR and average it with any +// existing responses for its field, elevation, azimuth, and ear. +static void AverageHrirMagnitude(const uint points, const uint n, const double *hrir, const double f, double *mag) +{ + uint m = 1 + (n / 2), i; + std::vector<complex_d> h(n); + std::vector<double> r(n); + + for(i = 0;i < points;i++) + h[i] = complex_d{hrir[i], 0.0}; + for(;i < n;i++) + h[i] = complex_d{0.0, 0.0}; + FftForward(n, h.data()); + MagnitudeResponse(n, h.data(), r.data()); + for(i = 0;i < m;i++) + mag[i] = Lerp(mag[i], r[i], f); +} + +// Process the list of sources in the data set definition. +int ProcessSources(const HeadModelT model, TokenReaderT *tr, HrirDataT *hData) +{ + uint channels = (hData->mChannelType == CT_STEREO) ? 2 : 1; + hData->mHrirsBase.resize(channels * hData->mIrCount * hData->mIrSize); + double *hrirs = hData->mHrirsBase.data(); + std::vector<double> hrir(hData->mIrPoints); + uint line, col, fi, ei, ai, ti; + int count; + + printf("Loading sources..."); + fflush(stdout); + count = 0; + while(TrIsOperator(tr, "[")) + { + double factor[2]{ 1.0, 1.0 }; + + TrIndication(tr, &line, &col); + TrReadOperator(tr, "["); + + if(TrIsOperator(tr, "*")) + { + SourceRefT src; + struct MYSOFA_EASY *sofa; + uint si; + + TrReadOperator(tr, "*"); + if(!TrReadOperator(tr, "]") || !TrReadOperator(tr, "=")) + return 0; + + TrIndication(tr, &line, &col); + if(!ReadSofaRef(tr, &src)) + return 0; + + if(hData->mChannelType == CT_STEREO) + { + char type[MAX_IDENT_LEN+1]; + ChannelTypeT channelType; + + if(!TrReadIdent(tr, MAX_IDENT_LEN, type)) + return 0; + + channelType = MatchChannelType(type); + + switch(channelType) + { + case CT_NONE: + TrErrorAt(tr, line, col, "Expected a channel type.\n"); + return 0; + case CT_MONO: + src.mChannel = 0; + break; + case CT_STEREO: + src.mChannel = 1; + break; + } + } + else + { + char type[MAX_IDENT_LEN+1]; + ChannelTypeT channelType; + + if(!TrReadIdent(tr, MAX_IDENT_LEN, type)) + return 0; + + channelType = MatchChannelType(type); + if(channelType != CT_MONO) + { + TrErrorAt(tr, line, col, "Expected a mono channel type.\n"); + return 0; + } + src.mChannel = 0; + } + + sofa = LoadSofaFile(&src, hData->mIrRate, hData->mIrPoints); + if(!sofa) return 0; + + for(si = 0;si < sofa->hrtf->M;si++) + { + printf("\rLoading sources... %d of %d", si+1, sofa->hrtf->M); + fflush(stdout); + + float aer[3] = { + sofa->hrtf->SourcePosition.values[3*si], + sofa->hrtf->SourcePosition.values[3*si + 1], + sofa->hrtf->SourcePosition.values[3*si + 2] + }; + mysofa_c2s(aer); + + if(std::fabs(aer[1]) >= 89.999f) + aer[0] = 0.0f; + else + aer[0] = std::fmod(360.0f - aer[0], 360.0f); + + for(fi = 0;fi < hData->mFdCount;fi++) + { + double delta = aer[2] - hData->mFds[fi].mDistance; + if(std::abs(delta) < 0.001) + break; + } + if(fi >= hData->mFdCount) + continue; + + double ef{(90.0 + aer[1]) * (hData->mFds[fi].mEvCount - 1) / 180.0}; + ei = (int)std::round(ef); + ef = (ef - ei) * 180.0f / (hData->mFds[fi].mEvCount - 1); + if(std::abs(ef) >= 0.1) + continue; + + double af{aer[0] * hData->mFds[fi].mEvs[ei].mAzCount / 360.0f}; + ai = (int)std::round(af); + af = (af - ai) * 360.0f / hData->mFds[fi].mEvs[ei].mAzCount; + ai = ai % hData->mFds[fi].mEvs[ei].mAzCount; + if(std::abs(af) >= 0.1) + continue; + + HrirAzT *azd = &hData->mFds[fi].mEvs[ei].mAzs[ai]; + + if(azd->mIrs[0] != nullptr) + { + TrErrorAt(tr, line, col, "Redefinition of source [ %d, %d, %d ].\n", fi, ei, ai); + return 0; + } + + ExtractSofaHrir(sofa, si, 0, src.mOffset, hData->mIrPoints, hrir.data()); + azd->mIrs[0] = &hrirs[hData->mIrSize * azd->mIndex]; + if(model == HM_DATASET) + azd->mDelays[0] = AverageHrirOnset(hData->mIrRate, hData->mIrPoints, hrir.data(), 1.0, azd->mDelays[0]); + AverageHrirMagnitude(hData->mIrPoints, hData->mFftSize, hrir.data(), 1.0, azd->mIrs[0]); + + if(src.mChannel == 1) + { + ExtractSofaHrir(sofa, si, 1, src.mOffset, hData->mIrPoints, hrir.data()); + azd->mIrs[1] = &hrirs[hData->mIrSize * (hData->mIrCount + azd->mIndex)]; + if(model == HM_DATASET) + azd->mDelays[1] = AverageHrirOnset(hData->mIrRate, hData->mIrPoints, hrir.data(), 1.0, azd->mDelays[1]); + AverageHrirMagnitude(hData->mIrPoints, hData->mFftSize, hrir.data(), 1.0, azd->mIrs[1]); + } + + // TODO: Since some SOFA files contain minimum phase HRIRs, + // it would be beneficial to check for per-measurement delays + // (when available) to reconstruct the HRTDs. + } + + continue; + } + + if(!ReadIndexTriplet(tr, hData, &fi, &ei, &ai)) + return 0; + if(!TrReadOperator(tr, "]")) + return 0; + HrirAzT *azd = &hData->mFds[fi].mEvs[ei].mAzs[ai]; + + if(azd->mIrs[0] != nullptr) + { + TrErrorAt(tr, line, col, "Redefinition of source.\n"); + return 0; + } + if(!TrReadOperator(tr, "=")) + return 0; + + for(;;) + { + SourceRefT src; + uint ti = 0; + + if(!ReadSourceRef(tr, &src)) + return 0; + + // TODO: Would be nice to display 'x of y files', but that would + // require preparing the source refs first to get a total count + // before loading them. + ++count; + printf("\rLoading sources... %d file%s", count, (count==1)?"":"s"); + fflush(stdout); + + if(!LoadSource(&src, hData->mIrRate, hData->mIrPoints, hrir.data())) + return 0; + + if(hData->mChannelType == CT_STEREO) + { + char ident[MAX_IDENT_LEN+1]; + + if(!TrReadIdent(tr, MAX_IDENT_LEN, ident)) + return 0; + ti = MatchTargetEar(ident); + if(static_cast<int>(ti) < 0) + { + TrErrorAt(tr, line, col, "Expected a target ear.\n"); + return 0; + } + } + azd->mIrs[ti] = &hrirs[hData->mIrSize * (ti * hData->mIrCount + azd->mIndex)]; + if(model == HM_DATASET) + azd->mDelays[ti] = AverageHrirOnset(hData->mIrRate, hData->mIrPoints, hrir.data(), 1.0 / factor[ti], azd->mDelays[ti]); + AverageHrirMagnitude(hData->mIrPoints, hData->mFftSize, hrir.data(), 1.0 / factor[ti], azd->mIrs[ti]); + factor[ti] += 1.0; + if(!TrIsOperator(tr, "+")) + break; + TrReadOperator(tr, "+"); + } + if(hData->mChannelType == CT_STEREO) + { + if(azd->mIrs[0] == nullptr) + { + TrErrorAt(tr, line, col, "Missing left ear source reference(s).\n"); + return 0; + } + else if(azd->mIrs[1] == nullptr) + { + TrErrorAt(tr, line, col, "Missing right ear source reference(s).\n"); + return 0; + } + } + } + printf("\n"); + for(fi = 0;fi < hData->mFdCount;fi++) + { + for(ei = 0;ei < hData->mFds[fi].mEvCount;ei++) + { + for(ai = 0;ai < hData->mFds[fi].mEvs[ei].mAzCount;ai++) + { + HrirAzT *azd = &hData->mFds[fi].mEvs[ei].mAzs[ai]; + if(azd->mIrs[0] != nullptr) + break; + } + if(ai < hData->mFds[fi].mEvs[ei].mAzCount) + break; + } + if(ei >= hData->mFds[fi].mEvCount) + { + TrError(tr, "Missing source references [ %d, *, * ].\n", fi); + return 0; + } + hData->mFds[fi].mEvStart = ei; + for(;ei < hData->mFds[fi].mEvCount;ei++) + { + for(ai = 0;ai < hData->mFds[fi].mEvs[ei].mAzCount;ai++) + { + HrirAzT *azd = &hData->mFds[fi].mEvs[ei].mAzs[ai]; + + if(azd->mIrs[0] == nullptr) + { + TrError(tr, "Missing source reference [ %d, %d, %d ].\n", fi, ei, ai); + return 0; + } + } + } + } + for(ti = 0;ti < channels;ti++) + { + for(fi = 0;fi < hData->mFdCount;fi++) + { + for(ei = 0;ei < hData->mFds[fi].mEvCount;ei++) + { + for(ai = 0;ai < hData->mFds[fi].mEvs[ei].mAzCount;ai++) + { + HrirAzT *azd = &hData->mFds[fi].mEvs[ei].mAzs[ai]; + + azd->mIrs[ti] = &hrirs[hData->mIrSize * (ti * hData->mIrCount + azd->mIndex)]; + } + } + } + } + if(!TrLoad(tr)) + { + mysofa_cache_release_all(); + return 1; + } + + TrError(tr, "Errant data at end of source list.\n"); + mysofa_cache_release_all(); + return 0; +} diff --git a/utils/makemhr/loaddef.h b/utils/makemhr/loaddef.h index e69de29b..5e8c5313 100644 --- a/utils/makemhr/loaddef.h +++ b/utils/makemhr/loaddef.h @@ -0,0 +1,30 @@ +#ifndef LOADDEF_H +#define LOADDEF_H + +#include <stdio.h> + +#include "makemhr.h" + + +// Constants for accessing the token reader's ring buffer. +#define TR_RING_BITS (16) +#define TR_RING_SIZE (1 << TR_RING_BITS) +#define TR_RING_MASK (TR_RING_SIZE - 1) + + +// Token reader state for parsing the data set definition. +struct TokenReaderT { + FILE *mFile; + const char *mName; + uint mLine; + uint mColumn; + char mRing[TR_RING_SIZE]; + size_t mIn; + size_t mOut; +}; + +void TrSetup(FILE *fp, const char *filename, TokenReaderT *tr); +int ProcessMetrics(TokenReaderT *tr, const uint fftSize, const uint truncSize, HrirDataT *hData); +int ProcessSources(const HeadModelT model, TokenReaderT *tr, HrirDataT *hData); + +#endif /* LOADDEF_H */ diff --git a/utils/makemhr/loadsofa.cpp b/utils/makemhr/loadsofa.cpp index e69de29b..a60526a8 100644 --- a/utils/makemhr/loadsofa.cpp +++ b/utils/makemhr/loadsofa.cpp @@ -0,0 +1,4 @@ + +#include "mysofa.h" + +#include "loadsofa.h" diff --git a/utils/makemhr/loadsofa.h b/utils/makemhr/loadsofa.h index e69de29b..9a73bd72 100644 --- a/utils/makemhr/loadsofa.h +++ b/utils/makemhr/loadsofa.h @@ -0,0 +1,6 @@ +#ifndef LOADSOFA_H +#define LOADSOFA_H + +#include "makemhr.h" + +#endif /* LOADSOFA_H */ diff --git a/utils/makemhr/makemhr.cpp b/utils/makemhr/makemhr.cpp index 27b1d69d..0c2c398c 100644 --- a/utils/makemhr/makemhr.cpp +++ b/utils/makemhr/makemhr.cpp @@ -92,6 +92,10 @@ #include "mysofa.h" +#include "makemhr.h" +#include "loaddef.h" +#include "loadsofa.h" + #include "win_main_utf8.h" namespace { @@ -108,68 +112,6 @@ using namespace std::placeholders; // The epsilon used to maintain signal stability. #define EPSILON (1e-9) -// Constants for accessing the token reader's ring buffer. -#define TR_RING_BITS (16) -#define TR_RING_SIZE (1 << TR_RING_BITS) -#define TR_RING_MASK (TR_RING_SIZE - 1) - -// The token reader's load interval in bytes. -#define TR_LOAD_SIZE (TR_RING_SIZE >> 2) - -// The maximum identifier length used when processing the data set -// definition. -#define MAX_IDENT_LEN (16) - -// The maximum path length used when processing filenames. -#define MAX_PATH_LEN (256) - -// The limits for the sample 'rate' metric in the data set definition and for -// resampling. -#define MIN_RATE (32000) -#define MAX_RATE (96000) - -// The limits for the HRIR 'points' metric in the data set definition. -#define MIN_POINTS (16) -#define MAX_POINTS (8192) - -// The limit to the number of 'distances' listed in the data set definition. -#define MAX_FD_COUNT (16) - -// The limits to the number of 'azimuths' listed in the data set definition. -#define MIN_EV_COUNT (5) -#define MAX_EV_COUNT (128) - -// The limits for each of the 'azimuths' listed in the data set definition. -#define MIN_AZ_COUNT (1) -#define MAX_AZ_COUNT (128) - -// The limits for the listener's head 'radius' in the data set definition. -#define MIN_RADIUS (0.05) -#define MAX_RADIUS (0.15) - -// The limits for the 'distance' from source to listener for each field in -// the definition file. -#define MIN_DISTANCE (0.05) -#define MAX_DISTANCE (2.50) - -// The maximum number of channels that can be addressed for a WAVE file -// source listed in the data set definition. -#define MAX_WAVE_CHANNELS (65535) - -// The limits to the byte size for a binary source listed in the definition -// file. -#define MIN_BIN_SIZE (2) -#define MAX_BIN_SIZE (4) - -// The minimum number of significant bits for binary sources listed in the -// data set definition. The maximum is calculated from the byte size. -#define MIN_BIN_BITS (16) - -// The limits to the number of significant bits for an ASCII source listed in -// the data set definition. -#define MIN_ASCII_BITS (16) -#define MAX_ASCII_BITS (32) - // The limits to the FFT window size override on the command line. #define MIN_FFTSIZE (65536) #define MAX_FFTSIZE (131072) @@ -199,21 +141,6 @@ using namespace std::placeholders; #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' -#define FOURCC_RIFX (0x58464952) // 'RIFX' -#define FOURCC_WAVE (0x45564157) // 'WAVE' -#define FOURCC_FMT (0x20746D66) // 'fmt ' -#define FOURCC_DATA (0x61746164) // 'data' -#define FOURCC_LIST (0x5453494C) // 'LIST' -#define FOURCC_WAVL (0x6C766177) // 'wavl' -#define FOURCC_SLNT (0x746E6C73) // 'slnt' - -// The supported wave formats. -#define WAVE_FORMAT_PCM (0x0001) -#define WAVE_FORMAT_IEEE_FLOAT (0x0003) -#define WAVE_FORMAT_EXTENSIBLE (0xFFFE) - // The maximum propagation delay value supported by OpenAL Soft. #define MAX_HRTD (63.0) @@ -221,19 +148,6 @@ using namespace std::placeholders; // response protocol 02. #define MHR_FORMAT ("MinPHR02") -// Sample and channel type enum values. -enum SampleTypeT { - ST_S16 = 0, - ST_S24 = 1 -}; - -// Certain iterations rely on these integer enum values. -enum ChannelTypeT { - CT_NONE = -1, - CT_MONO = 0, - CT_STEREO = 1 -}; - // Byte order for the serialization routines. enum ByteOrderT { BO_NONE, @@ -241,36 +155,6 @@ enum ByteOrderT { BO_BIG }; -// Source format for the references listed in the data set definition. -enum SourceFormatT { - SF_NONE, - SF_ASCII, // ASCII text file. - SF_BIN_LE, // Little-endian binary file. - SF_BIN_BE, // Big-endian binary file. - SF_WAVE, // RIFF/RIFX WAVE file. - SF_SOFA // Spatially Oriented Format for Accoustics (SOFA) file. -}; - -// Element types for the references listed in the data set definition. -enum ElementTypeT { - ET_NONE, - ET_INT, // Integer elements. - ET_FP // Floating-point elements. -}; - -// Head model used for calculating the impulse delays. -enum HeadModelT { - HM_NONE, - HM_DATASET, // Measure the onset from the dataset. - HM_SPHERE // Calculate the onset using a spherical head model. -}; - -/* Unsigned integer type. */ -using uint = unsigned int; - -/* Complex double type. */ -using complex_d = std::complex<double>; - /* Channel index enums. Mono uses LeftChannel only. */ enum ChannelIndex : uint { LeftChannel = 0u, @@ -278,538 +162,6 @@ enum ChannelIndex : uint { }; -// Token reader state for parsing the data set definition. -struct TokenReaderT { - FILE *mFile; - const char *mName; - uint mLine; - uint mColumn; - char mRing[TR_RING_SIZE]; - size_t mIn; - size_t mOut; -}; - -// Source reference state used when loading sources. -struct SourceRefT { - SourceFormatT mFormat; - ElementTypeT mType; - uint mSize; - int mBits; - uint mChannel; - double mAzimuth; - double mElevation; - double mRadius; - uint mSkip; - uint mOffset; - char mPath[MAX_PATH_LEN+1]; -}; - -// Structured HRIR storage for stereo azimuth pairs, elevations, and fields. -struct HrirAzT { - double mAzimuth{0.0}; - uint mIndex{0u}; - double mDelays[2]{0.0, 0.0}; - double *mIrs[2]{nullptr, nullptr}; -}; - -struct HrirEvT { - double mElevation{0.0}; - uint mIrCount{0u}; - uint mAzCount{0u}; - HrirAzT *mAzs{nullptr}; -}; - -struct HrirFdT { - double mDistance{0.0}; - uint mIrCount{0u}; - uint mEvCount{0u}; - uint mEvStart{0u}; - HrirEvT *mEvs{nullptr}; -}; - -// The HRIR metrics and data set used when loading, processing, and storing -// the resulting HRTF. -struct HrirDataT { - uint mIrRate{0u}; - SampleTypeT mSampleType{ST_S24}; - ChannelTypeT mChannelType{CT_NONE}; - uint mIrPoints{0u}; - uint mFftSize{0u}; - uint mIrSize{0u}; - double mRadius{0.0}; - uint mIrCount{0u}; - uint mFdCount{0u}; - - std::vector<double> mHrirsBase; - std::vector<HrirEvT> mEvsBase; - std::vector<HrirAzT> mAzsBase; - - std::vector<HrirFdT> mFds; -}; - -// The resampler metrics and FIR filter. -struct ResamplerT { - uint mP, mQ, mM, mL; - std::vector<double> mF; -}; - - -/***************************** - *** Token reader routines *** - *****************************/ - -/* Whitespace is not significant. It can process tokens as identifiers, numbers - * (integer and floating-point), strings, and operators. Strings must be - * encapsulated by double-quotes and cannot span multiple lines. - */ - -// Setup the reader on the given file. The filename can be NULL if no error -// output is desired. -static void TrSetup(FILE *fp, const char *filename, TokenReaderT *tr) -{ - const char *name = nullptr; - - if(filename) - { - const char *slash = strrchr(filename, '/'); - if(slash) - { - const char *bslash = strrchr(slash+1, '\\'); - if(bslash) name = bslash+1; - else name = slash+1; - } - else - { - const char *bslash = strrchr(filename, '\\'); - if(bslash) name = bslash+1; - else name = filename; - } - } - - tr->mFile = fp; - tr->mName = name; - tr->mLine = 1; - tr->mColumn = 1; - tr->mIn = 0; - tr->mOut = 0; -} - -// Prime the reader's ring buffer, and return a result indicating that there -// is text to process. -static int TrLoad(TokenReaderT *tr) -{ - size_t toLoad, in, count; - - toLoad = TR_RING_SIZE - (tr->mIn - tr->mOut); - if(toLoad >= TR_LOAD_SIZE && !feof(tr->mFile)) - { - // Load TR_LOAD_SIZE (or less if at the end of the file) per read. - toLoad = TR_LOAD_SIZE; - in = tr->mIn&TR_RING_MASK; - count = TR_RING_SIZE - in; - if(count < toLoad) - { - tr->mIn += fread(&tr->mRing[in], 1, count, tr->mFile); - tr->mIn += fread(&tr->mRing[0], 1, toLoad-count, tr->mFile); - } - else - tr->mIn += fread(&tr->mRing[in], 1, toLoad, tr->mFile); - - if(tr->mOut >= TR_RING_SIZE) - { - tr->mOut -= TR_RING_SIZE; - tr->mIn -= TR_RING_SIZE; - } - } - if(tr->mIn > tr->mOut) - return 1; - return 0; -} - -// Error display routine. Only displays when the base name is not NULL. -static void TrErrorVA(const TokenReaderT *tr, uint line, uint column, const char *format, va_list argPtr) -{ - if(!tr->mName) - return; - fprintf(stderr, "\nError (%s:%u:%u): ", tr->mName, line, column); - vfprintf(stderr, format, argPtr); -} - -// Used to display an error at a saved line/column. -static void TrErrorAt(const TokenReaderT *tr, uint line, uint column, const char *format, ...) -{ - va_list argPtr; - - va_start(argPtr, format); - TrErrorVA(tr, line, column, format, argPtr); - va_end(argPtr); -} - -// Used to display an error at the current line/column. -static void TrError(const TokenReaderT *tr, const char *format, ...) -{ - va_list argPtr; - - va_start(argPtr, format); - TrErrorVA(tr, tr->mLine, tr->mColumn, format, argPtr); - va_end(argPtr); -} - -// Skips to the next line. -static void TrSkipLine(TokenReaderT *tr) -{ - char ch; - - while(TrLoad(tr)) - { - ch = tr->mRing[tr->mOut&TR_RING_MASK]; - tr->mOut++; - if(ch == '\n') - { - tr->mLine++; - tr->mColumn = 1; - break; - } - tr->mColumn ++; - } -} - -// Skips to the next token. -static int TrSkipWhitespace(TokenReaderT *tr) -{ - while(TrLoad(tr)) - { - char ch{tr->mRing[tr->mOut&TR_RING_MASK]}; - if(isspace(ch)) - { - tr->mOut++; - if(ch == '\n') - { - tr->mLine++; - tr->mColumn = 1; - } - else - tr->mColumn++; - } - else if(ch == '#') - TrSkipLine(tr); - else - return 1; - } - return 0; -} - -// Get the line and/or column of the next token (or the end of input). -static void TrIndication(TokenReaderT *tr, uint *line, uint *column) -{ - TrSkipWhitespace(tr); - if(line) *line = tr->mLine; - if(column) *column = tr->mColumn; -} - -// Checks to see if a token is (likely to be) an identifier. It does not -// display any errors and will not proceed to the next token. -static int TrIsIdent(TokenReaderT *tr) -{ - if(!TrSkipWhitespace(tr)) - return 0; - char ch{tr->mRing[tr->mOut&TR_RING_MASK]}; - return ch == '_' || isalpha(ch); -} - - -// Checks to see if a token is the given operator. It does not display any -// errors and will not proceed to the next token. -static int TrIsOperator(TokenReaderT *tr, const char *op) -{ - size_t out, len; - char ch; - - if(!TrSkipWhitespace(tr)) - return 0; - out = tr->mOut; - len = 0; - while(op[len] != '\0' && out < tr->mIn) - { - ch = tr->mRing[out&TR_RING_MASK]; - if(ch != op[len]) break; - len++; - out++; - } - if(op[len] == '\0') - return 1; - return 0; -} - -/* The TrRead*() routines obtain the value of a matching token type. They - * display type, form, and boundary errors and will proceed to the next - * token. - */ - -// Reads and validates an identifier token. -static int TrReadIdent(TokenReaderT *tr, const uint maxLen, char *ident) -{ - uint col, len; - char ch; - - col = tr->mColumn; - if(TrSkipWhitespace(tr)) - { - col = tr->mColumn; - ch = tr->mRing[tr->mOut&TR_RING_MASK]; - if(ch == '_' || isalpha(ch)) - { - len = 0; - do { - if(len < maxLen) - ident[len] = ch; - len++; - tr->mOut++; - if(!TrLoad(tr)) - break; - ch = tr->mRing[tr->mOut&TR_RING_MASK]; - } while(ch == '_' || isdigit(ch) || isalpha(ch)); - - tr->mColumn += len; - if(len < maxLen) - { - ident[len] = '\0'; - return 1; - } - TrErrorAt(tr, tr->mLine, col, "Identifier is too long.\n"); - return 0; - } - } - TrErrorAt(tr, tr->mLine, col, "Expected an identifier.\n"); - return 0; -} - -// Reads and validates (including bounds) an integer token. -static int TrReadInt(TokenReaderT *tr, const int loBound, const int hiBound, int *value) -{ - uint col, digis, len; - char ch, temp[64+1]; - - col = tr->mColumn; - if(TrSkipWhitespace(tr)) - { - col = tr->mColumn; - len = 0; - ch = tr->mRing[tr->mOut&TR_RING_MASK]; - if(ch == '+' || ch == '-') - { - temp[len] = ch; - len++; - tr->mOut++; - } - digis = 0; - while(TrLoad(tr)) - { - ch = tr->mRing[tr->mOut&TR_RING_MASK]; - if(!isdigit(ch)) break; - if(len < 64) - temp[len] = ch; - len++; - digis++; - tr->mOut++; - } - tr->mColumn += len; - if(digis > 0 && ch != '.' && !isalpha(ch)) - { - if(len > 64) - { - TrErrorAt(tr, tr->mLine, col, "Integer is too long."); - return 0; - } - temp[len] = '\0'; - *value = strtol(temp, nullptr, 10); - if(*value < loBound || *value > hiBound) - { - TrErrorAt(tr, tr->mLine, col, "Expected a value from %d to %d.\n", loBound, hiBound); - return 0; - } - return 1; - } - } - TrErrorAt(tr, tr->mLine, col, "Expected an integer.\n"); - return 0; -} - -// Reads and validates (including bounds) a float token. -static int TrReadFloat(TokenReaderT *tr, const double loBound, const double hiBound, double *value) -{ - uint col, digis, len; - char ch, temp[64+1]; - - col = tr->mColumn; - if(TrSkipWhitespace(tr)) - { - col = tr->mColumn; - len = 0; - ch = tr->mRing[tr->mOut&TR_RING_MASK]; - if(ch == '+' || ch == '-') - { - temp[len] = ch; - len++; - tr->mOut++; - } - - digis = 0; - while(TrLoad(tr)) - { - ch = tr->mRing[tr->mOut&TR_RING_MASK]; - if(!isdigit(ch)) break; - if(len < 64) - temp[len] = ch; - len++; - digis++; - tr->mOut++; - } - if(ch == '.') - { - if(len < 64) - temp[len] = ch; - len++; - tr->mOut++; - } - while(TrLoad(tr)) - { - ch = tr->mRing[tr->mOut&TR_RING_MASK]; - if(!isdigit(ch)) break; - if(len < 64) - temp[len] = ch; - len++; - digis++; - tr->mOut++; - } - if(digis > 0) - { - if(ch == 'E' || ch == 'e') - { - if(len < 64) - temp[len] = ch; - len++; - digis = 0; - tr->mOut++; - if(ch == '+' || ch == '-') - { - if(len < 64) - temp[len] = ch; - len++; - tr->mOut++; - } - while(TrLoad(tr)) - { - ch = tr->mRing[tr->mOut&TR_RING_MASK]; - if(!isdigit(ch)) break; - if(len < 64) - temp[len] = ch; - len++; - digis++; - tr->mOut++; - } - } - tr->mColumn += len; - if(digis > 0 && ch != '.' && !isalpha(ch)) - { - if(len > 64) - { - TrErrorAt(tr, tr->mLine, col, "Float is too long."); - return 0; - } - temp[len] = '\0'; - *value = strtod(temp, nullptr); - if(*value < loBound || *value > hiBound) - { - TrErrorAt(tr, tr->mLine, col, "Expected a value from %f to %f.\n", loBound, hiBound); - return 0; - } - return 1; - } - } - else - tr->mColumn += len; - } - TrErrorAt(tr, tr->mLine, col, "Expected a float.\n"); - return 0; -} - -// Reads and validates a string token. -static int TrReadString(TokenReaderT *tr, const uint maxLen, char *text) -{ - uint col, len; - char ch; - - col = tr->mColumn; - if(TrSkipWhitespace(tr)) - { - col = tr->mColumn; - ch = tr->mRing[tr->mOut&TR_RING_MASK]; - if(ch == '\"') - { - tr->mOut++; - len = 0; - while(TrLoad(tr)) - { - ch = tr->mRing[tr->mOut&TR_RING_MASK]; - tr->mOut++; - if(ch == '\"') - break; - if(ch == '\n') - { - TrErrorAt(tr, tr->mLine, col, "Unterminated string at end of line.\n"); - return 0; - } - if(len < maxLen) - text[len] = ch; - len++; - } - if(ch != '\"') - { - tr->mColumn += 1 + len; - TrErrorAt(tr, tr->mLine, col, "Unterminated string at end of input.\n"); - return 0; - } - tr->mColumn += 2 + len; - if(len > maxLen) - { - TrErrorAt(tr, tr->mLine, col, "String is too long.\n"); - return 0; - } - text[len] = '\0'; - return 1; - } - } - TrErrorAt(tr, tr->mLine, col, "Expected a string.\n"); - return 0; -} - -// Reads and validates the given operator. -static int TrReadOperator(TokenReaderT *tr, const char *op) -{ - uint col, len; - char ch; - - col = tr->mColumn; - if(TrSkipWhitespace(tr)) - { - col = tr->mColumn; - len = 0; - while(op[len] != '\0' && TrLoad(tr)) - { - ch = tr->mRing[tr->mOut&TR_RING_MASK]; - if(ch != op[len]) break; - len++; - tr->mOut++; - } - tr->mColumn += len; - if(op[len] == '\0') - return 1; - } - TrErrorAt(tr, tr->mLine, col, "Expected '%s' operator.\n", op); - return 0; -} - /* Performs a string substitution. Any case-insensitive occurrences of the * pattern string are replaced with the replacement string. The result is * truncated if necessary. @@ -863,12 +215,6 @@ static double Clamp(const double val, const double lower, const double upper) return std::min(std::max(val, lower), upper); } -// Performs linear interpolation. -static double Lerp(const double a, const double b, const double f) -{ - return a + f * (b - a); -} - static inline uint dither_rng(uint *seed) { *seed = *seed * 96314165 + 907633515; @@ -940,14 +286,14 @@ static void FftSummation(const int n, const double s, complex_d *cplx) } // Performs a forward FFT. -static void FftForward(const uint n, complex_d *inout) +void FftForward(const uint n, complex_d *inout) { FftArrange(n, inout); FftSummation(n, 1.0, inout); } // Performs an inverse FFT. -static void FftInverse(const uint n, complex_d *inout) +void FftInverse(const uint n, complex_d *inout) { FftArrange(n, inout); FftSummation(n, -1.0, inout); @@ -984,7 +330,7 @@ static void Hilbert(const uint n, complex_d *inout) * minimum phase reconstruction. The mirrored half of the response is also * discarded. */ -static void MagnitudeResponse(const uint n, const complex_d *in, double *out) +void MagnitudeResponse(const uint n, const complex_d *in, double *out) { const uint m = 1 + (n / 2); uint i; @@ -1202,7 +548,7 @@ static double SincFilter(const int l, const double b, const double gain, const d // Calculate the resampling metrics and build the Kaiser-windowed sinc filter // that's used to cut frequencies above the destination nyquist. -static void ResamplerSetup(ResamplerT *rs, const uint srcRate, const uint dstRate) +void ResamplerSetup(ResamplerT *rs, const uint srcRate, const uint dstRate) { double cutoff, width, beta; uint gcd, l; @@ -1238,7 +584,7 @@ static void ResamplerSetup(ResamplerT *rs, const uint srcRate, const uint dstRat // Perform the upsample-filter-downsample resampling operation using a // polyphase filter implementation. -static void ResamplerRun(ResamplerT *rs, const uint inN, const double *in, const uint outN, double *out) +void ResamplerRun(ResamplerT *rs, const uint inN, const double *in, const uint outN, double *out) { const uint p = rs->mP, q = rs->mQ, m = rs->mM, l = rs->mL; std::vector<double> workspace; @@ -1285,606 +631,6 @@ static void ResamplerRun(ResamplerT *rs, const uint inN, const double *in, const } } -/************************* - *** File source input *** - *************************/ - -// Read a binary value of the specified byte order and byte size from a file, -// storing it as a 32-bit unsigned integer. -static int ReadBin4(FILE *fp, const char *filename, const ByteOrderT order, const uint bytes, uint32_t *out) -{ - uint8_t in[4]; - uint32_t accum; - uint i; - - if(fread(in, 1, bytes, fp) != bytes) - { - fprintf(stderr, "\nError: Bad read from file '%s'.\n", filename); - return 0; - } - accum = 0; - switch(order) - { - case BO_LITTLE: - for(i = 0;i < bytes;i++) - accum = (accum<<8) | in[bytes - i - 1]; - break; - case BO_BIG: - for(i = 0;i < bytes;i++) - accum = (accum<<8) | in[i]; - break; - default: - break; - } - *out = accum; - return 1; -} - -// Read a binary value of the specified byte order from a file, storing it as -// a 64-bit unsigned integer. -static int ReadBin8(FILE *fp, const char *filename, const ByteOrderT order, uint64_t *out) -{ - uint8_t in[8]; - uint64_t accum; - uint i; - - if(fread(in, 1, 8, fp) != 8) - { - fprintf(stderr, "\nError: Bad read from file '%s'.\n", filename); - return 0; - } - accum = 0ULL; - switch(order) - { - case BO_LITTLE: - for(i = 0;i < 8;i++) - accum = (accum<<8) | in[8 - i - 1]; - break; - case BO_BIG: - for(i = 0;i < 8;i++) - accum = (accum<<8) | in[i]; - break; - default: - break; - } - *out = accum; - return 1; -} - -/* Read a binary value of the specified type, byte order, and byte size from - * a file, converting it to a double. For integer types, the significant - * bits are used to normalize the result. The sign of bits determines - * whether they are padded toward the MSB (negative) or LSB (positive). - * Floating-point types are not normalized. - */ -static int ReadBinAsDouble(FILE *fp, const char *filename, const ByteOrderT order, const ElementTypeT type, const uint bytes, const int bits, double *out) -{ - union { - uint32_t ui; - int32_t i; - float f; - } v4; - union { - uint64_t ui; - double f; - } v8; - - *out = 0.0; - if(bytes > 4) - { - if(!ReadBin8(fp, filename, order, &v8.ui)) - return 0; - if(type == ET_FP) - *out = v8.f; - } - else - { - if(!ReadBin4(fp, filename, order, bytes, &v4.ui)) - return 0; - if(type == ET_FP) - *out = v4.f; - else - { - if(bits > 0) - v4.ui >>= (8*bytes) - (static_cast<uint>(bits)); - else - v4.ui &= (0xFFFFFFFF >> (32+bits)); - - if(v4.ui&static_cast<uint>(1<<(std::abs(bits)-1))) - v4.ui |= (0xFFFFFFFF << std::abs(bits)); - *out = v4.i / static_cast<double>(1<<(std::abs(bits)-1)); - } - } - return 1; -} - -/* Read an ascii value of the specified type from a file, converting it to a - * double. For integer types, the significant bits are used to normalize the - * result. The sign of the bits should always be positive. This also skips - * up to one separator character before the element itself. - */ -static int ReadAsciiAsDouble(TokenReaderT *tr, const char *filename, const ElementTypeT type, const uint bits, double *out) -{ - if(TrIsOperator(tr, ",")) - TrReadOperator(tr, ","); - else if(TrIsOperator(tr, ":")) - TrReadOperator(tr, ":"); - else if(TrIsOperator(tr, ";")) - TrReadOperator(tr, ";"); - else if(TrIsOperator(tr, "|")) - TrReadOperator(tr, "|"); - - if(type == ET_FP) - { - if(!TrReadFloat(tr, -std::numeric_limits<double>::infinity(), - std::numeric_limits<double>::infinity(), out)) - { - fprintf(stderr, "\nError: Bad read from file '%s'.\n", filename); - return 0; - } - } - else - { - int v; - if(!TrReadInt(tr, -(1<<(bits-1)), (1<<(bits-1))-1, &v)) - { - fprintf(stderr, "\nError: Bad read from file '%s'.\n", filename); - return 0; - } - *out = v / static_cast<double>((1<<(bits-1))-1); - } - return 1; -} - -// Read the RIFF/RIFX WAVE format chunk from a file, validating it against -// the source parameters and data set metrics. -static int ReadWaveFormat(FILE *fp, const ByteOrderT order, const uint hrirRate, SourceRefT *src) -{ - uint32_t fourCC, chunkSize; - uint32_t format, channels, rate, dummy, block, size, bits; - - chunkSize = 0; - do { - if(chunkSize > 0) - fseek(fp, static_cast<long>(chunkSize), SEEK_CUR); - if(!ReadBin4(fp, src->mPath, BO_LITTLE, 4, &fourCC) || - !ReadBin4(fp, src->mPath, order, 4, &chunkSize)) - return 0; - } while(fourCC != FOURCC_FMT); - if(!ReadBin4(fp, src->mPath, order, 2, &format) || - !ReadBin4(fp, src->mPath, order, 2, &channels) || - !ReadBin4(fp, src->mPath, order, 4, &rate) || - !ReadBin4(fp, src->mPath, order, 4, &dummy) || - !ReadBin4(fp, src->mPath, order, 2, &block)) - return 0; - block /= channels; - if(chunkSize > 14) - { - if(!ReadBin4(fp, src->mPath, order, 2, &size)) - return 0; - size /= 8; - if(block > size) - size = block; - } - else - size = block; - if(format == WAVE_FORMAT_EXTENSIBLE) - { - fseek(fp, 2, SEEK_CUR); - if(!ReadBin4(fp, src->mPath, order, 2, &bits)) - return 0; - if(bits == 0) - bits = 8 * size; - fseek(fp, 4, SEEK_CUR); - if(!ReadBin4(fp, src->mPath, order, 2, &format)) - return 0; - fseek(fp, static_cast<long>(chunkSize - 26), SEEK_CUR); - } - else - { - bits = 8 * size; - if(chunkSize > 14) - fseek(fp, static_cast<long>(chunkSize - 16), SEEK_CUR); - else - fseek(fp, static_cast<long>(chunkSize - 14), SEEK_CUR); - } - if(format != WAVE_FORMAT_PCM && format != WAVE_FORMAT_IEEE_FLOAT) - { - fprintf(stderr, "\nError: Unsupported WAVE format in file '%s'.\n", src->mPath); - return 0; - } - if(src->mChannel >= channels) - { - fprintf(stderr, "\nError: Missing source channel in WAVE file '%s'.\n", src->mPath); - return 0; - } - if(rate != hrirRate) - { - fprintf(stderr, "\nError: Mismatched source sample rate in WAVE file '%s'.\n", src->mPath); - return 0; - } - if(format == WAVE_FORMAT_PCM) - { - if(size < 2 || size > 4) - { - fprintf(stderr, "\nError: Unsupported sample size in WAVE file '%s'.\n", src->mPath); - return 0; - } - if(bits < 16 || bits > (8*size)) - { - fprintf(stderr, "\nError: Bad significant bits in WAVE file '%s'.\n", src->mPath); - return 0; - } - src->mType = ET_INT; - } - else - { - if(size != 4 && size != 8) - { - fprintf(stderr, "\nError: Unsupported sample size in WAVE file '%s'.\n", src->mPath); - return 0; - } - src->mType = ET_FP; - } - src->mSize = size; - src->mBits = static_cast<int>(bits); - src->mSkip = channels; - return 1; -} - -// Read a RIFF/RIFX WAVE data chunk, converting all elements to doubles. -static int ReadWaveData(FILE *fp, const SourceRefT *src, const ByteOrderT order, const uint n, double *hrir) -{ - int pre, post, skip; - uint i; - - pre = static_cast<int>(src->mSize * src->mChannel); - post = static_cast<int>(src->mSize * (src->mSkip - src->mChannel - 1)); - skip = 0; - for(i = 0;i < n;i++) - { - skip += pre; - if(skip > 0) - fseek(fp, skip, SEEK_CUR); - if(!ReadBinAsDouble(fp, src->mPath, order, src->mType, src->mSize, src->mBits, &hrir[i])) - return 0; - skip = post; - } - if(skip > 0) - fseek(fp, skip, SEEK_CUR); - return 1; -} - -// Read the RIFF/RIFX WAVE list or data chunk, converting all elements to -// doubles. -static int ReadWaveList(FILE *fp, const SourceRefT *src, const ByteOrderT order, const uint n, double *hrir) -{ - uint32_t fourCC, chunkSize, listSize, count; - uint block, skip, offset, i; - double lastSample; - - for(;;) - { - if(!ReadBin4(fp, src->mPath, BO_LITTLE, 4, &fourCC) || - !ReadBin4(fp, src->mPath, order, 4, &chunkSize)) - return 0; - - if(fourCC == FOURCC_DATA) - { - block = src->mSize * src->mSkip; - count = chunkSize / block; - if(count < (src->mOffset + n)) - { - fprintf(stderr, "\nError: Bad read from file '%s'.\n", src->mPath); - return 0; - } - fseek(fp, static_cast<long>(src->mOffset * block), SEEK_CUR); - if(!ReadWaveData(fp, src, order, n, &hrir[0])) - return 0; - return 1; - } - else if(fourCC == FOURCC_LIST) - { - if(!ReadBin4(fp, src->mPath, BO_LITTLE, 4, &fourCC)) - return 0; - chunkSize -= 4; - if(fourCC == FOURCC_WAVL) - break; - } - if(chunkSize > 0) - fseek(fp, static_cast<long>(chunkSize), SEEK_CUR); - } - listSize = chunkSize; - block = src->mSize * src->mSkip; - skip = src->mOffset; - offset = 0; - lastSample = 0.0; - while(offset < n && listSize > 8) - { - if(!ReadBin4(fp, src->mPath, BO_LITTLE, 4, &fourCC) || - !ReadBin4(fp, src->mPath, order, 4, &chunkSize)) - return 0; - listSize -= 8 + chunkSize; - if(fourCC == FOURCC_DATA) - { - count = chunkSize / block; - if(count > skip) - { - fseek(fp, static_cast<long>(skip * block), SEEK_CUR); - chunkSize -= skip * block; - count -= skip; - skip = 0; - if(count > (n - offset)) - count = n - offset; - if(!ReadWaveData(fp, src, order, count, &hrir[offset])) - return 0; - chunkSize -= count * block; - offset += count; - lastSample = hrir[offset - 1]; - } - else - { - skip -= count; - count = 0; - } - } - else if(fourCC == FOURCC_SLNT) - { - if(!ReadBin4(fp, src->mPath, order, 4, &count)) - return 0; - chunkSize -= 4; - if(count > skip) - { - count -= skip; - skip = 0; - if(count > (n - offset)) - count = n - offset; - for(i = 0; i < count; i ++) - hrir[offset + i] = lastSample; - offset += count; - } - else - { - skip -= count; - count = 0; - } - } - if(chunkSize > 0) - fseek(fp, static_cast<long>(chunkSize), SEEK_CUR); - } - if(offset < n) - { - fprintf(stderr, "\nError: Bad read from file '%s'.\n", src->mPath); - return 0; - } - return 1; -} - -// Load a source HRIR from an ASCII text file containing a list of elements -// separated by whitespace or common list operators (',', ';', ':', '|'). -static int LoadAsciiSource(FILE *fp, const SourceRefT *src, const uint n, double *hrir) -{ - TokenReaderT tr; - uint i, j; - double dummy; - - TrSetup(fp, nullptr, &tr); - for(i = 0;i < src->mOffset;i++) - { - if(!ReadAsciiAsDouble(&tr, src->mPath, src->mType, static_cast<uint>(src->mBits), &dummy)) - return 0; - } - for(i = 0;i < n;i++) - { - if(!ReadAsciiAsDouble(&tr, src->mPath, src->mType, static_cast<uint>(src->mBits), &hrir[i])) - return 0; - for(j = 0;j < src->mSkip;j++) - { - if(!ReadAsciiAsDouble(&tr, src->mPath, src->mType, static_cast<uint>(src->mBits), &dummy)) - return 0; - } - } - return 1; -} - -// Load a source HRIR from a binary file. -static int LoadBinarySource(FILE *fp, const SourceRefT *src, const ByteOrderT order, const uint n, double *hrir) -{ - uint i; - - fseek(fp, static_cast<long>(src->mOffset), SEEK_SET); - for(i = 0;i < n;i++) - { - if(!ReadBinAsDouble(fp, src->mPath, order, src->mType, src->mSize, src->mBits, &hrir[i])) - return 0; - if(src->mSkip > 0) - fseek(fp, static_cast<long>(src->mSkip), SEEK_CUR); - } - return 1; -} - -// Load a source HRIR from a RIFF/RIFX WAVE file. -static int LoadWaveSource(FILE *fp, SourceRefT *src, const uint hrirRate, const uint n, double *hrir) -{ - uint32_t fourCC, dummy; - ByteOrderT order; - - if(!ReadBin4(fp, src->mPath, BO_LITTLE, 4, &fourCC) || - !ReadBin4(fp, src->mPath, BO_LITTLE, 4, &dummy)) - return 0; - if(fourCC == FOURCC_RIFF) - order = BO_LITTLE; - else if(fourCC == FOURCC_RIFX) - order = BO_BIG; - else - { - fprintf(stderr, "\nError: No RIFF/RIFX chunk in file '%s'.\n", src->mPath); - return 0; - } - - if(!ReadBin4(fp, src->mPath, BO_LITTLE, 4, &fourCC)) - return 0; - if(fourCC != FOURCC_WAVE) - { - fprintf(stderr, "\nError: Not a RIFF/RIFX WAVE file '%s'.\n", src->mPath); - return 0; - } - if(!ReadWaveFormat(fp, order, hrirRate, src)) - return 0; - if(!ReadWaveList(fp, src, order, n, hrir)) - return 0; - return 1; -} - -// Load a Spatially Oriented Format for Accoustics (SOFA) file. -static struct MYSOFA_EASY* LoadSofaFile(SourceRefT *src, const uint hrirRate, const uint n) -{ - struct MYSOFA_EASY *sofa{mysofa_cache_lookup(src->mPath, (float)hrirRate)}; - if(sofa) return sofa; - - sofa = static_cast<MYSOFA_EASY*>(calloc(1, sizeof(*sofa))); - if(sofa == nullptr) - { - fprintf(stderr, "\nError: Out of memory.\n"); - return nullptr; - } - sofa->lookup = nullptr; - sofa->neighborhood = nullptr; - - int err; - sofa->hrtf = mysofa_load(src->mPath, &err); - if(!sofa->hrtf) - { - mysofa_close(sofa); - fprintf(stderr, "\nError: Could not load source file '%s'.\n", src->mPath); - return nullptr; - } - err = mysofa_check(sofa->hrtf); - if(err != MYSOFA_OK) -/* NOTE: Some valid SOFA files are failing this check. - { - mysofa_close(sofa); - fprintf(stderr, "\nError: Malformed source file '%s'.\n", src->mPath); - return nullptr; - }*/ - fprintf(stderr, "\nWarning: Supposedly malformed source file '%s'.\n", src->mPath); - if((src->mOffset + n) > sofa->hrtf->N) - { - mysofa_close(sofa); - fprintf(stderr, "\nError: Not enough samples in SOFA file '%s'.\n", src->mPath); - return nullptr; - } - if(src->mChannel >= sofa->hrtf->R) - { - mysofa_close(sofa); - fprintf(stderr, "\nError: Missing source receiver in SOFA file '%s'.\n", src->mPath); - return nullptr; - } - mysofa_tocartesian(sofa->hrtf); - sofa->lookup = mysofa_lookup_init(sofa->hrtf); - if(sofa->lookup == nullptr) - { - mysofa_close(sofa); - fprintf(stderr, "\nError: Out of memory.\n"); - return nullptr; - } - return mysofa_cache_store(sofa, src->mPath, (float)hrirRate); -} - -// Copies the HRIR data from a particular SOFA measurement. -static void ExtractSofaHrir(const struct MYSOFA_EASY *sofa, const uint index, const uint channel, const uint offset, const uint n, double *hrir) -{ - for(uint i{0u};i < n;i++) - hrir[i] = sofa->hrtf->DataIR.values[(index*sofa->hrtf->R + channel)*sofa->hrtf->N + offset + i]; -} - -// Load a source HRIR from a Spatially Oriented Format for Accoustics (SOFA) -// file. -static int LoadSofaSource(SourceRefT *src, const uint hrirRate, const uint n, double *hrir) -{ - struct MYSOFA_EASY *sofa; - float target[3]; - int nearest; - float *coords; - - sofa = LoadSofaFile(src, hrirRate, n); - if(sofa == nullptr) - return 0; - - /* NOTE: At some point it may be benficial or necessary to consider the - various coordinate systems, listener/source orientations, and - direciontal vectors defined in the SOFA file. - */ - target[0] = src->mAzimuth; - target[1] = src->mElevation; - target[2] = src->mRadius; - mysofa_s2c(target); - - nearest = mysofa_lookup(sofa->lookup, target); - if(nearest < 0) - { - fprintf(stderr, "\nError: Lookup failed in source file '%s'.\n", src->mPath); - return 0; - } - - coords = &sofa->hrtf->SourcePosition.values[3 * nearest]; - if(std::fabs(coords[0] - target[0]) > 0.001 || std::fabs(coords[1] - target[1]) > 0.001 || std::fabs(coords[2] - target[2]) > 0.001) - { - fprintf(stderr, "\nError: No impulse response at coordinates (%.3fr, %.1fev, %.1faz) in file '%s'.\n", src->mRadius, src->mElevation, src->mAzimuth, src->mPath); - target[0] = coords[0]; - target[1] = coords[1]; - target[2] = coords[2]; - mysofa_c2s(target); - fprintf(stderr, " Nearest candidate at (%.3fr, %.1fev, %.1faz).\n", target[2], target[1], target[0]); - return 0; - } - - ExtractSofaHrir(sofa, nearest, src->mChannel, src->mOffset, n, hrir); - - return 1; -} - -// Load a source HRIR from a supported file type. -static int LoadSource(SourceRefT *src, const uint hrirRate, const uint n, double *hrir) -{ - FILE *fp{nullptr}; - if(src->mFormat != SF_SOFA) - { - if(src->mFormat == SF_ASCII) - fp = fopen(src->mPath, "r"); - else - fp = fopen(src->mPath, "rb"); - if(fp == nullptr) - { - fprintf(stderr, "\nError: Could not open source file '%s'.\n", src->mPath); - return 0; - } - } - int result; - switch(src->mFormat) - { - case SF_ASCII: - result = LoadAsciiSource(fp, src, n, hrir); - break; - case SF_BIN_LE: - result = LoadBinarySource(fp, src, BO_LITTLE, n, hrir); - break; - case SF_BIN_BE: - result = LoadBinarySource(fp, src, BO_BIG, n, hrir); - break; - case SF_WAVE: - result = LoadWaveSource(fp, src, hrirRate, n, hrir); - break; - case SF_SOFA: - result = LoadSofaSource(src, hrirRate, n, hrir); - break; - default: - result = 0; - } - if(fp) fclose(fp); - return result; -} - /*************************** *** File storage output *** @@ -2029,49 +775,6 @@ static int StoreMhr(const HrirDataT *hData, const char *filename) *** HRTF processing *** ***********************/ -// Calculate the onset time of an HRIR and average it with any existing -// timing for its field, elevation, azimuth, and ear. -static double AverageHrirOnset(const uint rate, const uint n, const double *hrir, const double f, const double onset) -{ - std::vector<double> upsampled(10 * n); - { - ResamplerT rs; - ResamplerSetup(&rs, rate, 10 * rate); - ResamplerRun(&rs, n, hrir, 10 * n, upsampled.data()); - } - - double mag{0.0}; - for(uint i{0u};i < 10*n;i++) - mag = std::max(std::abs(upsampled[i]), mag); - - mag *= 0.15; - uint i{0u}; - for(;i < 10*n;i++) - { - if(std::abs(upsampled[i]) >= mag) - break; - } - return Lerp(onset, static_cast<double>(i) / (10*rate), f); -} - -// Calculate the magnitude response of an HRIR and average it with any -// existing responses for its field, elevation, azimuth, and ear. -static void AverageHrirMagnitude(const uint points, const uint n, const double *hrir, const double f, double *mag) -{ - uint m = 1 + (n / 2), i; - std::vector<complex_d> h(n); - std::vector<double> r(n); - - for(i = 0;i < points;i++) - h[i] = complex_d{hrir[i], 0.0}; - for(;i < n;i++) - h[i] = complex_d{0.0, 0.0}; - FftForward(n, h.data()); - MagnitudeResponse(n, h.data(), r.data()); - for(i = 0;i < m;i++) - mag[i] = Lerp(mag[i], r[i], f); -} - /* Balances the maximum HRIR magnitudes of multi-field data sets by * independently normalizing each field in relation to the overall maximum. * This is done to ignore distance attenuation. @@ -2781,7 +1484,7 @@ static void CalculateHrtds(const HeadModelT model, const double radius, HrirData } // Allocate and configure dynamic HRIR structures. -static int PrepareHrirData(const uint fdCount, const double distances[MAX_FD_COUNT], const uint evCounts[MAX_FD_COUNT], const uint azCounts[MAX_FD_COUNT * MAX_EV_COUNT], HrirDataT *hData) +int PrepareHrirData(const uint fdCount, const double distances[MAX_FD_COUNT], const uint evCounts[MAX_FD_COUNT], const uint azCounts[MAX_FD_COUNT * MAX_EV_COUNT], HrirDataT *hData) { uint evTotal = 0, azTotal = 0, fi, ei, ai; @@ -2832,775 +1535,6 @@ static int PrepareHrirData(const uint fdCount, const double distances[MAX_FD_COU return 1; } -// Match the channel type from a given identifier. -static ChannelTypeT MatchChannelType(const char *ident) -{ - if(strcasecmp(ident, "mono") == 0) - return CT_MONO; - if(strcasecmp(ident, "stereo") == 0) - return CT_STEREO; - return CT_NONE; -} - -// Process the data set definition to read and validate the data set metrics. -static int ProcessMetrics(TokenReaderT *tr, const uint fftSize, const uint truncSize, HrirDataT *hData) -{ - int hasRate = 0, hasType = 0, hasPoints = 0, hasRadius = 0; - int hasDistance = 0, hasAzimuths = 0; - char ident[MAX_IDENT_LEN+1]; - uint line, col; - double fpVal; - uint points; - int intVal; - double distances[MAX_FD_COUNT]; - uint fdCount = 0; - uint evCounts[MAX_FD_COUNT]; - std::vector<uint> azCounts(MAX_FD_COUNT * MAX_EV_COUNT); - - TrIndication(tr, &line, &col); - while(TrIsIdent(tr)) - { - TrIndication(tr, &line, &col); - if(!TrReadIdent(tr, MAX_IDENT_LEN, ident)) - return 0; - if(strcasecmp(ident, "rate") == 0) - { - if(hasRate) - { - TrErrorAt(tr, line, col, "Redefinition of 'rate'.\n"); - return 0; - } - if(!TrReadOperator(tr, "=")) - return 0; - if(!TrReadInt(tr, MIN_RATE, MAX_RATE, &intVal)) - return 0; - hData->mIrRate = static_cast<uint>(intVal); - hasRate = 1; - } - else if(strcasecmp(ident, "type") == 0) - { - char type[MAX_IDENT_LEN+1]; - - if(hasType) - { - TrErrorAt(tr, line, col, "Redefinition of 'type'.\n"); - return 0; - } - if(!TrReadOperator(tr, "=")) - return 0; - - if(!TrReadIdent(tr, MAX_IDENT_LEN, type)) - return 0; - hData->mChannelType = MatchChannelType(type); - if(hData->mChannelType == CT_NONE) - { - TrErrorAt(tr, line, col, "Expected a channel type.\n"); - return 0; - } - hasType = 1; - } - else if(strcasecmp(ident, "points") == 0) - { - if(hasPoints) - { - TrErrorAt(tr, line, col, "Redefinition of 'points'.\n"); - return 0; - } - if(!TrReadOperator(tr, "=")) - return 0; - TrIndication(tr, &line, &col); - if(!TrReadInt(tr, MIN_POINTS, MAX_POINTS, &intVal)) - return 0; - points = static_cast<uint>(intVal); - if(fftSize > 0 && points > fftSize) - { - TrErrorAt(tr, line, col, "Value exceeds the overridden FFT size.\n"); - return 0; - } - if(points < truncSize) - { - TrErrorAt(tr, line, col, "Value is below the truncation size.\n"); - return 0; - } - hData->mIrPoints = points; - if(fftSize <= 0) - { - hData->mFftSize = DEFAULT_FFTSIZE; - hData->mIrSize = 1 + (DEFAULT_FFTSIZE / 2); - } - else - { - hData->mFftSize = fftSize; - hData->mIrSize = 1 + (fftSize / 2); - if(points > hData->mIrSize) - hData->mIrSize = points; - } - hasPoints = 1; - } - else if(strcasecmp(ident, "radius") == 0) - { - if(hasRadius) - { - TrErrorAt(tr, line, col, "Redefinition of 'radius'.\n"); - return 0; - } - if(!TrReadOperator(tr, "=")) - return 0; - if(!TrReadFloat(tr, MIN_RADIUS, MAX_RADIUS, &fpVal)) - return 0; - hData->mRadius = fpVal; - hasRadius = 1; - } - else if(strcasecmp(ident, "distance") == 0) - { - uint count = 0; - - if(hasDistance) - { - TrErrorAt(tr, line, col, "Redefinition of 'distance'.\n"); - return 0; - } - if(!TrReadOperator(tr, "=")) - return 0; - - for(;;) - { - if(!TrReadFloat(tr, MIN_DISTANCE, MAX_DISTANCE, &fpVal)) - return 0; - if(count > 0 && fpVal <= distances[count - 1]) - { - TrError(tr, "Distances are not ascending.\n"); - return 0; - } - distances[count++] = fpVal; - if(!TrIsOperator(tr, ",")) - break; - if(count >= MAX_FD_COUNT) - { - TrError(tr, "Exceeded the maximum of %d fields.\n", MAX_FD_COUNT); - return 0; - } - TrReadOperator(tr, ","); - } - if(fdCount != 0 && count != fdCount) - { - TrError(tr, "Did not match the specified number of %d fields.\n", fdCount); - return 0; - } - fdCount = count; - hasDistance = 1; - } - else if(strcasecmp(ident, "azimuths") == 0) - { - uint count = 0; - - if(hasAzimuths) - { - TrErrorAt(tr, line, col, "Redefinition of 'azimuths'.\n"); - return 0; - } - if(!TrReadOperator(tr, "=")) - return 0; - - evCounts[0] = 0; - for(;;) - { - if(!TrReadInt(tr, MIN_AZ_COUNT, MAX_AZ_COUNT, &intVal)) - return 0; - azCounts[(count * MAX_EV_COUNT) + evCounts[count]++] = static_cast<uint>(intVal); - if(TrIsOperator(tr, ",")) - { - if(evCounts[count] >= MAX_EV_COUNT) - { - TrError(tr, "Exceeded the maximum of %d elevations.\n", MAX_EV_COUNT); - return 0; - } - TrReadOperator(tr, ","); - } - else - { - if(evCounts[count] < MIN_EV_COUNT) - { - TrErrorAt(tr, line, col, "Did not reach the minimum of %d azimuth counts.\n", MIN_EV_COUNT); - return 0; - } - if(azCounts[count * MAX_EV_COUNT] != 1 || azCounts[(count * MAX_EV_COUNT) + evCounts[count] - 1] != 1) - { - TrError(tr, "Poles are not singular for field %d.\n", count - 1); - return 0; - } - count++; - if(!TrIsOperator(tr, ";")) - break; - - if(count >= MAX_FD_COUNT) - { - TrError(tr, "Exceeded the maximum number of %d fields.\n", MAX_FD_COUNT); - return 0; - } - evCounts[count] = 0; - TrReadOperator(tr, ";"); - } - } - if(fdCount != 0 && count != fdCount) - { - TrError(tr, "Did not match the specified number of %d fields.\n", fdCount); - return 0; - } - fdCount = count; - hasAzimuths = 1; - } - else - { - TrErrorAt(tr, line, col, "Expected a metric name.\n"); - return 0; - } - TrSkipWhitespace(tr); - } - if(!(hasRate && hasPoints && hasRadius && hasDistance && hasAzimuths)) - { - TrErrorAt(tr, line, col, "Expected a metric name.\n"); - return 0; - } - if(distances[0] < hData->mRadius) - { - TrError(tr, "Distance cannot start below head radius.\n"); - return 0; - } - if(hData->mChannelType == CT_NONE) - hData->mChannelType = CT_MONO; - if(!PrepareHrirData(fdCount, distances, evCounts, azCounts.data(), hData)) - { - fprintf(stderr, "Error: Out of memory.\n"); - exit(-1); - } - return 1; -} - -// Parse an index triplet from the data set definition. -static int ReadIndexTriplet(TokenReaderT *tr, const HrirDataT *hData, uint *fi, uint *ei, uint *ai) -{ - int intVal; - - if(hData->mFdCount > 1) - { - if(!TrReadInt(tr, 0, static_cast<int>(hData->mFdCount) - 1, &intVal)) - return 0; - *fi = static_cast<uint>(intVal); - if(!TrReadOperator(tr, ",")) - return 0; - } - else - { - *fi = 0; - } - if(!TrReadInt(tr, 0, static_cast<int>(hData->mFds[*fi].mEvCount) - 1, &intVal)) - return 0; - *ei = static_cast<uint>(intVal); - if(!TrReadOperator(tr, ",")) - return 0; - if(!TrReadInt(tr, 0, static_cast<int>(hData->mFds[*fi].mEvs[*ei].mAzCount) - 1, &intVal)) - return 0; - *ai = static_cast<uint>(intVal); - return 1; -} - -// Match the source format from a given identifier. -static SourceFormatT MatchSourceFormat(const char *ident) -{ - if(strcasecmp(ident, "ascii") == 0) - return SF_ASCII; - if(strcasecmp(ident, "bin_le") == 0) - return SF_BIN_LE; - if(strcasecmp(ident, "bin_be") == 0) - return SF_BIN_BE; - if(strcasecmp(ident, "wave") == 0) - return SF_WAVE; - if(strcasecmp(ident, "sofa") == 0) - return SF_SOFA; - return SF_NONE; -} - -// Match the source element type from a given identifier. -static ElementTypeT MatchElementType(const char *ident) -{ - if(strcasecmp(ident, "int") == 0) - return ET_INT; - if(strcasecmp(ident, "fp") == 0) - return ET_FP; - return ET_NONE; -} - -// Parse and validate a source reference from the data set definition. -static int ReadSourceRef(TokenReaderT *tr, SourceRefT *src) -{ - char ident[MAX_IDENT_LEN+1]; - uint line, col; - double fpVal; - int intVal; - - TrIndication(tr, &line, &col); - if(!TrReadIdent(tr, MAX_IDENT_LEN, ident)) - return 0; - src->mFormat = MatchSourceFormat(ident); - if(src->mFormat == SF_NONE) - { - TrErrorAt(tr, line, col, "Expected a source format.\n"); - return 0; - } - if(!TrReadOperator(tr, "(")) - return 0; - if(src->mFormat == SF_SOFA) - { - if(!TrReadFloat(tr, MIN_DISTANCE, MAX_DISTANCE, &fpVal)) - return 0; - src->mRadius = fpVal; - if(!TrReadOperator(tr, ",")) - return 0; - if(!TrReadFloat(tr, -90.0, 90.0, &fpVal)) - return 0; - src->mElevation = fpVal; - if(!TrReadOperator(tr, ",")) - return 0; - if(!TrReadFloat(tr, -360.0, 360.0, &fpVal)) - return 0; - src->mAzimuth = fpVal; - if(!TrReadOperator(tr, ":")) - return 0; - if(!TrReadInt(tr, 0, MAX_WAVE_CHANNELS, &intVal)) - return 0; - src->mType = ET_NONE; - src->mSize = 0; - src->mBits = 0; - src->mChannel = (uint)intVal; - src->mSkip = 0; - } - else if(src->mFormat == SF_WAVE) - { - if(!TrReadInt(tr, 0, MAX_WAVE_CHANNELS, &intVal)) - return 0; - src->mType = ET_NONE; - src->mSize = 0; - src->mBits = 0; - src->mChannel = static_cast<uint>(intVal); - src->mSkip = 0; - } - else - { - TrIndication(tr, &line, &col); - if(!TrReadIdent(tr, MAX_IDENT_LEN, ident)) - return 0; - src->mType = MatchElementType(ident); - if(src->mType == ET_NONE) - { - TrErrorAt(tr, line, col, "Expected a source element type.\n"); - return 0; - } - if(src->mFormat == SF_BIN_LE || src->mFormat == SF_BIN_BE) - { - if(!TrReadOperator(tr, ",")) - return 0; - if(src->mType == ET_INT) - { - if(!TrReadInt(tr, MIN_BIN_SIZE, MAX_BIN_SIZE, &intVal)) - return 0; - src->mSize = static_cast<uint>(intVal); - if(!TrIsOperator(tr, ",")) - src->mBits = static_cast<int>(8*src->mSize); - else - { - TrReadOperator(tr, ","); - TrIndication(tr, &line, &col); - if(!TrReadInt(tr, -2147483647-1, 2147483647, &intVal)) - return 0; - if(std::abs(intVal) < MIN_BIN_BITS || static_cast<uint>(std::abs(intVal)) > (8*src->mSize)) - { - TrErrorAt(tr, line, col, "Expected a value of (+/-) %d to %d.\n", MIN_BIN_BITS, 8*src->mSize); - return 0; - } - src->mBits = intVal; - } - } - else - { - TrIndication(tr, &line, &col); - if(!TrReadInt(tr, -2147483647-1, 2147483647, &intVal)) - return 0; - if(intVal != 4 && intVal != 8) - { - TrErrorAt(tr, line, col, "Expected a value of 4 or 8.\n"); - return 0; - } - src->mSize = static_cast<uint>(intVal); - src->mBits = 0; - } - } - else if(src->mFormat == SF_ASCII && src->mType == ET_INT) - { - if(!TrReadOperator(tr, ",")) - return 0; - if(!TrReadInt(tr, MIN_ASCII_BITS, MAX_ASCII_BITS, &intVal)) - return 0; - src->mSize = 0; - src->mBits = intVal; - } - else - { - src->mSize = 0; - src->mBits = 0; - } - - if(!TrIsOperator(tr, ";")) - src->mSkip = 0; - else - { - TrReadOperator(tr, ";"); - if(!TrReadInt(tr, 0, 0x7FFFFFFF, &intVal)) - return 0; - src->mSkip = static_cast<uint>(intVal); - } - } - if(!TrReadOperator(tr, ")")) - return 0; - if(TrIsOperator(tr, "@")) - { - TrReadOperator(tr, "@"); - if(!TrReadInt(tr, 0, 0x7FFFFFFF, &intVal)) - return 0; - src->mOffset = static_cast<uint>(intVal); - } - else - src->mOffset = 0; - if(!TrReadOperator(tr, ":")) - return 0; - if(!TrReadString(tr, MAX_PATH_LEN, src->mPath)) - return 0; - return 1; -} - -// Parse and validate a SOFA source reference from the data set definition. -static int ReadSofaRef(TokenReaderT *tr, SourceRefT *src) -{ - char ident[MAX_IDENT_LEN+1]; - uint line, col; - int intVal; - - TrIndication(tr, &line, &col); - if(!TrReadIdent(tr, MAX_IDENT_LEN, ident)) - return 0; - src->mFormat = MatchSourceFormat(ident); - if(src->mFormat != SF_SOFA) - { - TrErrorAt(tr, line, col, "Expected the SOFA source format.\n"); - return 0; - } - - src->mType = ET_NONE; - src->mSize = 0; - src->mBits = 0; - src->mChannel = 0; - src->mSkip = 0; - - if(TrIsOperator(tr, "@")) - { - TrReadOperator(tr, "@"); - if(!TrReadInt(tr, 0, 0x7FFFFFFF, &intVal)) - return 0; - src->mOffset = (uint)intVal; - } - else - src->mOffset = 0; - if(!TrReadOperator(tr, ":")) - return 0; - if(!TrReadString(tr, MAX_PATH_LEN, src->mPath)) - return 0; - return 1; -} - -// Match the target ear (index) from a given identifier. -static int MatchTargetEar(const char *ident) -{ - if(strcasecmp(ident, "left") == 0) - return 0; - if(strcasecmp(ident, "right") == 0) - return 1; - return -1; -} - -// Process the list of sources in the data set definition. -static int ProcessSources(const HeadModelT model, TokenReaderT *tr, HrirDataT *hData) -{ - uint channels = (hData->mChannelType == CT_STEREO) ? 2 : 1; - hData->mHrirsBase.resize(channels * hData->mIrCount * hData->mIrSize); - double *hrirs = hData->mHrirsBase.data(); - std::vector<double> hrir(hData->mIrPoints); - uint line, col, fi, ei, ai, ti; - int count; - - printf("Loading sources..."); - fflush(stdout); - count = 0; - while(TrIsOperator(tr, "[")) - { - double factor[2]{ 1.0, 1.0 }; - - TrIndication(tr, &line, &col); - TrReadOperator(tr, "["); - - if(TrIsOperator(tr, "*")) - { - SourceRefT src; - struct MYSOFA_EASY *sofa; - uint si; - - TrReadOperator(tr, "*"); - if(!TrReadOperator(tr, "]") || !TrReadOperator(tr, "=")) - return 0; - - TrIndication(tr, &line, &col); - if(!ReadSofaRef(tr, &src)) - return 0; - - if(hData->mChannelType == CT_STEREO) - { - char type[MAX_IDENT_LEN+1]; - ChannelTypeT channelType; - - if(!TrReadIdent(tr, MAX_IDENT_LEN, type)) - return 0; - - channelType = MatchChannelType(type); - - switch(channelType) - { - case CT_NONE: - TrErrorAt(tr, line, col, "Expected a channel type.\n"); - return 0; - case CT_MONO: - src.mChannel = 0; - break; - case CT_STEREO: - src.mChannel = 1; - break; - } - } - else - { - char type[MAX_IDENT_LEN+1]; - ChannelTypeT channelType; - - if(!TrReadIdent(tr, MAX_IDENT_LEN, type)) - return 0; - - channelType = MatchChannelType(type); - if(channelType != CT_MONO) - { - TrErrorAt(tr, line, col, "Expected a mono channel type.\n"); - return 0; - } - src.mChannel = 0; - } - - sofa = LoadSofaFile(&src, hData->mIrRate, hData->mIrPoints); - if(!sofa) return 0; - - for(si = 0;si < sofa->hrtf->M;si++) - { - printf("\rLoading sources... %d of %d", si+1, sofa->hrtf->M); - fflush(stdout); - - float aer[3] = { - sofa->hrtf->SourcePosition.values[3*si], - sofa->hrtf->SourcePosition.values[3*si + 1], - sofa->hrtf->SourcePosition.values[3*si + 2] - }; - mysofa_c2s(aer); - - if(std::fabs(aer[1]) >= 89.999f) - aer[0] = 0.0f; - else - aer[0] = std::fmod(360.0f - aer[0], 360.0f); - - for(fi = 0;fi < hData->mFdCount;fi++) - { - double delta = aer[2] - hData->mFds[fi].mDistance; - if(std::abs(delta) < 0.001) - break; - } - if(fi >= hData->mFdCount) - continue; - - double ef{(90.0 + aer[1]) * (hData->mFds[fi].mEvCount - 1) / 180.0}; - ei = (int)std::round(ef); - ef = (ef - ei) * 180.0f / (hData->mFds[fi].mEvCount - 1); - if(std::abs(ef) >= 0.1) - continue; - - double af{aer[0] * hData->mFds[fi].mEvs[ei].mAzCount / 360.0f}; - ai = (int)std::round(af); - af = (af - ai) * 360.0f / hData->mFds[fi].mEvs[ei].mAzCount; - ai = ai % hData->mFds[fi].mEvs[ei].mAzCount; - if(std::abs(af) >= 0.1) - continue; - - HrirAzT *azd = &hData->mFds[fi].mEvs[ei].mAzs[ai]; - - if(azd->mIrs[0] != nullptr) - { - TrErrorAt(tr, line, col, "Redefinition of source [ %d, %d, %d ].\n", fi, ei, ai); - return 0; - } - - ExtractSofaHrir(sofa, si, 0, src.mOffset, hData->mIrPoints, hrir.data()); - azd->mIrs[0] = &hrirs[hData->mIrSize * azd->mIndex]; - if(model == HM_DATASET) - azd->mDelays[0] = AverageHrirOnset(hData->mIrRate, hData->mIrPoints, hrir.data(), 1.0, azd->mDelays[0]); - AverageHrirMagnitude(hData->mIrPoints, hData->mFftSize, hrir.data(), 1.0, azd->mIrs[0]); - - if(src.mChannel == 1) - { - ExtractSofaHrir(sofa, si, 1, src.mOffset, hData->mIrPoints, hrir.data()); - azd->mIrs[1] = &hrirs[hData->mIrSize * (hData->mIrCount + azd->mIndex)]; - if(model == HM_DATASET) - azd->mDelays[1] = AverageHrirOnset(hData->mIrRate, hData->mIrPoints, hrir.data(), 1.0, azd->mDelays[1]); - AverageHrirMagnitude(hData->mIrPoints, hData->mFftSize, hrir.data(), 1.0, azd->mIrs[1]); - } - - // TODO: Since some SOFA files contain minimum phase HRIRs, - // it would be beneficial to check for per-measurement delays - // (when available) to reconstruct the HRTDs. - } - - continue; - } - - if(!ReadIndexTriplet(tr, hData, &fi, &ei, &ai)) - return 0; - if(!TrReadOperator(tr, "]")) - return 0; - HrirAzT *azd = &hData->mFds[fi].mEvs[ei].mAzs[ai]; - - if(azd->mIrs[0] != nullptr) - { - TrErrorAt(tr, line, col, "Redefinition of source.\n"); - return 0; - } - if(!TrReadOperator(tr, "=")) - return 0; - - for(;;) - { - SourceRefT src; - uint ti = 0; - - if(!ReadSourceRef(tr, &src)) - return 0; - - // TODO: Would be nice to display 'x of y files', but that would - // require preparing the source refs first to get a total count - // before loading them. - ++count; - printf("\rLoading sources... %d file%s", count, (count==1)?"":"s"); - fflush(stdout); - - if(!LoadSource(&src, hData->mIrRate, hData->mIrPoints, hrir.data())) - return 0; - - if(hData->mChannelType == CT_STEREO) - { - char ident[MAX_IDENT_LEN+1]; - - if(!TrReadIdent(tr, MAX_IDENT_LEN, ident)) - return 0; - ti = MatchTargetEar(ident); - if(static_cast<int>(ti) < 0) - { - TrErrorAt(tr, line, col, "Expected a target ear.\n"); - return 0; - } - } - azd->mIrs[ti] = &hrirs[hData->mIrSize * (ti * hData->mIrCount + azd->mIndex)]; - if(model == HM_DATASET) - azd->mDelays[ti] = AverageHrirOnset(hData->mIrRate, hData->mIrPoints, hrir.data(), 1.0 / factor[ti], azd->mDelays[ti]); - AverageHrirMagnitude(hData->mIrPoints, hData->mFftSize, hrir.data(), 1.0 / factor[ti], azd->mIrs[ti]); - factor[ti] += 1.0; - if(!TrIsOperator(tr, "+")) - break; - TrReadOperator(tr, "+"); - } - if(hData->mChannelType == CT_STEREO) - { - if(azd->mIrs[0] == nullptr) - { - TrErrorAt(tr, line, col, "Missing left ear source reference(s).\n"); - return 0; - } - else if(azd->mIrs[1] == nullptr) - { - TrErrorAt(tr, line, col, "Missing right ear source reference(s).\n"); - return 0; - } - } - } - printf("\n"); - for(fi = 0;fi < hData->mFdCount;fi++) - { - for(ei = 0;ei < hData->mFds[fi].mEvCount;ei++) - { - for(ai = 0;ai < hData->mFds[fi].mEvs[ei].mAzCount;ai++) - { - HrirAzT *azd = &hData->mFds[fi].mEvs[ei].mAzs[ai]; - if(azd->mIrs[0] != nullptr) - break; - } - if(ai < hData->mFds[fi].mEvs[ei].mAzCount) - break; - } - if(ei >= hData->mFds[fi].mEvCount) - { - TrError(tr, "Missing source references [ %d, *, * ].\n", fi); - return 0; - } - hData->mFds[fi].mEvStart = ei; - for(;ei < hData->mFds[fi].mEvCount;ei++) - { - for(ai = 0;ai < hData->mFds[fi].mEvs[ei].mAzCount;ai++) - { - HrirAzT *azd = &hData->mFds[fi].mEvs[ei].mAzs[ai]; - - if(azd->mIrs[0] == nullptr) - { - TrError(tr, "Missing source reference [ %d, %d, %d ].\n", fi, ei, ai); - return 0; - } - } - } - } - for(ti = 0;ti < channels;ti++) - { - for(fi = 0;fi < hData->mFdCount;fi++) - { - for(ei = 0;ei < hData->mFds[fi].mEvCount;ei++) - { - for(ai = 0;ai < hData->mFds[fi].mEvs[ei].mAzCount;ai++) - { - HrirAzT *azd = &hData->mFds[fi].mEvs[ei].mAzs[ai]; - - azd->mIrs[ti] = &hrirs[hData->mIrSize * (ti * hData->mIrCount + azd->mIndex)]; - } - } - } - } - if(!TrLoad(tr)) - { - mysofa_cache_release_all(); - return 1; - } - - TrError(tr, "Errant data at end of source list.\n"); - mysofa_cache_release_all(); - return 0; -} /* Parse the data set definition and process the source data, storing the * resulting data set as desired. If the input name is NULL it will read @@ -3730,7 +1664,7 @@ int main(int argc, char *argv[]) outName = "./oalsoft_hrtf_%r.mhr"; outRate = 0; - fftSize = 0; + fftSize = DEFAULT_FFTSIZE; equalize = DEFAULT_EQUALIZE; surface = DEFAULT_SURFACE; limit = DEFAULT_LIMIT; diff --git a/utils/makemhr/makemhr.h b/utils/makemhr/makemhr.h new file mode 100644 index 00000000..d60aa300 --- /dev/null +++ b/utils/makemhr/makemhr.h @@ -0,0 +1,128 @@ +#ifndef MAKEMHR_H +#define MAKEMHR_H + +#include <vector> +#include <complex> + + +// The maximum path length used when processing filenames. +#define MAX_PATH_LEN (256) + +// The limit to the number of 'distances' listed in the data set definition. +#define MAX_FD_COUNT (16) + +// The limits to the number of 'azimuths' listed in the data set definition. +#define MIN_EV_COUNT (5) +#define MAX_EV_COUNT (128) + +// The limits for each of the 'azimuths' listed in the data set definition. +#define MIN_AZ_COUNT (1) +#define MAX_AZ_COUNT (128) + +// The limits for the 'distance' from source to listener for each field in +// the definition file. +#define MIN_DISTANCE (0.05) +#define MAX_DISTANCE (2.50) + +// The limits for the sample 'rate' metric in the data set definition and for +// resampling. +#define MIN_RATE (32000) +#define MAX_RATE (96000) + +// The limits for the HRIR 'points' metric in the data set definition. +#define MIN_POINTS (16) +#define MAX_POINTS (8192) + + +using uint = unsigned int; + +/* Complex double type. */ +using complex_d = std::complex<double>; + + +// Head model used for calculating the impulse delays. +enum HeadModelT { + HM_NONE, + HM_DATASET, // Measure the onset from the dataset. + HM_SPHERE // Calculate the onset using a spherical head model. +}; + + +// Sample and channel type enum values. +enum SampleTypeT { + ST_S16 = 0, + ST_S24 = 1 +}; + +// Certain iterations rely on these integer enum values. +enum ChannelTypeT { + CT_NONE = -1, + CT_MONO = 0, + CT_STEREO = 1 +}; + +// Structured HRIR storage for stereo azimuth pairs, elevations, and fields. +struct HrirAzT { + double mAzimuth{0.0}; + uint mIndex{0u}; + double mDelays[2]{0.0, 0.0}; + double *mIrs[2]{nullptr, nullptr}; +}; + +struct HrirEvT { + double mElevation{0.0}; + uint mIrCount{0u}; + uint mAzCount{0u}; + HrirAzT *mAzs{nullptr}; +}; + +struct HrirFdT { + double mDistance{0.0}; + uint mIrCount{0u}; + uint mEvCount{0u}; + uint mEvStart{0u}; + HrirEvT *mEvs{nullptr}; +}; + +// The HRIR metrics and data set used when loading, processing, and storing +// the resulting HRTF. +struct HrirDataT { + uint mIrRate{0u}; + SampleTypeT mSampleType{ST_S24}; + ChannelTypeT mChannelType{CT_NONE}; + uint mIrPoints{0u}; + uint mFftSize{0u}; + uint mIrSize{0u}; + double mRadius{0.0}; + uint mIrCount{0u}; + uint mFdCount{0u}; + + std::vector<double> mHrirsBase; + std::vector<HrirEvT> mEvsBase; + std::vector<HrirAzT> mAzsBase; + + std::vector<HrirFdT> mFds; +}; + + +int PrepareHrirData(const uint fdCount, const double distances[MAX_FD_COUNT], const uint evCounts[MAX_FD_COUNT], const uint azCounts[MAX_FD_COUNT * MAX_EV_COUNT], HrirDataT *hData); +void MagnitudeResponse(const uint n, const complex_d *in, double *out); +void FftForward(const uint n, complex_d *inout); +void FftInverse(const uint n, complex_d *inout); + + +// The resampler metrics and FIR filter. +struct ResamplerT { + uint mP, mQ, mM, mL; + std::vector<double> mF; +}; + +void ResamplerSetup(ResamplerT *rs, const uint srcRate, const uint dstRate); +void ResamplerRun(ResamplerT *rs, const uint inN, const double *in, const uint outN, double *out); + + +// Performs linear interpolation. +inline double Lerp(const double a, const double b, const double f) +{ return a + f * (b - a); } + +#endif /* MAKEMHR_H */ |