libAACenc/src/intensity.cpp

   1
   2 /* -----------------------------------------------------------------------------------------------------------
   3 Software License for The Fraunhofer FDK AAC Codec Library for Android
   4
   5 © Copyright  1995 - 2013 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
 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
 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
 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       FDK_ASSERT(50 >= 49);
 329       /* max width of scalefactorband is 96; width's are always even */
 330       /* inv_n is scaled with factor 2 to compensate fMultDiv2() in subsequent loops */
 331       inv_n = GetInvInt((sfbOffset[sfb + sfboffs + 1] - sfbOffset[sfb + sfboffs])>>1);
 332
 333       if (inv_n > FL2FXCONST_DBL(0.0f)) {
 334         INT s,sL,sR;
 335
 336         /* correlation := Pearson's product-moment coefficient */
 337         /* compute correlation between channels and check if it is over threshold */
 338         ml       = FL2FXCONST_DBL(0.0f);
 339         mr       = FL2FXCONST_DBL(0.0f);
 340         prod_lr  = FL2FXCONST_DBL(0.0f);
 341         square_l = FL2FXCONST_DBL(0.0f);
 342         square_r = FL2FXCONST_DBL(0.0f);
 343
 344         sL = calcSfbMaxScale(mdctSpectrumLeft,sfbOffset[sfb+sfboffs],sfbOffset[sfb+sfboffs+1]);
 345         sR = calcSfbMaxScale(mdctSpectrumRight,sfbOffset[sfb+sfboffs],sfbOffset[sfb+sfboffs+1]);
 346         s = fixMin(sL,sR);
 347
 348         for (j = sfbOffset[sfb + sfboffs]; j < sfbOffset[sfb + sfboffs + 1]; j++) {
 349           ml += fMultDiv2((mdctSpectrumLeft[j]  << s),inv_n);             // scaled with mdctScale - s + inv_n
 350           mr += fMultDiv2((mdctSpectrumRight[j] << s),inv_n);             // scaled with mdctScale - s + inv_n
 351         }
 352         ml = fMultDiv2(ml,inv_n);                                         // scaled with mdctScale - s + inv_n
 353         mr = fMultDiv2(mr,inv_n);                                         // scaled with mdctScale - s + inv_n
 354
 355         for (j = sfbOffset[sfb + sfboffs]; j < sfbOffset[sfb + sfboffs + 1]; j++) {
 356           tmp_l = fMultDiv2((mdctSpectrumLeft[j]  << s),inv_n) - ml;      // scaled with mdctScale - s + inv_n
 357           tmp_r = fMultDiv2((mdctSpectrumRight[j] << s),inv_n) - mr;      // scaled with mdctScale - s + inv_n
 358
 359           prod_lr  += fMultDiv2(tmp_l,tmp_r);                             // scaled with 2*(mdctScale - s + inv_n) + 1
 360           square_l += fPow2Div2(tmp_l);                                   // scaled with 2*(mdctScale - s + inv_n) + 1
 361           square_r += fPow2Div2(tmp_r);                                   // scaled with 2*(mdctScale - s + inv_n) + 1
 362         }
 363         prod_lr  = prod_lr  << 1;                                         // scaled with 2*(mdctScale - s + inv_n)
 364         square_l = square_l << 1;                                         // scaled with 2*(mdctScale - s + inv_n)
 365         square_r = square_r << 1;                                         // scaled with 2*(mdctScale - s + inv_n)
 366
 367         if (square_l > FL2FXCONST_DBL(0.0f) && square_r > FL2FXCONST_DBL(0.0f)) {
 368           INT channelCorrSF = 0;
 369
 370           /* local scaling of square_l and square_r is compensated after sqrt calculation */
 371           sL  = fixMax(0,(CntLeadingZeros(square_l)-1));
 372           sR  = fixMax(0,(CntLeadingZeros(square_r)-1));
 373           s   = ((sL + sR)>>1)<<1;
 374           sL  = fixMin(sL,s);
 375           sR  = s-sL;
 376           tmp = fMult(square_l<<sL,square_r<<sR);
 377           tmp = sqrtFixp(tmp);
 378
 379           FDK_ASSERT(tmp > FL2FXCONST_DBL(0.0f));
 380
 381           /* numerator and denominator have the same scaling */
 382           if (prod_lr < FL2FXCONST_DBL(0.0f) ) {
 383             channelCorr[sfb + sfboffs] = -(fDivNorm(-prod_lr,tmp,&channelCorrSF));
 384
 385           }
 386           else {
 387             channelCorr[sfb + sfboffs] =  (fDivNorm( prod_lr,tmp,&channelCorrSF));
 388           }
 389           channelCorrSF = fixMin(fixMax(( channelCorrSF + ((sL+sR)>>1)),-(DFRACT_BITS-1)),DFRACT_BITS-1);
 390
 391           if (channelCorrSF < 0) {
 392             channelCorr[sfb + sfboffs] = channelCorr[sfb + sfboffs] >> (-channelCorrSF);
 393           }
 394           else {
 395             /* avoid overflows due to limited computational accuracy */
 396             if ( fAbs(channelCorr[sfb + sfboffs]) > (((FIXP_DBL)MAXVAL_DBL)>>channelCorrSF) ) {
 397               if (channelCorr[sfb + sfboffs] < FL2FXCONST_DBL(0.0f))
 398                 channelCorr[sfb + sfboffs] = -(FIXP_DBL) MAXVAL_DBL;
 399               else
 400                 channelCorr[sfb + sfboffs] =  (FIXP_DBL) MAXVAL_DBL;
 401             }
 402             else {
 403               channelCorr[sfb + sfboffs] = channelCorr[sfb + sfboffs] << channelCorrSF;
 404             }
 405           }
 406         }
 407       }
 408
 409       /* for post processing: hrrErr is the error in terms of (too little) correlation
 410        * weighted with the loudness of the SFB; SFBs with small hrrErr can be merged */
 411       if (hrrErr[sfb + sfboffs] == FL2FXCONST_DBL(1.0/8.0)) {
 412         continue;
 413       }
 414
 415       hrrErr[sfb + sfboffs] = fMultDiv2((FL2FXCONST_DBL(0.25f)-(channelCorr[sfb + sfboffs]>>2)),normSfbLoudness[sfb + sfboffs]);
 416
 417       /* set IS mask/vector to 1, if correlation is high enough */
 418       if (fAbs(channelCorr[sfb + sfboffs]) >= isParams->corr_thresh) {
 419         isMask[sfb + sfboffs] = 1;
 420       }
 421     }
 422   }
 423 }
 424
 425
 426 /*****************************************************************************
 427
 428     functionname: FDKaacEnc_finalizeIntensityDecision
 429
 430     description:  Finalizes intensity decision
 431
 432     input:        isParams          scale: none
 433                   hrrErr            scale: none
 434                   realIsScale       scale: LD_DATA_SHIFT + REAL_SCALE_SF
 435                   normSfbLoudness   scale: none
 436
 437     output:       isMask            scale: none
 438
 439     returns:      none
 440
 441 *****************************************************************************/
 442 static void
 443 FDKaacEnc_finalizeIntensityDecision(const FIXP_DBL *hrrErr,
 444                                     INT            *isMask,
 445                                     const FIXP_DBL *realIsScale,
 446                                     const FIXP_DBL *normSfbLoudness,
 447                                     const INTENSITY_PARAMETERS *isParams,
 448                                     const INT       sfbCnt,
 449                                     const INT       sfbPerGroup,
 450                                     const INT       maxSfbPerGroup)
 451 {
 452   INT sfb,sfboffs, j;
 453   FIXP_DBL isScaleLast = FL2FXCONST_DBL(0.0f);
 454   INT isStartValueFound = 0;
 455
 456   for (sfboffs = 0; sfboffs < sfbCnt; sfboffs += sfbPerGroup) {
 457     INT startIsSfb = 0;
 458     INT inIsBlock = 0;
 459     INT currentIsSfbCount = 0;
 460     FIXP_DBL overallHrrError = FL2FXCONST_DBL(0.0f);
 461     FIXP_DBL isRegionLoudness = FL2FXCONST_DBL(0.0f);
 462
 463     for (sfb = 0; sfb < maxSfbPerGroup; sfb++) {
 464       if (isMask[sfboffs + sfb] == 1) {
 465         if (currentIsSfbCount == 0) {
 466           startIsSfb = sfboffs + sfb;
 467         }
 468         if (isStartValueFound==0) {
 469           isScaleLast = realIsScale[sfboffs + sfb];
 470           isStartValueFound = 1;
 471         }
 472         inIsBlock = 1;
 473         currentIsSfbCount++;
 474         overallHrrError  += hrrErr[sfboffs + sfb] >> (MAX_SFB_PER_GROUP_SF-3);
 475         isRegionLoudness += normSfbLoudness[sfboffs + sfb] >> MAX_SFB_PER_GROUP_SF;
 476       }
 477       else {
 478         /* based on correlation, IS should not be used
 479          * -> use it anyway, if overall error is below threshold
 480          *    and if local error does not exceed threshold
 481          * otherwise: check if there are enough IS SFBs
 482          */
 483         if (inIsBlock) {
 484           overallHrrError  += hrrErr[sfboffs + sfb] >> (MAX_SFB_PER_GROUP_SF-3);
 485           isRegionLoudness += normSfbLoudness[sfboffs + sfb] >> MAX_SFB_PER_GROUP_SF;
 486
 487           if ( (hrrErr[sfboffs + sfb] < (isParams->local_error_thresh>>3)) && (overallHrrError < (isParams->total_error_thresh>>MAX_SFB_PER_GROUP_SF)) ) {
 488             currentIsSfbCount++;
 489             /* overwrite correlation based decision */
 490             isMask[sfboffs + sfb] = 1;
 491           } else {
 492             inIsBlock = 0;
 493           }
 494         }
 495       }
 496       /* check for large direction deviation */
 497       if (inIsBlock) {
 498         if( fAbs(isScaleLast-realIsScale[sfboffs + sfb]) < (isParams->direction_deviation_thresh>>(REAL_SCALE_SF+LD_DATA_SHIFT-IS_DIRECTION_DEVIATION_THRESH_SF)) ) {
 499           isScaleLast = realIsScale[sfboffs + sfb];
 500         }
 501         else{
 502           isMask[sfboffs + sfb] = 0;
 503           inIsBlock = 0;
 504           currentIsSfbCount--;
 505         }
 506       }
 507
 508       if (currentIsSfbCount > 0 && (!inIsBlock || sfb == maxSfbPerGroup - 1)) {
 509         /* not enough SFBs -> do not use IS */
 510         if (currentIsSfbCount < isParams->min_is_sfbs || (isRegionLoudness < isParams->is_region_min_loudness>>MAX_SFB_PER_GROUP_SF)) {
 511           for(j = startIsSfb; j <= sfboffs + sfb; j++) {
 512             isMask[j] = 0;
 513           }
 514           isScaleLast = FL2FXCONST_DBL(0.0f);
 515           isStartValueFound = 0;
 516           for (j=0; j < startIsSfb; j++) {
 517             if (isMask[j]!=0) {
 518               isScaleLast = realIsScale[j];
 519               isStartValueFound = 1;
 520             }
 521           }
 522         }
 523         currentIsSfbCount = 0;
 524         overallHrrError = FL2FXCONST_DBL(0.0f);
 525         isRegionLoudness = FL2FXCONST_DBL(0.0f);
 526       }
 527     }
 528   }
 529 }
 530
 531
 532 /*****************************************************************************
 533
 534     functionname: FDKaacEnc_IntensityStereoProcessing
 535
 536     description:  Intensity stereo processing tool
 537
 538     input:        sfbEnergyLeft
 539                   sfbEnergyRight
 540                   mdctSpectrumLeft
 541                   mdctSpectrumRight
 542                   sfbThresholdLeft
 543                   sfbThresholdRight
 544                   sfbSpreadEnLeft
 545                   sfbSpreadEnRight
 546                   sfbEnergyLdDataLeft
 547                   sfbEnergyLdDataRight
 548
 549     output:       isBook
 550                   isScale
 551                   pnsData->pnsFlag
 552                   msDigest                 zeroed from start to sfbCnt
 553                   msMask                   zeroed from start to sfbCnt
 554                   mdctSpectrumRight        zeroed where isBook!=0
 555                   sfbEnergyRight           zeroed where isBook!=0
 556                   sfbSpreadEnRight       zeroed where isBook!=0
 557                   sfbThresholdRight        zeroed where isBook!=0
 558                   sfbEnergyLdDataRight     FL2FXCONST_DBL(-1.0) where isBook!=0
 559                   sfbThresholdLdDataRight  FL2FXCONST_DBL(-0.515625f) where isBook!=0
 560
 561     returns:      none
 562
 563 *****************************************************************************/
 564 void FDKaacEnc_IntensityStereoProcessing(
 565         FIXP_DBL                  *sfbEnergyLeft,
 566         FIXP_DBL                  *sfbEnergyRight,
 567         FIXP_DBL                  *mdctSpectrumLeft,
 568         FIXP_DBL                  *mdctSpectrumRight,
 569         FIXP_DBL                  *sfbThresholdLeft,
 570         FIXP_DBL                  *sfbThresholdRight,
 571         FIXP_DBL                  *sfbThresholdLdDataRight,
 572         FIXP_DBL                  *sfbSpreadEnLeft,
 573         FIXP_DBL                  *sfbSpreadEnRight,
 574         FIXP_DBL                  *sfbEnergyLdDataLeft,
 575         FIXP_DBL                  *sfbEnergyLdDataRight,
 576         INT                       *msDigest,
 577         INT                       *msMask,
 578         const INT                  sfbCnt,
 579         const INT                  sfbPerGroup,
 580         const INT                  maxSfbPerGroup,
 581         const INT                 *sfbOffset,
 582         const INT                  allowIS,
 583         INT                       *isBook,
 584         INT                       *isScale,
 585         PNS_DATA         *RESTRICT pnsData[2]
 586         )
 587 {
 588   INT sfb,sfboffs, j;
 589   FIXP_DBL scale;
 590   FIXP_DBL lr;
 591   FIXP_DBL hrrErr[MAX_GROUPED_SFB];
 592   FIXP_DBL normSfbLoudness[MAX_GROUPED_SFB];
 593   FIXP_DBL realIsScale[MAX_GROUPED_SFB];
 594   INTENSITY_PARAMETERS isParams;
 595   INT isMask[MAX_GROUPED_SFB];
 596
 597   FDKmemclear((void*)isBook,sfbCnt*sizeof(INT));
 598   FDKmemclear((void*)isMask,sfbCnt*sizeof(INT));
 599   FDKmemclear((void*)realIsScale,sfbCnt*sizeof(FIXP_DBL));
 600   FDKmemclear((void*)isScale,sfbCnt*sizeof(INT));
 601   FDKmemclear((void*)hrrErr,sfbCnt*sizeof(FIXP_DBL));
 602
 603   if (!allowIS)
 604     return;
 605
 606   FDKaacEnc_initIsParams(&isParams);
 607
 608   /* compute / set the following values per SFB:
 609    * - left/right ratio between channels
 610    * - normalized loudness
 611    *   + loudness == average of energy in channels to 0.25
 612    *   + normalization: division by sum of all SFB loudnesses
 613    * - isMask (is set to 0 if channels are the same or one is 0)
 614    */
 615    FDKaacEnc_prepareIntensityDecision(sfbEnergyLeft,
 616                                       sfbEnergyRight,
 617                                       sfbEnergyLdDataLeft,
 618                                       sfbEnergyLdDataRight,
 619                                       mdctSpectrumLeft,
 620                                       mdctSpectrumRight,
 621                                       &isParams,
 622                                       hrrErr,
 623                                       isMask,
 624                                       realIsScale,
 625                                       normSfbLoudness,
 626                                       sfbCnt,
 627                                       sfbPerGroup,
 628                                       maxSfbPerGroup,
 629                                       sfbOffset);
 630
 631   FDKaacEnc_finalizeIntensityDecision(hrrErr,
 632                                       isMask,
 633                                       realIsScale,
 634                                       normSfbLoudness,
 635                                       &isParams,
 636                                       sfbCnt,
 637                                       sfbPerGroup,
 638                                       maxSfbPerGroup);
 639
 640   for (sfb=0; sfb<sfbCnt; sfb+=sfbPerGroup) {
 641     for (sfboffs=0; sfboffs<maxSfbPerGroup; sfboffs++) {
 642       INT sL, sR;
 643       FIXP_DBL inv_n;
 644
 645       msMask[sfb+sfboffs] = 0;
 646       if (isMask[sfb+sfboffs] == 0) {
 647         continue;
 648       }
 649
 650       if (   (sfbEnergyLeft[sfb+sfboffs] < sfbThresholdLeft[sfb+sfboffs])
 651           &&(fMult(FL2FXCONST_DBL(1.0f/1.5f),sfbEnergyRight[sfb+sfboffs]) > sfbThresholdRight[sfb+sfboffs]) ) {
 652         continue;
 653       }
 654       /* NEW: if there is a big-enough IS region, switch off PNS */
 655       if (pnsData[0]) {
 656         if(pnsData[0]->pnsFlag[sfb+sfboffs]) {
 657           pnsData[0]->pnsFlag[sfb+sfboffs] = 0;
 658         }
 659         if(pnsData[1]->pnsFlag[sfb+sfboffs]) {
 660           pnsData[1]->pnsFlag[sfb+sfboffs] = 0;
 661         }
 662       }
 663
 664       inv_n = GetInvInt((sfbOffset[sfb + sfboffs + 1] - sfbOffset[sfb + sfboffs])>>1);  // scaled with 2 to compensate fMultDiv2() in subsequent loop
 665       sL = calcSfbMaxScale(mdctSpectrumLeft,sfbOffset[sfb+sfboffs],sfbOffset[sfb+sfboffs+1]);
 666       sR = calcSfbMaxScale(mdctSpectrumRight,sfbOffset[sfb+sfboffs],sfbOffset[sfb+sfboffs+1]);
 667
 668       lr = FL2FXCONST_DBL(0.0f);
 669       for (j=sfbOffset[sfb+sfboffs]; j<sfbOffset[sfb+sfboffs+1]; j++)
 670         lr += fMultDiv2(fMultDiv2(mdctSpectrumLeft[j]<<sL,mdctSpectrumRight[j]<<sR),inv_n);
 671       lr = lr<<1;
 672
 673       if (lr < FL2FXCONST_DBL(0.0f)) {
 674         /* This means OUT OF phase intensity stereo, cf. standard */
 675         INT s0, s1, s2;
 676         FIXP_DBL tmp, d, ed = FL2FXCONST_DBL(0.0f);
 677
 678         s0 = fixMin(sL,sR);
 679         for (j=sfbOffset[sfb+sfboffs]; j<sfbOffset[sfb+sfboffs+1]; j++) {
 680           d = ((mdctSpectrumLeft[j]<<s0)>>1) - ((mdctSpectrumRight[j]<<s0)>>1);
 681           ed += fMultDiv2(d,d)>>(MDCT_SPEC_SF-1);
 682         }
 683         msMask[sfb+sfboffs] = 1;
 684         tmp = fDivNorm(sfbEnergyLeft[sfb+sfboffs],ed,&s1);
 685         s2 = (s1) + (2*s0) - 2 - MDCT_SPEC_SF;
 686         if (s2 & 1) {
 687           tmp = tmp>>1;
 688           s2 = s2+1;
 689         }
 690         s2 = (s2>>1) + 1;  // +1 compensate fMultDiv2() in subsequent loop
 691         s2 = fixMin(fixMax(s2,-(DFRACT_BITS-1)),(DFRACT_BITS-1));
 692         scale = sqrtFixp(tmp);
 693         if (s2 < 0) {
 694           s2 = -s2;
 695           for (j=sfbOffset[sfb+sfboffs]; j<sfbOffset[sfb+sfboffs+1]; j++) {
 696             mdctSpectrumLeft[j] = (fMultDiv2(mdctSpectrumLeft[j],scale) - fMultDiv2(mdctSpectrumRight[j],scale)) >> s2;
 697             mdctSpectrumRight[j] = FL2FXCONST_DBL(0.0f);
 698           }
 699         }
 700         else {
 701           for (j=sfbOffset[sfb+sfboffs]; j<sfbOffset[sfb+sfboffs+1]; j++) {
 702             mdctSpectrumLeft[j] = (fMultDiv2(mdctSpectrumLeft[j],scale) - fMultDiv2(mdctSpectrumRight[j],scale)) << s2;
 703             mdctSpectrumRight[j] = FL2FXCONST_DBL(0.0f);
 704           }
 705         }
 706       }
 707       else {
 708         /* This means IN phase intensity stereo, cf. standard */
 709         INT s0,s1,s2;
 710         FIXP_DBL tmp, s, es = FL2FXCONST_DBL(0.0f);
 711
 712         s0 = fixMin(sL,sR);
 713         for (j=sfbOffset[sfb+sfboffs]; j<sfbOffset[sfb+sfboffs+1]; j++) {
 714           s   = ((mdctSpectrumLeft[j]<<s0)>>1) + ((mdctSpectrumRight[j]<<s0)>>1);
 715           es += fMultDiv2(s,s)>>(MDCT_SPEC_SF-1);     // scaled 2*(mdctScale - s0 + 1) + MDCT_SPEC_SF
 716         }
 717         msMask[sfb+sfboffs] = 0;
 718         tmp = fDivNorm(sfbEnergyLeft[sfb+sfboffs],es,&s1);
 719         s2 = (s1) + (2*s0) - 2 - MDCT_SPEC_SF;
 720         if (s2 & 1) {
 721           tmp = tmp>>1;
 722           s2 = s2 + 1;
 723         }
 724         s2 = (s2>>1) + 1; // +1 compensate fMultDiv2() in subsequent loop
 725         s2 = fixMin(fixMax(s2,-(DFRACT_BITS-1)),(DFRACT_BITS-1));
 726         scale = sqrtFixp(tmp);
 727         if (s2 < 0) {
 728           s2 = -s2;
 729           for (j=sfbOffset[sfb+sfboffs]; j<sfbOffset[sfb+sfboffs+1]; j++) {
 730             mdctSpectrumLeft[j] = (fMultDiv2(mdctSpectrumLeft[j],scale) + fMultDiv2(mdctSpectrumRight[j],scale)) >> s2;
 731             mdctSpectrumRight[j] = FL2FXCONST_DBL(0.0f);
 732           }
 733         }
 734         else {
 735           for (j=sfbOffset[sfb+sfboffs]; j<sfbOffset[sfb+sfboffs+1]; j++) {
 736             mdctSpectrumLeft[j] = (fMultDiv2(mdctSpectrumLeft[j],scale) + fMultDiv2(mdctSpectrumRight[j],scale)) << s2;
 737             mdctSpectrumRight[j] = FL2FXCONST_DBL(0.0f);
 738           }
 739         }
 740       }
 741
 742       isBook[sfb+sfboffs] = CODE_BOOK_IS_IN_PHASE_NO;
 743
 744       if ( realIsScale[sfb+sfboffs] < FL2FXCONST_DBL(0.0f) ) {
 745         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;
 746       }
 747       else {
 748         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));
 749       }
 750
 751       sfbEnergyRight[sfb+sfboffs] = FL2FXCONST_DBL(0.0f);
 752       sfbEnergyLdDataRight[sfb+sfboffs] = FL2FXCONST_DBL(-1.0f);
 753       sfbThresholdRight[sfb+sfboffs] = FL2FXCONST_DBL(0.0f);
 754       sfbThresholdLdDataRight[sfb+sfboffs] = FL2FXCONST_DBL(-0.515625f);
 755       sfbSpreadEnRight[sfb+sfboffs] = FL2FXCONST_DBL(0.0f);
 756
 757       *msDigest = MS_SOME;
 758     }
 759   }
 760 }
 761