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