Coverage Report

Created: 2026-02-26 07:03

next uncovered line (L), next uncovered region (R), next uncovered branch (B)
/src/aac/libFDK/src/mdct.cpp
Line
Count
Source
1
/* -----------------------------------------------------------------------------
2
Software License for The Fraunhofer FDK AAC Codec Library for Android
3
4
© Copyright  1995 - 2019 Fraunhofer-Gesellschaft zur Förderung der angewandten
5
Forschung e.V. All rights reserved.
6
7
 1.    INTRODUCTION
8
The Fraunhofer FDK AAC Codec Library for Android ("FDK AAC Codec") is software
9
that implements the MPEG Advanced Audio Coding ("AAC") encoding and decoding
10
scheme for digital audio. This FDK AAC Codec software is intended to be used on
11
a wide variety of Android devices.
12
13
AAC's HE-AAC and HE-AAC v2 versions are regarded as today's most efficient
14
general perceptual audio codecs. AAC-ELD is considered the best-performing
15
full-bandwidth communications codec by independent studies and is widely
16
deployed. AAC has been standardized by ISO and IEC as part of the MPEG
17
specifications.
18
19
Patent licenses for necessary patent claims for the FDK AAC Codec (including
20
those of Fraunhofer) may be obtained through Via Licensing
21
(www.vialicensing.com) or through the respective patent owners individually for
22
the purpose of encoding or decoding bit streams in products that are compliant
23
with the ISO/IEC MPEG audio standards. Please note that most manufacturers of
24
Android devices already license these patent claims through Via Licensing or
25
directly from the patent owners, and therefore FDK AAC Codec software may
26
already be covered under those patent licenses when it is used for those
27
licensed purposes only.
28
29
Commercially-licensed AAC software libraries, including floating-point versions
30
with enhanced sound quality, are also available from Fraunhofer. Users are
31
encouraged to check the Fraunhofer website for additional applications
32
information and documentation.
33
34
2.    COPYRIGHT LICENSE
35
36
Redistribution and use in source and binary forms, with or without modification,
37
are permitted without payment of copyright license fees provided that you
38
satisfy the following conditions:
39
40
You must retain the complete text of this software license in redistributions of
41
the FDK AAC Codec or your modifications thereto in source code form.
42
43
You must retain the complete text of this software license in the documentation
44
and/or other materials provided with redistributions of the FDK AAC Codec or
45
your modifications thereto in binary form. You must make available free of
46
charge copies of the complete source code of the FDK AAC Codec and your
47
modifications thereto to recipients of copies in binary form.
48
49
The name of Fraunhofer may not be used to endorse or promote products derived
50
from this library without prior written permission.
51
52
You may not charge copyright license fees for anyone to use, copy or distribute
53
the FDK AAC Codec software or your modifications thereto.
54
55
Your modified versions of the FDK AAC Codec must carry prominent notices stating
56
that you changed the software and the date of any change. For modified versions
57
of the FDK AAC Codec, the term "Fraunhofer FDK AAC Codec Library for Android"
58
must be replaced by the term "Third-Party Modified Version of the Fraunhofer FDK
59
AAC Codec Library for Android."
60
61
3.    NO PATENT LICENSE
62
63
NO EXPRESS OR IMPLIED LICENSES TO ANY PATENT CLAIMS, including without
64
limitation the patents of Fraunhofer, ARE GRANTED BY THIS SOFTWARE LICENSE.
65
Fraunhofer provides no warranty of patent non-infringement with respect to this
66
software.
67
68
You may use this FDK AAC Codec software or modifications thereto only for
69
purposes that are authorized by appropriate patent licenses.
70
71
4.    DISCLAIMER
72
73
This FDK AAC Codec software is provided by Fraunhofer on behalf of the copyright
74
holders and contributors "AS IS" and WITHOUT ANY EXPRESS OR IMPLIED WARRANTIES,
75
including but not limited to the implied warranties of merchantability and
76
fitness for a particular purpose. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR
77
CONTRIBUTORS BE LIABLE for any direct, indirect, incidental, special, exemplary,
78
or consequential damages, including but not limited to procurement of substitute
79
goods or services; loss of use, data, or profits, or business interruption,
80
however caused and on any theory of liability, whether in contract, strict
81
liability, or tort (including negligence), arising in any way out of the use of
82
this software, even if advised of the possibility of such damage.
83
84
5.    CONTACT INFORMATION
85
86
Fraunhofer Institute for Integrated Circuits IIS
87
Attention: Audio and Multimedia Departments - FDK AAC LL
88
Am Wolfsmantel 33
89
91058 Erlangen, Germany
90
91
www.iis.fraunhofer.de/amm
92
amm-info@iis.fraunhofer.de
93
----------------------------------------------------------------------------- */
94
95
/******************* Library for basic calculation routines ********************
96
97
   Author(s):   Josef Hoepfl, Manuel Jander, Youliy Ninov, Daniel Hagel
98
99
   Description: MDCT/MDST routines
100
101
*******************************************************************************/
102
103
#include "mdct.h"
104
105
#include "FDK_tools_rom.h"
106
#include "dct.h"
107
#include "fixpoint_math.h"
108
109
881k
void mdct_init(H_MDCT hMdct, FIXP_DBL *overlap, INT overlapBufferSize) {
110
881k
  hMdct->overlap.freq = overlap;
111
  // FDKmemclear(overlap, overlapBufferSize*sizeof(FIXP_DBL));
112
881k
  hMdct->prev_fr = 0;
113
881k
  hMdct->prev_nr = 0;
114
881k
  hMdct->prev_tl = 0;
115
881k
  hMdct->ov_size = overlapBufferSize;
116
881k
  hMdct->prevAliasSymmetry = 0;
117
881k
  hMdct->prevPrevAliasSymmetry = 0;
118
881k
  hMdct->pFacZir = NULL;
119
881k
  hMdct->pAsymOvlp = NULL;
120
881k
}
121
122
/*
123
This program implements the forward MDCT transform on an input block of data.
124
The input block is in a form (A,B,C,D) where A,B,C and D are the respective
125
1/4th segments of the block. The program takes the input block and folds it in
126
the form:
127
(-D-Cr,A-Br). This block is twice shorter and here the 'r' suffix denotes
128
flipping of the sequence (reversing the order of the samples). While folding the
129
input block in the above mentioned shorter block the program windows the data.
130
Because the two operations (windowing and folding) are not implemented
131
sequentially, but together the program's structure is not easy to understand.
132
Once the output (already windowed) block (-D-Cr,A-Br) is ready it is passed to
133
the DCT IV for processing.
134
*/
135
INT mdct_block(H_MDCT hMdct, const INT_PCM *RESTRICT timeData,
136
               const INT noInSamples, FIXP_DBL *RESTRICT mdctData,
137
               const INT nSpec, const INT tl, const FIXP_WTP *pRightWindowPart,
138
0
               const INT fr, SHORT *pMdctData_e) {
139
0
  int i, n;
140
  /* tl: transform length
141
     fl: left window slope length
142
     nl: left window slope offset
143
     fr: right window slope length
144
     nr: right window slope offset
145
     See FDK_tools/doc/intern/mdct.tex for more detail. */
146
0
  int fl, nl, nr;
147
0
  const FIXP_WTP *wls, *wrs;
148
149
0
  wrs = pRightWindowPart;
150
151
  /* Detect FRprevious / FL mismatches and override parameters accordingly */
152
0
  if (hMdct->prev_fr ==
153
0
      0) { /* At start just initialize and pass parameters as they are */
154
0
    hMdct->prev_fr = fr;
155
0
    hMdct->prev_wrs = wrs;
156
0
    hMdct->prev_tl = tl;
157
0
  }
158
159
  /* Derive NR */
160
0
  nr = (tl - fr) >> 1;
161
162
  /* Skip input samples if tl is smaller than block size */
163
0
  timeData += (noInSamples - tl) >> 1;
164
165
  /* windowing */
166
0
  for (n = 0; n < nSpec; n++) {
167
    /*
168
     * MDCT scale:
169
     * + 1: fMultDiv2() in windowing.
170
     * + 1: Because of factor 1/2 in Princen-Bradley compliant windowed TDAC.
171
     */
172
0
    INT mdctData_e = 1 + 1;
173
174
    /* Derive left parameters */
175
0
    wls = hMdct->prev_wrs;
176
0
    fl = hMdct->prev_fr;
177
0
    nl = (tl - fl) >> 1;
178
179
    /* Here we implement a simplified version of what happens after the this
180
    piece of code (see the comments below). We implement the folding of A and B
181
    segments to (A-Br) but A is zero, because in this part of the MDCT sequence
182
    the window coefficients with which A must be multiplied are zero.    */
183
0
    for (i = 0; i < nl; i++) {
184
#if SAMPLE_BITS == DFRACT_BITS /* SPC_BITS and DFRACT_BITS should be equal. */
185
      mdctData[(tl / 2) + i] = -((FIXP_DBL)timeData[tl - i - 1] >> (1));
186
#else
187
0
      mdctData[(tl / 2) + i] = -(FIXP_DBL)timeData[tl - i - 1]
188
0
                               << (DFRACT_BITS - SAMPLE_BITS - 1); /* 0(A)-Br */
189
0
#endif
190
0
    }
191
192
    /* Implements the folding and windowing of the left part of the sequence,
193
    that is segments A and B. The A segment is multiplied by the respective left
194
    window coefficient and placed in a temporary variable.
195
196
    tmp0 = fMultDiv2((FIXP_PCM)timeData[i+nl], pLeftWindowPart[i].v.im);
197
198
    After this the B segment taken in reverse order is multiplied by the left
199
    window and subtracted from the previously derived temporary variable, so
200
    that finally we implement the A-Br operation. This output is written to the
201
    right part of the MDCT output : (-D-Cr,A-Br).
202
203
    mdctData[(tl/2)+i+nl] = fMultSubDiv2(tmp0, (FIXP_PCM)timeData[tl-nl-i-1],
204
    pLeftWindowPart[i].v.re);//A*window-Br*window
205
206
    The (A-Br) data is written to the output buffer (mdctData) without being
207
    flipped.     */
208
0
    for (i = 0; i < fl / 2; i++) {
209
0
      FIXP_DBL tmp0;
210
0
      tmp0 = fMultDiv2((FIXP_PCM)timeData[i + nl], wls[i].v.im); /* a*window */
211
0
      mdctData[(tl / 2) + i + nl] =
212
0
          fMultSubDiv2(tmp0, (FIXP_PCM)timeData[tl - nl - i - 1],
213
0
                       wls[i].v.re); /* A*window-Br*window */
214
0
    }
215
216
    /* Right window slope offset */
217
    /* Here we implement a simplified version of what happens after the this
218
    piece of code (see the comments below). We implement the folding of C and D
219
    segments to (-D-Cr) but D is zero, because in this part of the MDCT sequence
220
    the window coefficients with which D must be multiplied are zero.    */
221
0
    for (i = 0; i < nr; i++) {
222
#if SAMPLE_BITS == \
223
    DFRACT_BITS /* This should be SPC_BITS instead of DFRACT_BITS. */
224
      mdctData[(tl / 2) - 1 - i] = -((FIXP_DBL)timeData[tl + i] >> (1));
225
#else
226
0
      mdctData[(tl / 2) - 1 - i] =
227
0
          -(FIXP_DBL)timeData[tl + i]
228
0
          << (DFRACT_BITS - SAMPLE_BITS - 1); /* -C flipped at placing */
229
0
#endif
230
0
    }
231
232
    /* Implements the folding and windowing of the right part of the sequence,
233
    that is, segments C and D. The C segment is multiplied by the respective
234
    right window coefficient and placed in a temporary variable.
235
236
    tmp1 = fMultDiv2((FIXP_PCM)timeData[tl+nr+i], pRightWindowPart[i].v.re);
237
238
    After this the D segment taken in reverse order is multiplied by the right
239
    window and added from the previously derived temporary variable, so that we
240
    get (C+Dr) operation. This output is negated to get (-C-Dr) and written to
241
    the left part of the MDCT output while being reversed (flipped) at the same
242
    time, so that from (-C-Dr) we get (-D-Cr)=> (-D-Cr,A-Br).
243
244
    mdctData[(tl/2)-nr-i-1] = -fMultAddDiv2(tmp1,
245
    (FIXP_PCM)timeData[(tl*2)-nr-i-1], pRightWindowPart[i].v.im);*/
246
0
    for (i = 0; i < fr / 2; i++) {
247
0
      FIXP_DBL tmp1;
248
0
      tmp1 = fMultDiv2((FIXP_PCM)timeData[tl + nr + i],
249
0
                       wrs[i].v.re); /* C*window */
250
0
      mdctData[(tl / 2) - nr - i - 1] =
251
0
          -fMultAddDiv2(tmp1, (FIXP_PCM)timeData[(tl * 2) - nr - i - 1],
252
0
                        wrs[i].v.im); /* -(C*window+Dr*window) and flip before
253
                                         placing -> -Cr - D */
254
0
    }
255
256
    /* We pass the shortened folded data (-D-Cr,A-Br) to the MDCT function */
257
0
    dct_IV(mdctData, tl, &mdctData_e);
258
259
0
    pMdctData_e[n] = (SHORT)mdctData_e;
260
261
0
    timeData += tl;
262
0
    mdctData += tl;
263
264
0
    hMdct->prev_wrs = wrs;
265
0
    hMdct->prev_fr = fr;
266
0
    hMdct->prev_tl = tl;
267
0
  }
268
269
0
  return nSpec * tl;
270
0
}
271
272
647k
void imdct_gain(FIXP_DBL *pGain_m, int *pGain_e, int tl) {
273
647k
  FIXP_DBL gain_m = *pGain_m;
274
647k
  int gain_e = *pGain_e;
275
647k
  int log2_tl;
276
277
647k
  gain_e += -MDCT_OUTPUT_GAIN - MDCT_OUT_HEADROOM + 1;
278
647k
  if (tl == 0) {
279
    /* Dont regard the 2/N factor from the IDCT. It is compensated for somewhere
280
     * else. */
281
23.6k
    *pGain_e = gain_e;
282
23.6k
    return;
283
23.6k
  }
284
285
623k
  log2_tl = DFRACT_BITS - 1 - fNormz((FIXP_DBL)tl);
286
623k
  gain_e += -log2_tl;
287
288
  /* Detect non-radix 2 transform length and add amplitude compensation factor
289
     which cannot be included into the exponent above */
290
623k
  switch ((tl) >> (log2_tl - 2)) {
291
164k
    case 0x7: /* 10 ms, 1/tl = 1.0/(FDKpow(2.0, -log2_tl) *
292
                 0.53333333333333333333) */
293
164k
      if (gain_m == (FIXP_DBL)0) {
294
164k
        gain_m = FL2FXCONST_DBL(0.53333333333333333333f);
295
164k
      } else {
296
0
        gain_m = fMult(gain_m, FL2FXCONST_DBL(0.53333333333333333333f));
297
0
      }
298
164k
      break;
299
177k
    case 0x6: /* 3/4 of radix 2, 1/tl = 1.0/(FDKpow(2.0, -log2_tl) * 2.0/3.0) */
300
177k
      if (gain_m == (FIXP_DBL)0) {
301
151k
        gain_m = FL2FXCONST_DBL(2.0 / 3.0f);
302
151k
      } else {
303
26.4k
        gain_m = fMult(gain_m, FL2FXCONST_DBL(2.0 / 3.0f));
304
26.4k
      }
305
177k
      break;
306
111
    case 0x5: /* 0.8 of radix 2 (e.g. tl 160), 1/tl = 1.0/(FDKpow(2.0, -log2_tl)
307
               * 0.8/1.5) */
308
111
      if (gain_m == (FIXP_DBL)0) {
309
111
        gain_m = FL2FXCONST_DBL(0.53333333333333333333f);
310
111
      } else {
311
0
        gain_m = fMult(gain_m, FL2FXCONST_DBL(0.53333333333333333333f));
312
0
      }
313
111
      break;
314
281k
    case 0x4:
315
      /* radix 2, nothing to do. */
316
281k
      break;
317
0
    default:
318
      /* unsupported */
319
0
      FDK_ASSERT(0);
320
0
      break;
321
623k
  }
322
323
623k
  *pGain_m = gain_m;
324
623k
  *pGain_e = gain_e;
325
623k
}
326
327
58.3k
INT imdct_drain(H_MDCT hMdct, FIXP_DBL *output, INT nrSamplesRoom) {
328
58.3k
  int buffered_samples = 0;
329
330
58.3k
  if (nrSamplesRoom > 0) {
331
25.1k
    buffered_samples = hMdct->ov_offset;
332
333
25.1k
    FDK_ASSERT(buffered_samples <= nrSamplesRoom);
334
335
25.1k
    if (buffered_samples > 0) {
336
3.34k
      FDKmemcpy(output, hMdct->overlap.time,
337
3.34k
                buffered_samples * sizeof(FIXP_DBL));
338
3.34k
      hMdct->ov_offset = 0;
339
3.34k
    }
340
25.1k
  }
341
58.3k
  return buffered_samples;
342
58.3k
}
343
344
44.2k
INT imdct_copy_ov_and_nr(H_MDCT hMdct, FIXP_DBL *pTimeData, INT nrSamples) {
345
44.2k
  FIXP_DBL *pOvl;
346
44.2k
  int nt, nf, i;
347
348
44.2k
  nt = fMin(hMdct->ov_offset, nrSamples);
349
44.2k
  nrSamples -= nt;
350
44.2k
  nf = fMin(hMdct->prev_nr, nrSamples);
351
44.2k
  FDKmemcpy(pTimeData, hMdct->overlap.time, nt * sizeof(FIXP_DBL));
352
44.2k
  pTimeData += nt;
353
354
44.2k
  pOvl = hMdct->overlap.freq + hMdct->ov_size - 1;
355
44.2k
  if (hMdct->prevPrevAliasSymmetry == 0) {
356
534k
    for (i = 0; i < nf; i++) {
357
489k
      FIXP_DBL x = -(*pOvl--);
358
489k
      *pTimeData = IMDCT_SCALE_DBL(x);
359
489k
      pTimeData++;
360
489k
    }
361
44.2k
  } else {
362
0
    for (i = 0; i < nf; i++) {
363
0
      FIXP_DBL x = (*pOvl--);
364
0
      *pTimeData = IMDCT_SCALE_DBL(x);
365
0
      pTimeData++;
366
0
    }
367
0
  }
368
369
44.2k
  return (nt + nf);
370
44.2k
}
371
372
void imdct_adapt_parameters(H_MDCT hMdct, int *pfl, int *pnl, int tl,
373
105k
                            const FIXP_WTP *wls, int noOutSamples) {
374
105k
  int fl = *pfl, nl = *pnl;
375
105k
  int window_diff, use_current = 0, use_previous = 0;
376
105k
  if (hMdct->prev_tl == 0) {
377
31.4k
    hMdct->prev_wrs = wls;
378
31.4k
    hMdct->prev_fr = fl;
379
31.4k
    hMdct->prev_nr = (noOutSamples - fl) >> 1;
380
31.4k
    hMdct->prev_tl = noOutSamples;
381
31.4k
    hMdct->ov_offset = 0;
382
31.4k
    use_current = 1;
383
31.4k
  }
384
385
105k
  window_diff = (hMdct->prev_fr - fl) >> 1;
386
387
  /* check if the previous window slope can be adjusted to match the current
388
   * window slope */
389
105k
  if (hMdct->prev_nr + window_diff > 0) {
390
47.0k
    use_current = 1;
391
47.0k
  }
392
  /* check if the current window slope can be adjusted to match the previous
393
   * window slope */
394
105k
  if (nl - window_diff > 0) {
395
48.6k
    use_previous = 1;
396
48.6k
  }
397
398
  /* if both is possible choose the larger of both window slope lengths */
399
105k
  if (use_current && use_previous) {
400
6.98k
    if (fl < hMdct->prev_fr) {
401
327
      use_current = 0;
402
327
    }
403
6.98k
  }
404
  /*
405
   * If the previous transform block is big enough, enlarge previous window
406
   * overlap, if not, then shrink current window overlap.
407
   */
408
105k
  if (use_current) {
409
63.8k
    hMdct->prev_nr += window_diff;
410
63.8k
    hMdct->prev_fr = fl;
411
63.8k
    hMdct->prev_wrs = wls;
412
63.8k
  } else {
413
42.0k
    nl -= window_diff;
414
42.0k
    fl = hMdct->prev_fr;
415
42.0k
  }
416
417
105k
  *pfl = fl;
418
105k
  *pnl = nl;
419
105k
}
420
421
/*
422
This program implements the inverse modulated lapped transform, a generalized
423
version of the inverse MDCT transform. Setting none of the MLT_*_ALIAS_FLAG
424
flags computes the IMDCT, setting all of them computes the IMDST. Other
425
combinations of these flags compute type III transforms used by the RSVD60
426
multichannel tool for transitions between MDCT/MDST. The following description
427
relates to the IMDCT only.
428
429
If we pass the data block (A,B,C,D,E,F) to the FORWARD MDCT it will produce two
430
outputs. The first one will be over the (A,B,C,D) part =>(-D-Cr,A-Br) and the
431
second one will be over the (C,D,E,F) part => (-F-Er,C-Dr), since there is a
432
overlap between consequtive passes of the algorithm. This overlap is over the
433
(C,D) segments. The two outputs will be given sequentially to the DCT IV
434
algorithm. At the INVERSE MDCT side we get two consecutive outputs from the IDCT
435
IV algorithm, namely the same blocks: (-D-Cr,A-Br) and (-F-Er,C-Dr). The first
436
of them lands in the Overlap buffer and the second is in the working one, which,
437
one algorithm pass later will substitute the one residing in the overlap
438
register. The IMDCT algorithm has to produce the C and D segments from the two
439
buffers. In order to do this we take the left part of the overlap
440
buffer(-D-Cr,A-Br), namely (-D-Cr) and add it appropriately to the right part of
441
the working buffer (-F-Er,C-Dr), namely (C-Dr), so that we get first the C
442
segment and later the D segment. We do this in the following way: From the right
443
part of the working buffer(C-Dr) we subtract the flipped left part of the
444
overlap buffer(-D-Cr):
445
446
Result = (C-Dr) - flipped(-D-Cr) = C -Dr + Dr + C = 2C
447
We divide by two and get the C segment. What we did is adding the right part of
448
the first frame to the left part of the second one.   While applying these
449
operation we multiply the respective segments with the appropriate window
450
functions.
451
452
In order to get the D segment we do the following:
453
From the negated second part of the working buffer(C-Dr) we subtract the flipped
454
first part of the overlap buffer (-D-Cr):
455
456
Result= - (C -Dr) - flipped(-D-Cr)= -C +Dr +Dr +C = 2Dr.
457
After dividing by two and flipping we get the D segment.What we did is adding
458
the right part of the first frame to the left part of the second one.   While
459
applying these operation we multiply the respective segments with the
460
appropriate window functions.
461
462
Once we have obtained the C and D segments the overlap buffer is emptied and the
463
current buffer is sent in it, so that the E and F segments are available for
464
decoding in the next algorithm pass.*/
465
INT imlt_block(H_MDCT hMdct, FIXP_DBL *output, FIXP_DBL *spectrum,
466
               const SHORT scalefactor[], const INT nSpec,
467
               const INT noOutSamples, const INT tl, const FIXP_WTP *wls,
468
               INT fl, const FIXP_WTP *wrs, const INT fr, FIXP_DBL gain,
469
337k
               int flags) {
470
337k
  FIXP_DBL *pOvl;
471
337k
  FIXP_DBL *pOut0 = output, *pOut1;
472
337k
  INT nl, nr;
473
337k
  int w, i, nrSamples = 0, specShiftScale, transform_gain_e = 0;
474
337k
  int currAliasSymmetry = (flags & MLT_FLAG_CURR_ALIAS_SYMMETRY);
475
476
  /* Derive NR and NL */
477
337k
  nr = (tl - fr) >> 1;
478
337k
  nl = (tl - fl) >> 1;
479
480
  /* Include 2/N IMDCT gain into gain factor and exponent. */
481
337k
  imdct_gain(&gain, &transform_gain_e, tl);
482
483
  /* Detect FRprevious / FL mismatches and override parameters accordingly */
484
337k
  if (hMdct->prev_fr != fl) {
485
93.0k
    imdct_adapt_parameters(hMdct, &fl, &nl, tl, wls, noOutSamples);
486
93.0k
  }
487
488
337k
  pOvl = hMdct->overlap.freq + hMdct->ov_size - 1;
489
490
337k
  if (noOutSamples > nrSamples) {
491
    /* Purge buffered output. */
492
39.5M
    for (i = 0; i < hMdct->ov_offset; i++) {
493
39.2M
      *pOut0 = hMdct->overlap.time[i];
494
39.2M
      pOut0++;
495
39.2M
    }
496
335k
    nrSamples = hMdct->ov_offset;
497
335k
    hMdct->ov_offset = 0;
498
335k
  }
499
500
1.31M
  for (w = 0; w < nSpec; w++) {
501
973k
    FIXP_DBL *pSpec, *pCurr;
502
973k
    const FIXP_WTP *pWindow;
503
504
    /* Detect FRprevious / FL mismatches and override parameters accordingly */
505
973k
    if (hMdct->prev_fr != fl) {
506
0
      imdct_adapt_parameters(hMdct, &fl, &nl, tl, wls, noOutSamples);
507
0
    }
508
509
973k
    specShiftScale = transform_gain_e;
510
511
    /* Setup window pointers */
512
973k
    pWindow = hMdct->prev_wrs;
513
514
    /* Current spectrum */
515
973k
    pSpec = spectrum + w * tl;
516
517
    /* DCT IV of current spectrum. */
518
973k
    if (currAliasSymmetry == 0) {
519
973k
      if (hMdct->prevAliasSymmetry == 0) {
520
973k
        dct_IV(pSpec, tl, &specShiftScale);
521
973k
      } else {
522
0
        FIXP_DBL _tmp[1024 + ALIGNMENT_DEFAULT / sizeof(FIXP_DBL)];
523
0
        FIXP_DBL *tmp = (FIXP_DBL *)ALIGN_PTR(_tmp);
524
0
        C_ALLOC_ALIGNED_REGISTER(tmp, sizeof(_tmp));
525
0
        dct_III(pSpec, tmp, tl, &specShiftScale);
526
0
        C_ALLOC_ALIGNED_UNREGISTER(tmp);
527
0
      }
528
973k
    } else {
529
0
      if (hMdct->prevAliasSymmetry == 0) {
530
0
        FIXP_DBL _tmp[1024 + ALIGNMENT_DEFAULT / sizeof(FIXP_DBL)];
531
0
        FIXP_DBL *tmp = (FIXP_DBL *)ALIGN_PTR(_tmp);
532
0
        C_ALLOC_ALIGNED_REGISTER(tmp, sizeof(_tmp));
533
0
        dst_III(pSpec, tmp, tl, &specShiftScale);
534
0
        C_ALLOC_ALIGNED_UNREGISTER(tmp);
535
0
      } else {
536
0
        dst_IV(pSpec, tl, &specShiftScale);
537
0
      }
538
0
    }
539
540
    /* Optional scaling of time domain - no yet windowed - of current spectrum
541
     */
542
    /* and de-scale current spectrum signal (time domain, no yet windowed) */
543
973k
    if (gain != (FIXP_DBL)0) {
544
130M
      for (i = 0; i < tl; i++) {
545
129M
        pSpec[i] = fMult(pSpec[i], gain);
546
129M
      }
547
586k
    }
548
549
973k
    {
550
973k
      int loc_scale =
551
973k
          fixmin_I(scalefactor[w] + specShiftScale, (INT)DFRACT_BITS - 1);
552
973k
      DWORD_ALIGNED(pSpec);
553
973k
      scaleValuesSaturate(pSpec, tl, loc_scale);
554
973k
    }
555
556
973k
    if (noOutSamples <= nrSamples) {
557
      /* Divert output first half to overlap buffer if we already got enough
558
       * output samples. */
559
274k
      pOut0 = hMdct->overlap.time + hMdct->ov_offset;
560
274k
      hMdct->ov_offset += hMdct->prev_nr + fl / 2;
561
698k
    } else {
562
      /* Account output samples */
563
698k
      nrSamples += hMdct->prev_nr + fl / 2;
564
698k
    }
565
566
    /* NR output samples 0 .. NR. -overlap[TL/2..TL/2-NR] */
567
973k
    if ((hMdct->pFacZir != 0) && (hMdct->prev_nr == fl / 2)) {
568
      /* In the case of ACELP -> TCX20 -> FD short add FAC ZIR on nr signal part
569
       */
570
35.2k
      for (i = 0; i < hMdct->prev_nr; i++) {
571
34.6k
        FIXP_DBL x = -(*pOvl--);
572
34.6k
        *pOut0 = fAddSaturate(x, IMDCT_SCALE_DBL(hMdct->pFacZir[i]));
573
34.6k
        pOut0++;
574
34.6k
      }
575
635
      hMdct->pFacZir = NULL;
576
972k
    } else {
577
      /* Here we implement a simplified version of what happens after the this
578
      piece of code (see the comments below). We implement the folding of C and
579
      D segments from (-D-Cr) but D is zero, because in this part of the MDCT
580
      sequence the window coefficients with which D must be multiplied are zero.
581
      "pOut0" writes sequentially the C block from left to right.   */
582
972k
      if (hMdct->prevPrevAliasSymmetry == 0) {
583
27.5M
        for (i = 0; i < hMdct->prev_nr; i++) {
584
26.5M
          FIXP_DBL x = -(*pOvl--);
585
26.5M
          *pOut0 = IMDCT_SCALE_DBL(x);
586
26.5M
          pOut0++;
587
26.5M
        }
588
972k
      } else {
589
0
        for (i = 0; i < hMdct->prev_nr; i++) {
590
0
          FIXP_DBL x = *pOvl--;
591
0
          *pOut0 = IMDCT_SCALE_DBL(x);
592
0
          pOut0++;
593
0
        }
594
0
      }
595
972k
    }
596
597
973k
    if (noOutSamples <= nrSamples) {
598
      /* Divert output second half to overlap buffer if we already got enough
599
       * output samples. */
600
375k
      pOut1 = hMdct->overlap.time + hMdct->ov_offset + fl / 2 - 1;
601
375k
      hMdct->ov_offset += fl / 2 + nl;
602
598k
    } else {
603
598k
      pOut1 = pOut0 + (fl - 1);
604
598k
      nrSamples += fl / 2 + nl;
605
598k
    }
606
607
    /* output samples before window crossing point NR .. TL/2.
608
     * -overlap[TL/2-NR..TL/2-NR-FL/2] + current[NR..TL/2] */
609
    /* output samples after window crossing point TL/2 .. TL/2+FL/2.
610
     * -overlap[0..FL/2] - current[TL/2..FL/2] */
611
973k
    pCurr = pSpec + tl - fl / 2;
612
973k
    DWORD_ALIGNED(pCurr);
613
973k
    C_ALLOC_ALIGNED_REGISTER(pWindow, fl);
614
973k
    DWORD_ALIGNED(pWindow);
615
973k
    C_ALLOC_ALIGNED_UNREGISTER(pWindow);
616
617
973k
    if (hMdct->prevPrevAliasSymmetry == 0) {
618
973k
      if (hMdct->prevAliasSymmetry == 0) {
619
973k
        if (!hMdct->pAsymOvlp) {
620
123M
          for (i = 0; i < fl / 2; i++) {
621
122M
            FIXP_DBL x0, x1;
622
122M
            cplxMultDiv2(&x1, &x0, *pCurr++, -*pOvl--, pWindow[i]);
623
122M
            *pOut0 = IMDCT_SCALE_DBL_LSH1(x0);
624
122M
            *pOut1 = IMDCT_SCALE_DBL_LSH1(-x1);
625
122M
            pOut0++;
626
122M
            pOut1--;
627
122M
          }
628
973k
        } else {
629
0
          FIXP_DBL *pAsymOvl = hMdct->pAsymOvlp + fl / 2 - 1;
630
0
          for (i = 0; i < fl / 2; i++) {
631
0
            FIXP_DBL x0, x1;
632
0
            x1 = -fMultDiv2(*pCurr, pWindow[i].v.re) +
633
0
                 fMultDiv2(*pAsymOvl, pWindow[i].v.im);
634
0
            x0 = fMultDiv2(*pCurr, pWindow[i].v.im) -
635
0
                 fMultDiv2(*pOvl, pWindow[i].v.re);
636
0
            pCurr++;
637
0
            pOvl--;
638
0
            pAsymOvl--;
639
0
            *pOut0++ = IMDCT_SCALE_DBL_LSH1(x0);
640
0
            *pOut1-- = IMDCT_SCALE_DBL_LSH1(x1);
641
0
          }
642
0
          hMdct->pAsymOvlp = NULL;
643
0
        }
644
973k
      } else { /* prevAliasingSymmetry == 1 */
645
0
        for (i = 0; i < fl / 2; i++) {
646
0
          FIXP_DBL x0, x1;
647
0
          cplxMultDiv2(&x1, &x0, *pCurr++, -*pOvl--, pWindow[i]);
648
0
          *pOut0 = IMDCT_SCALE_DBL_LSH1(x0);
649
0
          *pOut1 = IMDCT_SCALE_DBL_LSH1(x1);
650
0
          pOut0++;
651
0
          pOut1--;
652
0
        }
653
0
      }
654
973k
    } else { /* prevPrevAliasingSymmetry == 1 */
655
0
      if (hMdct->prevAliasSymmetry == 0) {
656
0
        for (i = 0; i < fl / 2; i++) {
657
0
          FIXP_DBL x0, x1;
658
0
          cplxMultDiv2(&x1, &x0, *pCurr++, *pOvl--, pWindow[i]);
659
0
          *pOut0 = IMDCT_SCALE_DBL_LSH1(x0);
660
0
          *pOut1 = IMDCT_SCALE_DBL_LSH1(-x1);
661
0
          pOut0++;
662
0
          pOut1--;
663
0
        }
664
0
      } else { /* prevAliasingSymmetry == 1 */
665
0
        for (i = 0; i < fl / 2; i++) {
666
0
          FIXP_DBL x0, x1;
667
0
          cplxMultDiv2(&x1, &x0, *pCurr++, *pOvl--, pWindow[i]);
668
0
          *pOut0 = IMDCT_SCALE_DBL_LSH1(x0);
669
0
          *pOut1 = IMDCT_SCALE_DBL_LSH1(x1);
670
0
          pOut0++;
671
0
          pOut1--;
672
0
        }
673
0
      }
674
0
    }
675
676
973k
    if (hMdct->pFacZir != 0) {
677
      /* add FAC ZIR of previous ACELP -> mdct transition */
678
3.03k
      FIXP_DBL *pOut = pOut0 - fl / 2;
679
3.03k
      FDK_ASSERT(fl / 2 <= 128);
680
364k
      for (i = 0; i < fl / 2; i++) {
681
361k
        pOut[i] = fAddSaturate(pOut[i], IMDCT_SCALE_DBL(hMdct->pFacZir[i]));
682
361k
      }
683
3.03k
      hMdct->pFacZir = NULL;
684
3.03k
    }
685
973k
    pOut0 += (fl / 2) + nl;
686
687
    /* NL output samples TL/2+FL/2..TL. - current[FL/2..0] */
688
973k
    pOut1 += (fl / 2) + 1;
689
973k
    pCurr = pSpec + tl - fl / 2 - 1;
690
    /* Here we implement a simplified version of what happens above the this
691
    piece of code (see the comments above). We implement the folding of C and D
692
    segments from (C-Dr) but C is zero, because in this part of the MDCT
693
    sequence the window coefficients with which C must be multiplied are zero.
694
    "pOut1" writes sequentially the D block from left to right.   */
695
973k
    if (hMdct->prevAliasSymmetry == 0) {
696
26.4M
      for (i = 0; i < nl; i++) {
697
25.4M
        FIXP_DBL x = -(*pCurr--);
698
25.4M
        *pOut1++ = IMDCT_SCALE_DBL(x);
699
25.4M
      }
700
973k
    } else {
701
0
      for (i = 0; i < nl; i++) {
702
0
        FIXP_DBL x = *pCurr--;
703
0
        *pOut1++ = IMDCT_SCALE_DBL(x);
704
0
      }
705
0
    }
706
707
    /* Set overlap source pointer for next window pOvl = pSpec + tl/2 - 1; */
708
973k
    pOvl = pSpec + tl / 2 - 1;
709
710
    /* Previous window values. */
711
973k
    hMdct->prev_nr = nr;
712
973k
    hMdct->prev_fr = fr;
713
973k
    hMdct->prev_tl = tl;
714
973k
    hMdct->prev_wrs = wrs;
715
716
    /* Previous aliasing symmetry */
717
973k
    hMdct->prevPrevAliasSymmetry = hMdct->prevAliasSymmetry;
718
973k
    hMdct->prevAliasSymmetry = currAliasSymmetry;
719
973k
  }
720
721
  /* Save overlap */
722
723
337k
  pOvl = hMdct->overlap.freq + hMdct->ov_size - tl / 2;
724
337k
  FDKmemcpy(pOvl, &spectrum[(nSpec - 1) * tl], (tl / 2) * sizeof(FIXP_DBL));
725
726
337k
  return nrSamples;
727
337k
}