lib/res0.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-2000             *
   9  * by Monty <monty@xiph.org> and The XIPHOPHORUS Company            *
  10  * http://www.xiph.org/                                             *
  11  *                                                                  *
  12  ********************************************************************
  13
  14  function: residue backend 0 implementation
  15  last mod: $Id: res0.c,v 1.11 2000/04/06 15:55:41 xiphmont Exp $
  16
  17  ********************************************************************/
  18
  19 /* Slow, slow, slow, simpleminded and did I mention it was slow?  The
  20    encode/decode loops are coded for clarity and performance is not
  21    yet even a nagging little idea lurking in the shadows.  Oh and BTW,
  22    it's slow. */
  23
  24 #include <stdlib.h>
  25 #include <string.h>
  26 #include <math.h>
  27 #include <stdio.h>
  28 #include "vorbis/codec.h"
  29 #include "bitwise.h"
  30 #include "registry.h"
  31 #include "scales.h"
  32 #include "bookinternal.h"
  33 #include "misc.h"
  34 #include "os.h"
  35
  36 typedef struct {
  37   vorbis_info_residue0 *info;
  38
  39   int         parts;
  40   codebook   *phrasebook;
  41
  42   codebook ***partbooks;
  43   int        *partstages;
  44
  45   int         partvals;
  46   int       **decodemap;
  47 } vorbis_look_residue0;
  48
  49 void free_info(vorbis_info_residue *i){
  50   if(i){
  51     memset(i,0,sizeof(vorbis_info_residue0));
  52     free(i);
  53   }
  54 }
  55
  56 void free_look(vorbis_look_residue *i){
  57   int j;
  58   if(i){
  59     vorbis_look_residue0 *look=(vorbis_look_residue0 *)i;
  60     for(j=0;j<look->parts;j++)
  61       if(look->partbooks[j])free(look->partbooks[j]);
  62     free(look->partbooks);
  63     for(j=0;j<look->partvals;j++)
  64       free(look->decodemap[j]);
  65     free(look->decodemap);
  66     if(look->partstages)free(look->partstages);
  67     memset(i,0,sizeof(vorbis_look_residue0));
  68     free(i);
  69   }
  70 }
  71
  72 void pack(vorbis_info_residue *vr,oggpack_buffer *opb){
  73   vorbis_info_residue0 *info=(vorbis_info_residue0 *)vr;
  74   int j,acc=0;
  75   _oggpack_write(opb,info->begin,24);
  76   _oggpack_write(opb,info->end,24);
  77
  78   _oggpack_write(opb,info->grouping-1,24);  /* residue vectors to group and
  79                                              code with a partitioned book */
  80   _oggpack_write(opb,info->partitions-1,6); /* possible partition choices */
  81   _oggpack_write(opb,info->groupbook,8);  /* group huffman book */
  82   for(j=0;j<info->partitions;j++){
  83     _oggpack_write(opb,info->secondstages[j],4); /* zero *is* a valid choice */
  84     acc+=info->secondstages[j];
  85   }
  86   for(j=0;j<acc;j++)
  87     _oggpack_write(opb,info->booklist[j],8);
  88
  89 }
  90
  91 /* vorbis_info is for range checking */
  92 vorbis_info_residue *unpack(vorbis_info *vi,oggpack_buffer *opb){
  93   int j,acc=0;
  94   vorbis_info_residue0 *info=calloc(1,sizeof(vorbis_info_residue0));
  95
  96   info->begin=_oggpack_read(opb,24);
  97   info->end=_oggpack_read(opb,24);
  98   info->grouping=_oggpack_read(opb,24)+1;
  99   info->partitions=_oggpack_read(opb,6)+1;
 100   info->groupbook=_oggpack_read(opb,8);
 101   for(j=0;j<info->partitions;j++)
 102     acc+=info->secondstages[j]=_oggpack_read(opb,4);
 103   for(j=0;j<acc;j++)
 104     info->booklist[j]=_oggpack_read(opb,8);
 105
 106   if(info->groupbook>=vi->books)goto errout;
 107   for(j=0;j<acc;j++)
 108     if(info->booklist[j]>=vi->books)goto errout;
 109
 110   return(info);
 111  errout:
 112   free_info(info);
 113   return(NULL);
 114 }
 115
 116 vorbis_look_residue *look (vorbis_dsp_state *vd,vorbis_info_mode *vm,
 117                           vorbis_info_residue *vr){
 118   vorbis_info_residue0 *info=(vorbis_info_residue0 *)vr;
 119   vorbis_look_residue0 *look=calloc(1,sizeof(vorbis_look_residue0));
 120   int j,k,acc=0;
 121   int dim;
 122   look->info=info;
 123
 124   look->parts=info->partitions;
 125   look->phrasebook=vd->fullbooks+info->groupbook;
 126   dim=look->phrasebook->dim;
 127
 128   look->partbooks=calloc(look->parts,sizeof(codebook **));
 129   look->partstages=calloc(look->parts,sizeof(int));
 130
 131   for(j=0;j<look->parts;j++){
 132     int stages=info->secondstages[j];
 133     if(stages){
 134       look->partbooks[j]=malloc(stages*sizeof(codebook *));
 135       for(k=0;k<stages;k++)
 136         look->partbooks[j][k]=vd->fullbooks+info->booklist[acc++];
 137     }
 138     look->partstages[j]=stages;
 139   }
 140
 141   look->partvals=pow(look->parts,dim);
 142   look->decodemap=malloc(look->partvals*sizeof(int *));
 143   for(j=0;j<look->partvals;j++){
 144     long val=j;
 145     long mult=look->partvals/look->parts;
 146     look->decodemap[j]=malloc(dim*sizeof(int));
 147     for(k=0;k<dim;k++){
 148       long deco=val/mult;
 149       val-=deco*mult;
 150       mult/=look->parts;
 151       look->decodemap[j][k]=deco;
 152     }
 153   }
 154
 155   return(look);
 156 }
 157
 158 /* returns the distance error from encoding with this book set */
 159 static double _testpart(double *vec,int n, int stages, codebook **books){
 160   int i,j;
 161
 162   double *work=alloca(n*sizeof(double)),acc=0.;
 163   memcpy(work,vec,n*sizeof(double));
 164
 165   if(stages==0){
 166     /* a mild hack.  We want partitions with samples values under
 167        fabs(.5) to be fully zeroed; if this case is met, we return an
 168        error of -1 (which cannot be beaten).  If the samples values
 169        don't meet this criteria, return the real error */
 170     for(i=0;i<n;i++)
 171       if(fabs(vec[i])>.5)break;
 172     if(i==n)return(-1.);
 173
 174     /* real (squared) error */
 175     for(i=0;i<n;i++)
 176       acc+=vec[i]*vec[i];
 177
 178   }else{
 179     for(j=0;j<stages;j++){
 180       acc=0.;
 181       for(i=0;i<n;i+=books[j]->dim)
 182         acc+=vorbis_book_vE(books[j],work+i);
 183     }
 184   }
 185
 186   return(acc);
 187 }
 188
 189 static int _testhack(double *vec,int n){
 190   int i;
 191   double acc=0.;
 192   double max=0.;
 193   for(i=0;i<n;i++)
 194     acc+=todB(fabs(vec[i]));
 195   acc=fromdB(acc/n);
 196   for(i=0;i<n;i++)
 197     max=(fabs(vec[i])>max?fabs(vec[i]):max);
 198
 199   if(max<.5)return(0);
 200   if(max<2.5 && acc<1.5)return(1);
 201   if(max<6.)return(2);
 202   return(3);
 203 }
 204
 205 static int _encodepart(oggpack_buffer *opb,double *vec, int n,
 206                        int stages, codebook **books){
 207   int i,j,bits=0;
 208
 209   double *work=alloca(n*sizeof(double));
 210   memcpy(work,vec,n*sizeof(double));
 211
 212   for(j=0;j<stages;j++)
 213     for(i=0;i<n;i+=books[j]->dim)
 214       bits+=vorbis_book_encodevE(books[j],work+i,opb);
 215
 216   return(bits);
 217 }
 218
 219 static int _decodepart(oggpack_buffer *opb,double *work,double *vec, int n,
 220                        int stages, codebook **books){
 221   int i,j;
 222
 223   memset(work,0,n*sizeof(double));
 224   for(j=0;j<stages;j++)
 225     for(i=0;i<n;i+=books[j]->dim)
 226       vorbis_book_decodev(books[j],work+i,opb);
 227
 228   for(i=0;i<n;i++)
 229     vec[i]*=work[i];
 230
 231   return(0);
 232 }
 233
 234 int forward(vorbis_block *vb,vorbis_look_residue *vl,
 235             double **in,int ch){
 236   long i,j,k,l;
 237   vorbis_look_residue0 *look=(vorbis_look_residue0 *)vl;
 238   vorbis_info_residue0 *info=look->info;
 239
 240   /* move all this setup out later */
 241   int samples_per_partition=info->grouping;
 242   int possible_partitions=info->partitions;
 243   int partitions_per_word=look->phrasebook->dim;
 244   int n=info->end-info->begin;
 245   long phrasebits=0,resbitsT=0;
 246   long *resbits=alloca(sizeof(long)*possible_partitions);
 247   long *resvals=alloca(sizeof(long)*possible_partitions);
 248
 249   int partvals=n/samples_per_partition;
 250   int partwords=(partvals+partitions_per_word-1)/partitions_per_word;
 251   long **partword=_vorbis_block_alloc(vb,ch*sizeof(long *));
 252   partvals=partwords*partitions_per_word;
 253
 254   /* we find/encode the patition type for each partition of each
 255      channel.  We'll go back and do the interleaved encoding in a
 256      bit.  For now, clarity */
 257
 258   if(ch)_analysis_output("a_res",vb->sequence,in[0],n);
 259   memset(resbits,0,sizeof(long)*possible_partitions);
 260   memset(resvals,0,sizeof(long)*possible_partitions);
 261
 262   for(i=0;i<ch;i++){
 263     partword[i]=_vorbis_block_alloc(vb,n/samples_per_partition*sizeof(long));
 264     memset(partword[i],0,n/samples_per_partition*sizeof(long));
 265   }
 266
 267   for(i=info->begin,l=0;i<info->end;i+=samples_per_partition,l++){
 268     for(j=0;j<ch;j++){
 269
 270       /* find the best encoding for the partition using each possible
 271          book.  Use the book that had lowest error (we arrange the
 272          books to make optimal choice very obvious and not even think
 273          about bits) */
 274       partword[j][l]=_testhack(in[j]+i,samples_per_partition);
 275
 276 #if 0
 277       double best=_testpart(in[j]+i,samples_per_partition,
 278                             look->partstages[0],look->partbooks[0]);
 279       for(k=1;k<info->partitions;k++){
 280         double this=_testpart(in[j]+i,samples_per_partition,
 281                               look->partstages[k],look->partbooks[k]);
 282         if(this<best){
 283           best=this;
 284           partword[j][l]=k;
 285         }
 286       }
 287 #endif
 288     }
 289   }
 290
 291   /* we code the partition words for each channel, then the residual
 292      words for a partition per channel until we've written all the
 293      partitions for that partition word.  Then write the next parition
 294      channel words... */
 295
 296   for(i=info->begin,l=0;i<info->end;){
 297     /* first we encode a partition codeword for each channel */
 298     for(j=0;j<ch;j++){
 299       long val=partword[j][l];
 300       for(k=1;k<partitions_per_word;k++)
 301         val= val*possible_partitions+partword[j][l+k];
 302       phrasebits+=vorbis_book_encode(look->phrasebook,val,&vb->opb);
 303     }
 304     /* now we encode interleaved residual values for the partitions */
 305     for(k=0;k<partitions_per_word;k++,l++,i+=samples_per_partition)
 306       for(j=0;j<ch;j++){
 307         resbits[partword[j][l]]+=
 308           _encodepart(&vb->opb,in[j]+i,samples_per_partition,
 309                       look->partstages[partword[j][l]],
 310                       look->partbooks[partword[j][l]]);
 311         resvals[partword[j][l]]+=samples_per_partition;
 312       }
 313
 314   }
 315
 316   for(i=0;i<possible_partitions;i++)resbitsT+=resbits[i];
 317   fprintf(stderr,"Encoded %ld res vectors in %ld phrasing and %ld res bits\n\t",
 318           ch*(info->end-info->begin),phrasebits,resbitsT);
 319   for(i=0;i<possible_partitions;i++)
 320     fprintf(stderr,"%ld(%ld):%ld ",i,resvals[i],resbits[i]);
 321   fprintf(stderr,"\n");
 322
 323   return(0);
 324 }
 325
 326 int inverse(vorbis_block *vb,vorbis_look_residue *vl,double **in,int ch){
 327   long i,j,k,l;
 328   vorbis_look_residue0 *look=(vorbis_look_residue0 *)vl;
 329   vorbis_info_residue0 *info=look->info;
 330
 331   /* move all this setup out later */
 332   int samples_per_partition=info->grouping;
 333   int partitions_per_word=look->phrasebook->dim;
 334   int n=info->end-info->begin;
 335
 336   int partvals=n/samples_per_partition;
 337   int partwords=(partvals+partitions_per_word-1)/partitions_per_word;
 338   int **partword=alloca(ch*sizeof(long *));
 339   double *work=alloca(sizeof(double)*samples_per_partition);
 340   partvals=partwords*partitions_per_word;
 341
 342   for(i=info->begin,l=0;i<info->end;){
 343     /* fetch the partition word for each channel */
 344     for(j=0;j<ch;j++){
 345       int temp=vorbis_book_decode(look->phrasebook,&vb->opb);
 346       partword[j]=look->decodemap[temp];
 347       if(partword[j]==NULL)exit(1);
 348     }
 349
 350     /* now we decode interleaved residual values for the partitions */
 351     for(k=0;k<partitions_per_word;k++,l++,i+=samples_per_partition)
 352       for(j=0;j<ch;j++){
 353         int part=partword[j][k];
 354         _decodepart(&vb->opb,work,in[j]+i,samples_per_partition,
 355                     look->partstages[part],
 356                     look->partbooks[part]);
 357       }
 358   }
 359
 360   if(ch)_analysis_output("s_res",vb->sequence,in[0],n);
 361
 362   return(0);
 363 }
 364
 365 vorbis_func_residue residue0_exportbundle={
 366   &pack,
 367   &unpack,
 368   &look,
 369   &free_info,
 370   &free_look,
 371   &forward,
 372   &inverse
 373 };