lib/block.c

   1 /********************************************************************
   2  *                                                                  *
   3  * THIS FILE IS PART OF THE Ogg Vorbis SOFTWARE CODEC SOURCE CODE.  *
   4  * USE, DISTRIBUTION AND REPRODUCTION OF THIS SOURCE IS GOVERNED BY *
   5  * THE GNU PUBLIC LICENSE 2, WHICH IS INCLUDED WITH THIS SOURCE.    *
   6  * PLEASE READ THESE TERMS DISTRIBUTING.                            *
   7  *                                                                  *
   8  * THE OggSQUISH SOURCE CODE IS (C) COPYRIGHT 1994-1999             *
   9  * by 1999 Monty <monty@xiph.org> and The XIPHOPHORUS Company       *
  10  * http://www.xiph.org/                                             *
  11  *                                                                  *
  12  ********************************************************************
  13
  14  function: PCM data vector blocking, windowing and dis/reassembly
  15  author: Monty <xiphmont@mit.edu>
  16  modifications by: Monty
  17  last modification date: Oct 2 1999
  18
  19  Handle windowing, overlap-add, etc of the PCM vectors.  This is made
  20  more amusing by Vorbis' current two allowed block sizes.
  21
  22  Vorbis manipulates the dynamic range of the incoming PCM data
  23  envelope to minimise time-domain energy leakage from percussive and
  24  plosive waveforms being quantized in the MDCT domain.
  25
  26  ********************************************************************/
  27
  28 #include <stdlib.h>
  29 #include <string.h>
  30 #include "codec.h"
  31 #include "window.h"
  32 #include "envelope.h"
  33 #include "mdct.h"
  34 #include "lpc.h"
  35
  36 /* pcm accumulator examples (not exhaustive):
  37
  38  <-------------- lW ---------------->
  39                    <--------------- W ---------------->
  40 :            .....|.....       _______________         |
  41 :        .'''     |     '''_---      |       |\        |
  42 :.....'''         |_____--- '''......|       | \_______|
  43 :.................|__________________|_______|__|______|
  44                   |<------ Sl ------>|      > Sr <     |endW
  45                   |beginSl           |endSl  |  |endSr
  46                   |beginW            |endlW  |beginSr
  47
  48
  49                       |< lW >|
  50                    <--------------- W ---------------->
  51                   |   |  ..  ______________            |
  52                   |   | '  `/        |     ---_        |
  53                   |___.'___/`.       |         ---_____|
  54                   |_______|__|_______|_________________|
  55                   |      >|Sl|<      |<------ Sr ----->|endW
  56                   |       |  |endSl  |beginSr          |endSr
  57                   |beginW |  |endlW
  58                   mult[0] |beginSl                     mult[n]
  59
  60  <-------------- lW ----------------->
  61                           |<--W-->|
  62 :            ..............  ___  |   |
  63 :        .'''             |`/   \ |   |
  64 :.....'''                 |/`....\|...|
  65 :.........................|___|___|___|
  66                           |Sl |Sr |endW
  67                           |   |   |endSr
  68                           |   |beginSr
  69                           |   |endSl
  70                           |beginSl
  71                           |beginW
  72 */
  73
  74 /* block abstraction setup *********************************************/
  75
  76 int vorbis_block_init(vorbis_dsp_state *v, vorbis_block *vb){
  77   int i;
  78   memset(vb,0,sizeof(vorbis_block));
  79   vb->vd=v;
  80
  81   vb->pcm_storage=v->block_size[1];
  82   vb->pcm_channels=v->pcm_channels;
  83   vb->mult_storage=v->block_size[1]/v->samples_per_envelope_step;
  84   vb->mult_channels=v->envelope_channels;
  85   vb->floor_channels=v->vi.floorch;
  86   vb->floor_storage=v->vi.floororder;
  87
  88   vb->pcm=malloc(vb->pcm_channels*sizeof(double *));
  89   for(i=0;i<vb->pcm_channels;i++)
  90     vb->pcm[i]=malloc(vb->pcm_storage*sizeof(double));
  91
  92   vb->mult=malloc(vb->mult_channels*sizeof(double *));
  93   for(i=0;i<vb->mult_channels;i++)
  94     vb->mult[i]=malloc(vb->mult_storage*sizeof(double));
  95
  96   vb->lsp=malloc(vb->floor_channels*sizeof(double *));
  97   vb->lpc=malloc(vb->floor_channels*sizeof(double *));
  98   vb->amp=malloc(vb->floor_channels*sizeof(double));
  99   for(i=0;i<vb->floor_channels;i++){
 100     vb->lsp[i]=malloc(vb->floor_storage*sizeof(double));
 101     vb->lpc[i]=malloc(vb->floor_storage*sizeof(double));
 102   }
 103
 104   return(0);
 105 }
 106
 107 int vorbis_block_clear(vorbis_block *vb){
 108   int i;
 109   if(vb->pcm){
 110     for(i=0;i<vb->pcm_channels;i++)
 111       free(vb->pcm[i]);
 112     free(vb->pcm);
 113   }
 114   if(vb->mult){
 115     for(i=0;i<vb->mult_channels;i++)
 116       free(vb->mult[i]);
 117     free(vb->mult);
 118   }
 119   memset(vb,0,sizeof(vorbis_block));
 120   return(0);
 121 }
 122
 123 /* Analysis side code, but directly related to blocking.  Thus it's
 124    here and not in analysis.c (which is for analysis transforms only).
 125    The init is here because some of it is shared */
 126
 127 static int _vds_shared_init(vorbis_dsp_state *v,vorbis_info *vi){
 128   memset(v,0,sizeof(vorbis_dsp_state));
 129
 130   memcpy(&v->vi,vi,sizeof(vorbis_info));
 131   _ve_envelope_init(&v->ve,vi->envelopesa);
 132   mdct_init(&v->vm[0],vi->smallblock);
 133   mdct_init(&v->vm[1],vi->largeblock);
 134
 135   v->samples_per_envelope_step=vi->envelopesa;
 136   v->block_size[0]=vi->smallblock;
 137   v->block_size[1]=vi->largeblock;
 138
 139   v->window[0][0][0]=_vorbis_window(v->block_size[0],
 140                                    v->block_size[0]/2,v->block_size[0]/2);
 141   v->window[0][0][1]=v->window[0][0][0];
 142   v->window[0][1][0]=v->window[0][0][0];
 143   v->window[0][1][1]=v->window[0][0][0];
 144
 145   v->window[1][0][0]=_vorbis_window(v->block_size[1],
 146                                    v->block_size[0]/2,v->block_size[0]/2);
 147   v->window[1][0][1]=_vorbis_window(v->block_size[1],
 148                                    v->block_size[0]/2,v->block_size[1]/2);
 149   v->window[1][1][0]=_vorbis_window(v->block_size[1],
 150                                    v->block_size[1]/2,v->block_size[0]/2);
 151   v->window[1][1][1]=_vorbis_window(v->block_size[1],
 152                                     v->block_size[1]/2,v->block_size[1]/2);
 153
 154   /* initialize the storage vectors to a decent size greater than the
 155      minimum */
 156
 157   v->pcm_storage=8192; /* we'll assume later that we have
 158                           a minimum of twice the blocksize of
 159                           accumulated samples in analysis */
 160   v->pcm_channels=v->vi.channels=vi->channels;
 161   v->pcm=malloc(vi->channels*sizeof(double *));
 162   v->pcmret=malloc(vi->channels*sizeof(double *));
 163   {
 164     int i;
 165     for(i=0;i<vi->channels;i++)
 166       v->pcm[i]=calloc(v->pcm_storage,sizeof(double));
 167   }
 168
 169   /* Initialize the envelope multiplier storage */
 170
 171   if(vi->envelopech){
 172     v->envelope_storage=v->pcm_storage/v->samples_per_envelope_step;
 173     v->envelope_channels=vi->envelopech;
 174     v->multipliers=calloc(v->envelope_channels,sizeof(double *));
 175     {
 176       int i;
 177       for(i=0;i<v->envelope_channels;i++){
 178         v->multipliers[i]=calloc(v->envelope_storage,sizeof(double));
 179       }
 180     }
 181   }
 182
 183   /* all 1 (large block) or 0 (small block) */
 184   /* explicitly set for the sake of clarity */
 185   v->lW=0; /* previous window size */
 186   v->W=0;  /* current window size */
 187
 188   /* all vector indexes; multiples of samples_per_envelope_step */
 189   v->centerW=v->block_size[1]/2;
 190
 191   v->pcm_current=v->centerW;
 192   v->envelope_current=v->centerW/v->samples_per_envelope_step;
 193   return(0);
 194 }
 195
 196 /* arbitrary settings and spec-mandated numbers get filled in here */
 197 int vorbis_analysis_init(vorbis_dsp_state *v,vorbis_info *vi){
 198   _vds_shared_init(v,vi);
 199
 200   /* Yes, wasteful to have four lookups.  This will get collapsed once
 201      things crystallize */
 202   lpc_init(&v->vl[0],vi->smallblock/2,vi->smallblock/2,
 203            vi->floororder,vi->flooroctaves,1);
 204   lpc_init(&v->vl[1],vi->largeblock/2,vi->largeblock/2,
 205            vi->floororder,vi->flooroctaves,1);
 206
 207   lpc_init(&v->vbal[0],vi->smallblock/2,256,
 208            vi->balanceorder,vi->balanceoctaves,1);
 209   lpc_init(&v->vbal[1],vi->largeblock/2,256,
 210            vi->balanceorder,vi->balanceoctaves,1);
 211
 212   return(0);
 213 }
 214
 215 void vorbis_dsp_clear(vorbis_dsp_state *v){
 216   int i,j,k;
 217   if(v){
 218
 219     if(v->window[0][0][0])free(v->window[0][0][0]);
 220     for(j=0;j<2;j++)
 221       for(k=0;k<2;k++)
 222         if(v->window[1][j][k])free(v->window[1][j][k]);
 223     if(v->pcm){
 224       for(i=0;i<v->pcm_channels;i++)
 225         if(v->pcm[i])free(v->pcm[i]);
 226       free(v->pcm);
 227       if(v->pcmret)free(v->pcmret);
 228     }
 229     if(v->multipliers){
 230       for(i=0;i<v->envelope_channels;i++)
 231         if(v->multipliers[i])free(v->multipliers[i]);
 232       free(v->multipliers);
 233     }
 234     _ve_envelope_clear(&v->ve);
 235     mdct_clear(&v->vm[0]);
 236     mdct_clear(&v->vm[1]);
 237     lpc_clear(&v->vl[0]);
 238     lpc_clear(&v->vl[1]);
 239     lpc_clear(&v->vbal[0]);
 240     lpc_clear(&v->vbal[1]);
 241     memset(v,0,sizeof(vorbis_dsp_state));
 242   }
 243 }
 244
 245 double **vorbis_analysis_buffer(vorbis_dsp_state *v, int vals){
 246   int i;
 247
 248   /* Do we have enough storage space for the requested buffer? If not,
 249      expand the PCM (and envelope) storage */
 250
 251   if(v->pcm_current+vals>=v->pcm_storage){
 252     v->pcm_storage=v->pcm_current+vals*2;
 253     v->envelope_storage=v->pcm_storage/v->samples_per_envelope_step;
 254
 255     for(i=0;i<v->pcm_channels;i++){
 256       v->pcm[i]=realloc(v->pcm[i],v->pcm_storage*sizeof(double));
 257     }
 258     for(i=0;i<v->envelope_channels;i++){
 259       v->multipliers[i]=realloc(v->multipliers[i],
 260                                 v->envelope_storage*sizeof(double));
 261     }
 262   }
 263
 264   for(i=0;i<v->pcm_channels;i++)
 265     v->pcmret[i]=v->pcm[i]+v->pcm_current;
 266
 267   return(v->pcmret);
 268 }
 269
 270 /* call with val<=0 to set eof */
 271
 272 int vorbis_analysis_wrote(vorbis_dsp_state *v, int vals){
 273   if(vals<=0){
 274     /* We're encoding the end of the stream.  Just make sure we have
 275        [at least] a full block of zeroes at the end. */
 276
 277     int i;
 278     vorbis_analysis_buffer(v,v->block_size[1]*2);
 279     v->eofflag=v->pcm_current;
 280     v->pcm_current+=v->block_size[1]*2;
 281     for(i=0;i<v->pcm_channels;i++)
 282       memset(v->pcm[i]+v->eofflag,0,
 283              (v->pcm_current-v->eofflag)*sizeof(double));
 284   }else{
 285
 286     if(v->pcm_current+vals>v->pcm_storage)
 287       return(-1);
 288
 289     v->pcm_current+=vals;
 290   }
 291   return(0);
 292 }
 293
 294 /* do the deltas, envelope shaping, pre-echo and determine the size of
 295    the next block on which to continue analysis */
 296 int vorbis_analysis_blockout(vorbis_dsp_state *v,vorbis_block *vb){
 297   int i,j;
 298   long beginW=v->centerW-v->block_size[v->W]/2,centerNext;
 299   long beginM=beginW/v->samples_per_envelope_step;
 300
 301   /* check to see if we're done... */
 302   if(v->eofflag==-1)return(0);
 303
 304   /* if we have any unfilled envelope blocks for which we have PCM
 305      data, fill them up in before proceeding. */
 306
 307   if(v->pcm_current/v->samples_per_envelope_step>v->envelope_current){
 308     /* This generates the multipliers, but does not sparsify the vector.
 309        That's done by block before coding */
 310     _ve_envelope_multipliers(v);
 311   }
 312
 313   /* By our invariant, we have lW, W and centerW set.  Search for
 314      the next boundary so we can determine nW (the next window size)
 315      which lets us compute the shape of the current block's window */
 316
 317   /* overconserve for now; any block with a non-placeholder multiplier
 318      should be minimal size. We can be greedy and only look at nW size */
 319
 320   if(v->vi.smallblock<v->vi.largeblock){
 321
 322     if(v->W)
 323       /* this is a long window; we start the search forward of centerW
 324          because that's the fastest we could react anyway */
 325       i=v->centerW+v->block_size[1]/4-v->block_size[0]/4;
 326     else
 327       /* short window.  Search from centerW */
 328       i=v->centerW;
 329     i/=v->samples_per_envelope_step;
 330
 331     for(;i<v->envelope_current;i++){
 332       for(j=0;j<v->envelope_channels;j++)
 333         if(v->multipliers[j][i-1]*v->vi.preecho_thresh<
 334            v->multipliers[j][i])break;
 335       if(j<v->envelope_channels)break;
 336     }
 337
 338     if(i<v->envelope_current){
 339       /* Ooo, we hit a multiplier. Is it beyond the boundary to make the
 340          upcoming block large ? */
 341       long largebound;
 342       if(v->W)
 343         largebound=v->centerW+v->block_size[1];
 344       else
 345         largebound=v->centerW+v->block_size[0]/4+v->block_size[1]*3/4;
 346       largebound/=v->samples_per_envelope_step;
 347
 348       if(i>=largebound)
 349         v->nW=1;
 350       else
 351         v->nW=0;
 352
 353     }else{
 354       /* Assume maximum; if the block is incomplete given current
 355          buffered data, this will be detected below */
 356       v->nW=1;
 357     }
 358   }else
 359     v->nW=1;
 360     v->nW=1;
 361
 362   /* Do we actually have enough data *now* for the next block? The
 363      reason to check is that if we had no multipliers, that could
 364      simply been due to running out of data.  In that case, we don;t
 365      know the size of the next block for sure and we need that now to
 366      figure out the window shape of this block */
 367
 368   centerNext=v->centerW+v->block_size[v->W]/4+v->block_size[v->nW]/4;
 369
 370   {
 371     long blockbound=centerNext+v->block_size[v->nW]/2;
 372     if(v->pcm_current<blockbound)return(0); /* not enough data yet */
 373   }
 374
 375   /* fill in the block */
 376   vb->lW=v->lW;
 377   vb->W=v->W;
 378   vb->nW=v->nW;
 379   vb->vd=v;
 380
 381   vb->pcmend=v->block_size[v->W];
 382   vb->multend=vb->pcmend / v->samples_per_envelope_step;
 383
 384   if(vb->floor_channels!=v->vi.floorch ||
 385      vb->floor_storage!=v->vi.floororder ||
 386      v->pcm_channels!=vb->pcm_channels ||
 387      v->block_size[1]!=vb->pcm_storage ||
 388      v->envelope_channels!=vb->mult_channels){
 389
 390     /* Storage not initialized or initilized for some other codec
 391        instance with different settings */
 392
 393     vorbis_block_clear(vb);
 394     vorbis_block_init(v,vb);
 395   }
 396
 397   /* copy the vectors */
 398   for(i=0;i<v->pcm_channels;i++)
 399     memcpy(vb->pcm[i],v->pcm[i]+beginW,v->block_size[v->W]*sizeof(double));
 400   for(i=0;i<v->envelope_channels;i++)
 401     memcpy(vb->mult[i],v->multipliers[i]+beginM,v->block_size[v->W]/
 402            v->samples_per_envelope_step*sizeof(double));
 403
 404   vb->frameno=v->frame;
 405
 406   /* handle eof detection: eof==0 means that we've not yet received EOF
 407                            eof>0  marks the last 'real' sample in pcm[]
 408                            eof<0  'no more to do'; doesn't get here */
 409
 410   if(v->eofflag){
 411     if(v->centerW>=v->eofflag){
 412       v->eofflag=-1;
 413       vb->eofflag=1;
 414     }
 415   }
 416
 417   /* advance storage vectors and clean up */
 418   {
 419     int new_centerNext=v->block_size[1]/2;
 420     int movementW=centerNext-new_centerNext;
 421     int movementM=movementW/v->samples_per_envelope_step;
 422
 423     /* the multipliers and pcm stay synced up because the blocksizes
 424        must be multiples of samples_per_envelope_step (minimum
 425        multiple is 2) */
 426
 427     for(i=0;i<v->pcm_channels;i++)
 428       memmove(v->pcm[i],v->pcm[i]+movementW,
 429               (v->pcm_current-movementW)*sizeof(double));
 430
 431     for(i=0;i<v->envelope_channels;i++){
 432       memmove(v->multipliers[i],v->multipliers[i]+movementM,
 433               (v->envelope_current-movementM)*sizeof(double));
 434     }
 435
 436     v->pcm_current-=movementW;
 437     v->envelope_current-=movementM;
 438
 439     v->lW=v->W;
 440     v->W=v->nW;
 441     v->centerW=new_centerNext;
 442
 443     v->frame++;
 444     v->samples+=movementW;
 445
 446     if(v->eofflag)
 447       v->eofflag-=movementW;
 448   }
 449
 450   /* done */
 451   return(1);
 452 }
 453
 454 int vorbis_synthesis_init(vorbis_dsp_state *v,vorbis_info *vi){
 455   int temp=vi->envelopech;
 456   vi->envelopech=0; /* we don't need multiplier buffering in syn */
 457   _vds_shared_init(v,vi);
 458   vi->envelopech=temp;
 459
 460   /* Yes, wasteful to have four lookups.  This will get collapsed once
 461      things crystallize */
 462   lpc_init(&v->vl[0],vi->smallblock/2,vi->smallblock/2,
 463            vi->floororder,vi->flooroctaves,0);
 464   lpc_init(&v->vl[1],vi->largeblock/2,vi->largeblock/2,
 465            vi->floororder,vi->flooroctaves,0);
 466   lpc_init(&v->vbal[0],vi->smallblock/2,256,
 467            vi->balanceorder,vi->balanceoctaves,0);
 468   lpc_init(&v->vbal[1],vi->largeblock/2,256,
 469            vi->balanceorder,vi->balanceoctaves,0);
 470
 471
 472   /* Adjust centerW to allow an easier mechanism for determining output */
 473   v->pcm_returned=v->centerW;
 474   v->centerW-= v->block_size[v->W]/4+v->block_size[v->lW]/4;
 475   return(0);
 476 }
 477
 478 /* Unike in analysis, the window is only partially applied.  Envelope
 479    is previously applied (the whole envelope, if any, is shipped in
 480    each block) */
 481
 482 int vorbis_synthesis_blockin(vorbis_dsp_state *v,vorbis_block *vb){
 483
 484   /* Shift out any PCM that we returned previously */
 485
 486   if(v->pcm_returned  && v->centerW>v->block_size[1]/2){
 487
 488     /* don't shift too much; we need to have a minimum PCM buffer of
 489        1/2 long block */
 490
 491     int shift=v->centerW-v->block_size[1]/2;
 492     shift=(v->pcm_returned<shift?v->pcm_returned:shift);
 493
 494     v->pcm_current-=shift;
 495     v->centerW-=shift;
 496     v->pcm_returned-=shift;
 497
 498     if(shift){
 499       int i;
 500       for(i=0;i<v->pcm_channels;i++)
 501         memmove(v->pcm[i],v->pcm[i]+shift,
 502                 v->pcm_current*sizeof(double));
 503     }
 504   }
 505
 506   {
 507     int sizeW=v->block_size[vb->W];
 508     int centerW=v->centerW+v->block_size[vb->lW]/4+sizeW/4;
 509     int beginW=centerW-sizeW/2;
 510     int endW=beginW+sizeW;
 511     int beginSl;
 512     int endSl;
 513     int i,j;
 514     double *window;
 515
 516     /* Do we have enough PCM storage for the block? */
 517     if(endW>v->pcm_storage){
 518       /* expand the PCM storage */
 519
 520       v->pcm_storage=endW+v->block_size[1];
 521
 522       for(i=0;i<v->pcm_channels;i++)
 523         v->pcm[i]=realloc(v->pcm[i],v->pcm_storage*sizeof(double));
 524     }
 525
 526     /* Overlap/add */
 527     switch(vb->W){
 528     case 0:
 529       beginSl=0;
 530       endSl=v->block_size[0]/2;
 531       break;
 532     case 1:
 533       beginSl=v->block_size[1]/4-v->block_size[vb->lW]/4;
 534       endSl=beginSl+v->block_size[vb->lW]/2;
 535       break;
 536     }
 537
 538     window=v->window[vb->W][0][vb->lW]+v->block_size[vb->W]/2;
 539
 540     for(j=0;j<v->pcm_channels;j++){
 541       double *pcm=v->pcm[j]+beginW;
 542
 543       /* the add section */
 544       for(i=beginSl;i<endSl;i++)
 545         pcm[i]=pcm[i]*window[i]+vb->pcm[j][i];
 546       /* the remaining section */
 547       for(;i<sizeW;i++)
 548         pcm[i]=vb->pcm[j][i];
 549     }
 550
 551     /* Update, cleanup */
 552
 553     v->centerW=centerW;
 554     v->pcm_current=endW;
 555
 556     if(vb->eofflag)v->eofflag=1;
 557   }
 558   return(0);
 559 }
 560
 561 int vorbis_synthesis_pcmout(vorbis_dsp_state *v,double ***pcm){
 562   if(v->pcm_returned<v->centerW){
 563     int i;
 564     for(i=0;i<v->pcm_channels;i++)
 565       v->pcmret[i]=v->pcm[i]+v->pcm_returned;
 566     *pcm=v->pcmret;
 567     return(v->centerW-v->pcm_returned);
 568   }
 569   return(0);
 570 }
 571
 572 int vorbis_synthesis_read(vorbis_dsp_state *v,int bytes){
 573   if(bytes && v->pcm_returned+bytes>v->centerW)return(-1);
 574   v->pcm_returned+=bytes;
 575   return(0);
 576 }
 577