lib/mapping0.c

   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-2009             *
   9  * by the Xiph.Org Foundation http://www.xiph.org/                  *
  10  *                                                                  *
  11  ********************************************************************
  12
  13  function: channel mapping 0 implementation
  14  last mod: $Id$
  15
  16  ********************************************************************/
  17
  18 #include <stdlib.h>
  19 #include <stdio.h>
  20 #include <string.h>
  21 #include <math.h>
  22 #include <ogg/ogg.h>
  23 #include "vorbis/codec.h"
  24 #include "codec_internal.h"
  25 #include "codebook.h"
  26 #include "window.h"
  27 #include "registry.h"
  28 #include "psy.h"
  29 #include "misc.h"
  30
  31 /* simplistic, wasteful way of doing this (unique lookup for each
  32    mode/submapping); there should be a central repository for
  33    identical lookups.  That will require minor work, so I'm putting it
  34    off as low priority.
  35
  36    Why a lookup for each backend in a given mode?  Because the
  37    blocksize is set by the mode, and low backend lookups may require
  38    parameters from other areas of the mode/mapping */
  39
  40 static void mapping0_free_info(vorbis_info_mapping *i){
  41   vorbis_info_mapping0 *info=(vorbis_info_mapping0 *)i;
  42   if(info){
  43     memset(info,0,sizeof(*info));
  44     _ogg_free(info);
  45   }
  46 }
  47
  48 static int ilog(unsigned int v){
  49   int ret=0;
  50   if(v)--v;
  51   while(v){
  52     ret++;
  53     v>>=1;
  54   }
  55   return(ret);
  56 }
  57
  58 static void mapping0_pack(vorbis_info *vi,vorbis_info_mapping *vm,
  59                           oggpack_buffer *opb){
  60   int i;
  61   vorbis_info_mapping0 *info=(vorbis_info_mapping0 *)vm;
  62
  63   /* another 'we meant to do it this way' hack...  up to beta 4, we
  64      packed 4 binary zeros here to signify one submapping in use.  We
  65      now redefine that to mean four bitflags that indicate use of
  66      deeper features; bit0:submappings, bit1:coupling,
  67      bit2,3:reserved. This is backward compatable with all actual uses
  68      of the beta code. */
  69
  70   if(info->submaps>1){
  71     oggpack_write(opb,1,1);
  72     oggpack_write(opb,info->submaps-1,4);
  73   }else
  74     oggpack_write(opb,0,1);
  75
  76   if(info->coupling_steps>0){
  77     oggpack_write(opb,1,1);
  78     oggpack_write(opb,info->coupling_steps-1,8);
  79
  80     for(i=0;i<info->coupling_steps;i++){
  81       oggpack_write(opb,info->coupling_mag[i],ilog(vi->channels));
  82       oggpack_write(opb,info->coupling_ang[i],ilog(vi->channels));
  83     }
  84   }else
  85     oggpack_write(opb,0,1);
  86
  87   oggpack_write(opb,0,2); /* 2,3:reserved */
  88
  89   /* we don't write the channel submappings if we only have one... */
  90   if(info->submaps>1){
  91     for(i=0;i<vi->channels;i++)
  92       oggpack_write(opb,info->chmuxlist[i],4);
  93   }
  94   for(i=0;i<info->submaps;i++){
  95     oggpack_write(opb,0,8); /* time submap unused */
  96     oggpack_write(opb,info->floorsubmap[i],8);
  97     oggpack_write(opb,info->residuesubmap[i],8);
  98   }
  99 }
 100
 101 /* also responsible for range checking */
 102 static vorbis_info_mapping *mapping0_unpack(vorbis_info *vi,oggpack_buffer *opb){
 103   int i,b;
 104   vorbis_info_mapping0 *info=_ogg_calloc(1,sizeof(*info));
 105   codec_setup_info     *ci=vi->codec_setup;
 106   memset(info,0,sizeof(*info));
 107
 108   b=oggpack_read(opb,1);
 109   if(b<0)goto err_out;
 110   if(b){
 111     info->submaps=oggpack_read(opb,4)+1;
 112     if(info->submaps<=0)goto err_out;
 113   }else
 114     info->submaps=1;
 115
 116   b=oggpack_read(opb,1);
 117   if(b<0)goto err_out;
 118   if(b){
 119     info->coupling_steps=oggpack_read(opb,8)+1;
 120     if(info->coupling_steps<=0)goto err_out;
 121     for(i=0;i<info->coupling_steps;i++){
 122       int testM=info->coupling_mag[i]=oggpack_read(opb,ilog(vi->channels));
 123       int testA=info->coupling_ang[i]=oggpack_read(opb,ilog(vi->channels));
 124
 125       if(testM<0 ||
 126          testA<0 ||
 127          testM==testA ||
 128          testM>=vi->channels ||
 129          testA>=vi->channels) goto err_out;
 130     }
 131
 132   }
 133
 134   if(oggpack_read(opb,2)!=0)goto err_out; /* 2,3:reserved */
 135
 136   if(info->submaps>1){
 137     for(i=0;i<vi->channels;i++){
 138       info->chmuxlist[i]=oggpack_read(opb,4);
 139       if(info->chmuxlist[i]>=info->submaps || info->chmuxlist[i]<0)goto err_out;
 140     }
 141   }
 142   for(i=0;i<info->submaps;i++){
 143     oggpack_read(opb,8); /* time submap unused */
 144     info->floorsubmap[i]=oggpack_read(opb,8);
 145     if(info->floorsubmap[i]>=ci->floors || info->floorsubmap[i]<0)goto err_out;
 146     info->residuesubmap[i]=oggpack_read(opb,8);
 147     if(info->residuesubmap[i]>=ci->residues || info->residuesubmap[i]<0)goto err_out;
 148   }
 149
 150   return info;
 151
 152  err_out:
 153   mapping0_free_info(info);
 154   return(NULL);
 155 }
 156
 157 #include "os.h"
 158 #include "lpc.h"
 159 #include "lsp.h"
 160 #include "envelope.h"
 161 #include "mdct.h"
 162 #include "psy.h"
 163 #include "scales.h"
 164
 165 #if 0
 166 static long seq=0;
 167 static ogg_int64_t total=0;
 168 static float FLOOR1_fromdB_LOOKUP[256]={
 169   1.0649863e-07F, 1.1341951e-07F, 1.2079015e-07F, 1.2863978e-07F,
 170   1.3699951e-07F, 1.4590251e-07F, 1.5538408e-07F, 1.6548181e-07F,
 171   1.7623575e-07F, 1.8768855e-07F, 1.9988561e-07F, 2.128753e-07F,
 172   2.2670913e-07F, 2.4144197e-07F, 2.5713223e-07F, 2.7384213e-07F,
 173   2.9163793e-07F, 3.1059021e-07F, 3.3077411e-07F, 3.5226968e-07F,
 174   3.7516214e-07F, 3.9954229e-07F, 4.2550680e-07F, 4.5315863e-07F,
 175   4.8260743e-07F, 5.1396998e-07F, 5.4737065e-07F, 5.8294187e-07F,
 176   6.2082472e-07F, 6.6116941e-07F, 7.0413592e-07F, 7.4989464e-07F,
 177   7.9862701e-07F, 8.5052630e-07F, 9.0579828e-07F, 9.6466216e-07F,
 178   1.0273513e-06F, 1.0941144e-06F, 1.1652161e-06F, 1.2409384e-06F,
 179   1.3215816e-06F, 1.4074654e-06F, 1.4989305e-06F, 1.5963394e-06F,
 180   1.7000785e-06F, 1.8105592e-06F, 1.9282195e-06F, 2.0535261e-06F,
 181   2.1869758e-06F, 2.3290978e-06F, 2.4804557e-06F, 2.6416497e-06F,
 182   2.8133190e-06F, 2.9961443e-06F, 3.1908506e-06F, 3.3982101e-06F,
 183   3.6190449e-06F, 3.8542308e-06F, 4.1047004e-06F, 4.3714470e-06F,
 184   4.6555282e-06F, 4.9580707e-06F, 5.2802740e-06F, 5.6234160e-06F,
 185   5.9888572e-06F, 6.3780469e-06F, 6.7925283e-06F, 7.2339451e-06F,
 186   7.7040476e-06F, 8.2047000e-06F, 8.7378876e-06F, 9.3057248e-06F,
 187   9.9104632e-06F, 1.0554501e-05F, 1.1240392e-05F, 1.1970856e-05F,
 188   1.2748789e-05F, 1.3577278e-05F, 1.4459606e-05F, 1.5399272e-05F,
 189   1.6400004e-05F, 1.7465768e-05F, 1.8600792e-05F, 1.9809576e-05F,
 190   2.1096914e-05F, 2.2467911e-05F, 2.3928002e-05F, 2.5482978e-05F,
 191   2.7139006e-05F, 2.8902651e-05F, 3.0780908e-05F, 3.2781225e-05F,
 192   3.4911534e-05F, 3.7180282e-05F, 3.9596466e-05F, 4.2169667e-05F,
 193   4.4910090e-05F, 4.7828601e-05F, 5.0936773e-05F, 5.4246931e-05F,
 194   5.7772202e-05F, 6.1526565e-05F, 6.5524908e-05F, 6.9783085e-05F,
 195   7.4317983e-05F, 7.9147585e-05F, 8.4291040e-05F, 8.9768747e-05F,
 196   9.5602426e-05F, 0.00010181521F, 0.00010843174F, 0.00011547824F,
 197   0.00012298267F, 0.00013097477F, 0.00013948625F, 0.00014855085F,
 198   0.00015820453F, 0.00016848555F, 0.00017943469F, 0.00019109536F,
 199   0.00020351382F, 0.00021673929F, 0.00023082423F, 0.00024582449F,
 200   0.00026179955F, 0.00027881276F, 0.00029693158F, 0.00031622787F,
 201   0.00033677814F, 0.00035866388F, 0.00038197188F, 0.00040679456F,
 202   0.00043323036F, 0.00046138411F, 0.00049136745F, 0.00052329927F,
 203   0.00055730621F, 0.00059352311F, 0.00063209358F, 0.00067317058F,
 204   0.00071691700F, 0.00076350630F, 0.00081312324F, 0.00086596457F,
 205   0.00092223983F, 0.00098217216F, 0.0010459992F, 0.0011139742F,
 206   0.0011863665F, 0.0012634633F, 0.0013455702F, 0.0014330129F,
 207   0.0015261382F, 0.0016253153F, 0.0017309374F, 0.0018434235F,
 208   0.0019632195F, 0.0020908006F, 0.0022266726F, 0.0023713743F,
 209   0.0025254795F, 0.0026895994F, 0.0028643847F, 0.0030505286F,
 210   0.0032487691F, 0.0034598925F, 0.0036847358F, 0.0039241906F,
 211   0.0041792066F, 0.0044507950F, 0.0047400328F, 0.0050480668F,
 212   0.0053761186F, 0.0057254891F, 0.0060975636F, 0.0064938176F,
 213   0.0069158225F, 0.0073652516F, 0.0078438871F, 0.0083536271F,
 214   0.0088964928F, 0.009474637F, 0.010090352F, 0.010746080F,
 215   0.011444421F, 0.012188144F, 0.012980198F, 0.013823725F,
 216   0.014722068F, 0.015678791F, 0.016697687F, 0.017782797F,
 217   0.018938423F, 0.020169149F, 0.021479854F, 0.022875735F,
 218   0.024362330F, 0.025945531F, 0.027631618F, 0.029427276F,
 219   0.031339626F, 0.033376252F, 0.035545228F, 0.037855157F,
 220   0.040315199F, 0.042935108F, 0.045725273F, 0.048696758F,
 221   0.051861348F, 0.055231591F, 0.058820850F, 0.062643361F,
 222   0.066714279F, 0.071049749F, 0.075666962F, 0.080584227F,
 223   0.085821044F, 0.091398179F, 0.097337747F, 0.10366330F,
 224   0.11039993F, 0.11757434F, 0.12521498F, 0.13335215F,
 225   0.14201813F, 0.15124727F, 0.16107617F, 0.17154380F,
 226   0.18269168F, 0.19456402F, 0.20720788F, 0.22067342F,
 227   0.23501402F, 0.25028656F, 0.26655159F, 0.28387361F,
 228   0.30232132F, 0.32196786F, 0.34289114F, 0.36517414F,
 229   0.38890521F, 0.41417847F, 0.44109412F, 0.46975890F,
 230   0.50028648F, 0.53279791F, 0.56742212F, 0.60429640F,
 231   0.64356699F, 0.68538959F, 0.72993007F, 0.77736504F,
 232   0.82788260F, 0.88168307F, 0.9389798F, 1.F,
 233 };
 234
 235 #endif
 236
 237
 238 static int mapping0_forward(vorbis_block *vb){
 239   vorbis_dsp_state      *vd=vb->vd;
 240   vorbis_info           *vi=vd->vi;
 241   codec_setup_info      *ci=vi->codec_setup;
 242   private_state         *b=vb->vd->backend_state;
 243   vorbis_block_internal *vbi=(vorbis_block_internal *)vb->internal;
 244   int                    n=vb->pcmend;
 245   int i,j,k;
 246
 247   int    *nonzero    = alloca(sizeof(*nonzero)*vi->channels);
 248   float  **gmdct     = _vorbis_block_alloc(vb,vi->channels*sizeof(*gmdct));
 249   int    **iwork      = _vorbis_block_alloc(vb,vi->channels*sizeof(*iwork));
 250   int ***floor_posts = _vorbis_block_alloc(vb,vi->channels*sizeof(*floor_posts));
 251
 252   float global_ampmax=vbi->ampmax;
 253   float *local_ampmax=alloca(sizeof(*local_ampmax)*vi->channels);
 254   int blocktype=vbi->blocktype;
 255
 256   int modenumber=vb->W;
 257   vorbis_info_mapping0 *info=ci->map_param[modenumber];
 258   vorbis_look_psy *psy_look=
 259     b->psy+blocktype+(vb->W?2:0);
 260
 261   vb->mode=modenumber;
 262
 263   for(i=0;i<vi->channels;i++){
 264     float scale=4.f/n;
 265     float scale_dB;
 266
 267     float *pcm     =vb->pcm[i];
 268     float *logfft  =pcm;
 269
 270     gmdct[i]=_vorbis_block_alloc(vb,n/2*sizeof(**gmdct));
 271
 272     scale_dB=todB(&scale) + .345; /* + .345 is a hack; the original
 273                                      todB estimation used on IEEE 754
 274                                      compliant machines had a bug that
 275                                      returned dB values about a third
 276                                      of a decibel too high.  The bug
 277                                      was harmless because tunings
 278                                      implicitly took that into
 279                                      account.  However, fixing the bug
 280                                      in the estimator requires
 281                                      changing all the tunings as well.
 282                                      For now, it's easier to sync
 283                                      things back up here, and
 284                                      recalibrate the tunings in the
 285                                      next major model upgrade. */
 286
 287 #if 0
 288     if(vi->channels==2){
 289       if(i==0)
 290         _analysis_output("pcmL",seq,pcm,n,0,0,total-n/2);
 291       else
 292         _analysis_output("pcmR",seq,pcm,n,0,0,total-n/2);
 293     }else{
 294       _analysis_output("pcm",seq,pcm,n,0,0,total-n/2);
 295     }
 296 #endif
 297
 298     /* window the PCM data */
 299     _vorbis_apply_window(pcm,b->window,ci->blocksizes,vb->lW,vb->W,vb->nW);
 300
 301 #if 0
 302     if(vi->channels==2){
 303       if(i==0)
 304         _analysis_output("windowedL",seq,pcm,n,0,0,total-n/2);
 305       else
 306         _analysis_output("windowedR",seq,pcm,n,0,0,total-n/2);
 307     }else{
 308       _analysis_output("windowed",seq,pcm,n,0,0,total-n/2);
 309     }
 310 #endif
 311
 312     /* transform the PCM data */
 313     /* only MDCT right now.... */
 314     mdct_forward(b->transform[vb->W][0],pcm,gmdct[i]);
 315
 316     /* FFT yields more accurate tonal estimation (not phase sensitive) */
 317     drft_forward(&b->fft_look[vb->W],pcm);
 318     logfft[0]=scale_dB+todB(pcm)  + .345; /* + .345 is a hack; the
 319                                      original todB estimation used on
 320                                      IEEE 754 compliant machines had a
 321                                      bug that returned dB values about
 322                                      a third of a decibel too high.
 323                                      The bug was harmless because
 324                                      tunings implicitly took that into
 325                                      account.  However, fixing the bug
 326                                      in the estimator requires
 327                                      changing all the tunings as well.
 328                                      For now, it's easier to sync
 329                                      things back up here, and
 330                                      recalibrate the tunings in the
 331                                      next major model upgrade. */
 332     local_ampmax[i]=logfft[0];
 333     for(j=1;j<n-1;j+=2){
 334       float temp=pcm[j]*pcm[j]+pcm[j+1]*pcm[j+1];
 335       temp=logfft[(j+1)>>1]=scale_dB+.5f*todB(&temp)  + .345; /* +
 336                                      .345 is a hack; the original todB
 337                                      estimation used on IEEE 754
 338                                      compliant machines had a bug that
 339                                      returned dB values about a third
 340                                      of a decibel too high.  The bug
 341                                      was harmless because tunings
 342                                      implicitly took that into
 343                                      account.  However, fixing the bug
 344                                      in the estimator requires
 345                                      changing all the tunings as well.
 346                                      For now, it's easier to sync
 347                                      things back up here, and
 348                                      recalibrate the tunings in the
 349                                      next major model upgrade. */
 350       if(temp>local_ampmax[i])local_ampmax[i]=temp;
 351     }
 352
 353     if(local_ampmax[i]>0.f)local_ampmax[i]=0.f;
 354     if(local_ampmax[i]>global_ampmax)global_ampmax=local_ampmax[i];
 355
 356 #if 0
 357     if(vi->channels==2){
 358       if(i==0){
 359         _analysis_output("fftL",seq,logfft,n/2,1,0,0);
 360       }else{
 361         _analysis_output("fftR",seq,logfft,n/2,1,0,0);
 362       }
 363     }else{
 364       _analysis_output("fft",seq,logfft,n/2,1,0,0);
 365     }
 366 #endif
 367
 368   }
 369
 370   {
 371     float   *noise        = _vorbis_block_alloc(vb,n/2*sizeof(*noise));
 372     float   *tone         = _vorbis_block_alloc(vb,n/2*sizeof(*tone));
 373
 374     for(i=0;i<vi->channels;i++){
 375       /* the encoder setup assumes that all the modes used by any
 376          specific bitrate tweaking use the same floor */
 377
 378       int submap=info->chmuxlist[i];
 379
 380       /* the following makes things clearer to *me* anyway */
 381       float *mdct    =gmdct[i];
 382       float *logfft  =vb->pcm[i];
 383
 384       float *logmdct =logfft+n/2;
 385       float *logmask =logfft;
 386
 387       vb->mode=modenumber;
 388
 389       floor_posts[i]=_vorbis_block_alloc(vb,PACKETBLOBS*sizeof(**floor_posts));
 390       memset(floor_posts[i],0,sizeof(**floor_posts)*PACKETBLOBS);
 391
 392       for(j=0;j<n/2;j++)
 393         logmdct[j]=todB(mdct+j)  + .345; /* + .345 is a hack; the original
 394                                      todB estimation used on IEEE 754
 395                                      compliant machines had a bug that
 396                                      returned dB values about a third
 397                                      of a decibel too high.  The bug
 398                                      was harmless because tunings
 399                                      implicitly took that into
 400                                      account.  However, fixing the bug
 401                                      in the estimator requires
 402                                      changing all the tunings as well.
 403                                      For now, it's easier to sync
 404                                      things back up here, and
 405                                      recalibrate the tunings in the
 406                                      next major model upgrade. */
 407
 408 #if 0
 409       if(vi->channels==2){
 410         if(i==0)
 411           _analysis_output("mdctL",seq,logmdct,n/2,1,0,0);
 412         else
 413           _analysis_output("mdctR",seq,logmdct,n/2,1,0,0);
 414       }else{
 415         _analysis_output("mdct",seq,logmdct,n/2,1,0,0);
 416       }
 417 #endif
 418
 419       /* first step; noise masking.  Not only does 'noise masking'
 420          give us curves from which we can decide how much resolution
 421          to give noise parts of the spectrum, it also implicitly hands
 422          us a tonality estimate (the larger the value in the
 423          'noise_depth' vector, the more tonal that area is) */
 424
 425       _vp_noisemask(psy_look,
 426                     logmdct,
 427                     noise); /* noise does not have by-frequency offset
 428                                bias applied yet */
 429 #if 0
 430       if(vi->channels==2){
 431         if(i==0)
 432           _analysis_output("noiseL",seq,noise,n/2,1,0,0);
 433         else
 434           _analysis_output("noiseR",seq,noise,n/2,1,0,0);
 435       }else{
 436         _analysis_output("noise",seq,noise,n/2,1,0,0);
 437       }
 438 #endif
 439
 440       /* second step: 'all the other crap'; all the stuff that isn't
 441          computed/fit for bitrate management goes in the second psy
 442          vector.  This includes tone masking, peak limiting and ATH */
 443
 444       _vp_tonemask(psy_look,
 445                    logfft,
 446                    tone,
 447                    global_ampmax,
 448                    local_ampmax[i]);
 449
 450 #if 0
 451       if(vi->channels==2){
 452         if(i==0)
 453           _analysis_output("toneL",seq,tone,n/2,1,0,0);
 454         else
 455           _analysis_output("toneR",seq,tone,n/2,1,0,0);
 456       }else{
 457         _analysis_output("tone",seq,tone,n/2,1,0,0);
 458       }
 459 #endif
 460
 461       /* third step; we offset the noise vectors, overlay tone
 462          masking.  We then do a floor1-specific line fit.  If we're
 463          performing bitrate management, the line fit is performed
 464          multiple times for up/down tweakage on demand. */
 465
 466 #if 0
 467       {
 468       float aotuv[psy_look->n];
 469 #endif
 470
 471         _vp_offset_and_mix(psy_look,
 472                            noise,
 473                            tone,
 474                            1,
 475                            logmask,
 476                            mdct,
 477                            logmdct);
 478
 479 #if 0
 480         if(vi->channels==2){
 481           if(i==0)
 482             _analysis_output("aotuvM1_L",seq,aotuv,psy_look->n,1,1,0);
 483           else
 484             _analysis_output("aotuvM1_R",seq,aotuv,psy_look->n,1,1,0);
 485         }else{
 486           _analysis_output("aotuvM1",seq,aotuv,psy_look->n,1,1,0);
 487         }
 488       }
 489 #endif
 490
 491
 492 #if 0
 493       if(vi->channels==2){
 494         if(i==0)
 495           _analysis_output("mask1L",seq,logmask,n/2,1,0,0);
 496         else
 497           _analysis_output("mask1R",seq,logmask,n/2,1,0,0);
 498       }else{
 499         _analysis_output("mask1",seq,logmask,n/2,1,0,0);
 500       }
 501 #endif
 502
 503       /* this algorithm is hardwired to floor 1 for now; abort out if
 504          we're *not* floor1.  This won't happen unless someone has
 505          broken the encode setup lib.  Guard it anyway. */
 506       if(ci->floor_type[info->floorsubmap[submap]]!=1)return(-1);
 507
 508       floor_posts[i][PACKETBLOBS/2]=
 509         floor1_fit(vb,b->flr[info->floorsubmap[submap]],
 510                    logmdct,
 511                    logmask);
 512
 513       /* are we managing bitrate?  If so, perform two more fits for
 514          later rate tweaking (fits represent hi/lo) */
 515       if(vorbis_bitrate_managed(vb) && floor_posts[i][PACKETBLOBS/2]){
 516         /* higher rate by way of lower noise curve */
 517
 518         _vp_offset_and_mix(psy_look,
 519                            noise,
 520                            tone,
 521                            2,
 522                            logmask,
 523                            mdct,
 524                            logmdct);
 525
 526 #if 0
 527         if(vi->channels==2){
 528           if(i==0)
 529             _analysis_output("mask2L",seq,logmask,n/2,1,0,0);
 530           else
 531             _analysis_output("mask2R",seq,logmask,n/2,1,0,0);
 532         }else{
 533           _analysis_output("mask2",seq,logmask,n/2,1,0,0);
 534         }
 535 #endif
 536
 537         floor_posts[i][PACKETBLOBS-1]=
 538           floor1_fit(vb,b->flr[info->floorsubmap[submap]],
 539                      logmdct,
 540                      logmask);
 541
 542         /* lower rate by way of higher noise curve */
 543         _vp_offset_and_mix(psy_look,
 544                            noise,
 545                            tone,
 546                            0,
 547                            logmask,
 548                            mdct,
 549                            logmdct);
 550
 551 #if 0
 552         if(vi->channels==2){
 553           if(i==0)
 554             _analysis_output("mask0L",seq,logmask,n/2,1,0,0);
 555           else
 556             _analysis_output("mask0R",seq,logmask,n/2,1,0,0);
 557         }else{
 558           _analysis_output("mask0",seq,logmask,n/2,1,0,0);
 559         }
 560 #endif
 561
 562         floor_posts[i][0]=
 563           floor1_fit(vb,b->flr[info->floorsubmap[submap]],
 564                      logmdct,
 565                      logmask);
 566
 567         /* we also interpolate a range of intermediate curves for
 568            intermediate rates */
 569         for(k=1;k<PACKETBLOBS/2;k++)
 570           floor_posts[i][k]=
 571             floor1_interpolate_fit(vb,b->flr[info->floorsubmap[submap]],
 572                                    floor_posts[i][0],
 573                                    floor_posts[i][PACKETBLOBS/2],
 574                                    k*65536/(PACKETBLOBS/2));
 575         for(k=PACKETBLOBS/2+1;k<PACKETBLOBS-1;k++)
 576           floor_posts[i][k]=
 577             floor1_interpolate_fit(vb,b->flr[info->floorsubmap[submap]],
 578                                    floor_posts[i][PACKETBLOBS/2],
 579                                    floor_posts[i][PACKETBLOBS-1],
 580                                    (k-PACKETBLOBS/2)*65536/(PACKETBLOBS/2));
 581       }
 582     }
 583   }
 584   vbi->ampmax=global_ampmax;
 585
 586   /*
 587     the next phases are performed once for vbr-only and PACKETBLOB
 588     times for bitrate managed modes.
 589
 590     1) encode actual mode being used
 591     2) encode the floor for each channel, compute coded mask curve/res
 592     3) normalize and couple.
 593     4) encode residue
 594     5) save packet bytes to the packetblob vector
 595
 596   */
 597
 598   /* iterate over the many masking curve fits we've created */
 599
 600   {
 601     int **couple_bundle=alloca(sizeof(*couple_bundle)*vi->channels);
 602     int *zerobundle=alloca(sizeof(*zerobundle)*vi->channels);
 603     int **sortindex=alloca(sizeof(*sortindex)*vi->channels);
 604     float **mag_memo=NULL;
 605     int **mag_sort=NULL;
 606
 607     if(info->coupling_steps){
 608       mag_memo=_vp_quantize_couple_memo(vb,
 609                                         &ci->psy_g_param,
 610                                         psy_look,
 611                                         info,
 612                                         gmdct);
 613
 614       mag_sort=_vp_quantize_couple_sort(vb,
 615                                         psy_look,
 616                                         info,
 617                                         mag_memo);
 618
 619       hf_reduction(&ci->psy_g_param,
 620                    psy_look,
 621                    info,
 622                    mag_memo);
 623     }
 624
 625     memset(sortindex,0,sizeof(*sortindex)*vi->channels);
 626     if(psy_look->vi->normal_channel_p){
 627       for(i=0;i<vi->channels;i++){
 628         float *mdct    =gmdct[i];
 629         sortindex[i]=alloca(sizeof(**sortindex)*n/2);
 630         iwork[i]=_vorbis_block_alloc(vb,n/2*sizeof(**iwork));
 631         _vp_noise_normalize_sort(psy_look,mdct,sortindex[i]);
 632       }
 633     }
 634
 635     for(k=(vorbis_bitrate_managed(vb)?0:PACKETBLOBS/2);
 636         k<=(vorbis_bitrate_managed(vb)?PACKETBLOBS-1:PACKETBLOBS/2);
 637         k++){
 638       oggpack_buffer *opb=vbi->packetblob[k];
 639
 640       /* start out our new packet blob with packet type and mode */
 641       /* Encode the packet type */
 642       oggpack_write(opb,0,1);
 643       /* Encode the modenumber */
 644       /* Encode frame mode, pre,post windowsize, then dispatch */
 645       oggpack_write(opb,modenumber,b->modebits);
 646       if(vb->W){
 647         oggpack_write(opb,vb->lW,1);
 648         oggpack_write(opb,vb->nW,1);
 649       }
 650
 651       /* encode floor, compute masking curve, sep out residue */
 652       for(i=0;i<vi->channels;i++){
 653         int submap=info->chmuxlist[i];
 654         float *mdct    =gmdct[i];
 655         float *res     =vb->pcm[i];
 656         int   *ilogmask=iwork[i];
 657
 658         nonzero[i]=floor1_encode(opb,vb,b->flr[info->floorsubmap[submap]],
 659                                  floor_posts[i][k],
 660                                  ilogmask);
 661 #if 0
 662         {
 663           char buf[80];
 664           sprintf(buf,"maskI%c%d",i?'R':'L',k);
 665           float work[n/2];
 666           for(j=0;j<n/2;j++)
 667             work[j]=FLOOR1_fromdB_LOOKUP[iwork[j]];
 668           _analysis_output(buf,seq,work,n/2,1,1,0);
 669         }
 670 #endif
 671         _vp_remove_floor(psy_look,
 672                          mdct,
 673                          ilogmask,
 674                          res,
 675                          ci->psy_g_param.sliding_lowpass[vb->W][k]);
 676
 677         _vp_noise_normalize(psy_look,res,res+n/2,sortindex[i]);
 678
 679
 680 #if 0
 681         {
 682           char buf[80];
 683           float work[n/2];
 684           for(j=0;j<n/2;j++)
 685             work[j]=FLOOR1_fromdB_LOOKUP[iwork[j]]*(res+n/2)[j];
 686           sprintf(buf,"resI%c%d",i?'R':'L',k);
 687           _analysis_output(buf,seq,work,n/2,1,1,0);
 688
 689         }
 690 #endif
 691       }
 692
 693       /* our iteration is now based on masking curve, not prequant and
 694          coupling.  Only one prequant/coupling step */
 695
 696       /* quantize/couple */
 697       /* incomplete implementation that assumes the tree is all depth
 698          one, or no tree at all */
 699       if(info->coupling_steps){
 700         _vp_couple(k,
 701                    &ci->psy_g_param,
 702                    psy_look,
 703                    info,
 704                    vb->pcm,
 705                    mag_memo,
 706                    mag_sort,
 707                    iwork,
 708                    nonzero,
 709                    ci->psy_g_param.sliding_lowpass[vb->W][k]);
 710       }
 711       for(j=0;j<vi->channels;j++){
 712         int *ires = iwork[j];
 713         float *res = vb->pcm[j]+n/2;
 714         for(i=0;i<n/2;i++)
 715           ires[i] = (int)rint(res[i]);
 716       }
 717
 718       /* classify and encode by submap */
 719       for(i=0;i<info->submaps;i++){
 720         int ch_in_bundle=0;
 721         long **classifications;
 722         int resnum=info->residuesubmap[i];
 723
 724         for(j=0;j<vi->channels;j++){
 725           if(info->chmuxlist[j]==i){
 726             zerobundle[ch_in_bundle]=0;
 727             if(nonzero[j])zerobundle[ch_in_bundle]=1;
 728             couple_bundle[ch_in_bundle++]=iwork[j];
 729           }
 730         }
 731
 732         classifications=_residue_P[ci->residue_type[resnum]]->
 733           class(vb,b->residue[resnum],couple_bundle,zerobundle,ch_in_bundle);
 734
 735         ch_in_bundle=0;
 736         for(j=0;j<vi->channels;j++)
 737           if(info->chmuxlist[j]==i)
 738             couple_bundle[ch_in_bundle++]=iwork[j];
 739
 740         _residue_P[ci->residue_type[resnum]]->
 741           forward(opb,vb,b->residue[resnum],
 742                   couple_bundle,zerobundle,ch_in_bundle,classifications);
 743       }
 744
 745       /* ok, done encoding.  Next protopacket. */
 746     }
 747
 748   }
 749
 750 #if 0
 751   seq++;
 752   total+=ci->blocksizes[vb->W]/4+ci->blocksizes[vb->nW]/4;
 753 #endif
 754   return(0);
 755 }
 756
 757 static int mapping0_inverse(vorbis_block *vb,vorbis_info_mapping *l){
 758   vorbis_dsp_state     *vd=vb->vd;
 759   vorbis_info          *vi=vd->vi;
 760   codec_setup_info     *ci=vi->codec_setup;
 761   private_state        *b=vd->backend_state;
 762   vorbis_info_mapping0 *info=(vorbis_info_mapping0 *)l;
 763
 764   int                   i,j;
 765   long                  n=vb->pcmend=ci->blocksizes[vb->W];
 766
 767   float **pcmbundle=alloca(sizeof(*pcmbundle)*vi->channels);
 768   int    *zerobundle=alloca(sizeof(*zerobundle)*vi->channels);
 769
 770   int   *nonzero  =alloca(sizeof(*nonzero)*vi->channels);
 771   void **floormemo=alloca(sizeof(*floormemo)*vi->channels);
 772
 773   /* recover the spectral envelope; store it in the PCM vector for now */
 774   for(i=0;i<vi->channels;i++){
 775     int submap=info->chmuxlist[i];
 776     floormemo[i]=_floor_P[ci->floor_type[info->floorsubmap[submap]]]->
 777       inverse1(vb,b->flr[info->floorsubmap[submap]]);
 778     if(floormemo[i])
 779       nonzero[i]=1;
 780     else
 781       nonzero[i]=0;
 782     memset(vb->pcm[i],0,sizeof(*vb->pcm[i])*n/2);
 783   }
 784
 785   /* channel coupling can 'dirty' the nonzero listing */
 786   for(i=0;i<info->coupling_steps;i++){
 787     if(nonzero[info->coupling_mag[i]] ||
 788        nonzero[info->coupling_ang[i]]){
 789       nonzero[info->coupling_mag[i]]=1;
 790       nonzero[info->coupling_ang[i]]=1;
 791     }
 792   }
 793
 794   /* recover the residue into our working vectors */
 795   for(i=0;i<info->submaps;i++){
 796     int ch_in_bundle=0;
 797     for(j=0;j<vi->channels;j++){
 798       if(info->chmuxlist[j]==i){
 799         if(nonzero[j])
 800           zerobundle[ch_in_bundle]=1;
 801         else
 802           zerobundle[ch_in_bundle]=0;
 803         pcmbundle[ch_in_bundle++]=vb->pcm[j];
 804       }
 805     }
 806
 807     _residue_P[ci->residue_type[info->residuesubmap[i]]]->
 808       inverse(vb,b->residue[info->residuesubmap[i]],
 809               pcmbundle,zerobundle,ch_in_bundle);
 810   }
 811
 812   /* channel coupling */
 813   for(i=info->coupling_steps-1;i>=0;i--){
 814     float *pcmM=vb->pcm[info->coupling_mag[i]];
 815     float *pcmA=vb->pcm[info->coupling_ang[i]];
 816
 817     for(j=0;j<n/2;j++){
 818       float mag=pcmM[j];
 819       float ang=pcmA[j];
 820
 821       if(mag>0)
 822         if(ang>0){
 823           pcmM[j]=mag;
 824           pcmA[j]=mag-ang;
 825         }else{
 826           pcmA[j]=mag;
 827           pcmM[j]=mag+ang;
 828         }
 829       else
 830         if(ang>0){
 831           pcmM[j]=mag;
 832           pcmA[j]=mag+ang;
 833         }else{
 834           pcmA[j]=mag;
 835           pcmM[j]=mag-ang;
 836         }
 837     }
 838   }
 839
 840   /* compute and apply spectral envelope */
 841   for(i=0;i<vi->channels;i++){
 842     float *pcm=vb->pcm[i];
 843     int submap=info->chmuxlist[i];
 844     _floor_P[ci->floor_type[info->floorsubmap[submap]]]->
 845       inverse2(vb,b->flr[info->floorsubmap[submap]],
 846                floormemo[i],pcm);
 847   }
 848
 849   /* transform the PCM data; takes PCM vector, vb; modifies PCM vector */
 850   /* only MDCT right now.... */
 851   for(i=0;i<vi->channels;i++){
 852     float *pcm=vb->pcm[i];
 853     mdct_backward(b->transform[vb->W][0],pcm,pcm);
 854   }
 855
 856   /* all done! */
 857   return(0);
 858 }
 859
 860 /* export hooks */
 861 const vorbis_func_mapping mapping0_exportbundle={
 862   &mapping0_pack,
 863   &mapping0_unpack,
 864   &mapping0_free_info,
 865   &mapping0_forward,
 866   &mapping0_inverse
 867 };