tor-browser

The Tor Browser
git clone https://git.dasho.dev/tor-browser.git
Log | Files | Refs | README | LICENSE

vorbis_floor1.c (30922B)


      1 /********************************************************************
      2 *                                                                  *
      3 * THIS FILE IS PART OF THE OggVorbis SOFTWARE CODEC SOURCE CODE.   *
      4 * USE, DISTRIBUTION AND REPRODUCTION OF THIS LIBRARY SOURCE IS     *
      5 * GOVERNED BY A BSD-STYLE SOURCE LICENSE INCLUDED WITH THIS SOURCE *
      6 * IN 'COPYING'. PLEASE READ THESE TERMS BEFORE DISTRIBUTING.       *
      7 *                                                                  *
      8 * THE OggVorbis SOURCE CODE IS (C) COPYRIGHT 1994-2015             *
      9 * by the Xiph.Org Foundation https://xiph.org/                     *
     10 *                                                                  *
     11 ********************************************************************
     12 
     13 function: floor backend 1 implementation
     14 
     15 ********************************************************************/
     16 
     17 #include <stdlib.h>
     18 #include <string.h>
     19 #include <math.h>
     20 #include <ogg/ogg.h>
     21 #include "vorbis/codec.h"
     22 #include "codec_internal.h"
     23 #include "registry.h"
     24 #include "codebook.h"
     25 #include "misc.h"
     26 #include "scales.h"
     27 
     28 #include <stdio.h>
     29 
     30 #define floor1_rangedB 140 /* floor 1 fixed at -140dB to 0dB range */
     31 
     32 typedef struct lsfit_acc{
     33  int x0;
     34  int x1;
     35 
     36  int xa;
     37  int ya;
     38  int x2a;
     39  int y2a;
     40  int xya;
     41  int an;
     42 
     43  int xb;
     44  int yb;
     45  int x2b;
     46  int y2b;
     47  int xyb;
     48  int bn;
     49 } lsfit_acc;
     50 
     51 /***********************************************/
     52 
     53 static void floor1_free_info(vorbis_info_floor *i){
     54  vorbis_info_floor1 *info=(vorbis_info_floor1 *)i;
     55  if(info){
     56    memset(info,0,sizeof(*info));
     57    _ogg_free(info);
     58  }
     59 }
     60 
     61 static void floor1_free_look(vorbis_look_floor *i){
     62  vorbis_look_floor1 *look=(vorbis_look_floor1 *)i;
     63  if(look){
     64    /*fprintf(stderr,"floor 1 bit usage %f:%f (%f total)\n",
     65            (float)look->phrasebits/look->frames,
     66            (float)look->postbits/look->frames,
     67            (float)(look->postbits+look->phrasebits)/look->frames);*/
     68 
     69    memset(look,0,sizeof(*look));
     70    _ogg_free(look);
     71  }
     72 }
     73 
     74 static void floor1_pack (vorbis_info_floor *i,oggpack_buffer *opb){
     75  vorbis_info_floor1 *info=(vorbis_info_floor1 *)i;
     76  int j,k;
     77  int count=0;
     78  int rangebits;
     79  int maxposit=info->postlist[1];
     80  int maxclass=-1;
     81 
     82  /* save out partitions */
     83  oggpack_write(opb,info->partitions,5); /* only 0 to 31 legal */
     84  for(j=0;j<info->partitions;j++){
     85    oggpack_write(opb,info->partitionclass[j],4); /* only 0 to 15 legal */
     86    if(maxclass<info->partitionclass[j])maxclass=info->partitionclass[j];
     87  }
     88 
     89  /* save out partition classes */
     90  for(j=0;j<maxclass+1;j++){
     91    oggpack_write(opb,info->class_dim[j]-1,3); /* 1 to 8 */
     92    oggpack_write(opb,info->class_subs[j],2); /* 0 to 3 */
     93    if(info->class_subs[j])oggpack_write(opb,info->class_book[j],8);
     94    for(k=0;k<(1<<info->class_subs[j]);k++)
     95      oggpack_write(opb,info->class_subbook[j][k]+1,8);
     96  }
     97 
     98  /* save out the post list */
     99  oggpack_write(opb,info->mult-1,2);     /* only 1,2,3,4 legal now */
    100  /* maxposit cannot legally be less than 1; this is encode-side, we
    101     can assume our setup is OK */
    102  oggpack_write(opb,ov_ilog(maxposit-1),4);
    103  rangebits=ov_ilog(maxposit-1);
    104 
    105  for(j=0,k=0;j<info->partitions;j++){
    106    count+=info->class_dim[info->partitionclass[j]];
    107    for(;k<count;k++)
    108      oggpack_write(opb,info->postlist[k+2],rangebits);
    109  }
    110 }
    111 
    112 static int icomp(const void *a,const void *b){
    113  return(**(int **)a-**(int **)b);
    114 }
    115 
    116 static vorbis_info_floor *floor1_unpack (vorbis_info *vi,oggpack_buffer *opb){
    117  codec_setup_info     *ci=vi->codec_setup;
    118  int j,k,count=0,maxclass=-1,rangebits;
    119 
    120  vorbis_info_floor1 *info=_ogg_calloc(1,sizeof(*info));
    121  /* read partitions */
    122  info->partitions=oggpack_read(opb,5); /* only 0 to 31 legal */
    123  for(j=0;j<info->partitions;j++){
    124    info->partitionclass[j]=oggpack_read(opb,4); /* only 0 to 15 legal */
    125    if(info->partitionclass[j]<0)goto err_out;
    126    if(maxclass<info->partitionclass[j])maxclass=info->partitionclass[j];
    127  }
    128 
    129  /* read partition classes */
    130  for(j=0;j<maxclass+1;j++){
    131    info->class_dim[j]=oggpack_read(opb,3)+1; /* 1 to 8 */
    132    info->class_subs[j]=oggpack_read(opb,2); /* 0,1,2,3 bits */
    133    if(info->class_subs[j]<0)
    134      goto err_out;
    135    if(info->class_subs[j])info->class_book[j]=oggpack_read(opb,8);
    136    if(info->class_book[j]<0 || info->class_book[j]>=ci->books)
    137      goto err_out;
    138    for(k=0;k<(1<<info->class_subs[j]);k++){
    139      info->class_subbook[j][k]=oggpack_read(opb,8)-1;
    140      if(info->class_subbook[j][k]<-1 || info->class_subbook[j][k]>=ci->books)
    141        goto err_out;
    142    }
    143  }
    144 
    145  /* read the post list */
    146  info->mult=oggpack_read(opb,2)+1;     /* only 1,2,3,4 legal now */
    147  rangebits=oggpack_read(opb,4);
    148  if(rangebits<0)goto err_out;
    149 
    150  for(j=0,k=0;j<info->partitions;j++){
    151    count+=info->class_dim[info->partitionclass[j]];
    152    if(count>VIF_POSIT) goto err_out;
    153    for(;k<count;k++){
    154      int t=info->postlist[k+2]=oggpack_read(opb,rangebits);
    155      if(t<0 || t>=(1<<rangebits))
    156        goto err_out;
    157    }
    158  }
    159  info->postlist[0]=0;
    160  info->postlist[1]=1<<rangebits;
    161 
    162  /* don't allow repeated values in post list as they'd result in
    163     zero-length segments */
    164  {
    165    int *sortpointer[VIF_POSIT+2];
    166    for(j=0;j<count+2;j++)sortpointer[j]=info->postlist+j;
    167    qsort(sortpointer,count+2,sizeof(*sortpointer),icomp);
    168 
    169    for(j=1;j<count+2;j++)
    170      if(*sortpointer[j-1]==*sortpointer[j])goto err_out;
    171  }
    172 
    173  return(info);
    174 
    175 err_out:
    176  floor1_free_info(info);
    177  return(NULL);
    178 }
    179 
    180 static vorbis_look_floor *floor1_look(vorbis_dsp_state *vd,
    181                                      vorbis_info_floor *in){
    182 
    183  int *sortpointer[VIF_POSIT+2];
    184  vorbis_info_floor1 *info=(vorbis_info_floor1 *)in;
    185  vorbis_look_floor1 *look=_ogg_calloc(1,sizeof(*look));
    186  int i,j,n=0;
    187 
    188  (void)vd;
    189 
    190  look->vi=info;
    191  look->n=info->postlist[1];
    192 
    193  /* we drop each position value in-between already decoded values,
    194     and use linear interpolation to predict each new value past the
    195     edges.  The positions are read in the order of the position
    196     list... we precompute the bounding positions in the lookup.  Of
    197     course, the neighbors can change (if a position is declined), but
    198     this is an initial mapping */
    199 
    200  for(i=0;i<info->partitions;i++)n+=info->class_dim[info->partitionclass[i]];
    201  n+=2;
    202  look->posts=n;
    203 
    204  /* also store a sorted position index */
    205  for(i=0;i<n;i++)sortpointer[i]=info->postlist+i;
    206  qsort(sortpointer,n,sizeof(*sortpointer),icomp);
    207 
    208  /* points from sort order back to range number */
    209  for(i=0;i<n;i++)look->forward_index[i]=sortpointer[i]-info->postlist;
    210  /* points from range order to sorted position */
    211  for(i=0;i<n;i++)look->reverse_index[look->forward_index[i]]=i;
    212  /* we actually need the post values too */
    213  for(i=0;i<n;i++)look->sorted_index[i]=info->postlist[look->forward_index[i]];
    214 
    215  /* quantize values to multiplier spec */
    216  switch(info->mult){
    217  case 1: /* 1024 -> 256 */
    218    look->quant_q=256;
    219    break;
    220  case 2: /* 1024 -> 128 */
    221    look->quant_q=128;
    222    break;
    223  case 3: /* 1024 -> 86 */
    224    look->quant_q=86;
    225    break;
    226  case 4: /* 1024 -> 64 */
    227    look->quant_q=64;
    228    break;
    229  }
    230 
    231  /* discover our neighbors for decode where we don't use fit flags
    232     (that would push the neighbors outward) */
    233  for(i=0;i<n-2;i++){
    234    int lo=0;
    235    int hi=1;
    236    int lx=0;
    237    int hx=look->n;
    238    int currentx=info->postlist[i+2];
    239    for(j=0;j<i+2;j++){
    240      int x=info->postlist[j];
    241      if(x>lx && x<currentx){
    242        lo=j;
    243        lx=x;
    244      }
    245      if(x<hx && x>currentx){
    246        hi=j;
    247        hx=x;
    248      }
    249    }
    250    look->loneighbor[i]=lo;
    251    look->hineighbor[i]=hi;
    252  }
    253 
    254  return(look);
    255 }
    256 
    257 static int render_point(int x0,int x1,int y0,int y1,int x){
    258  y0&=0x7fff; /* mask off flag */
    259  y1&=0x7fff;
    260 
    261  {
    262    int dy=y1-y0;
    263    int adx=x1-x0;
    264    int ady=abs(dy);
    265    int err=ady*(x-x0);
    266 
    267    int off=err/adx;
    268    if(dy<0)return(y0-off);
    269    return(y0+off);
    270  }
    271 }
    272 
    273 static int vorbis_dBquant(const float *x){
    274  int i= *x*7.3142857f+1023.5f;
    275  if(i>1023)return(1023);
    276  if(i<0)return(0);
    277  return i;
    278 }
    279 
    280 static const float FLOOR1_fromdB_LOOKUP[256]={
    281  1.0649863e-07F, 1.1341951e-07F, 1.2079015e-07F, 1.2863978e-07F,
    282  1.3699951e-07F, 1.4590251e-07F, 1.5538408e-07F, 1.6548181e-07F,
    283  1.7623575e-07F, 1.8768855e-07F, 1.9988561e-07F, 2.128753e-07F,
    284  2.2670913e-07F, 2.4144197e-07F, 2.5713223e-07F, 2.7384213e-07F,
    285  2.9163793e-07F, 3.1059021e-07F, 3.3077411e-07F, 3.5226968e-07F,
    286  3.7516214e-07F, 3.9954229e-07F, 4.2550680e-07F, 4.5315863e-07F,
    287  4.8260743e-07F, 5.1396998e-07F, 5.4737065e-07F, 5.8294187e-07F,
    288  6.2082472e-07F, 6.6116941e-07F, 7.0413592e-07F, 7.4989464e-07F,
    289  7.9862701e-07F, 8.5052630e-07F, 9.0579828e-07F, 9.6466216e-07F,
    290  1.0273513e-06F, 1.0941144e-06F, 1.1652161e-06F, 1.2409384e-06F,
    291  1.3215816e-06F, 1.4074654e-06F, 1.4989305e-06F, 1.5963394e-06F,
    292  1.7000785e-06F, 1.8105592e-06F, 1.9282195e-06F, 2.0535261e-06F,
    293  2.1869758e-06F, 2.3290978e-06F, 2.4804557e-06F, 2.6416497e-06F,
    294  2.8133190e-06F, 2.9961443e-06F, 3.1908506e-06F, 3.3982101e-06F,
    295  3.6190449e-06F, 3.8542308e-06F, 4.1047004e-06F, 4.3714470e-06F,
    296  4.6555282e-06F, 4.9580707e-06F, 5.2802740e-06F, 5.6234160e-06F,
    297  5.9888572e-06F, 6.3780469e-06F, 6.7925283e-06F, 7.2339451e-06F,
    298  7.7040476e-06F, 8.2047000e-06F, 8.7378876e-06F, 9.3057248e-06F,
    299  9.9104632e-06F, 1.0554501e-05F, 1.1240392e-05F, 1.1970856e-05F,
    300  1.2748789e-05F, 1.3577278e-05F, 1.4459606e-05F, 1.5399272e-05F,
    301  1.6400004e-05F, 1.7465768e-05F, 1.8600792e-05F, 1.9809576e-05F,
    302  2.1096914e-05F, 2.2467911e-05F, 2.3928002e-05F, 2.5482978e-05F,
    303  2.7139006e-05F, 2.8902651e-05F, 3.0780908e-05F, 3.2781225e-05F,
    304  3.4911534e-05F, 3.7180282e-05F, 3.9596466e-05F, 4.2169667e-05F,
    305  4.4910090e-05F, 4.7828601e-05F, 5.0936773e-05F, 5.4246931e-05F,
    306  5.7772202e-05F, 6.1526565e-05F, 6.5524908e-05F, 6.9783085e-05F,
    307  7.4317983e-05F, 7.9147585e-05F, 8.4291040e-05F, 8.9768747e-05F,
    308  9.5602426e-05F, 0.00010181521F, 0.00010843174F, 0.00011547824F,
    309  0.00012298267F, 0.00013097477F, 0.00013948625F, 0.00014855085F,
    310  0.00015820453F, 0.00016848555F, 0.00017943469F, 0.00019109536F,
    311  0.00020351382F, 0.00021673929F, 0.00023082423F, 0.00024582449F,
    312  0.00026179955F, 0.00027881276F, 0.00029693158F, 0.00031622787F,
    313  0.00033677814F, 0.00035866388F, 0.00038197188F, 0.00040679456F,
    314  0.00043323036F, 0.00046138411F, 0.00049136745F, 0.00052329927F,
    315  0.00055730621F, 0.00059352311F, 0.00063209358F, 0.00067317058F,
    316  0.00071691700F, 0.00076350630F, 0.00081312324F, 0.00086596457F,
    317  0.00092223983F, 0.00098217216F, 0.0010459992F, 0.0011139742F,
    318  0.0011863665F, 0.0012634633F, 0.0013455702F, 0.0014330129F,
    319  0.0015261382F, 0.0016253153F, 0.0017309374F, 0.0018434235F,
    320  0.0019632195F, 0.0020908006F, 0.0022266726F, 0.0023713743F,
    321  0.0025254795F, 0.0026895994F, 0.0028643847F, 0.0030505286F,
    322  0.0032487691F, 0.0034598925F, 0.0036847358F, 0.0039241906F,
    323  0.0041792066F, 0.0044507950F, 0.0047400328F, 0.0050480668F,
    324  0.0053761186F, 0.0057254891F, 0.0060975636F, 0.0064938176F,
    325  0.0069158225F, 0.0073652516F, 0.0078438871F, 0.0083536271F,
    326  0.0088964928F, 0.009474637F, 0.010090352F, 0.010746080F,
    327  0.011444421F, 0.012188144F, 0.012980198F, 0.013823725F,
    328  0.014722068F, 0.015678791F, 0.016697687F, 0.017782797F,
    329  0.018938423F, 0.020169149F, 0.021479854F, 0.022875735F,
    330  0.024362330F, 0.025945531F, 0.027631618F, 0.029427276F,
    331  0.031339626F, 0.033376252F, 0.035545228F, 0.037855157F,
    332  0.040315199F, 0.042935108F, 0.045725273F, 0.048696758F,
    333  0.051861348F, 0.055231591F, 0.058820850F, 0.062643361F,
    334  0.066714279F, 0.071049749F, 0.075666962F, 0.080584227F,
    335  0.085821044F, 0.091398179F, 0.097337747F, 0.10366330F,
    336  0.11039993F, 0.11757434F, 0.12521498F, 0.13335215F,
    337  0.14201813F, 0.15124727F, 0.16107617F, 0.17154380F,
    338  0.18269168F, 0.19456402F, 0.20720788F, 0.22067342F,
    339  0.23501402F, 0.25028656F, 0.26655159F, 0.28387361F,
    340  0.30232132F, 0.32196786F, 0.34289114F, 0.36517414F,
    341  0.38890521F, 0.41417847F, 0.44109412F, 0.46975890F,
    342  0.50028648F, 0.53279791F, 0.56742212F, 0.60429640F,
    343  0.64356699F, 0.68538959F, 0.72993007F, 0.77736504F,
    344  0.82788260F, 0.88168307F, 0.9389798F, 1.F,
    345 };
    346 
    347 static void render_line(int n, int x0,int x1,int y0,int y1,float *d){
    348  int dy=y1-y0;
    349  int adx=x1-x0;
    350  int ady=abs(dy);
    351  int base=dy/adx;
    352  int sy=(dy<0?base-1:base+1);
    353  int x=x0;
    354  int y=y0;
    355  int err=0;
    356 
    357  ady-=abs(base*adx);
    358 
    359  if(n>x1)n=x1;
    360 
    361  if(x<n)
    362    d[x]*=FLOOR1_fromdB_LOOKUP[y];
    363 
    364  while(++x<n){
    365    err=err+ady;
    366    if(err>=adx){
    367      err-=adx;
    368      y+=sy;
    369    }else{
    370      y+=base;
    371    }
    372    d[x]*=FLOOR1_fromdB_LOOKUP[y];
    373  }
    374 }
    375 
    376 static void render_line0(int n, int x0,int x1,int y0,int y1,int *d){
    377  int dy=y1-y0;
    378  int adx=x1-x0;
    379  int ady=abs(dy);
    380  int base=dy/adx;
    381  int sy=(dy<0?base-1:base+1);
    382  int x=x0;
    383  int y=y0;
    384  int err=0;
    385 
    386  ady-=abs(base*adx);
    387 
    388  if(n>x1)n=x1;
    389 
    390  if(x<n)
    391    d[x]=y;
    392 
    393  while(++x<n){
    394    err=err+ady;
    395    if(err>=adx){
    396      err-=adx;
    397      y+=sy;
    398    }else{
    399      y+=base;
    400    }
    401    d[x]=y;
    402  }
    403 }
    404 
    405 /* the floor has already been filtered to only include relevant sections */
    406 static int accumulate_fit(const float *flr,const float *mdct,
    407                          int x0, int x1,lsfit_acc *a,
    408                          int n,vorbis_info_floor1 *info){
    409  long i;
    410 
    411  int xa=0,ya=0,x2a=0,y2a=0,xya=0,na=0, xb=0,yb=0,x2b=0,y2b=0,xyb=0,nb=0;
    412 
    413  memset(a,0,sizeof(*a));
    414  a->x0=x0;
    415  a->x1=x1;
    416  if(x1>=n)x1=n-1;
    417 
    418  for(i=x0;i<=x1;i++){
    419    int quantized=vorbis_dBquant(flr+i);
    420    if(quantized){
    421      if(mdct[i]+info->twofitatten>=flr[i]){
    422        xa  += i;
    423        ya  += quantized;
    424        x2a += i*i;
    425        y2a += quantized*quantized;
    426        xya += i*quantized;
    427        na++;
    428      }else{
    429        xb  += i;
    430        yb  += quantized;
    431        x2b += i*i;
    432        y2b += quantized*quantized;
    433        xyb += i*quantized;
    434        nb++;
    435      }
    436    }
    437  }
    438 
    439  a->xa=xa;
    440  a->ya=ya;
    441  a->x2a=x2a;
    442  a->y2a=y2a;
    443  a->xya=xya;
    444  a->an=na;
    445 
    446  a->xb=xb;
    447  a->yb=yb;
    448  a->x2b=x2b;
    449  a->y2b=y2b;
    450  a->xyb=xyb;
    451  a->bn=nb;
    452 
    453  return(na);
    454 }
    455 
    456 static int fit_line(lsfit_acc *a,int fits,int *y0,int *y1,
    457                    vorbis_info_floor1 *info){
    458  double xb=0,yb=0,x2b=0,y2b=0,xyb=0,bn=0;
    459  int i;
    460  int x0=a[0].x0;
    461  int x1=a[fits-1].x1;
    462 
    463  for(i=0;i<fits;i++){
    464    double weight = (a[i].bn+a[i].an)*info->twofitweight/(a[i].an+1)+1.;
    465 
    466    xb+=a[i].xb + a[i].xa * weight;
    467    yb+=a[i].yb + a[i].ya * weight;
    468    x2b+=a[i].x2b + a[i].x2a * weight;
    469    y2b+=a[i].y2b + a[i].y2a * weight;
    470    xyb+=a[i].xyb + a[i].xya * weight;
    471    bn+=a[i].bn + a[i].an * weight;
    472  }
    473 
    474  if(*y0>=0){
    475    xb+=   x0;
    476    yb+=  *y0;
    477    x2b+=  x0 *  x0;
    478    y2b+= *y0 * *y0;
    479    xyb+= *y0 *  x0;
    480    bn++;
    481  }
    482 
    483  if(*y1>=0){
    484    xb+=   x1;
    485    yb+=  *y1;
    486    x2b+=  x1 *  x1;
    487    y2b+= *y1 * *y1;
    488    xyb+= *y1 *  x1;
    489    bn++;
    490  }
    491 
    492  {
    493    double denom=(bn*x2b-xb*xb);
    494 
    495    if(denom>0.){
    496      double a=(yb*x2b-xyb*xb)/denom;
    497      double b=(bn*xyb-xb*yb)/denom;
    498      *y0=rint(a+b*x0);
    499      *y1=rint(a+b*x1);
    500 
    501      /* limit to our range! */
    502      if(*y0>1023)*y0=1023;
    503      if(*y1>1023)*y1=1023;
    504      if(*y0<0)*y0=0;
    505      if(*y1<0)*y1=0;
    506 
    507      return 0;
    508    }else{
    509      *y0=0;
    510      *y1=0;
    511      return 1;
    512    }
    513  }
    514 }
    515 
    516 static int inspect_error(int x0,int x1,int y0,int y1,const float *mask,
    517                         const float *mdct,
    518                         vorbis_info_floor1 *info){
    519  int dy=y1-y0;
    520  int adx=x1-x0;
    521  int ady=abs(dy);
    522  int base=dy/adx;
    523  int sy=(dy<0?base-1:base+1);
    524  int x=x0;
    525  int y=y0;
    526  int err=0;
    527  int val=vorbis_dBquant(mask+x);
    528  int mse=0;
    529  int n=0;
    530 
    531  ady-=abs(base*adx);
    532 
    533  mse=(y-val);
    534  mse*=mse;
    535  n++;
    536  if(mdct[x]+info->twofitatten>=mask[x]){
    537    if(y+info->maxover<val)return(1);
    538    if(y-info->maxunder>val)return(1);
    539  }
    540 
    541  while(++x<x1){
    542    err=err+ady;
    543    if(err>=adx){
    544      err-=adx;
    545      y+=sy;
    546    }else{
    547      y+=base;
    548    }
    549 
    550    val=vorbis_dBquant(mask+x);
    551    mse+=((y-val)*(y-val));
    552    n++;
    553    if(mdct[x]+info->twofitatten>=mask[x]){
    554      if(val){
    555        if(y+info->maxover<val)return(1);
    556        if(y-info->maxunder>val)return(1);
    557      }
    558    }
    559  }
    560 
    561  if(info->maxover*info->maxover/n>info->maxerr)return(0);
    562  if(info->maxunder*info->maxunder/n>info->maxerr)return(0);
    563  if(mse/n>info->maxerr)return(1);
    564  return(0);
    565 }
    566 
    567 static int post_Y(int *A,int *B,int pos){
    568  if(A[pos]<0)
    569    return B[pos];
    570  if(B[pos]<0)
    571    return A[pos];
    572 
    573  return (A[pos]+B[pos])>>1;
    574 }
    575 
    576 int *floor1_fit(vorbis_block *vb,vorbis_look_floor1 *look,
    577                          const float *logmdct,   /* in */
    578                          const float *logmask){
    579  long i,j;
    580  vorbis_info_floor1 *info=look->vi;
    581  long n=look->n;
    582  long posts=look->posts;
    583  long nonzero=0;
    584  lsfit_acc fits[VIF_POSIT+1];
    585  int fit_valueA[VIF_POSIT+2]; /* index by range list position */
    586  int fit_valueB[VIF_POSIT+2]; /* index by range list position */
    587 
    588  int loneighbor[VIF_POSIT+2]; /* sorted index of range list position (+2) */
    589  int hineighbor[VIF_POSIT+2];
    590  int *output=NULL;
    591  int memo[VIF_POSIT+2];
    592 
    593  for(i=0;i<posts;i++)fit_valueA[i]=-200; /* mark all unused */
    594  for(i=0;i<posts;i++)fit_valueB[i]=-200; /* mark all unused */
    595  for(i=0;i<posts;i++)loneighbor[i]=0; /* 0 for the implicit 0 post */
    596  for(i=0;i<posts;i++)hineighbor[i]=1; /* 1 for the implicit post at n */
    597  for(i=0;i<posts;i++)memo[i]=-1;      /* no neighbor yet */
    598 
    599  /* quantize the relevant floor points and collect them into line fit
    600     structures (one per minimal division) at the same time */
    601  if(posts==0){
    602    nonzero+=accumulate_fit(logmask,logmdct,0,n,fits,n,info);
    603  }else{
    604    for(i=0;i<posts-1;i++)
    605      nonzero+=accumulate_fit(logmask,logmdct,look->sorted_index[i],
    606                              look->sorted_index[i+1],fits+i,
    607                              n,info);
    608  }
    609 
    610  if(nonzero){
    611    /* start by fitting the implicit base case.... */
    612    int y0=-200;
    613    int y1=-200;
    614    fit_line(fits,posts-1,&y0,&y1,info);
    615 
    616    fit_valueA[0]=y0;
    617    fit_valueB[0]=y0;
    618    fit_valueB[1]=y1;
    619    fit_valueA[1]=y1;
    620 
    621    /* Non degenerate case */
    622    /* start progressive splitting.  This is a greedy, non-optimal
    623       algorithm, but simple and close enough to the best
    624       answer. */
    625    for(i=2;i<posts;i++){
    626      int sortpos=look->reverse_index[i];
    627      int ln=loneighbor[sortpos];
    628      int hn=hineighbor[sortpos];
    629 
    630      /* eliminate repeat searches of a particular range with a memo */
    631      if(memo[ln]!=hn){
    632        /* haven't performed this error search yet */
    633        int lsortpos=look->reverse_index[ln];
    634        int hsortpos=look->reverse_index[hn];
    635        memo[ln]=hn;
    636 
    637        {
    638          /* A note: we want to bound/minimize *local*, not global, error */
    639          int lx=info->postlist[ln];
    640          int hx=info->postlist[hn];
    641          int ly=post_Y(fit_valueA,fit_valueB,ln);
    642          int hy=post_Y(fit_valueA,fit_valueB,hn);
    643 
    644          if(ly==-1 || hy==-1){
    645            exit(1);
    646          }
    647 
    648          if(inspect_error(lx,hx,ly,hy,logmask,logmdct,info)){
    649            /* outside error bounds/begin search area.  Split it. */
    650            int ly0=-200;
    651            int ly1=-200;
    652            int hy0=-200;
    653            int hy1=-200;
    654            int ret0=fit_line(fits+lsortpos,sortpos-lsortpos,&ly0,&ly1,info);
    655            int ret1=fit_line(fits+sortpos,hsortpos-sortpos,&hy0,&hy1,info);
    656 
    657            if(ret0){
    658              ly0=ly;
    659              ly1=hy0;
    660            }
    661            if(ret1){
    662              hy0=ly1;
    663              hy1=hy;
    664            }
    665 
    666            if(ret0 && ret1){
    667              fit_valueA[i]=-200;
    668              fit_valueB[i]=-200;
    669            }else{
    670              /* store new edge values */
    671              fit_valueB[ln]=ly0;
    672              if(ln==0)fit_valueA[ln]=ly0;
    673              fit_valueA[i]=ly1;
    674              fit_valueB[i]=hy0;
    675              fit_valueA[hn]=hy1;
    676              if(hn==1)fit_valueB[hn]=hy1;
    677 
    678              if(ly1>=0 || hy0>=0){
    679                /* store new neighbor values */
    680                for(j=sortpos-1;j>=0;j--)
    681                  if(hineighbor[j]==hn)
    682                    hineighbor[j]=i;
    683                  else
    684                    break;
    685                for(j=sortpos+1;j<posts;j++)
    686                  if(loneighbor[j]==ln)
    687                    loneighbor[j]=i;
    688                  else
    689                    break;
    690              }
    691            }
    692          }else{
    693            fit_valueA[i]=-200;
    694            fit_valueB[i]=-200;
    695          }
    696        }
    697      }
    698    }
    699 
    700    output=_vorbis_block_alloc(vb,sizeof(*output)*posts);
    701 
    702    output[0]=post_Y(fit_valueA,fit_valueB,0);
    703    output[1]=post_Y(fit_valueA,fit_valueB,1);
    704 
    705    /* fill in posts marked as not using a fit; we will zero
    706       back out to 'unused' when encoding them so long as curve
    707       interpolation doesn't force them into use */
    708    for(i=2;i<posts;i++){
    709      int ln=look->loneighbor[i-2];
    710      int hn=look->hineighbor[i-2];
    711      int x0=info->postlist[ln];
    712      int x1=info->postlist[hn];
    713      int y0=output[ln];
    714      int y1=output[hn];
    715 
    716      int predicted=render_point(x0,x1,y0,y1,info->postlist[i]);
    717      int vx=post_Y(fit_valueA,fit_valueB,i);
    718 
    719      if(vx>=0 && predicted!=vx){
    720        output[i]=vx;
    721      }else{
    722        output[i]= predicted|0x8000;
    723      }
    724    }
    725  }
    726 
    727  return(output);
    728 
    729 }
    730 
    731 int *floor1_interpolate_fit(vorbis_block *vb,vorbis_look_floor1 *look,
    732                          int *A,int *B,
    733                          int del){
    734 
    735  long i;
    736  long posts=look->posts;
    737  int *output=NULL;
    738 
    739  if(A && B){
    740    output=_vorbis_block_alloc(vb,sizeof(*output)*posts);
    741 
    742    /* overly simpleminded--- look again post 1.2 */
    743    for(i=0;i<posts;i++){
    744      output[i]=((65536-del)*(A[i]&0x7fff)+del*(B[i]&0x7fff)+32768)>>16;
    745      if(A[i]&0x8000 && B[i]&0x8000)output[i]|=0x8000;
    746    }
    747  }
    748 
    749  return(output);
    750 }
    751 
    752 
    753 int floor1_encode(oggpack_buffer *opb,vorbis_block *vb,
    754                  vorbis_look_floor1 *look,
    755                  int *post,int *ilogmask){
    756 
    757  long i,j;
    758  vorbis_info_floor1 *info=look->vi;
    759  long posts=look->posts;
    760  codec_setup_info *ci=vb->vd->vi->codec_setup;
    761  int out[VIF_POSIT+2];
    762  static_codebook **sbooks=ci->book_param;
    763  codebook *books=ci->fullbooks;
    764 
    765  /* quantize values to multiplier spec */
    766  if(post){
    767    for(i=0;i<posts;i++){
    768      int val=post[i]&0x7fff;
    769      switch(info->mult){
    770      case 1: /* 1024 -> 256 */
    771        val>>=2;
    772        break;
    773      case 2: /* 1024 -> 128 */
    774        val>>=3;
    775        break;
    776      case 3: /* 1024 -> 86 */
    777        val/=12;
    778        break;
    779      case 4: /* 1024 -> 64 */
    780        val>>=4;
    781        break;
    782      }
    783      post[i]=val | (post[i]&0x8000);
    784    }
    785 
    786    out[0]=post[0];
    787    out[1]=post[1];
    788 
    789    /* find prediction values for each post and subtract them */
    790    for(i=2;i<posts;i++){
    791      int ln=look->loneighbor[i-2];
    792      int hn=look->hineighbor[i-2];
    793      int x0=info->postlist[ln];
    794      int x1=info->postlist[hn];
    795      int y0=post[ln];
    796      int y1=post[hn];
    797 
    798      int predicted=render_point(x0,x1,y0,y1,info->postlist[i]);
    799 
    800      if((post[i]&0x8000) || (predicted==post[i])){
    801        post[i]=predicted|0x8000; /* in case there was roundoff jitter
    802                                     in interpolation */
    803        out[i]=0;
    804      }else{
    805        int headroom=(look->quant_q-predicted<predicted?
    806                      look->quant_q-predicted:predicted);
    807 
    808        int val=post[i]-predicted;
    809 
    810        /* at this point the 'deviation' value is in the range +/- max
    811           range, but the real, unique range can always be mapped to
    812           only [0-maxrange).  So we want to wrap the deviation into
    813           this limited range, but do it in the way that least screws
    814           an essentially gaussian probability distribution. */
    815 
    816        if(val<0)
    817          if(val<-headroom)
    818            val=headroom-val-1;
    819          else
    820            val=-1-(val<<1);
    821        else
    822          if(val>=headroom)
    823            val= val+headroom;
    824          else
    825            val<<=1;
    826 
    827        out[i]=val;
    828        post[ln]&=0x7fff;
    829        post[hn]&=0x7fff;
    830      }
    831    }
    832 
    833    /* we have everything we need. pack it out */
    834    /* mark nontrivial floor */
    835    oggpack_write(opb,1,1);
    836 
    837    /* beginning/end post */
    838    look->frames++;
    839    look->postbits+=ov_ilog(look->quant_q-1)*2;
    840    oggpack_write(opb,out[0],ov_ilog(look->quant_q-1));
    841    oggpack_write(opb,out[1],ov_ilog(look->quant_q-1));
    842 
    843 
    844    /* partition by partition */
    845    for(i=0,j=2;i<info->partitions;i++){
    846      int class=info->partitionclass[i];
    847      int cdim=info->class_dim[class];
    848      int csubbits=info->class_subs[class];
    849      int csub=1<<csubbits;
    850      int bookas[8]={0,0,0,0,0,0,0,0};
    851      int cval=0;
    852      int cshift=0;
    853      int k,l;
    854 
    855      /* generate the partition's first stage cascade value */
    856      if(csubbits){
    857        int maxval[8]={0,0,0,0,0,0,0,0}; /* gcc's static analysis
    858                                            issues a warning without
    859                                            initialization */
    860        for(k=0;k<csub;k++){
    861          int booknum=info->class_subbook[class][k];
    862          if(booknum<0){
    863            maxval[k]=1;
    864          }else{
    865            maxval[k]=sbooks[info->class_subbook[class][k]]->entries;
    866          }
    867        }
    868        for(k=0;k<cdim;k++){
    869          for(l=0;l<csub;l++){
    870            int val=out[j+k];
    871            if(val<maxval[l]){
    872              bookas[k]=l;
    873              break;
    874            }
    875          }
    876          cval|= bookas[k]<<cshift;
    877          cshift+=csubbits;
    878        }
    879        /* write it */
    880        look->phrasebits+=
    881          vorbis_book_encode(books+info->class_book[class],cval,opb);
    882 
    883 #ifdef TRAIN_FLOOR1
    884        {
    885          FILE *of;
    886          char buffer[80];
    887          sprintf(buffer,"line_%dx%ld_class%d.vqd",
    888                  vb->pcmend/2,posts-2,class);
    889          of=fopen(buffer,"a");
    890          fprintf(of,"%d\n",cval);
    891          fclose(of);
    892        }
    893 #endif
    894      }
    895 
    896      /* write post values */
    897      for(k=0;k<cdim;k++){
    898        int book=info->class_subbook[class][bookas[k]];
    899        if(book>=0){
    900          /* hack to allow training with 'bad' books */
    901          if(out[j+k]<(books+book)->entries)
    902            look->postbits+=vorbis_book_encode(books+book,
    903                                               out[j+k],opb);
    904          /*else
    905            fprintf(stderr,"+!");*/
    906 
    907 #ifdef TRAIN_FLOOR1
    908          {
    909            FILE *of;
    910            char buffer[80];
    911            sprintf(buffer,"line_%dx%ld_%dsub%d.vqd",
    912                    vb->pcmend/2,posts-2,class,bookas[k]);
    913            of=fopen(buffer,"a");
    914            fprintf(of,"%d\n",out[j+k]);
    915            fclose(of);
    916          }
    917 #endif
    918        }
    919      }
    920      j+=cdim;
    921    }
    922 
    923    {
    924      /* generate quantized floor equivalent to what we'd unpack in decode */
    925      /* render the lines */
    926      int hx=0;
    927      int lx=0;
    928      int ly=post[0]*info->mult;
    929      int n=ci->blocksizes[vb->W]/2;
    930 
    931      for(j=1;j<look->posts;j++){
    932        int current=look->forward_index[j];
    933        int hy=post[current]&0x7fff;
    934        if(hy==post[current]){
    935 
    936          hy*=info->mult;
    937          hx=info->postlist[current];
    938 
    939          render_line0(n,lx,hx,ly,hy,ilogmask);
    940 
    941          lx=hx;
    942          ly=hy;
    943        }
    944      }
    945      for(j=hx;j<vb->pcmend/2;j++)ilogmask[j]=ly; /* be certain */
    946      return(1);
    947    }
    948  }else{
    949    oggpack_write(opb,0,1);
    950    memset(ilogmask,0,vb->pcmend/2*sizeof(*ilogmask));
    951    return(0);
    952  }
    953 }
    954 
    955 static void *floor1_inverse1(vorbis_block *vb,vorbis_look_floor *in){
    956  vorbis_look_floor1 *look=(vorbis_look_floor1 *)in;
    957  vorbis_info_floor1 *info=look->vi;
    958  codec_setup_info   *ci=vb->vd->vi->codec_setup;
    959 
    960  int i,j,k;
    961  codebook *books=ci->fullbooks;
    962 
    963  /* unpack wrapped/predicted values from stream */
    964  if(oggpack_read(&vb->opb,1)==1){
    965    int *fit_value=_vorbis_block_alloc(vb,(look->posts)*sizeof(*fit_value));
    966 
    967    fit_value[0]=oggpack_read(&vb->opb,ov_ilog(look->quant_q-1));
    968    fit_value[1]=oggpack_read(&vb->opb,ov_ilog(look->quant_q-1));
    969 
    970    /* partition by partition */
    971    for(i=0,j=2;i<info->partitions;i++){
    972      int class=info->partitionclass[i];
    973      int cdim=info->class_dim[class];
    974      int csubbits=info->class_subs[class];
    975      int csub=1<<csubbits;
    976      int cval=0;
    977 
    978      /* decode the partition's first stage cascade value */
    979      if(csubbits){
    980        cval=vorbis_book_decode(books+info->class_book[class],&vb->opb);
    981 
    982        if(cval==-1)goto eop;
    983      }
    984 
    985      for(k=0;k<cdim;k++){
    986        int book=info->class_subbook[class][cval&(csub-1)];
    987        cval>>=csubbits;
    988        if(book>=0){
    989          if((fit_value[j+k]=vorbis_book_decode(books+book,&vb->opb))==-1)
    990            goto eop;
    991        }else{
    992          fit_value[j+k]=0;
    993        }
    994      }
    995      j+=cdim;
    996    }
    997 
    998    /* unwrap positive values and reconsitute via linear interpolation */
    999    for(i=2;i<look->posts;i++){
   1000      int predicted=render_point(info->postlist[look->loneighbor[i-2]],
   1001                                 info->postlist[look->hineighbor[i-2]],
   1002                                 fit_value[look->loneighbor[i-2]],
   1003                                 fit_value[look->hineighbor[i-2]],
   1004                                 info->postlist[i]);
   1005      int hiroom=look->quant_q-predicted;
   1006      int loroom=predicted;
   1007      int room=(hiroom<loroom?hiroom:loroom)<<1;
   1008      int val=fit_value[i];
   1009 
   1010      if(val){
   1011        if(val>=room){
   1012          if(hiroom>loroom){
   1013            val = val-loroom;
   1014          }else{
   1015            val = -1-(val-hiroom);
   1016          }
   1017        }else{
   1018          if(val&1){
   1019            val= -((val+1)>>1);
   1020          }else{
   1021            val>>=1;
   1022          }
   1023        }
   1024 
   1025        fit_value[i]=(val+predicted)&0x7fff;
   1026        fit_value[look->loneighbor[i-2]]&=0x7fff;
   1027        fit_value[look->hineighbor[i-2]]&=0x7fff;
   1028 
   1029      }else{
   1030        fit_value[i]=predicted|0x8000;
   1031      }
   1032 
   1033    }
   1034 
   1035    return(fit_value);
   1036  }
   1037 eop:
   1038  return(NULL);
   1039 }
   1040 
   1041 static int floor1_inverse2(vorbis_block *vb,vorbis_look_floor *in,void *memo,
   1042                          float *out){
   1043  vorbis_look_floor1 *look=(vorbis_look_floor1 *)in;
   1044  vorbis_info_floor1 *info=look->vi;
   1045 
   1046  codec_setup_info   *ci=vb->vd->vi->codec_setup;
   1047  int                  n=ci->blocksizes[vb->W]/2;
   1048  int j;
   1049 
   1050  if(memo){
   1051    /* render the lines */
   1052    int *fit_value=(int *)memo;
   1053    int hx=0;
   1054    int lx=0;
   1055    int ly=fit_value[0]*info->mult;
   1056    /* guard lookup against out-of-range values */
   1057    ly=(ly<0?0:ly>255?255:ly);
   1058 
   1059    for(j=1;j<look->posts;j++){
   1060      int current=look->forward_index[j];
   1061      int hy=fit_value[current]&0x7fff;
   1062      if(hy==fit_value[current]){
   1063 
   1064        hx=info->postlist[current];
   1065        hy*=info->mult;
   1066        /* guard lookup against out-of-range values */
   1067        hy=(hy<0?0:hy>255?255:hy);
   1068 
   1069        render_line(n,lx,hx,ly,hy,out);
   1070 
   1071        lx=hx;
   1072        ly=hy;
   1073      }
   1074    }
   1075    for(j=hx;j<n;j++)out[j]*=FLOOR1_fromdB_LOOKUP[ly]; /* be certain */
   1076    return(1);
   1077  }
   1078  memset(out,0,sizeof(*out)*n);
   1079  return(0);
   1080 }
   1081 
   1082 /* export hooks */
   1083 const vorbis_func_floor floor1_exportbundle={
   1084  &floor1_pack,&floor1_unpack,&floor1_look,&floor1_free_info,
   1085  &floor1_free_look,&floor1_inverse1,&floor1_inverse2
   1086 };