Fade between HRTF coefficients, to reduce noise from sudden changes

1 files changed, 138 insertions, 2 deletions

diff --git a/Alc/hrtf.c b/Alc/hrtf.c
index 4e03e198..6edc37e7 100644
--- a/Alc/hrtf.c
+++ b/Alc/hrtf.c
@@ -80,6 +80,37 @@ static void CalcAzIndices(ALuint evidx, ALfloat az, ALuint *azidx, ALfloat *azmu
     *azmu = az - (ALuint)az;
 }
 
+// Calculates the normalized HRTF transition factor (delta) from the changes
+// in gain and listener to source angle between updates.  The result is a
+// normalized delta factor than can be used to calculate moving HRIR stepping
+// values.
+ALfloat CalcHrtfDelta(ALfloat oldGain, ALfloat newGain, const ALfloat olddir[3], const ALfloat newdir[3])
+{
+    ALfloat gainChange, angleChange, delta;
+
+    // Calculate the normalized dB gain change.
+    gainChange = aluFabs((log10(__max(newGain, 0.0001f)) -
+                          log10(__max(oldGain, 0.0001f))) / log10(0.0001f));
+
+    // Calculate the normalized listener to source angle change when there is
+    // enough gain to notice it.
+    angleChange = 0.0f;
+    if(gainChange > 0.0001f || newGain > 0.0001f)
+    {
+        // No angle change when the directions are equal or degenerate (when
+        // both have zero length).
+        if(newdir[0]-olddir[0] || newdir[1]-olddir[1] || newdir[2]-olddir[2])
+            angleChange = aluAcos(olddir[0]*newdir[0] +
+                                  olddir[1]*newdir[1] +
+                                  olddir[2]*newdir[2]) / M_PI;
+
+    }
+
+    // Use the largest of the two changes for the delta factor.
+    delta = __max(gainChange, angleChange) * 2.0f;
+    return __min(delta, 1.0f);
+}
+
 // Calculates static HRIR coefficients and delays for the given polar
 // elevation and azimuth in radians.  Linear interpolation is used to
 // increase the apparent resolution of the HRIR dataset.  The coefficients
@@ -144,10 +175,115 @@ void GetLerpedHrtfCoeffs(ALfloat elevation, ALfloat azimuth, ALfloat gain, ALflo
     // Calculate the HRIR delays using linear interpolation.
     delays[0] = (ALuint)(lerp(lerp(Hrtf.delays[lidx[0]], Hrtf.delays[lidx[1]], mu[0]),
                               lerp(Hrtf.delays[lidx[2]], Hrtf.delays[lidx[3]], mu[1]),
-                              mu[2]) + 0.5f);
+                              mu[2]) * 65536.0f);
+    delays[1] = (ALuint)(lerp(lerp(Hrtf.delays[ridx[0]], Hrtf.delays[ridx[1]], mu[0]),
+                              lerp(Hrtf.delays[ridx[2]], Hrtf.delays[ridx[3]], mu[1]),
+                              mu[2]) * 65536.0f);
+}
+
+// Calculates the moving HRIR target coefficients, target delays, and
+// stepping values for the given polar elevation and azimuth in radians.
+// Linear interpolation is used to increase the apparent resolution of the
+// HRIR dataset.  The coefficients are also normalized and attenuated by the
+// specified gain.  Stepping resolution and count is determined using the
+// given delta factor between 0.0 and 1.0.
+ALint GetMovingHrtfCoeffs(ALfloat elevation, ALfloat azimuth, ALfloat gain, ALfloat delta, ALint counter, ALfloat (*coeffs)[2], ALuint *delays, ALfloat (*coeffStep)[2], ALint *delayStep)
+{
+    ALfloat step;
+    ALuint evidx[2], azidx[2];
+    ALfloat mu[3];
+    ALuint lidx[4], ridx[4];
+    ALuint i;
+    ALfloat left, right;
+
+    // Claculate elevation indices and interpolation factor.
+    CalcEvIndices(elevation, evidx, &mu[2]);
+
+    // Calculate azimuth indices and interpolation factor for the first
+    // elevation.
+    CalcAzIndices(evidx[0], azimuth, azidx, &mu[0]);
+
+    // Calculate the first set of linear HRIR indices for left and right
+    // channels.
+    lidx[0] = evOffset[evidx[0]] + azidx[0];
+    lidx[1] = evOffset[evidx[0]] + azidx[1];
+    ridx[0] = evOffset[evidx[0]] + ((azCount[evidx[0]]-azidx[0]) % azCount[evidx[0]]);
+    ridx[1] = evOffset[evidx[0]] + ((azCount[evidx[0]]-azidx[1]) % azCount[evidx[0]]);
+
+    // Calculate azimuth indices and interpolation factor for the second
+    // elevation.
+    CalcAzIndices(evidx[1], azimuth, azidx, &mu[1]);
+
+    // Calculate the second set of linear HRIR indices for left and right
+    // channels.
+    lidx[2] = evOffset[evidx[1]] + azidx[0];
+    lidx[3] = evOffset[evidx[1]] + azidx[1];
+    ridx[2] = evOffset[evidx[1]] + ((azCount[evidx[1]]-azidx[0]) % azCount[evidx[1]]);
+    ridx[3] = evOffset[evidx[1]] + ((azCount[evidx[1]]-azidx[1]) % azCount[evidx[1]]);
+
+    // Calculate the stepping parameters.
+    delta = __max(0.015f * delta * Hrtf.sampleRate, 1.0f);
+    step = 1.0f / delta;
+
+    // Calculate the normalized and attenuated target HRIR coefficients using
+    // linear interpolation when there is enough gain to warrant it.  Zero
+    // the target coefficients if gain is too low.  Then calculate the
+    // coefficient stepping values using the target and previous running
+    // coefficients.
+    if(gain > 0.0001f)
+    {
+        ALdouble scale = gain * (1.0/32767.0);
+        for(i = 0;i < HRIR_LENGTH;i++)
+        {
+            left = coeffs[i][0] - (coeffStep[i][0] * counter);
+            right = coeffs[i][1] - (coeffStep[i][1] * counter);
+
+            coeffs[i][0] = lerp(lerp(Hrtf.coeffs[lidx[0]][i], Hrtf.coeffs[lidx[1]][i], mu[0]),
+                                lerp(Hrtf.coeffs[lidx[2]][i], Hrtf.coeffs[lidx[3]][i], mu[1]),
+                                mu[2]) * scale;
+            coeffs[i][1] = lerp(lerp(Hrtf.coeffs[ridx[0]][i], Hrtf.coeffs[ridx[1]][i], mu[0]),
+                                lerp(Hrtf.coeffs[ridx[2]][i], Hrtf.coeffs[ridx[3]][i], mu[1]),
+                                mu[2]) * scale;
+
+            coeffStep[i][0] = step * (coeffs[i][0] - left);
+            coeffStep[i][1] = step * (coeffs[i][1] - right);
+        }
+    }
+    else
+    {
+        for(i = 0;i < HRIR_LENGTH;i++)
+        {
+            left = coeffs[i][0] - (coeffStep[i][0] * counter);
+            right = coeffs[i][1] - (coeffStep[i][1] * counter);
+
+            coeffs[i][0] = 0.0f;
+            coeffs[i][1] = 0.0f;
+
+            coeffStep[i][0] = step * -left;
+            coeffStep[i][1] = step * -right;
+        }
+    }
+
+    // Calculate the HRIR delays using linear interpolation.  Then calculate
+    // the delay stepping values using the target and previous running
+    // delays.
+    left = delays[0] - (delayStep[0] * counter);
+    right = delays[1] - (delayStep[1] * counter);
+
+    delays[0] = (ALuint)(lerp(lerp(Hrtf.delays[lidx[0]], Hrtf.delays[lidx[1]], mu[0]),
+                              lerp(Hrtf.delays[lidx[2]], Hrtf.delays[lidx[3]], mu[1]),
+                              mu[2]) * 65536.0f);
     delays[1] = (ALuint)(lerp(lerp(Hrtf.delays[ridx[0]], Hrtf.delays[ridx[1]], mu[0]),
                               lerp(Hrtf.delays[ridx[2]], Hrtf.delays[ridx[3]], mu[1]),
-                              mu[2]) + 0.5f);
+                              mu[2]) * 65536.0f);
+
+    delayStep[0] = (ALint)(step * (delays[0] - left));
+    delayStep[1] = (ALint)(step * (delays[1] - right));
+
+    // The stepping count is the number of samples necessary for the HRIR to
+    // complete its transition.  The mixer will only apply stepping for this
+    // many samples.
+    return (ALuint)delta;
 }
 
 ALCboolean IsHrtfCompatible(ALCdevice *device)