1 
2 /* -----------------------------------------------------------------------------------------------------------
3 Software License for The Fraunhofer FDK AAC Codec Library for Android
4 
5 � Copyright  1995 - 2015 Fraunhofer-Gesellschaft zur F�rderung der angewandten Forschung e.V.
6   All rights reserved.
7 
8  1.    INTRODUCTION
9 The Fraunhofer FDK AAC Codec Library for Android ("FDK AAC Codec") is software that implements
10 the MPEG Advanced Audio Coding ("AAC") encoding and decoding scheme for digital audio.
11 This FDK AAC Codec software is intended to be used on a wide variety of Android devices.
12 
13 AAC's HE-AAC and HE-AAC v2 versions are regarded as today's most efficient general perceptual
14 audio codecs. AAC-ELD is considered the best-performing full-bandwidth communications codec by
15 independent studies and is widely deployed. AAC has been standardized by ISO and IEC as part
16 of the MPEG specifications.
17 
18 Patent licenses for necessary patent claims for the FDK AAC Codec (including those of Fraunhofer)
19 may be obtained through Via Licensing (www.vialicensing.com) or through the respective patent owners
20 individually for the purpose of encoding or decoding bit streams in products that are compliant with
21 the ISO/IEC MPEG audio standards. Please note that most manufacturers of Android devices already license
22 these patent claims through Via Licensing or directly from the patent owners, and therefore FDK AAC Codec
23 software may already be covered under those patent licenses when it is used for those licensed purposes only.
24 
25 Commercially-licensed AAC software libraries, including floating-point versions with enhanced sound quality,
26 are also available from Fraunhofer. Users are encouraged to check the Fraunhofer website for additional
27 applications information and documentation.
28 
29 2.    COPYRIGHT LICENSE
30 
31 Redistribution and use in source and binary forms, with or without modification, are permitted without
32 payment of copyright license fees provided that you satisfy the following conditions:
33 
34 You must retain the complete text of this software license in redistributions of the FDK AAC Codec or
35 your modifications thereto in source code form.
36 
37 You must retain the complete text of this software license in the documentation and/or other materials
38 provided with redistributions of the FDK AAC Codec or your modifications thereto in binary form.
39 You must make available free of charge copies of the complete source code of the FDK AAC Codec and your
40 modifications thereto to recipients of copies in binary form.
41 
42 The name of Fraunhofer may not be used to endorse or promote products derived from this library without
43 prior written permission.
44 
45 You may not charge copyright license fees for anyone to use, copy or distribute the FDK AAC Codec
46 software or your modifications thereto.
47 
48 Your modified versions of the FDK AAC Codec must carry prominent notices stating that you changed the software
49 and the date of any change. For modified versions of the FDK AAC Codec, the term
50 "Fraunhofer FDK AAC Codec Library for Android" must be replaced by the term
51 "Third-Party Modified Version of the Fraunhofer FDK AAC Codec Library for Android."
52 
53 3.    NO PATENT LICENSE
54 
55 NO EXPRESS OR IMPLIED LICENSES TO ANY PATENT CLAIMS, including without limitation the patents of Fraunhofer,
56 ARE GRANTED BY THIS SOFTWARE LICENSE. Fraunhofer provides no warranty of patent non-infringement with
57 respect to this software.
58 
59 You may use this FDK AAC Codec software or modifications thereto only for purposes that are authorized
60 by appropriate patent licenses.
61 
62 4.    DISCLAIMER
63 
64 This FDK AAC Codec software is provided by Fraunhofer on behalf of the copyright holders and contributors
65 "AS IS" and WITHOUT ANY EXPRESS OR IMPLIED WARRANTIES, including but not limited to the implied warranties
66 of merchantability and fitness for a particular purpose. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR
67 CONTRIBUTORS BE LIABLE for any direct, indirect, incidental, special, exemplary, or consequential damages,
68 including but not limited to procurement of substitute goods or services; loss of use, data, or profits,
69 or business interruption, however caused and on any theory of liability, whether in contract, strict
70 liability, or tort (including negligence), arising in any way out of the use of this software, even if
71 advised of the possibility of such damage.
72 
73 5.    CONTACT INFORMATION
74 
75 Fraunhofer Institute for Integrated Circuits IIS
76 Attention: Audio and Multimedia Departments - FDK AAC LL
77 Am Wolfsmantel 33
78 91058 Erlangen, Germany
79 
80 www.iis.fraunhofer.de/amm
81 amm-info@iis.fraunhofer.de
82 ----------------------------------------------------------------------------------------------------------- */
83 
84 /******************************** MPEG Audio Encoder **************************
85 
86    Initial author:       A. Horndasch (code originally from lwr) / Josef Hoepfl (FDK)
87    contents/description: intensity stereo processing
88 
89 ******************************************************************************/
90 
91 #include "intensity.h"
92 #include "interface.h"
93 #include "psy_configuration.h"
94 #include "psy_const.h"
95 #include "qc_main.h"
96 #include "bit_cnt.h"
97 
98 /* only set an IS seed it left/right channel correlation is above IS_CORR_THRESH */
99 #define IS_CORR_THRESH                FL2FXCONST_DBL(0.95f)
100 
101 /* when expanding the IS region to more SFBs only accept an error that is
102  * not more than IS_TOTAL_ERROR_THRESH overall and
103  * not more than IS_LOCAL_ERROR_THRESH for the current SFB */
104 #define IS_TOTAL_ERROR_THRESH         FL2FXCONST_DBL(0.04f)
105 #define IS_LOCAL_ERROR_THRESH         FL2FXCONST_DBL(0.01f)
106 
107 /* the maximum allowed change of the intensity direction (unit: IS scale) - scaled with factor 0.25 - */
108 #define IS_DIRECTION_DEVIATION_THRESH_SF 2
109 #define IS_DIRECTION_DEVIATION_THRESH FL2FXCONST_DBL(2.0f/(1<<IS_DIRECTION_DEVIATION_THRESH_SF))
110 
111 /* IS regions need to have a minimal percentage of the overall loudness, e.g. 0.06 == 6% */
112 #define IS_REGION_MIN_LOUDNESS        FL2FXCONST_DBL(0.1f)
113 
114 /* only perform IS if IS_MIN_SFBS neighboring SFBs can be processed */
115 #define IS_MIN_SFBS                   6
116 
117 /* only do IS if
118  * if IS_LEFT_RIGHT_RATIO_THRESH < sfbEnergyLeft[sfb]/sfbEnergyRight[sfb] < 1 / IS_LEFT_RIGHT_RATIO_THRESH
119  * -> no IS if the panning angle is not far from the middle, MS will do */
120 /* this is equivalent to a scale of +/-1.02914634566 */
121 #define IS_LEFT_RIGHT_RATIO_THRESH    FL2FXCONST_DBL(0.7f)
122 
123 /* scalefactor of realScale */
124 #define REAL_SCALE_SF                    1
125 
126 /* scalefactor overallLoudness */
127 #define OVERALL_LOUDNESS_SF              6
128 
129 /* scalefactor for sum over max samples per goup */
130 #define MAX_SFB_PER_GROUP_SF             6
131 
132 /* scalefactor for sum of mdct spectrum */
133 #define MDCT_SPEC_SF                     6
134 
135 
136 typedef struct
137 {
138 
139   FIXP_DBL corr_thresh;                 /*!< Only set an IS seed it left/right channel correlation is above corr_thresh */
140 
141   FIXP_DBL total_error_thresh;          /*!< When expanding the IS region to more SFBs only accept an error that is
142                                              not more than 'total_error_thresh' overall. */
143 
144   FIXP_DBL local_error_thresh;          /*!< When expanding the IS region to more SFBs only accept an error that is
145                                              not more than 'local_error_thresh' for the current SFB. */
146 
147   FIXP_DBL direction_deviation_thresh;  /*!< The maximum allowed change of the intensity direction (unit: IS scale) */
148 
149   FIXP_DBL is_region_min_loudness;      /*!< IS regions need to have a minimal percentage of the overall loudness, e.g. 0.06 == 6% */
150 
151   INT      min_is_sfbs;                 /*!< Only perform IS if 'min_is_sfbs' neighboring SFBs can be processed */
152 
153   FIXP_DBL left_right_ratio_threshold;  /*!< No IS if the panning angle is not far from the middle, MS will do */
154 
155 } INTENSITY_PARAMETERS;
156 
157 
158 /*****************************************************************************
159 
160     functionname: calcSfbMaxScale
161 
162     description:  Calc max value in scalefactor band
163 
164     input:        *mdctSpectrum
165                    l1
166                    l2
167 
168     output:       none
169 
170     returns:      scalefactor
171 
172 *****************************************************************************/
173 static INT
calcSfbMaxScale(const FIXP_DBL * mdctSpectrum,const INT l1,const INT l2)174 calcSfbMaxScale(const FIXP_DBL *mdctSpectrum,
175                 const INT       l1,
176                 const INT       l2)
177 {
178   INT i;
179   INT sfbMaxScale;
180   FIXP_DBL maxSpc;
181 
182   maxSpc = FL2FXCONST_DBL(0.0);
183   for (i=l1; i<l2; i++) {
184     FIXP_DBL tmp = fixp_abs((FIXP_DBL)mdctSpectrum[i]);
185     maxSpc = fixMax(maxSpc, tmp);
186   }
187   sfbMaxScale = (maxSpc==FL2FXCONST_DBL(0.0)) ? (DFRACT_BITS-2) : CntLeadingZeros(maxSpc)-1;
188 
189   return sfbMaxScale;
190  }
191 
192 
193 /*****************************************************************************
194 
195     functionname: FDKaacEnc_initIsParams
196 
197     description:  Initialization of intensity parameters
198 
199     input:        isParams
200 
201     output:       isParams
202 
203     returns:      none
204 
205 *****************************************************************************/
206 static void
FDKaacEnc_initIsParams(INTENSITY_PARAMETERS * isParams)207 FDKaacEnc_initIsParams(INTENSITY_PARAMETERS *isParams)
208 {
209   isParams->corr_thresh                = IS_CORR_THRESH;
210   isParams->total_error_thresh         = IS_TOTAL_ERROR_THRESH;
211   isParams->local_error_thresh         = IS_LOCAL_ERROR_THRESH;
212   isParams->direction_deviation_thresh = IS_DIRECTION_DEVIATION_THRESH;
213   isParams->is_region_min_loudness     = IS_REGION_MIN_LOUDNESS;
214   isParams->min_is_sfbs                = IS_MIN_SFBS;
215   isParams->left_right_ratio_threshold = IS_LEFT_RIGHT_RATIO_THRESH;
216 }
217 
218 
219 /*****************************************************************************
220 
221     functionname: FDKaacEnc_prepareIntensityDecision
222 
223     description:  Prepares intensity decision
224 
225     input:        sfbEnergyLeft
226                   sfbEnergyRight
227                   sfbEnergyLdDataLeft
228                   sfbEnergyLdDataRight
229                   mdctSpectrumLeft
230                   sfbEnergyLdDataRight
231                   isParams
232 
233     output:       hrrErr            scale: none
234                   isMask            scale: none
235                   realScale         scale: LD_DATA_SHIFT + REAL_SCALE_SF
236                   normSfbLoudness   scale: none
237 
238     returns:      none
239 
240 *****************************************************************************/
241 static void
FDKaacEnc_prepareIntensityDecision(const FIXP_DBL * sfbEnergyLeft,const FIXP_DBL * sfbEnergyRight,const FIXP_DBL * sfbEnergyLdDataLeft,const FIXP_DBL * sfbEnergyLdDataRight,const FIXP_DBL * mdctSpectrumLeft,const FIXP_DBL * mdctSpectrumRight,const INTENSITY_PARAMETERS * isParams,FIXP_DBL * hrrErr,INT * isMask,FIXP_DBL * realScale,FIXP_DBL * normSfbLoudness,const INT sfbCnt,const INT sfbPerGroup,const INT maxSfbPerGroup,const INT * sfbOffset)242 FDKaacEnc_prepareIntensityDecision(const FIXP_DBL    *sfbEnergyLeft,
243                                    const FIXP_DBL    *sfbEnergyRight,
244                                    const FIXP_DBL    *sfbEnergyLdDataLeft,
245                                    const FIXP_DBL    *sfbEnergyLdDataRight,
246                                    const FIXP_DBL    *mdctSpectrumLeft,
247                                    const FIXP_DBL    *mdctSpectrumRight,
248                                    const INTENSITY_PARAMETERS *isParams,
249                                    FIXP_DBL    *hrrErr,
250                                    INT         *isMask,
251                                    FIXP_DBL    *realScale,
252                                    FIXP_DBL    *normSfbLoudness,
253                                    const INT    sfbCnt,
254                                    const INT    sfbPerGroup,
255                                    const INT    maxSfbPerGroup,
256                                    const INT   *sfbOffset)
257 {
258   INT j,sfb,sfboffs;
259   INT grpCounter;
260 
261   /* temporary variables to compute loudness */
262   FIXP_DBL overallLoudness[MAX_NO_OF_GROUPS];
263 
264   /* temporary variables to compute correlation */
265   FIXP_DBL channelCorr[MAX_GROUPED_SFB];
266   FIXP_DBL ml, mr;
267   FIXP_DBL prod_lr;
268   FIXP_DBL square_l, square_r;
269   FIXP_DBL tmp_l, tmp_r;
270   FIXP_DBL inv_n;
271 
272   FDKmemclear(channelCorr,     MAX_GROUPED_SFB*sizeof(FIXP_DBL));
273   FDKmemclear(normSfbLoudness, MAX_GROUPED_SFB*sizeof(FIXP_DBL));
274   FDKmemclear(overallLoudness, MAX_NO_OF_GROUPS*sizeof(FIXP_DBL));
275   FDKmemclear(realScale,       MAX_GROUPED_SFB*sizeof(FIXP_DBL));
276 
277   for (grpCounter = 0, sfboffs = 0; sfboffs < sfbCnt; sfboffs += sfbPerGroup, grpCounter++) {
278     overallLoudness[grpCounter] = FL2FXCONST_DBL(0.0f);
279     for (sfb = 0; sfb < maxSfbPerGroup; sfb++) {
280       INT sL,sR,s;
281       FIXP_DBL isValue = sfbEnergyLdDataLeft[sfb+sfboffs]-sfbEnergyLdDataRight[sfb+sfboffs];
282 
283       /* delimitate intensity scale value to representable range */
284       realScale[sfb + sfboffs] = fixMin(FL2FXCONST_DBL(60.f/(1<<(REAL_SCALE_SF+LD_DATA_SHIFT))), fixMax(FL2FXCONST_DBL(-60.f/(1<<(REAL_SCALE_SF+LD_DATA_SHIFT))), isValue));
285 
286       sL = fixMax(0,(CntLeadingZeros(sfbEnergyLeft[sfb + sfboffs])-1));
287       sR = fixMax(0,(CntLeadingZeros(sfbEnergyRight[sfb + sfboffs])-1));
288       s  = (fixMin(sL,sR)>>2)<<2;
289       normSfbLoudness[sfb + sfboffs] = sqrtFixp(sqrtFixp(((sfbEnergyLeft[sfb + sfboffs]<<s) >> 1) + ((sfbEnergyRight[sfb + sfboffs]<<s) >> 1))) >> (s>>2);
290 
291       overallLoudness[grpCounter] += normSfbLoudness[sfb + sfboffs] >> OVERALL_LOUDNESS_SF;
292       /* don't do intensity if
293        * - panning angle is too close to the middle or
294        * - one channel is non-existent or
295        * - if it is dual mono */
296       if(   (sfbEnergyLeft[sfb + sfboffs] >= fMult(isParams->left_right_ratio_threshold,sfbEnergyRight[sfb + sfboffs]))
297          && (fMult(isParams->left_right_ratio_threshold,sfbEnergyLeft[sfb + sfboffs]) <= sfbEnergyRight[sfb + sfboffs]) ) {
298 
299         /* this will prevent post processing from considering this SFB for merging */
300         hrrErr[sfb + sfboffs] = FL2FXCONST_DBL(1.0/8.0);
301       }
302     }
303   }
304 
305   for (grpCounter = 0, sfboffs = 0; sfboffs < sfbCnt; sfboffs += sfbPerGroup, grpCounter++) {
306     INT invOverallLoudnessSF;
307     FIXP_DBL invOverallLoudness;
308 
309     if (overallLoudness[grpCounter] == FL2FXCONST_DBL(0.0)) {
310       invOverallLoudness = FL2FXCONST_DBL(0.0);
311       invOverallLoudnessSF = 0;
312     }
313     else {
314       invOverallLoudness = fDivNorm((FIXP_DBL)MAXVAL_DBL, overallLoudness[grpCounter],&invOverallLoudnessSF);
315       invOverallLoudnessSF = invOverallLoudnessSF - OVERALL_LOUDNESS_SF + 1; /* +1: compensate fMultDiv2() in subsequent loop */
316     }
317     invOverallLoudnessSF = fixMin(fixMax(invOverallLoudnessSF,-(DFRACT_BITS-1)),DFRACT_BITS-1);
318 
319     for (sfb = 0; sfb < maxSfbPerGroup; sfb++) {
320       FIXP_DBL tmp;
321 
322       tmp = fMultDiv2((normSfbLoudness[sfb + sfboffs]>>OVERALL_LOUDNESS_SF)<<OVERALL_LOUDNESS_SF,invOverallLoudness);
323 
324       normSfbLoudness[sfb + sfboffs] = scaleValue(tmp, invOverallLoudnessSF);
325 
326       channelCorr[sfb + sfboffs] = FL2FXCONST_DBL(0.0f);
327 
328       /* max width of scalefactorband is 96; width's are always even */
329       /* inv_n is scaled with factor 2 to compensate fMultDiv2() in subsequent loops */
330       inv_n = GetInvInt((sfbOffset[sfb + sfboffs + 1] - sfbOffset[sfb + sfboffs])>>1);
331 
332       if (inv_n > FL2FXCONST_DBL(0.0f)) {
333         INT s,sL,sR;
334 
335         /* correlation := Pearson's product-moment coefficient */
336         /* compute correlation between channels and check if it is over threshold */
337         ml       = FL2FXCONST_DBL(0.0f);
338         mr       = FL2FXCONST_DBL(0.0f);
339         prod_lr  = FL2FXCONST_DBL(0.0f);
340         square_l = FL2FXCONST_DBL(0.0f);
341         square_r = FL2FXCONST_DBL(0.0f);
342 
343         sL = calcSfbMaxScale(mdctSpectrumLeft,sfbOffset[sfb+sfboffs],sfbOffset[sfb+sfboffs+1]);
344         sR = calcSfbMaxScale(mdctSpectrumRight,sfbOffset[sfb+sfboffs],sfbOffset[sfb+sfboffs+1]);
345         s = fixMin(sL,sR);
346 
347         for (j = sfbOffset[sfb + sfboffs]; j < sfbOffset[sfb + sfboffs + 1]; j++) {
348           ml += fMultDiv2((mdctSpectrumLeft[j]  << s),inv_n);             // scaled with mdctScale - s + inv_n
349           mr += fMultDiv2((mdctSpectrumRight[j] << s),inv_n);             // scaled with mdctScale - s + inv_n
350         }
351         ml = fMultDiv2(ml,inv_n);                                         // scaled with mdctScale - s + inv_n
352         mr = fMultDiv2(mr,inv_n);                                         // scaled with mdctScale - s + inv_n
353 
354         for (j = sfbOffset[sfb + sfboffs]; j < sfbOffset[sfb + sfboffs + 1]; j++) {
355           tmp_l = fMultDiv2((mdctSpectrumLeft[j]  << s),inv_n) - ml;      // scaled with mdctScale - s + inv_n
356           tmp_r = fMultDiv2((mdctSpectrumRight[j] << s),inv_n) - mr;      // scaled with mdctScale - s + inv_n
357 
358           prod_lr  += fMultDiv2(tmp_l,tmp_r);                             // scaled with 2*(mdctScale - s + inv_n) + 1
359           square_l += fPow2Div2(tmp_l);                                   // scaled with 2*(mdctScale - s + inv_n) + 1
360           square_r += fPow2Div2(tmp_r);                                   // scaled with 2*(mdctScale - s + inv_n) + 1
361         }
362         prod_lr  = prod_lr  << 1;                                         // scaled with 2*(mdctScale - s + inv_n)
363         square_l = square_l << 1;                                         // scaled with 2*(mdctScale - s + inv_n)
364         square_r = square_r << 1;                                         // scaled with 2*(mdctScale - s + inv_n)
365 
366         if (square_l > FL2FXCONST_DBL(0.0f) && square_r > FL2FXCONST_DBL(0.0f)) {
367           INT channelCorrSF = 0;
368 
369           /* local scaling of square_l and square_r is compensated after sqrt calculation */
370           sL  = fixMax(0,(CntLeadingZeros(square_l)-1));
371           sR  = fixMax(0,(CntLeadingZeros(square_r)-1));
372           s   = ((sL + sR)>>1)<<1;
373           sL  = fixMin(sL,s);
374           sR  = s-sL;
375           tmp = fMult(square_l<<sL,square_r<<sR);
376           tmp = sqrtFixp(tmp);
377 
378           FDK_ASSERT(tmp > FL2FXCONST_DBL(0.0f));
379 
380           /* numerator and denominator have the same scaling */
381           if (prod_lr < FL2FXCONST_DBL(0.0f) ) {
382             channelCorr[sfb + sfboffs] = -(fDivNorm(-prod_lr,tmp,&channelCorrSF));
383 
384           }
385           else {
386             channelCorr[sfb + sfboffs] =  (fDivNorm( prod_lr,tmp,&channelCorrSF));
387           }
388           channelCorrSF = fixMin(fixMax(( channelCorrSF + ((sL+sR)>>1)),-(DFRACT_BITS-1)),DFRACT_BITS-1);
389 
390           if (channelCorrSF < 0) {
391             channelCorr[sfb + sfboffs] = channelCorr[sfb + sfboffs] >> (-channelCorrSF);
392           }
393           else {
394             /* avoid overflows due to limited computational accuracy */
395             if ( fAbs(channelCorr[sfb + sfboffs]) > (((FIXP_DBL)MAXVAL_DBL)>>channelCorrSF) ) {
396               if (channelCorr[sfb + sfboffs] < FL2FXCONST_DBL(0.0f))
397                 channelCorr[sfb + sfboffs] = -(FIXP_DBL) MAXVAL_DBL;
398               else
399                 channelCorr[sfb + sfboffs] =  (FIXP_DBL) MAXVAL_DBL;
400             }
401             else {
402               channelCorr[sfb + sfboffs] = channelCorr[sfb + sfboffs] << channelCorrSF;
403             }
404           }
405         }
406       }
407 
408       /* for post processing: hrrErr is the error in terms of (too little) correlation
409        * weighted with the loudness of the SFB; SFBs with small hrrErr can be merged */
410       if (hrrErr[sfb + sfboffs] == FL2FXCONST_DBL(1.0/8.0)) {
411         continue;
412       }
413 
414       hrrErr[sfb + sfboffs] = fMultDiv2((FL2FXCONST_DBL(0.25f)-(channelCorr[sfb + sfboffs]>>2)),normSfbLoudness[sfb + sfboffs]);
415 
416       /* set IS mask/vector to 1, if correlation is high enough */
417       if (fAbs(channelCorr[sfb + sfboffs]) >= isParams->corr_thresh) {
418         isMask[sfb + sfboffs] = 1;
419       }
420     }
421   }
422 }
423 
424 
425 /*****************************************************************************
426 
427     functionname: FDKaacEnc_finalizeIntensityDecision
428 
429     description:  Finalizes intensity decision
430 
431     input:        isParams          scale: none
432                   hrrErr            scale: none
433                   realIsScale       scale: LD_DATA_SHIFT + REAL_SCALE_SF
434                   normSfbLoudness   scale: none
435 
436     output:       isMask            scale: none
437 
438     returns:      none
439 
440 *****************************************************************************/
441 static void
FDKaacEnc_finalizeIntensityDecision(const FIXP_DBL * hrrErr,INT * isMask,const FIXP_DBL * realIsScale,const FIXP_DBL * normSfbLoudness,const INTENSITY_PARAMETERS * isParams,const INT sfbCnt,const INT sfbPerGroup,const INT maxSfbPerGroup)442 FDKaacEnc_finalizeIntensityDecision(const FIXP_DBL *hrrErr,
443                                     INT            *isMask,
444                                     const FIXP_DBL *realIsScale,
445                                     const FIXP_DBL *normSfbLoudness,
446                                     const INTENSITY_PARAMETERS *isParams,
447                                     const INT       sfbCnt,
448                                     const INT       sfbPerGroup,
449                                     const INT       maxSfbPerGroup)
450 {
451   INT sfb,sfboffs, j;
452   FIXP_DBL isScaleLast = FL2FXCONST_DBL(0.0f);
453   INT isStartValueFound = 0;
454 
455   for (sfboffs = 0; sfboffs < sfbCnt; sfboffs += sfbPerGroup) {
456     INT startIsSfb = 0;
457     INT inIsBlock = 0;
458     INT currentIsSfbCount = 0;
459     FIXP_DBL overallHrrError = FL2FXCONST_DBL(0.0f);
460     FIXP_DBL isRegionLoudness = FL2FXCONST_DBL(0.0f);
461 
462     for (sfb = 0; sfb < maxSfbPerGroup; sfb++) {
463       if (isMask[sfboffs + sfb] == 1) {
464         if (currentIsSfbCount == 0) {
465           startIsSfb = sfboffs + sfb;
466         }
467         if (isStartValueFound==0) {
468           isScaleLast = realIsScale[sfboffs + sfb];
469           isStartValueFound = 1;
470         }
471         inIsBlock = 1;
472         currentIsSfbCount++;
473         overallHrrError  += hrrErr[sfboffs + sfb] >> (MAX_SFB_PER_GROUP_SF-3);
474         isRegionLoudness += normSfbLoudness[sfboffs + sfb] >> MAX_SFB_PER_GROUP_SF;
475       }
476       else {
477         /* based on correlation, IS should not be used
478          * -> use it anyway, if overall error is below threshold
479          *    and if local error does not exceed threshold
480          * otherwise: check if there are enough IS SFBs
481          */
482         if (inIsBlock) {
483           overallHrrError  += hrrErr[sfboffs + sfb] >> (MAX_SFB_PER_GROUP_SF-3);
484           isRegionLoudness += normSfbLoudness[sfboffs + sfb] >> MAX_SFB_PER_GROUP_SF;
485 
486           if ( (hrrErr[sfboffs + sfb] < (isParams->local_error_thresh>>3)) && (overallHrrError < (isParams->total_error_thresh>>MAX_SFB_PER_GROUP_SF)) ) {
487             currentIsSfbCount++;
488             /* overwrite correlation based decision */
489             isMask[sfboffs + sfb] = 1;
490           } else {
491             inIsBlock = 0;
492           }
493         }
494       }
495       /* check for large direction deviation */
496       if (inIsBlock) {
497         if( fAbs(isScaleLast-realIsScale[sfboffs + sfb]) < (isParams->direction_deviation_thresh>>(REAL_SCALE_SF+LD_DATA_SHIFT-IS_DIRECTION_DEVIATION_THRESH_SF)) ) {
498           isScaleLast = realIsScale[sfboffs + sfb];
499         }
500         else{
501           isMask[sfboffs + sfb] = 0;
502           inIsBlock = 0;
503           currentIsSfbCount--;
504         }
505       }
506 
507       if (currentIsSfbCount > 0 && (!inIsBlock || sfb == maxSfbPerGroup - 1)) {
508         /* not enough SFBs -> do not use IS */
509         if (currentIsSfbCount < isParams->min_is_sfbs || (isRegionLoudness < isParams->is_region_min_loudness>>MAX_SFB_PER_GROUP_SF)) {
510           for(j = startIsSfb; j <= sfboffs + sfb; j++) {
511             isMask[j] = 0;
512           }
513           isScaleLast = FL2FXCONST_DBL(0.0f);
514           isStartValueFound = 0;
515           for (j=0; j < startIsSfb; j++) {
516             if (isMask[j]!=0) {
517               isScaleLast = realIsScale[j];
518               isStartValueFound = 1;
519             }
520           }
521         }
522         currentIsSfbCount = 0;
523         overallHrrError = FL2FXCONST_DBL(0.0f);
524         isRegionLoudness = FL2FXCONST_DBL(0.0f);
525       }
526     }
527   }
528 }
529 
530 
531 /*****************************************************************************
532 
533     functionname: FDKaacEnc_IntensityStereoProcessing
534 
535     description:  Intensity stereo processing tool
536 
537     input:        sfbEnergyLeft
538                   sfbEnergyRight
539                   mdctSpectrumLeft
540                   mdctSpectrumRight
541                   sfbThresholdLeft
542                   sfbThresholdRight
543                   sfbSpreadEnLeft
544                   sfbSpreadEnRight
545                   sfbEnergyLdDataLeft
546                   sfbEnergyLdDataRight
547 
548     output:       isBook
549                   isScale
550                   pnsData->pnsFlag
551                   msDigest                 zeroed from start to sfbCnt
552                   msMask                   zeroed from start to sfbCnt
553                   mdctSpectrumRight        zeroed where isBook!=0
554                   sfbEnergyRight           zeroed where isBook!=0
555                   sfbSpreadEnRight       zeroed where isBook!=0
556                   sfbThresholdRight        zeroed where isBook!=0
557                   sfbEnergyLdDataRight     FL2FXCONST_DBL(-1.0) where isBook!=0
558                   sfbThresholdLdDataRight  FL2FXCONST_DBL(-0.515625f) where isBook!=0
559 
560     returns:      none
561 
562 *****************************************************************************/
FDKaacEnc_IntensityStereoProcessing(FIXP_DBL * sfbEnergyLeft,FIXP_DBL * sfbEnergyRight,FIXP_DBL * mdctSpectrumLeft,FIXP_DBL * mdctSpectrumRight,FIXP_DBL * sfbThresholdLeft,FIXP_DBL * sfbThresholdRight,FIXP_DBL * sfbThresholdLdDataRight,FIXP_DBL * sfbSpreadEnLeft,FIXP_DBL * sfbSpreadEnRight,FIXP_DBL * sfbEnergyLdDataLeft,FIXP_DBL * sfbEnergyLdDataRight,INT * msDigest,INT * msMask,const INT sfbCnt,const INT sfbPerGroup,const INT maxSfbPerGroup,const INT * sfbOffset,const INT allowIS,INT * isBook,INT * isScale,PNS_DATA * RESTRICT pnsData[2])563 void FDKaacEnc_IntensityStereoProcessing(
564         FIXP_DBL                  *sfbEnergyLeft,
565         FIXP_DBL                  *sfbEnergyRight,
566         FIXP_DBL                  *mdctSpectrumLeft,
567         FIXP_DBL                  *mdctSpectrumRight,
568         FIXP_DBL                  *sfbThresholdLeft,
569         FIXP_DBL                  *sfbThresholdRight,
570         FIXP_DBL                  *sfbThresholdLdDataRight,
571         FIXP_DBL                  *sfbSpreadEnLeft,
572         FIXP_DBL                  *sfbSpreadEnRight,
573         FIXP_DBL                  *sfbEnergyLdDataLeft,
574         FIXP_DBL                  *sfbEnergyLdDataRight,
575         INT                       *msDigest,
576         INT                       *msMask,
577         const INT                  sfbCnt,
578         const INT                  sfbPerGroup,
579         const INT                  maxSfbPerGroup,
580         const INT                 *sfbOffset,
581         const INT                  allowIS,
582         INT                       *isBook,
583         INT                       *isScale,
584         PNS_DATA         *RESTRICT pnsData[2]
585         )
586 {
587   INT sfb,sfboffs, j;
588   FIXP_DBL scale;
589   FIXP_DBL lr;
590   FIXP_DBL hrrErr[MAX_GROUPED_SFB];
591   FIXP_DBL normSfbLoudness[MAX_GROUPED_SFB];
592   FIXP_DBL realIsScale[MAX_GROUPED_SFB];
593   INTENSITY_PARAMETERS isParams;
594   INT isMask[MAX_GROUPED_SFB];
595 
596   FDKmemclear((void*)isBook,sfbCnt*sizeof(INT));
597   FDKmemclear((void*)isMask,sfbCnt*sizeof(INT));
598   FDKmemclear((void*)realIsScale,sfbCnt*sizeof(FIXP_DBL));
599   FDKmemclear((void*)isScale,sfbCnt*sizeof(INT));
600   FDKmemclear((void*)hrrErr,sfbCnt*sizeof(FIXP_DBL));
601 
602   if (!allowIS)
603     return;
604 
605   FDKaacEnc_initIsParams(&isParams);
606 
607   /* compute / set the following values per SFB:
608    * - left/right ratio between channels
609    * - normalized loudness
610    *   + loudness == average of energy in channels to 0.25
611    *   + normalization: division by sum of all SFB loudnesses
612    * - isMask (is set to 0 if channels are the same or one is 0)
613    */
614    FDKaacEnc_prepareIntensityDecision(sfbEnergyLeft,
615                                       sfbEnergyRight,
616                                       sfbEnergyLdDataLeft,
617                                       sfbEnergyLdDataRight,
618                                       mdctSpectrumLeft,
619                                       mdctSpectrumRight,
620                                       &isParams,
621                                       hrrErr,
622                                       isMask,
623                                       realIsScale,
624                                       normSfbLoudness,
625                                       sfbCnt,
626                                       sfbPerGroup,
627                                       maxSfbPerGroup,
628                                       sfbOffset);
629 
630   FDKaacEnc_finalizeIntensityDecision(hrrErr,
631                                       isMask,
632                                       realIsScale,
633                                       normSfbLoudness,
634                                       &isParams,
635                                       sfbCnt,
636                                       sfbPerGroup,
637                                       maxSfbPerGroup);
638 
639   for (sfb=0; sfb<sfbCnt; sfb+=sfbPerGroup) {
640     for (sfboffs=0; sfboffs<maxSfbPerGroup; sfboffs++) {
641       INT sL, sR;
642       FIXP_DBL inv_n;
643 
644       msMask[sfb+sfboffs] = 0;
645       if (isMask[sfb+sfboffs] == 0) {
646         continue;
647       }
648 
649       if (   (sfbEnergyLeft[sfb+sfboffs] < sfbThresholdLeft[sfb+sfboffs])
650           &&(fMult(FL2FXCONST_DBL(1.0f/1.5f),sfbEnergyRight[sfb+sfboffs]) > sfbThresholdRight[sfb+sfboffs]) ) {
651         continue;
652       }
653       /* NEW: if there is a big-enough IS region, switch off PNS */
654       if (pnsData[0]) {
655         if(pnsData[0]->pnsFlag[sfb+sfboffs]) {
656           pnsData[0]->pnsFlag[sfb+sfboffs] = 0;
657         }
658         if(pnsData[1]->pnsFlag[sfb+sfboffs]) {
659           pnsData[1]->pnsFlag[sfb+sfboffs] = 0;
660         }
661       }
662 
663       inv_n = GetInvInt((sfbOffset[sfb + sfboffs + 1] - sfbOffset[sfb + sfboffs])>>1);  // scaled with 2 to compensate fMultDiv2() in subsequent loop
664       sL = calcSfbMaxScale(mdctSpectrumLeft,sfbOffset[sfb+sfboffs],sfbOffset[sfb+sfboffs+1]);
665       sR = calcSfbMaxScale(mdctSpectrumRight,sfbOffset[sfb+sfboffs],sfbOffset[sfb+sfboffs+1]);
666 
667       lr = FL2FXCONST_DBL(0.0f);
668       for (j=sfbOffset[sfb+sfboffs]; j<sfbOffset[sfb+sfboffs+1]; j++)
669         lr += fMultDiv2(fMultDiv2(mdctSpectrumLeft[j]<<sL,mdctSpectrumRight[j]<<sR),inv_n);
670       lr = lr<<1;
671 
672       if (lr < FL2FXCONST_DBL(0.0f)) {
673         /* This means OUT OF phase intensity stereo, cf. standard */
674         INT s0, s1, s2;
675         FIXP_DBL tmp, d, ed = FL2FXCONST_DBL(0.0f);
676 
677         s0 = fixMin(sL,sR);
678         for (j=sfbOffset[sfb+sfboffs]; j<sfbOffset[sfb+sfboffs+1]; j++) {
679           d = ((mdctSpectrumLeft[j]<<s0)>>1) - ((mdctSpectrumRight[j]<<s0)>>1);
680           ed += fMultDiv2(d,d)>>(MDCT_SPEC_SF-1);
681         }
682         msMask[sfb+sfboffs] = 1;
683         tmp = fDivNorm(sfbEnergyLeft[sfb+sfboffs],ed,&s1);
684         s2 = (s1) + (2*s0) - 2 - MDCT_SPEC_SF;
685         if (s2 & 1) {
686           tmp = tmp>>1;
687           s2 = s2+1;
688         }
689         s2 = (s2>>1) + 1;  // +1 compensate fMultDiv2() in subsequent loop
690         s2 = fixMin(fixMax(s2,-(DFRACT_BITS-1)),(DFRACT_BITS-1));
691         scale = sqrtFixp(tmp);
692         if (s2 < 0) {
693           s2 = -s2;
694           for (j=sfbOffset[sfb+sfboffs]; j<sfbOffset[sfb+sfboffs+1]; j++) {
695             mdctSpectrumLeft[j] = (fMultDiv2(mdctSpectrumLeft[j],scale) - fMultDiv2(mdctSpectrumRight[j],scale)) >> s2;
696             mdctSpectrumRight[j] = FL2FXCONST_DBL(0.0f);
697           }
698         }
699         else {
700           for (j=sfbOffset[sfb+sfboffs]; j<sfbOffset[sfb+sfboffs+1]; j++) {
701             mdctSpectrumLeft[j] = (fMultDiv2(mdctSpectrumLeft[j],scale) - fMultDiv2(mdctSpectrumRight[j],scale)) << s2;
702             mdctSpectrumRight[j] = FL2FXCONST_DBL(0.0f);
703           }
704         }
705       }
706       else {
707         /* This means IN phase intensity stereo, cf. standard */
708         INT s0,s1,s2;
709         FIXP_DBL tmp, s, es = FL2FXCONST_DBL(0.0f);
710 
711         s0 = fixMin(sL,sR);
712         for (j=sfbOffset[sfb+sfboffs]; j<sfbOffset[sfb+sfboffs+1]; j++) {
713           s   = ((mdctSpectrumLeft[j]<<s0)>>1) + ((mdctSpectrumRight[j]<<s0)>>1);
714           es += fMultDiv2(s,s)>>(MDCT_SPEC_SF-1);     // scaled 2*(mdctScale - s0 + 1) + MDCT_SPEC_SF
715         }
716         msMask[sfb+sfboffs] = 0;
717         tmp = fDivNorm(sfbEnergyLeft[sfb+sfboffs],es,&s1);
718         s2 = (s1) + (2*s0) - 2 - MDCT_SPEC_SF;
719         if (s2 & 1) {
720           tmp = tmp>>1;
721           s2 = s2 + 1;
722         }
723         s2 = (s2>>1) + 1; // +1 compensate fMultDiv2() in subsequent loop
724         s2 = fixMin(fixMax(s2,-(DFRACT_BITS-1)),(DFRACT_BITS-1));
725         scale = sqrtFixp(tmp);
726         if (s2 < 0) {
727           s2 = -s2;
728           for (j=sfbOffset[sfb+sfboffs]; j<sfbOffset[sfb+sfboffs+1]; j++) {
729             mdctSpectrumLeft[j] = (fMultDiv2(mdctSpectrumLeft[j],scale) + fMultDiv2(mdctSpectrumRight[j],scale)) >> s2;
730             mdctSpectrumRight[j] = FL2FXCONST_DBL(0.0f);
731           }
732         }
733         else {
734           for (j=sfbOffset[sfb+sfboffs]; j<sfbOffset[sfb+sfboffs+1]; j++) {
735             mdctSpectrumLeft[j] = (fMultDiv2(mdctSpectrumLeft[j],scale) + fMultDiv2(mdctSpectrumRight[j],scale)) << s2;
736             mdctSpectrumRight[j] = FL2FXCONST_DBL(0.0f);
737           }
738         }
739       }
740 
741       isBook[sfb+sfboffs] = CODE_BOOK_IS_IN_PHASE_NO;
742 
743       if ( realIsScale[sfb+sfboffs] < FL2FXCONST_DBL(0.0f) ) {
744         isScale[sfb+sfboffs] = (INT)(((realIsScale[sfb+sfboffs]>>1)-FL2FXCONST_DBL(0.5f/(1<<(REAL_SCALE_SF+LD_DATA_SHIFT+1))))>>(DFRACT_BITS-1-REAL_SCALE_SF-LD_DATA_SHIFT-1)) + 1;
745       }
746       else {
747         isScale[sfb+sfboffs] = (INT)(((realIsScale[sfb+sfboffs]>>1)+FL2FXCONST_DBL(0.5f/(1<<(REAL_SCALE_SF+LD_DATA_SHIFT+1))))>>(DFRACT_BITS-1-REAL_SCALE_SF-LD_DATA_SHIFT-1));
748       }
749 
750       sfbEnergyRight[sfb+sfboffs] = FL2FXCONST_DBL(0.0f);
751       sfbEnergyLdDataRight[sfb+sfboffs] = FL2FXCONST_DBL(-1.0f);
752       sfbThresholdRight[sfb+sfboffs] = FL2FXCONST_DBL(0.0f);
753       sfbThresholdLdDataRight[sfb+sfboffs] = FL2FXCONST_DBL(-0.515625f);
754       sfbSpreadEnRight[sfb+sfboffs] = FL2FXCONST_DBL(0.0f);
755 
756       *msDigest = MS_SOME;
757     }
758   }
759 }
760 
761