"Fossies" - the Fresh Open Source Software Archive

Member "faac-1_30/libfaac/quantize.c" (16 Oct 2019, 12553 Bytes) of package /linux/misc/faac-1_30.tar.gz:


As a special service "Fossies" has tried to format the requested source page into HTML format using (guessed) C and C++ source code syntax highlighting (style: standard) with prefixed line numbers and code folding option. Alternatively you can here view or download the uninterpreted source code file. For more information about "quantize.c" see the Fossies "Dox" file reference documentation.

    1 /****************************************************************************
    2     Quantizer core functions
    3     quality setting, error distribution, etc.
    4 
    5     Copyright (C) 2017 Krzysztof Nikiel
    6 
    7     This program is free software: you can redistribute it and/or modify
    8     it under the terms of the GNU General Public License as published by
    9     the Free Software Foundation, either version 3 of the License, or
   10     (at your option) any later version.
   11 
   12     This program is distributed in the hope that it will be useful,
   13     but WITHOUT ANY WARRANTY; without even the implied warranty of
   14     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
   15     GNU General Public License for more details.
   16 
   17     You should have received a copy of the GNU General Public License
   18     along with this program.  If not, see <http://www.gnu.org/licenses/>.
   19 ****************************************************************************/
   20 
   21 #include <math.h>
   22 #include <stdio.h>
   23 #include "quantize.h"
   24 #include "huff2.h"
   25 
   26 #ifdef HAVE_IMMINTRIN_H
   27 # include <immintrin.h>
   28 #endif
   29 
   30 #ifdef __SSE2__
   31 # ifdef __GNUC__
   32 #  include <cpuid.h>
   33 # endif
   34 #endif
   35 
   36 #ifdef _MSC_VER
   37 # include <immintrin.h>
   38 # include <intrin.h>
   39 # define __SSE2__
   40 # define bit_SSE2 (1 << 26)
   41 #endif
   42 
   43 #ifdef __GNUC__
   44 #define GCC_VERSION (__GNUC__ * 10000 \
   45                      + __GNUC_MINOR__ * 100 \
   46                      + __GNUC_PATCHLEVEL__)
   47 #endif
   48 
   49 #define MAGIC_NUMBER  0.4054
   50 #define NOISEFLOOR 0.4
   51 
   52 // band sound masking
   53 static void bmask(CoderInfo *coderInfo, double *xr0, double *bandqual,
   54                   int gnum, double quality)
   55 {
   56   int sfb, start, end, cnt;
   57   int *cb_offset = coderInfo->sfb_offset;
   58   int last;
   59   double avgenrg;
   60   double powm = 0.4;
   61   double totenrg = 0.0;
   62   int gsize = coderInfo->groups.len[gnum];
   63   double *xr;
   64   int win;
   65   int enrgcnt = 0;
   66 
   67 
   68   for (sfb = 0; sfb < coderInfo->sfbn; sfb++)
   69   {
   70       start = coderInfo->sfb_offset[sfb];
   71       end = coderInfo->sfb_offset[sfb + 1];
   72 
   73       xr = xr0;
   74       for (win = 0; win < gsize; win++)
   75       {
   76           for (cnt = start; cnt < end; cnt++)
   77           {
   78               totenrg += xr[cnt] * xr[cnt];
   79               enrgcnt++;
   80           }
   81 
   82           xr += BLOCK_LEN_SHORT;
   83       }
   84   }
   85 
   86   if (totenrg < ((NOISEFLOOR * NOISEFLOOR) * (double)enrgcnt))
   87   {
   88       for (sfb = 0; sfb < coderInfo->sfbn; sfb++)
   89           bandqual[sfb] = 0.0;
   90 
   91       return;
   92   }
   93 
   94   for (sfb = 0; sfb < coderInfo->sfbn; sfb++)
   95   {
   96     double avge, maxe;
   97     double target;
   98 
   99     start = cb_offset[sfb];
  100     end = cb_offset[sfb + 1];
  101 
  102     avge = 0.0;
  103     maxe = 0.0;
  104     xr = xr0;
  105     for (win = 0; win < gsize; win++)
  106     {
  107         for (cnt = start; cnt < end; cnt++)
  108         {
  109             double e = xr[cnt]*xr[cnt];
  110             avge += e;
  111             if (maxe < e)
  112                 maxe = e;
  113         }
  114         xr += BLOCK_LEN_SHORT;
  115     }
  116     maxe *= gsize;
  117 
  118 #define NOISETONE 0.2
  119     if (coderInfo->block_type == ONLY_SHORT_WINDOW)
  120     {
  121         last = BLOCK_LEN_SHORT;
  122         avgenrg = totenrg / last;
  123         avgenrg *= end - start;
  124 
  125         target = NOISETONE * pow(avge/avgenrg, powm);
  126         target += (1.0 - NOISETONE) * 0.45 * pow(maxe/avgenrg, powm);
  127 
  128         target *= 1.5;
  129     }
  130     else
  131     {
  132         last = BLOCK_LEN_LONG;
  133         avgenrg = totenrg / last;
  134         avgenrg *= end - start;
  135 
  136         target = NOISETONE * pow(avge/avgenrg, powm);
  137         target += (1.0 - NOISETONE) * 0.45 * pow(maxe/avgenrg, powm);
  138     }
  139 
  140     target *= 10.0 / (1.0 + ((double)(start+end)/last));
  141 
  142     bandqual[sfb] = target * quality;
  143   }
  144 }
  145 
  146 enum {MAXSHORTBAND = 36};
  147 // use band quality levels to quantize a group of windows
  148 static void qlevel(CoderInfo *coderInfo,
  149                    const double *xr0,
  150                    const double *bandqual,
  151                    int gnum,
  152                    int pnslevel
  153                   )
  154 {
  155     int sb, cnt;
  156 #if !defined(__clang__) && defined(__GNUC__) && (GCC_VERSION >= 40600)
  157     /* 2^0.25 (1.50515 dB) step from AAC specs */
  158     static const double sfstep = 1.0 / log10(sqrt(sqrt(2.0)));
  159 #else
  160     static const double sfstep = 20 / 1.50515;
  161 #endif
  162     int gsize = coderInfo->groups.len[gnum];
  163     double pnsthr = 0.1 * pnslevel;
  164 #ifdef __SSE2__
  165     int cpuid[4];
  166     int sse2 = 0;
  167 
  168     cpuid[3] = 0;
  169 # ifdef __GNUC__
  170     __cpuid(1, cpuid[0], cpuid[1], cpuid[2], cpuid[3]);
  171 # endif
  172 # ifdef _MSC_VER
  173     __cpuid(cpuid, 1);
  174 # endif
  175     if (cpuid[3] & bit_SSE2)
  176         sse2 = 1;
  177 #endif
  178 
  179     for (sb = 0; sb < coderInfo->sfbn; sb++)
  180     {
  181       double sfacfix;
  182       int sfac;
  183       double rmsx;
  184       double etot;
  185       int xitab[8 * MAXSHORTBAND];
  186       int *xi;
  187       int start, end;
  188       const double *xr;
  189       int win;
  190 
  191       if (coderInfo->book[coderInfo->bandcnt] != HCB_NONE)
  192       {
  193           coderInfo->bandcnt++;
  194           continue;
  195       }
  196 
  197       start = coderInfo->sfb_offset[sb];
  198       end = coderInfo->sfb_offset[sb+1];
  199 
  200       etot = 0.0;
  201       xr = xr0;
  202       for (win = 0; win < gsize; win++)
  203       {
  204           for (cnt = start; cnt < end; cnt++)
  205           {
  206               double e = xr[cnt] * xr[cnt];
  207               etot += e;
  208           }
  209           xr += BLOCK_LEN_SHORT;
  210       }
  211       etot /= (double)gsize;
  212       rmsx = sqrt(etot / (end - start));
  213 
  214       if ((rmsx < NOISEFLOOR) || (!bandqual[sb]))
  215       {
  216           coderInfo->book[coderInfo->bandcnt++] = HCB_ZERO;
  217           continue;
  218       }
  219 
  220 #ifndef DRM
  221       if (bandqual[sb] < pnsthr)
  222       {
  223           coderInfo->book[coderInfo->bandcnt] = HCB_PNS;
  224           coderInfo->sf[coderInfo->bandcnt] +=
  225               lrint(log10(etot) * (0.5 * sfstep));
  226           coderInfo->bandcnt++;
  227           continue;
  228       }
  229 #endif
  230 
  231       sfac = lrint(log10(bandqual[sb] / rmsx) * sfstep);
  232       if ((SF_OFFSET - sfac) < 10)
  233           sfacfix = 0.0;
  234       else
  235           sfacfix = pow(10, sfac / sfstep);
  236 
  237       xr = xr0 + start;
  238       end -= start;
  239       xi = xitab;
  240       for (win = 0; win < gsize; win++)
  241       {
  242 #ifdef __SSE2__
  243           if (sse2)
  244           {
  245               for (cnt = 0; cnt < end; cnt += 4)
  246               {
  247                   __m128 x = {xr[cnt], xr[cnt + 1], xr[cnt + 2], xr[cnt + 3]};
  248 
  249                   x = _mm_max_ps(x, _mm_sub_ps((__m128){0, 0, 0, 0}, x));
  250                   x = _mm_mul_ps(x, (__m128){sfacfix, sfacfix, sfacfix, sfacfix});
  251                   x = _mm_mul_ps(x, _mm_sqrt_ps(x));
  252                   x = _mm_sqrt_ps(x);
  253                   x = _mm_add_ps(x, (__m128){MAGIC_NUMBER, MAGIC_NUMBER, MAGIC_NUMBER, MAGIC_NUMBER});
  254 
  255                   *(__m128i*)(xi + cnt) = _mm_cvttps_epi32(x);
  256               }
  257               for (cnt = 0; cnt < end; cnt++)
  258               {
  259                   if (xr[cnt] < 0)
  260                       xi[cnt] = -xi[cnt];
  261               }
  262               xi += cnt;
  263               xr += BLOCK_LEN_SHORT;
  264               continue;
  265           }
  266 #endif
  267 
  268           for (cnt = 0; cnt < end; cnt++)
  269           {
  270               double tmp = fabs(xr[cnt]);
  271 
  272               tmp *= sfacfix;
  273               tmp = sqrt(tmp * sqrt(tmp));
  274 
  275               xi[cnt] = (int)(tmp + MAGIC_NUMBER);
  276               if (xr[cnt] < 0)
  277                   xi[cnt] = -xi[cnt];
  278           }
  279           xi += cnt;
  280           xr += BLOCK_LEN_SHORT;
  281       }
  282       huffbook(coderInfo, xitab, gsize * end);
  283       coderInfo->sf[coderInfo->bandcnt++] += SF_OFFSET - sfac;
  284     }
  285 }
  286 
  287 int BlocQuant(CoderInfo *coder, double *xr, AACQuantCfg *aacquantCfg)
  288 {
  289     double bandlvl[MAX_SCFAC_BANDS];
  290     int cnt;
  291     double *gxr;
  292 
  293     coder->global_gain = 0;
  294 
  295     coder->bandcnt = 0;
  296     coder->datacnt = 0;
  297 #ifdef DRM
  298     coder->iLenReordSpData = 0; /* init length of reordered spectral data */
  299     coder->iLenLongestCW = 0; /* init length of longest codeword */
  300     coder->cur_cw = 0; /* init codeword counter */
  301 #endif
  302 
  303     {
  304         int lastis;
  305         int lastsf;
  306 
  307         gxr = xr;
  308         for (cnt = 0; cnt < coder->groups.n; cnt++)
  309         {
  310             bmask(coder, gxr, bandlvl, cnt,
  311                   (double)aacquantCfg->quality/DEFQUAL);
  312             qlevel(coder, gxr, bandlvl, cnt, aacquantCfg->pnslevel);
  313             gxr += coder->groups.len[cnt] * BLOCK_LEN_SHORT;
  314         }
  315 
  316         coder->global_gain = 0;
  317         for (cnt = 0; cnt < coder->bandcnt; cnt++)
  318         {
  319             int book = coder->book[cnt];
  320             if (!book)
  321                 continue;
  322             if ((book != HCB_INTENSITY) && (book != HCB_INTENSITY2))
  323             {
  324                 coder->global_gain = coder->sf[cnt];
  325                 break;
  326             }
  327         }
  328 
  329         lastsf = coder->global_gain;
  330         lastis = 0;
  331         // fixme: move SF range check to quantizer
  332         for (cnt = 0; cnt < coder->bandcnt; cnt++)
  333         {
  334             int book = coder->book[cnt];
  335             if ((book == HCB_INTENSITY) || (book == HCB_INTENSITY2))
  336             {
  337                 int diff = coder->sf[cnt] - lastis;
  338 
  339                 if (diff < -60)
  340                     diff = -60;
  341                 if (diff > 60)
  342                     diff = 60;
  343 
  344                 lastis += diff;
  345                 coder->sf[cnt] = lastis;
  346             }
  347             else if (book == HCB_ESC)
  348             {
  349                 int diff = coder->sf[cnt] - lastsf;
  350 
  351                 if (diff < -60)
  352                     diff = -60;
  353                 if (diff > 60)
  354                     diff = 60;
  355 
  356                 lastsf += diff;
  357                 coder->sf[cnt] = lastsf;
  358             }
  359         }
  360         return 1;
  361     }
  362     return 0;
  363 }
  364 
  365 void CalcBW(unsigned *bw, int rate, SR_INFO *sr, AACQuantCfg *aacquantCfg)
  366 {
  367     // find max short frame band
  368     int max = *bw * (BLOCK_LEN_SHORT << 1) / rate;
  369     int cnt;
  370     int l;
  371 
  372     l = 0;
  373     for (cnt = 0; cnt < sr->num_cb_short; cnt++)
  374     {
  375         if (l >= max)
  376             break;
  377         l += sr->cb_width_short[cnt];
  378     }
  379     aacquantCfg->max_cbs = cnt;
  380     if (aacquantCfg->pnslevel)
  381         *bw = (double)l * rate / (BLOCK_LEN_SHORT << 1);
  382 
  383     // find max long frame band
  384     max = *bw * (BLOCK_LEN_LONG << 1) / rate;
  385     l = 0;
  386     for (cnt = 0; cnt < sr->num_cb_long; cnt++)
  387     {
  388         if (l >= max)
  389             break;
  390         l += sr->cb_width_long[cnt];
  391     }
  392     aacquantCfg->max_cbl = cnt;
  393     aacquantCfg->max_l = l;
  394 
  395     *bw = (double)l * rate / (BLOCK_LEN_LONG << 1);
  396 }
  397 
  398 enum {MINSFB = 2};
  399 
  400 static void calce(double *xr, int *bands, double e[NSFB_SHORT], int maxsfb,
  401                   int maxl)
  402 {
  403     int sfb;
  404     int l;
  405 
  406     // mute lines above cutoff freq
  407     for (l = maxl; l < bands[maxsfb]; l++)
  408         xr[l] = 0.0;
  409 
  410     for (sfb = MINSFB; sfb < maxsfb; sfb++)
  411     {
  412         e[sfb] = 0;
  413         for (l = bands[sfb]; l < bands[sfb + 1]; l++)
  414             e[sfb] += xr[l] * xr[l];
  415     }
  416 }
  417 
  418 static void resete(double min[NSFB_SHORT], double max[NSFB_SHORT],
  419                    double e[NSFB_SHORT], int maxsfb)
  420 {
  421     int sfb;
  422     for (sfb = MINSFB; sfb < maxsfb; sfb++)
  423         min[sfb] = max[sfb] = e[sfb];
  424 }
  425 
  426 #define PRINTSTAT 0
  427 #if PRINTSTAT
  428 static int groups = 0;
  429 static int frames = 0;
  430 #endif
  431 void BlocGroup(double *xr, CoderInfo *coderInfo, AACQuantCfg *cfg)
  432 {
  433     int win, sfb;
  434     double e[NSFB_SHORT];
  435     double min[NSFB_SHORT];
  436     double max[NSFB_SHORT];
  437     const double thr = 3.0;
  438     int win0;
  439     int fastmin;
  440     int maxsfb, maxl;
  441 
  442     if (coderInfo->block_type != ONLY_SHORT_WINDOW)
  443     {
  444         coderInfo->groups.n = 1;
  445         coderInfo->groups.len[0] = 1;
  446         return;
  447     }
  448 
  449     maxl = cfg->max_l / 8;
  450     maxsfb = cfg->max_cbs;
  451     fastmin = ((maxsfb - MINSFB) * 3) >> 2;
  452 
  453 #ifdef DRM
  454     coderInfo->groups.n = 1;
  455     coderInfo->groups.len[0] = 8;
  456     return;
  457 #endif
  458 
  459 #if PRINTSTAT
  460     frames++;
  461 #endif
  462     calce(xr, coderInfo->sfb_offset, e, maxsfb, maxl);
  463     resete(min, max, e, maxsfb);
  464     win0 = 0;
  465     coderInfo->groups.n = 0;
  466     for (win = 1; win < MAX_SHORT_WINDOWS; win++)
  467     {
  468         int fast = 0;
  469 
  470         calce(xr + win * BLOCK_LEN_SHORT, coderInfo->sfb_offset, e, maxsfb, maxl);
  471         for (sfb = MINSFB; sfb < maxsfb; sfb++)
  472         {
  473             if (min[sfb] > e[sfb])
  474                 min[sfb] = e[sfb];
  475             if (max[sfb] < e[sfb])
  476                 max[sfb] = e[sfb];
  477 
  478             if (max[sfb] > thr * min[sfb])
  479                 fast++;
  480         }
  481         if (fast > fastmin)
  482         {
  483             coderInfo->groups.len[coderInfo->groups.n++] = win - win0;
  484             win0 = win;
  485             resete(min, max, e, maxsfb);
  486         }
  487     }
  488     coderInfo->groups.len[coderInfo->groups.n++] = win - win0;
  489 #if PRINTSTAT
  490     groups += coderInfo->groups.n;
  491 #endif
  492 }
  493 
  494 void BlocStat(void)
  495 {
  496 #if PRINTSTAT
  497     printf("frames:%d; groups:%d; g/f:%f\n", frames, groups, (double)groups/frames);
  498 #endif
  499 }