1// Copyright (C) 2022 The Qt Company Ltd.
2// SPDX-License-Identifier: LicenseRef-Qt-Commercial OR LGPL-3.0-only OR GPL-2.0-only OR GPL-3.0-only
3#include "qambisonicdecoder_p.h"
4
5#include "qambisonicdecoderdata_p.h"
6#include <cmath>
7#include <qdebug.h>
8
9QT_BEGIN_NAMESPACE
10
11// Ambisonic decoding is described in detail in https://ambisonics.dreamhosters.com/BLaH3.pdf.
12// We're using a phase matched band splitting filter to split the ambisonic signal into a low
13// and high frequency component and apply matrix conversions to those components individually
14// as described in the document.
15//
16// We are currently not using a near field compensation filter, something that could potentially
17// improve sound quality further.
18//
19// For mono and stereo decoding, we use a simpler algorithm to avoid artificially dampening signals
20// coming from the back, as we do not have any speakers in that direction and the calculations
21// through matlab would give us audible 'holes'.
22
23struct QAmbisonicDecoderData
24{
25 QAudioFormat::ChannelConfig config;
26 const float *lf[3];
27 const float *hf[3];
28 const float *reverb;
29};
30
31const float reverb_x_0[] =
32{
33 1.f, 0.f, // L
34 0.f, 1.f, // R
35 .7f, .7f, // C
36 1.f, 0.f, // Ls
37 0.f, 1.f, // Rs
38 1.f, 0.f, // Lb
39 0.f, 1.f, // Rb
40};
41
42const float reverb_x_1[] =
43{
44 1.f, 0.f, // L
45 0.f, 1.f, // R
46 .7f, .7f, // C
47 .0f, .0f, // LFE
48 1.f, 0.f, // Ls
49 0.f, 1.f, // Rs
50 1.f, 0.f, // Lb
51 0.f, 1.f, // Rb
52};
53
54static const QAmbisonicDecoderData decoderMap[] =
55{
56 { .config: QAudioFormat::ChannelConfigSurround5Dot0,
57 .lf: { decoderMatrix_5dot0_1_lf, decoderMatrix_5dot0_2_lf, decoderMatrix_5dot0_3_lf },
58 .hf: { decoderMatrix_5dot0_1_hf, decoderMatrix_5dot0_2_hf, decoderMatrix_5dot0_3_hf },
59 .reverb: reverb_x_0
60 },
61 { .config: QAudioFormat::ChannelConfigSurround5Dot1,
62 .lf: { decoderMatrix_5dot1_1_lf, decoderMatrix_5dot1_2_lf, decoderMatrix_5dot1_3_lf },
63 .hf: { decoderMatrix_5dot1_1_hf, decoderMatrix_5dot1_2_hf, decoderMatrix_5dot1_3_hf },
64 .reverb: reverb_x_1
65 },
66 { .config: QAudioFormat::ChannelConfigSurround7Dot0,
67 .lf: { decoderMatrix_7dot0_1_lf, decoderMatrix_7dot0_2_lf, decoderMatrix_7dot0_3_lf },
68 .hf: { decoderMatrix_7dot0_1_hf, decoderMatrix_7dot0_2_hf, decoderMatrix_7dot0_3_hf },
69 .reverb: reverb_x_0
70 },
71 { .config: QAudioFormat::ChannelConfigSurround7Dot1,
72 .lf: { decoderMatrix_7dot1_1_lf, decoderMatrix_7dot1_2_lf, decoderMatrix_7dot1_3_lf },
73 .hf: { decoderMatrix_7dot1_1_hf, decoderMatrix_7dot1_2_hf, decoderMatrix_7dot1_3_hf },
74 .reverb: reverb_x_1
75 }
76};
77
78// Implements a split second order IIR filter
79// The audio data is split into a phase synced low and high frequency part
80// This allows us to apply different factors to both parts for better sound
81// localization when converting from ambisonic formats
82//
83// Details are described in https://ambisonics.dreamhosters.com/BLaH3.pdf, Appendix A.2.
84class QAmbisonicDecoderFilter
85{
86public:
87 QAmbisonicDecoderFilter() = default;
88 void configure(float sampleRate, float cutoffFrequency = 380)
89 {
90 double k = tan(M_PI*cutoffFrequency/sampleRate);
91 a1 = float(2.*(k*k - 1.)/(k*k + 2*k + 1.));
92 a2 = float((k*k - 2*k + 1.)/(k*k + 2*k + 1.));
93
94 b0_lf = float(k*k/(k*k + 2*k + 1));
95 b1_lf = 2.f*b0_lf;
96
97 b0_hf = float(1./(k*k + 2*k + 1));
98 b1_hf = -2.f*b0_hf;
99 }
100
101 struct Output
102 {
103 float lf;
104 float hf;
105 };
106
107 Output next(float x)
108 {
109 float r_lf = x*b0_lf +
110 prevX[0]*b1_lf +
111 prevX[1]*b0_lf -
112 prevR_lf[0]*a1 -
113 prevR_lf[1]*a2;
114 float r_hf = x*b0_hf +
115 prevX[0]*b1_hf +
116 prevX[1]*b0_hf -
117 prevR_hf[0]*a1 -
118 prevR_hf[1]*a2;
119 prevX[1] = prevX[0];
120 prevX[0] = x;
121 prevR_lf[1] = prevR_lf[0];
122 prevR_lf[0] = r_lf;
123 prevR_hf[1] = prevR_hf[0];
124 prevR_hf[0] = r_hf;
125 return { .lf: r_lf, .hf: r_hf };
126 }
127
128private:
129 float a1 = 0.;
130 float a2 = 0.;
131
132 float b0_hf = 0.;
133 float b1_hf = 0.;
134
135 float b0_lf = 0.;
136 float b1_lf = 0.;
137
138 float prevX[2] = {};
139 float prevR_lf[2] = {};
140 float prevR_hf[2] = {};
141};
142
143
144QAmbisonicDecoder::QAmbisonicDecoder(AmbisonicLevel ambisonicLevel, const QAudioFormat &format)
145 : level(ambisonicLevel)
146{
147 Q_ASSERT(level > 0 && level <= 3);
148 inputChannels = (level+1)*(level+1);
149 outputChannels = format.channelCount();
150
151 channelConfig = format.channelConfig();
152 if (channelConfig == QAudioFormat::ChannelConfigUnknown)
153 channelConfig = format.defaultChannelConfigForChannelCount(channelCount: format.channelCount());
154
155 if (channelConfig == QAudioFormat::ChannelConfigMono ||
156 channelConfig == QAudioFormat::ChannelConfigStereo ||
157 channelConfig == QAudioFormat::ChannelConfig2Dot1 ||
158 channelConfig == QAudioFormat::ChannelConfig3Dot0 ||
159 channelConfig == QAudioFormat::ChannelConfig3Dot1) {
160 // these are non surround configs and handled manually to avoid
161 // audible holes for sounds coming from behing
162 //
163 // We use a simpler decoding process here, only taking first order
164 // ambisonics into account
165 //
166 // Left and right channels get 50% W and 50% X
167 // Center gets 50% W and 50% Y
168 // LFE gets 50% W
169 simpleDecoderFactors = new float[4*outputChannels];
170 float *r = new float[2*outputChannels]; // reverb output is in stereo
171 float *f = simpleDecoderFactors;
172 reverbFactors = r;
173 if (channelConfig & QAudioFormat::channelConfig(channels: QAudioFormat::FrontLeft)) {
174 f[0] = 0.5f; f[1] = 0.5f; f[2] = 0.; f[3] = 0.f;
175 f += 4;
176 r[0] = 1.; r[1] = 0.;
177 r += 2;
178 }
179 if (channelConfig & QAudioFormat::channelConfig(channels: QAudioFormat::FrontRight)) {
180 f[0] = 0.5f; f[1] = -0.5f; f[2] = 0.; f[3] = 0.f;
181 f += 4;
182 r[0] = 0.; r[1] = 1.;
183 r += 2;
184 }
185 if (channelConfig & QAudioFormat::channelConfig(channels: QAudioFormat::FrontCenter)) {
186 f[0] = 0.5f; f[1] = -0.f; f[2] = 0.; f[3] = 0.5f;
187 f += 4;
188 r[0] = .5; r[1] = .5;
189 r += 2;
190 }
191 if (channelConfig & QAudioFormat::channelConfig(channels: QAudioFormat::LFE)) {
192 f[0] = 0.5f; f[1] = -0.f; f[2] = 0.; f[3] = 0.0f;
193 f += 4;
194 r[0] = 0.; r[1] = 0.;
195 r += 2;
196 }
197 Q_UNUSED(f);
198 Q_UNUSED(r);
199 Q_ASSERT((f - simpleDecoderFactors) == 4*outputChannels);
200 Q_ASSERT((r - reverbFactors) == 2*outputChannels);
201
202 return;
203 }
204
205 for (const auto &d : decoderMap) {
206 if (d.config == channelConfig) {
207 decoderData = &d;
208 reverbFactors = decoderData->reverb;
209 break;
210 }
211 }
212 if (!decoderData) {
213 // can't handle this,
214 outputChannels = 0;
215 return;
216 }
217
218 filters = new QAmbisonicDecoderFilter[inputChannels];
219 for (int i = 0; i < inputChannels; ++i)
220 filters[i].configure(sampleRate: format.sampleRate());
221}
222
223QAmbisonicDecoder::~QAmbisonicDecoder()
224{
225 if (simpleDecoderFactors) {
226 delete simpleDecoderFactors;
227 delete reverbFactors;
228 }
229}
230
231void QAmbisonicDecoder::processBuffer(const float *input[], float *output, int nSamples)
232{
233 float *o = output;
234 memset(s: o, c: 0, n: nSamples*outputChannels*sizeof(float));
235
236 if (simpleDecoderFactors) {
237 for (int i = 0; i < nSamples; ++i) {
238 for (int j = 0; j < 4; ++j) {
239 for (int k = 0; k < outputChannels; ++k)
240 o[k] += simpleDecoderFactors[k*4 + j]*input[j][i];
241 }
242 o += outputChannels;
243 }
244 return;
245 }
246
247 const float *matrix_hi = decoderData->hf[level - 1];
248 const float *matrix_lo = decoderData->lf[level - 1];
249 for (int i = 0; i < nSamples; ++i) {
250 QAmbisonicDecoderFilter::Output buf[maxAmbisonicChannels];
251 for (int j = 0; j < inputChannels; ++j)
252 buf[j] = filters[j].next(x: input[j][i]);
253 for (int j = 0; j < inputChannels; ++j) {
254 for (int k = 0; k < outputChannels; ++k)
255 o[k] += matrix_lo[k*inputChannels + j]*buf[j].lf + matrix_hi[k*inputChannels + j]*buf[j].hf;
256 }
257 o += outputChannels;
258 }
259}
260
261void QAmbisonicDecoder::processBuffer(const float *input[], short *output, int nSamples)
262{
263 const float *reverb[] = { nullptr, nullptr };
264 return processBufferWithReverb(input, reverb, output, nSamples);
265}
266
267void QAmbisonicDecoder::processBufferWithReverb(const float *input[], const float *reverb[2], short *output, int nSamples)
268{
269 if (simpleDecoderFactors) {
270 for (int i = 0; i < nSamples; ++i) {
271 float o[4] = {};
272 for (int k = 0; k < outputChannels; ++k) {
273 for (int j = 0; j < 4; ++j)
274 o[k] += simpleDecoderFactors[k*4 + j]*input[j][i];
275 }
276 if (reverb[0]) {
277 for (int k = 0; k < outputChannels; ++k) {
278 o[k] += reverb[0][i]*reverbFactors[2*k] + reverb[1][i]*reverbFactors[2*k+1];
279 }
280 }
281
282 for (int k = 0; k < outputChannels; ++k)
283 output[k] = static_cast<short>(o[k]*32768.);
284 output += outputChannels;
285 }
286 return;
287 }
288
289 // qDebug() << "XXX" << inputChannels << outputChannels;
290 const float *matrix_hi = decoderData->hf[level - 1];
291 const float *matrix_lo = decoderData->lf[level - 1];
292 for (int i = 0; i < nSamples; ++i) {
293 QAmbisonicDecoderFilter::Output buf[maxAmbisonicChannels];
294 for (int j = 0; j < inputChannels; ++j)
295 buf[j] = filters[j].next(x: input[j][i]);
296 float o[32] = {}; // we can't support more than 32 channels from our API
297 for (int j = 0; j < inputChannels; ++j) {
298 for (int k = 0; k < outputChannels; ++k)
299 o[k] += matrix_lo[k*inputChannels + j]*buf[j].lf + matrix_hi[k*inputChannels + j]*buf[j].hf;
300 }
301 if (reverb[0]) {
302 for (int k = 0; k < outputChannels; ++k) {
303 o[k] += reverb[0][i]*reverbFactors[2*k] + reverb[1][i]*reverbFactors[2*k+1];
304 }
305 }
306 for (int k = 0; k < outputChannels; ++k)
307 output[k] = static_cast<short>(o[k]*32768.);
308 output += outputChannels;
309 }
310
311}
312
313QT_END_NAMESPACE
314
315

source code of qtmultimedia/src/spatialaudio/qambisonicdecoder.cpp