Coverage Report

Created: 2024-09-06 07:53

/src/fdk-aac/libFDK/src/mdct.cpp
Line
Count
Source (jump to first uncovered line)
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
0
void mdct_init(H_MDCT hMdct, FIXP_DBL *overlap, INT overlapBufferSize) {
110
0
  hMdct->overlap.freq = overlap;
111
  // FDKmemclear(overlap, overlapBufferSize*sizeof(FIXP_DBL));
112
0
  hMdct->prev_fr = 0;
113
0
  hMdct->prev_nr = 0;
114
0
  hMdct->prev_tl = 0;
115
0
  hMdct->ov_size = overlapBufferSize;
116
0
  hMdct->prevAliasSymmetry = 0;
117
0
  hMdct->prevPrevAliasSymmetry = 0;
118
0
  hMdct->pFacZir = NULL;
119
0
  hMdct->pAsymOvlp = NULL;
120
0
}
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
0
void imdct_gain(FIXP_DBL *pGain_m, int *pGain_e, int tl) {
273
0
  FIXP_DBL gain_m = *pGain_m;
274
0
  int gain_e = *pGain_e;
275
0
  int log2_tl;
276
277
0
  gain_e += -MDCT_OUTPUT_GAIN - MDCT_OUT_HEADROOM + 1;
278
0
  if (tl == 0) {
279
    /* Dont regard the 2/N factor from the IDCT. It is compensated for somewhere
280
     * else. */
281
0
    *pGain_e = gain_e;
282
0
    return;
283
0
  }
284
285
0
  log2_tl = DFRACT_BITS - 1 - fNormz((FIXP_DBL)tl);
286
0
  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
0
  switch ((tl) >> (log2_tl - 2)) {
291
0
    case 0x7: /* 10 ms, 1/tl = 1.0/(FDKpow(2.0, -log2_tl) *
292
                 0.53333333333333333333) */
293
0
      if (gain_m == (FIXP_DBL)0) {
294
0
        gain_m = FL2FXCONST_DBL(0.53333333333333333333f);
295
0
      } else {
296
0
        gain_m = fMult(gain_m, FL2FXCONST_DBL(0.53333333333333333333f));
297
0
      }
298
0
      break;
299
0
    case 0x6: /* 3/4 of radix 2, 1/tl = 1.0/(FDKpow(2.0, -log2_tl) * 2.0/3.0) */
300
0
      if (gain_m == (FIXP_DBL)0) {
301
0
        gain_m = FL2FXCONST_DBL(2.0 / 3.0f);
302
0
      } else {
303
0
        gain_m = fMult(gain_m, FL2FXCONST_DBL(2.0 / 3.0f));
304
0
      }
305
0
      break;
306
0
    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
0
      if (gain_m == (FIXP_DBL)0) {
309
0
        gain_m = FL2FXCONST_DBL(0.53333333333333333333f);
310
0
      } else {
311
0
        gain_m = fMult(gain_m, FL2FXCONST_DBL(0.53333333333333333333f));
312
0
      }
313
0
      break;
314
0
    case 0x4:
315
      /* radix 2, nothing to do. */
316
0
      break;
317
0
    default:
318
      /* unsupported */
319
0
      FDK_ASSERT(0);
320
0
      break;
321
0
  }
322
323
0
  *pGain_m = gain_m;
324
0
  *pGain_e = gain_e;
325
0
}
326
327
0
INT imdct_drain(H_MDCT hMdct, FIXP_DBL *output, INT nrSamplesRoom) {
328
0
  int buffered_samples = 0;
329
330
0
  if (nrSamplesRoom > 0) {
331
0
    buffered_samples = hMdct->ov_offset;
332
333
0
    FDK_ASSERT(buffered_samples <= nrSamplesRoom);
334
335
0
    if (buffered_samples > 0) {
336
0
      FDKmemcpy(output, hMdct->overlap.time,
337
0
                buffered_samples * sizeof(FIXP_DBL));
338
0
      hMdct->ov_offset = 0;
339
0
    }
340
0
  }
341
0
  return buffered_samples;
342
0
}
343
344
0
INT imdct_copy_ov_and_nr(H_MDCT hMdct, FIXP_DBL *pTimeData, INT nrSamples) {
345
0
  FIXP_DBL *pOvl;
346
0
  int nt, nf, i;
347
348
0
  nt = fMin(hMdct->ov_offset, nrSamples);
349
0
  nrSamples -= nt;
350
0
  nf = fMin(hMdct->prev_nr, nrSamples);
351
0
  FDKmemcpy(pTimeData, hMdct->overlap.time, nt * sizeof(FIXP_DBL));
352
0
  pTimeData += nt;
353
354
0
  pOvl = hMdct->overlap.freq + hMdct->ov_size - 1;
355
0
  if (hMdct->prevPrevAliasSymmetry == 0) {
356
0
    for (i = 0; i < nf; i++) {
357
0
      FIXP_DBL x = -(*pOvl--);
358
0
      *pTimeData = IMDCT_SCALE_DBL(x);
359
0
      pTimeData++;
360
0
    }
361
0
  } 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
0
  return (nt + nf);
370
0
}
371
372
void imdct_adapt_parameters(H_MDCT hMdct, int *pfl, int *pnl, int tl,
373
0
                            const FIXP_WTP *wls, int noOutSamples) {
374
0
  int fl = *pfl, nl = *pnl;
375
0
  int window_diff, use_current = 0, use_previous = 0;
376
0
  if (hMdct->prev_tl == 0) {
377
0
    hMdct->prev_wrs = wls;
378
0
    hMdct->prev_fr = fl;
379
0
    hMdct->prev_nr = (noOutSamples - fl) >> 1;
380
0
    hMdct->prev_tl = noOutSamples;
381
0
    hMdct->ov_offset = 0;
382
0
    use_current = 1;
383
0
  }
384
385
0
  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
0
  if (hMdct->prev_nr + window_diff > 0) {
390
0
    use_current = 1;
391
0
  }
392
  /* check if the current window slope can be adjusted to match the previous
393
   * window slope */
394
0
  if (nl - window_diff > 0) {
395
0
    use_previous = 1;
396
0
  }
397
398
  /* if both is possible choose the larger of both window slope lengths */
399
0
  if (use_current && use_previous) {
400
0
    if (fl < hMdct->prev_fr) {
401
0
      use_current = 0;
402
0
    }
403
0
  }
404
  /*
405
   * If the previous transform block is big enough, enlarge previous window
406
   * overlap, if not, then shrink current window overlap.
407
   */
408
0
  if (use_current) {
409
0
    hMdct->prev_nr += window_diff;
410
0
    hMdct->prev_fr = fl;
411
0
    hMdct->prev_wrs = wls;
412
0
  } else {
413
0
    nl -= window_diff;
414
0
    fl = hMdct->prev_fr;
415
0
  }
416
417
0
  *pfl = fl;
418
0
  *pnl = nl;
419
0
}
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
0
               int flags) {
470
0
  FIXP_DBL *pOvl;
471
0
  FIXP_DBL *pOut0 = output, *pOut1;
472
0
  INT nl, nr;
473
0
  int w, i, nrSamples = 0, specShiftScale, transform_gain_e = 0;
474
0
  int currAliasSymmetry = (flags & MLT_FLAG_CURR_ALIAS_SYMMETRY);
475
476
  /* Derive NR and NL */
477
0
  nr = (tl - fr) >> 1;
478
0
  nl = (tl - fl) >> 1;
479
480
  /* Include 2/N IMDCT gain into gain factor and exponent. */
481
0
  imdct_gain(&gain, &transform_gain_e, tl);
482
483
  /* Detect FRprevious / FL mismatches and override parameters accordingly */
484
0
  if (hMdct->prev_fr != fl) {
485
0
    imdct_adapt_parameters(hMdct, &fl, &nl, tl, wls, noOutSamples);
486
0
  }
487
488
0
  pOvl = hMdct->overlap.freq + hMdct->ov_size - 1;
489
490
0
  if (noOutSamples > nrSamples) {
491
    /* Purge buffered output. */
492
0
    for (i = 0; i < hMdct->ov_offset; i++) {
493
0
      *pOut0 = hMdct->overlap.time[i];
494
0
      pOut0++;
495
0
    }
496
0
    nrSamples = hMdct->ov_offset;
497
0
    hMdct->ov_offset = 0;
498
0
  }
499
500
0
  for (w = 0; w < nSpec; w++) {
501
0
    FIXP_DBL *pSpec, *pCurr;
502
0
    const FIXP_WTP *pWindow;
503
504
    /* Detect FRprevious / FL mismatches and override parameters accordingly */
505
0
    if (hMdct->prev_fr != fl) {
506
0
      imdct_adapt_parameters(hMdct, &fl, &nl, tl, wls, noOutSamples);
507
0
    }
508
509
0
    specShiftScale = transform_gain_e;
510
511
    /* Setup window pointers */
512
0
    pWindow = hMdct->prev_wrs;
513
514
    /* Current spectrum */
515
0
    pSpec = spectrum + w * tl;
516
517
    /* DCT IV of current spectrum. */
518
0
    if (currAliasSymmetry == 0) {
519
0
      if (hMdct->prevAliasSymmetry == 0) {
520
0
        dct_IV(pSpec, tl, &specShiftScale);
521
0
      } 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
0
    } 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
0
    if (gain != (FIXP_DBL)0) {
544
0
      for (i = 0; i < tl; i++) {
545
0
        pSpec[i] = fMult(pSpec[i], gain);
546
0
      }
547
0
    }
548
549
0
    {
550
0
      int loc_scale =
551
0
          fixmin_I(scalefactor[w] + specShiftScale, (INT)DFRACT_BITS - 1);
552
0
      DWORD_ALIGNED(pSpec);
553
0
      scaleValuesSaturate(pSpec, tl, loc_scale);
554
0
    }
555
556
0
    if (noOutSamples <= nrSamples) {
557
      /* Divert output first half to overlap buffer if we already got enough
558
       * output samples. */
559
0
      pOut0 = hMdct->overlap.time + hMdct->ov_offset;
560
0
      hMdct->ov_offset += hMdct->prev_nr + fl / 2;
561
0
    } else {
562
      /* Account output samples */
563
0
      nrSamples += hMdct->prev_nr + fl / 2;
564
0
    }
565
566
    /* NR output samples 0 .. NR. -overlap[TL/2..TL/2-NR] */
567
0
    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
0
      for (i = 0; i < hMdct->prev_nr; i++) {
571
0
        FIXP_DBL x = -(*pOvl--);
572
0
        *pOut0 = fAddSaturate(x, IMDCT_SCALE_DBL(hMdct->pFacZir[i]));
573
0
        pOut0++;
574
0
      }
575
0
      hMdct->pFacZir = NULL;
576
0
    } 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
0
      if (hMdct->prevPrevAliasSymmetry == 0) {
583
0
        for (i = 0; i < hMdct->prev_nr; i++) {
584
0
          FIXP_DBL x = -(*pOvl--);
585
0
          *pOut0 = IMDCT_SCALE_DBL(x);
586
0
          pOut0++;
587
0
        }
588
0
      } 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
0
    }
596
597
0
    if (noOutSamples <= nrSamples) {
598
      /* Divert output second half to overlap buffer if we already got enough
599
       * output samples. */
600
0
      pOut1 = hMdct->overlap.time + hMdct->ov_offset + fl / 2 - 1;
601
0
      hMdct->ov_offset += fl / 2 + nl;
602
0
    } else {
603
0
      pOut1 = pOut0 + (fl - 1);
604
0
      nrSamples += fl / 2 + nl;
605
0
    }
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
0
    pCurr = pSpec + tl - fl / 2;
612
0
    DWORD_ALIGNED(pCurr);
613
0
    C_ALLOC_ALIGNED_REGISTER(pWindow, fl);
614
0
    DWORD_ALIGNED(pWindow);
615
0
    C_ALLOC_ALIGNED_UNREGISTER(pWindow);
616
617
0
    if (hMdct->prevPrevAliasSymmetry == 0) {
618
0
      if (hMdct->prevAliasSymmetry == 0) {
619
0
        if (!hMdct->pAsymOvlp) {
620
0
          for (i = 0; i < fl / 2; i++) {
621
0
            FIXP_DBL x0, x1;
622
0
            cplxMultDiv2(&x1, &x0, *pCurr++, -*pOvl--, pWindow[i]);
623
0
            *pOut0 = IMDCT_SCALE_DBL_LSH1(x0);
624
0
            *pOut1 = IMDCT_SCALE_DBL_LSH1(-x1);
625
0
            pOut0++;
626
0
            pOut1--;
627
0
          }
628
0
        } 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
0
      } 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
0
    } 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
0
    if (hMdct->pFacZir != 0) {
677
      /* add FAC ZIR of previous ACELP -> mdct transition */
678
0
      FIXP_DBL *pOut = pOut0 - fl / 2;
679
0
      FDK_ASSERT(fl / 2 <= 128);
680
0
      for (i = 0; i < fl / 2; i++) {
681
0
        pOut[i] = fAddSaturate(pOut[i], IMDCT_SCALE_DBL(hMdct->pFacZir[i]));
682
0
      }
683
0
      hMdct->pFacZir = NULL;
684
0
    }
685
0
    pOut0 += (fl / 2) + nl;
686
687
    /* NL output samples TL/2+FL/2..TL. - current[FL/2..0] */
688
0
    pOut1 += (fl / 2) + 1;
689
0
    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
0
    if (hMdct->prevAliasSymmetry == 0) {
696
0
      for (i = 0; i < nl; i++) {
697
0
        FIXP_DBL x = -(*pCurr--);
698
0
        *pOut1++ = IMDCT_SCALE_DBL(x);
699
0
      }
700
0
    } 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
0
    pOvl = pSpec + tl / 2 - 1;
709
710
    /* Previous window values. */
711
0
    hMdct->prev_nr = nr;
712
0
    hMdct->prev_fr = fr;
713
0
    hMdct->prev_tl = tl;
714
0
    hMdct->prev_wrs = wrs;
715
716
    /* Previous aliasing symmetry */
717
0
    hMdct->prevPrevAliasSymmetry = hMdct->prevAliasSymmetry;
718
0
    hMdct->prevAliasSymmetry = currAliasSymmetry;
719
0
  }
720
721
  /* Save overlap */
722
723
0
  pOvl = hMdct->overlap.freq + hMdct->ov_size - tl / 2;
724
0
  FDKmemcpy(pOvl, &spectrum[(nSpec - 1) * tl], (tl / 2) * sizeof(FIXP_DBL));
725
726
0
  return nrSamples;
727
0
}