kcat/openal-soft
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Slight update to the UHJ coefficients

Browse Source 
The extended precision of the encoder's 0.6512*X term was guesswork, with no
real basis for it. Switch back to the original value until something better
actually comes up. Also updates the decoder to account for the change in the
encoder.
This commit is contained in:
Chris Robinson 2021-11-28 05:21:33 -08:00
parent bd254c5426
commit a73b64ce3c
3 changed files with 46 additions and 46 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

24
core/uhjfilter.cpp
View File

@ -28,7 +28,7 @@ const PhaseShifterT<UhjFilterBase::sFilterDelay*2> PShift{};
* *
* Left = (S + D)/2.0 * Left = (S + D)/2.0
* Right = (S - D)/2.0 * Right = (S - D)/2.0
* T = j(-0.1432*W + 0.6511746*X) - 0.7071068*Y * T = j(-0.1432*W + 0.6512*X) - 0.7071068*Y
* Q = 0.9772*Z * Q = 0.9772*Z
* *
* where j is a wide-band +90 degree phase shift. 3-channel UHJ excludes Q, * where j is a wide-band +90 degree phase shift. 3-channel UHJ excludes Q,
@ -95,9 +95,9 @@ void UhjEncoder::encode(float *LeftOut, float *RightOut, const FloatBufferLine *
* S = Left + Right * S = Left + Right
* D = Left - Right * D = Left - Right
* *
* W = 0.981530*S + 0.197484*j(0.828347*D + 0.767835*T) * W = 0.981532*S + 0.197484*j(0.828331*D + 0.767820*T)
* X = 0.418504*S - j(0.828347*D + 0.767835*T) * X = 0.418496*S - j(0.828331*D + 0.767820*T)
* Y = 0.795954*D - 0.676406*T + j(0.186626*S) * Y = 0.795968*D - 0.676392*T + j(0.186633*S)
* Z = 1.023332*Q * Z = 1.023332*Q
* *
* where j is a +90 degree phase shift. 3-channel UHJ excludes Q, while 2- * where j is a +90 degree phase shift. 3-channel UHJ excludes Q, while 2-
@ -130,19 +130,19 @@ void UhjDecoder::decode(const al::span<BufferLine> samples, const size_t offset,
float *RESTRICT xoutput{samples[1].data() + offset}; float *RESTRICT xoutput{samples[1].data() + offset};
float *RESTRICT youtput{samples[2].data() + offset}; float *RESTRICT youtput{samples[2].data() + offset};
/* Precompute j(0.828347*D + 0.767835*T) and store in xoutput. */ /* Precompute j(0.828331*D + 0.767820*T) and store in xoutput. */
auto tmpiter = std::copy(mDTHistory.cbegin(), mDTHistory.cend(), mTemp.begin()); auto tmpiter = std::copy(mDTHistory.cbegin(), mDTHistory.cend(), mTemp.begin());
std::transform(mD.cbegin(), mD.cbegin()+samplesToDo+sFilterDelay, mT.cbegin(), tmpiter, std::transform(mD.cbegin(), mD.cbegin()+samplesToDo+sFilterDelay, mT.cbegin(), tmpiter,
[](const float d, const float t) noexcept { return 0.828347f*d + 0.767835f*t; }); [](const float d, const float t) noexcept { return 0.828331f*d + 0.767820f*t; });
std::copy_n(mTemp.cbegin()+forwardSamples, mDTHistory.size(), mDTHistory.begin()); std::copy_n(mTemp.cbegin()+forwardSamples, mDTHistory.size(), mDTHistory.begin());
PShift.process({xoutput, samplesToDo}, mTemp.data()); PShift.process({xoutput, samplesToDo}, mTemp.data());
/* W = 0.981530*S + 0.197484*j(0.828347*D + 0.767835*T) */ /* W = 0.981532*S + 0.197484*j(0.828331*D + 0.767820*T) */
for(size_t i{0};i < samplesToDo;++i) for(size_t i{0};i < samplesToDo;++i)
woutput[i] = 0.981530f*mS[i] + 0.197484f*xoutput[i]; woutput[i] = 0.981532f*mS[i] + 0.197484f*xoutput[i];
/* X = 0.418504*S - j(0.828347*D + 0.767835*T) */ /* X = 0.418496*S - j(0.828331*D + 0.767820*T) */
for(size_t i{0};i < samplesToDo;++i) for(size_t i{0};i < samplesToDo;++i)
xoutput[i] = 0.418504f*mS[i] - xoutput[i]; xoutput[i] = 0.418496f*mS[i] - xoutput[i];
/* Precompute j*S and store in youtput. */ /* Precompute j*S and store in youtput. */
tmpiter = std::copy(mSHistory.cbegin(), mSHistory.cend(), mTemp.begin()); tmpiter = std::copy(mSHistory.cbegin(), mSHistory.cend(), mTemp.begin());
@ -150,9 +150,9 @@ void UhjDecoder::decode(const al::span<BufferLine> samples, const size_t offset,
std::copy_n(mTemp.cbegin()+forwardSamples, mSHistory.size(), mSHistory.begin()); std::copy_n(mTemp.cbegin()+forwardSamples, mSHistory.size(), mSHistory.begin());
PShift.process({youtput, samplesToDo}, mTemp.data()); PShift.process({youtput, samplesToDo}, mTemp.data());
/* Y = 0.795954*D - 0.676406*T + j(0.186626*S) */ /* Y = 0.795968*D - 0.676392*T + j(0.186633*S) */
for(size_t i{0};i < samplesToDo;++i) for(size_t i{0};i < samplesToDo;++i)
youtput[i] = 0.795954f*mD[i] - 0.676406f*mT[i] + 0.186626f*youtput[i]; youtput[i] = 0.795968f*mD[i] - 0.676392f*mT[i] + 0.186633f*youtput[i];
if(samples.size() > 3) if(samples.size() > 3)
{ {

66
utils/uhjdecoder.cpp
View File

@ -141,9 +141,9 @@ const PhaseShifterT<UhjDecoder::sFilterDelay*2> PShift{};
* S = Left + Right * S = Left + Right
* D = Left - Right * D = Left - Right
* *
* W = 0.981530*S + 0.197484*j(0.828347*D + 0.767835*T) * W = 0.981532*S + 0.197484*j(0.828331*D + 0.767820*T)
* X = 0.418504*S - j(0.828347*D + 0.767835*T) * X = 0.418496*S - j(0.828331*D + 0.767820*T)
* Y = 0.795954*D - 0.676406*T + j(0.186626*S) * Y = 0.795968*D - 0.676392*T + j(0.186633*S)
* Z = 1.023332*Q * Z = 1.023332*Q
* *
* where j is a +90 degree phase shift. 3-channel UHJ excludes Q, while 2- * where j is a +90 degree phase shift. 3-channel UHJ excludes Q, while 2-
@ -162,49 +162,49 @@ const PhaseShifterT<UhjDecoder::sFilterDelay*2> PShift{};
* Gerzon's original paper for deriving Sigma (S) or Delta (D) from the L and R * Gerzon's original paper for deriving Sigma (S) or Delta (D) from the L and R
* signals. As proof, taking Y for example: * signals. As proof, taking Y for example:
* *
* Y = 0.795954*D - 0.676406*T + j(0.186626*S) * Y = 0.795968*D - 0.676392*T + j(0.186633*S)
* *
* * Plug in the encoding parameters, using ? as a placeholder for whether S * * Plug in the encoding parameters, using ? as a placeholder for whether S
* and D should receive an extra 0.5 factor * and D should receive an extra 0.5 factor
* Y = 0.795954*(j(-0.3420201*W + 0.5098604*X) + 0.6554516*Y)*? - * Y = 0.795968*(j(-0.3420201*W + 0.5098604*X) + 0.6554516*Y)*? -
* 0.676406*(j(-0.1432*W + 0.6511746*X) - 0.7071068*Y) + * 0.676392*(j(-0.1432*W + 0.6512*X) - 0.7071068*Y) +
* 0.186626*j(0.9396926*W + 0.1855740*X)*? * 0.186633*j(0.9396926*W + 0.1855740*X)*?
* *
* * Move common factors in * * Move common factors in
* Y = (j(-0.3420201*0.795954*?*W + 0.5098604*0.795954*?*X) + 0.6554516*0.795954*?*Y) - * Y = (j(-0.3420201*0.795968*?*W + 0.5098604*0.795968*?*X) + 0.6554516*0.795968*?*Y) -
* (j(-0.1432*0.676406*W + 0.6511746*0.676406*X) - 0.7071068*0.676406*Y) + * (j(-0.1432*0.676392*W + 0.6512*0.676392*X) - 0.7071068*0.676392*Y) +
* j(0.9396926*0.186626*?*W + 0.1855740*0.186626*?*X) * j(0.9396926*0.186633*?*W + 0.1855740*0.186633*?*X)
* *
* * Clean up extraneous groupings * * Clean up extraneous groupings
* Y = j(-0.3420201*0.795954*?*W + 0.5098604*0.795954*?*X) + 0.6554516*0.795954*?*Y - * Y = j(-0.3420201*0.795968*?*W + 0.5098604*0.795968*?*X) + 0.6554516*0.795968*?*Y -
* j(-0.1432*0.676406*W + 0.6511746*0.676406*X) + 0.7071068*0.676406*Y + * j(-0.1432*0.676392*W + 0.6512*0.676392*X) + 0.7071068*0.676392*Y +
* j*(0.9396926*0.186626*?*W + 0.1855740*0.186626*?*X) * j*(0.9396926*0.186633*?*W + 0.1855740*0.186633*?*X)
* *
* * Move phase shifts together and combine them * * Move phase shifts together and combine them
* Y = j(-0.3420201*0.795954*?*W + 0.5098604*0.795954*?*X - -0.1432*0.676406*W - * Y = j(-0.3420201*0.795968*?*W + 0.5098604*0.795968*?*X - -0.1432*0.676392*W -
* 0.6511746*0.676406*X + 0.9396926*0.186626*?*W + 0.1855740*0.186626*?*X) + * 0.6512*0.676392*X + 0.9396926*0.186633*?*W + 0.1855740*0.186633*?*X) +
* 0.6554516*0.795954*?*Y + 0.7071068*0.676406*Y * 0.6554516*0.795968*?*Y + 0.7071068*0.676392*Y
* *
* * Reorder terms * * Reorder terms
* Y = j(-0.3420201*0.795954*?*W + 0.1432*0.676406*W + 0.9396926*0.186626*?*W + * Y = j(-0.3420201*0.795968*?*W + 0.1432*0.676392*W + 0.9396926*0.186633*?*W +
* 0.5098604*0.795954*?*X + -0.6511746*0.676406*X + 0.1855740*0.186626*?*X) + * 0.5098604*0.795968*?*X + -0.6512*0.676392*X + 0.1855740*0.186633*?*X) +
* 0.7071068*0.676406*Y + 0.6554516*0.795954*?*Y * 0.7071068*0.676392*Y + 0.6554516*0.795968*?*Y
* *
* * Move common factors out * * Move common factors out
* Y = j((-0.3420201*0.795954*? + 0.1432*0.676406 + 0.9396926*0.186626*?)*W + * Y = j((-0.3420201*0.795968*? + 0.1432*0.676392 + 0.9396926*0.186633*?)*W +
* ( 0.5098604*0.795954*? + -0.6511746*0.676406 + 0.1855740*0.186626*?)*X) + * ( 0.5098604*0.795968*? + -0.6512*0.676392 + 0.1855740*0.186633*?)*X) +
* (0.7071068*0.676406 + 0.6554516*0.795954*?)*Y * (0.7071068*0.676392 + 0.6554516*0.795968*?)*Y
* *
* * Result w/ 0.5 factor: * * Result w/ 0.5 factor:
* -0.3420201*0.795954*0.5 + 0.1432*0.676406 + 0.9396926*0.186626*0.5 = 0.04843*W * -0.3420201*0.795968*0.5 + 0.1432*0.676392 + 0.9396926*0.186633*0.5 = 0.04843*W
* 0.5098604*0.795954*0.5 + -0.6511746*0.676406 + 0.1855740*0.186626*0.5 = -0.22023*X * 0.5098604*0.795968*0.5 + -0.6512*0.676392 + 0.1855740*0.186633*0.5 = -0.22023*X
* 0.7071068*0.676406 + 0.6554516*0.795954*0.5 = 0.73915*Y * 0.7071068*0.676392 + 0.6554516*0.795968*0.5 = 0.73914*Y
* -> Y = j(0.04843*W + -0.22023*X) + 0.73915*Y * -> Y = j(0.04843*W + -0.22023*X) + 0.73914*Y
* *
* * Result w/o 0.5 factor: * * Result w/o 0.5 factor:
* -0.3420201*0.795954 + 0.1432*0.676406 + 0.9396926*0.186626 = 0.00000*W * -0.3420201*0.795968 + 0.1432*0.676392 + 0.9396926*0.186633 = 0.00000*W
* 0.5098604*0.795954 + -0.6511746*0.676406 + 0.1855740*0.186626 = 0.00000*X * 0.5098604*0.795968 + -0.6512*0.676392 + 0.1855740*0.186633 = 0.00000*X
* 0.7071068*0.676406 + 0.6554516*0.795954 = 1.00000*Y * 0.7071068*0.676392 + 0.6554516*0.795968 = 1.00000*Y
* -> Y = j(0.00000*W + 0.00000*X) + 1.00000*Y * -> Y = j(0.00000*W + 0.00000*X) + 1.00000*Y
* *
* Not halving produces a result matching the original input. * Not halving produces a result matching the original input.
@ -252,8 +252,8 @@ void UhjDecoder::decode(const float *RESTRICT InSamples, const size_t InChannels
for(size_t i{0};i < SamplesToDo;++i) for(size_t i{0};i < SamplesToDo;++i)
{ {
/* W = 0.981530*S + 0.197484*j(0.828347*D + 0.767835*T) */ /* W = 0.981532*S + 0.197484*j(0.828347*D + 0.767835*T) */
woutput[i] = 0.981530f*mS[i] + 0.197484f*xoutput[i]; woutput[i] = 0.981532f*mS[i] + 0.197484f*xoutput[i];
/* X = 0.418504*S - j(0.828347*D + 0.767835*T) */ /* X = 0.418504*S - j(0.828347*D + 0.767835*T) */
xoutput[i] = 0.418504f*mS[i] - xoutput[i]; xoutput[i] = 0.418504f*mS[i] - xoutput[i];
} }
@ -266,8 +266,8 @@ void UhjDecoder::decode(const float *RESTRICT InSamples, const size_t InChannels
for(size_t i{0};i < SamplesToDo;++i) for(size_t i{0};i < SamplesToDo;++i)
{ {
/* Y = 0.795954*D - 0.676406*T + j(0.186626*S) */ /* Y = 0.795968*D - 0.676392*T + j(0.186633*S) */
youtput[i] = 0.795954f*mD[i] - 0.676406f*mT[i] + 0.186626f*youtput[i]; youtput[i] = 0.795968f*mD[i] - 0.676392f*mT[i] + 0.186633f*youtput[i];
} }
if(OutSamples.size() > 3) if(OutSamples.size() > 3)

2
utils/uhjencoder.cpp
View File

@ -89,7 +89,7 @@ const PhaseShifterT<UhjEncoder::sFilterDelay*2> PShift{};
* *
* Left = (S + D)/2.0 * Left = (S + D)/2.0
* Right = (S - D)/2.0 * Right = (S - D)/2.0
* T = j(-0.1432*W + 0.6511746*X) - 0.7071068*Y * T = j(-0.1432*W + 0.6512*X) - 0.7071068*Y
* Q = 0.9772*Z * Q = 0.9772*Z
* *
* where j is a wide-band +90 degree phase shift. T is excluded from 2-channel * where j is a wide-band +90 degree phase shift. T is excluded from 2-channel