user/oggplay/Tremor/floor0.c

   1 /********************************************************************
   2  *                                                                  *
   3  * THIS FILE IS PART OF THE OggVorbis 'TREMOR' CODEC SOURCE CODE.   *
   4  *                                                                  *
   5  * USE, DISTRIBUTION AND REPRODUCTION OF THIS LIBRARY SOURCE IS     *
   6  * GOVERNED BY A BSD-STYLE SOURCE LICENSE INCLUDED WITH THIS SOURCE *
   7  * IN 'COPYING'. PLEASE READ THESE TERMS BEFORE DISTRIBUTING.       *
   8  *                                                                  *
   9  * THE OggVorbis 'TREMOR' SOURCE CODE IS (C) COPYRIGHT 1994-2002    *
  10  * BY THE Xiph.Org FOUNDATION http://www.xiph.org/                  *
  11  *                                                                  *
  12  ********************************************************************
  13
  14  function: floor backend 0 implementation
  15
  16  ********************************************************************/
  17
  18 #include <stdlib.h>
  19 #include <string.h>
  20 #include <math.h>
  21 #include "ogg.h"
  22 #include "ivorbiscodec.h"
  23 #include "codec_internal.h"
  24 #include "registry.h"
  25 #include "codebook.h"
  26 #include "misc.h"
  27 #include "os.h"
  28
  29 #define LSP_FRACBITS 14
  30
  31 typedef struct {
  32   long n;
  33   int ln;
  34   int  m;
  35   int *linearmap;
  36
  37   vorbis_info_floor0 *vi;
  38   ogg_int32_t *lsp_look;
  39
  40 } vorbis_look_floor0;
  41
  42 /*************** LSP decode ********************/
  43
  44 #include "lsp_lookup.h"
  45
  46 /* interpolated 1./sqrt(p) where .5 <= a < 1. (.100000... to .111111...) in
  47    16.16 format
  48    returns in m.8 format */
  49
  50 static long ADJUST_SQRT2[2]={8192,5792};
  51 static inline ogg_int32_t vorbis_invsqlook_i(long a,long e){
  52   long i=(a&0x7fff)>>(INVSQ_LOOKUP_I_SHIFT-1);
  53   long d=a&INVSQ_LOOKUP_I_MASK;                              /*  0.10 */
  54   long val=INVSQ_LOOKUP_I[i]-                                /*  1.16 */
  55     ((INVSQ_LOOKUP_IDel[i]*d)>>INVSQ_LOOKUP_I_SHIFT);        /* result 1.16 */
  56   val*=ADJUST_SQRT2[e&1];
  57   e=(e>>1)+21;
  58   return(val>>e);
  59 }
  60
  61 /* interpolated lookup based fromdB function, domain -140dB to 0dB only */
  62 /* a is in n.12 format */
  63 static inline ogg_int32_t vorbis_fromdBlook_i(long a){
  64   int i=(-a)>>(12-FROMdB2_SHIFT);
  65   if(i<0) return 0x7fffffff;
  66   if(i>=(FROMdB_LOOKUP_SZ<<FROMdB_SHIFT))return 0;
  67
  68   return FROMdB_LOOKUP[i>>FROMdB_SHIFT] * FROMdB2_LOOKUP[i&FROMdB2_MASK];
  69 }
  70
  71 /* interpolated lookup based cos function, domain 0 to PI only */
  72 /* a is in 0.16 format, where 0==0, 2^^16-1==PI, return 0.14 */
  73 static inline ogg_int32_t vorbis_coslook_i(long a){
  74   int i=a>>COS_LOOKUP_I_SHIFT;
  75   int d=a&COS_LOOKUP_I_MASK;
  76   return COS_LOOKUP_I[i]- ((d*(COS_LOOKUP_I[i]-COS_LOOKUP_I[i+1]))>>
  77                            COS_LOOKUP_I_SHIFT);
  78 }
  79
  80 /* interpolated lookup based cos function */
  81 /* a is in 0.16 format, where 0==0, 2^^16==PI, return .LSP_FRACBITS */
  82 static inline ogg_int32_t vorbis_coslook2_i(long a){
  83   a=a&0x1ffff;
  84
  85   if(a>0x10000)a=0x20000-a;
  86   {
  87     int i=a>>COS_LOOKUP_I_SHIFT;
  88     int d=a&COS_LOOKUP_I_MASK;
  89     a=((COS_LOOKUP_I[i]<<COS_LOOKUP_I_SHIFT)-
  90        d*(COS_LOOKUP_I[i]-COS_LOOKUP_I[i+1]))>>
  91       (COS_LOOKUP_I_SHIFT-LSP_FRACBITS+14);
  92   }
  93
  94   return(a);
  95 }
  96
  97 static int barklook[28]={
  98   0,100,200,301,          405,516,635,766,
  99   912,1077,1263,1476,     1720,2003,2333,2721,
 100   3184,3742,4428,5285,    6376,7791,9662,12181,
 101   15624,20397,27087,36554
 102 };
 103
 104 /* used in init only; interpolate the long way */
 105 static inline ogg_int32_t toBARK(int n){
 106   int i;
 107   for(i=0;i<27;i++)
 108     if(n>=barklook[i] && n<barklook[i+1])break;
 109
 110   if(i==27){
 111     return 27<<15;
 112   }else{
 113     int gap=barklook[i+1]-barklook[i];
 114     int del=n-barklook[i];
 115
 116     return((i<<15)+((del<<15)/gap));
 117   }
 118 }
 119
 120 static int MLOOP_1[64]={
 121    0,10,11,11, 12,12,12,12, 13,13,13,13, 13,13,13,13,
 122   14,14,14,14, 14,14,14,14, 14,14,14,14, 14,14,14,14,
 123   15,15,15,15, 15,15,15,15, 15,15,15,15, 15,15,15,15,
 124   15,15,15,15, 15,15,15,15, 15,15,15,15, 15,15,15,15,
 125 };
 126
 127 static int MLOOP_2[64]={
 128   0,4,5,5, 6,6,6,6, 7,7,7,7, 7,7,7,7,
 129   8,8,8,8, 8,8,8,8, 8,8,8,8, 8,8,8,8,
 130   9,9,9,9, 9,9,9,9, 9,9,9,9, 9,9,9,9,
 131   9,9,9,9, 9,9,9,9, 9,9,9,9, 9,9,9,9,
 132 };
 133
 134 static int MLOOP_3[8]={0,1,2,2,3,3,3,3};
 135
 136 void vorbis_lsp_to_curve(ogg_int32_t *curve,int *map,int n,int ln,
 137                          ogg_int32_t *lsp,int m,
 138                          ogg_int32_t amp,
 139                          ogg_int32_t ampoffset,
 140                          ogg_int32_t *icos){
 141
 142   /* 0 <= m < 256 */
 143
 144   /* set up for using all int later */
 145   int i;
 146   int ampoffseti=ampoffset*4096;
 147   int ampi=amp;
 148   ogg_int32_t *ilsp=(ogg_int32_t *)alloca(m*sizeof(*ilsp));
 149   ogg_int32_t invsq=0x517cc2;
 150   /* lsp is in 8.24, range 0 to PI; coslook wants it in .16 0 to 1*/
 151   for(i=0;i<m;i++){
 152     ogg_int32_t val=MULT32(lsp[i],invsq);
 153
 154     /* safeguard against a malicious stream */
 155     if(val<0 || (val>>COS_LOOKUP_I_SHIFT)>=COS_LOOKUP_I_SZ){
 156       memset(curve,0,sizeof(*curve)*n);
 157       return;
 158     }
 159
 160     ilsp[i]=vorbis_coslook_i(val);
 161   }
 162
 163   i=0;
 164   while(i<n){
 165     int j,k=map[i];
 166     ogg_uint32_t pi=46341; /* 2**-.5 in 0.16 */
 167     ogg_uint32_t qi=46341;
 168     ogg_int32_t qexp=0,shift;
 169     ogg_int32_t wi=icos[k];
 170
 171 #ifdef _V_LSP_MATH_ASM
 172     lsp_loop_asm(&qi,&pi,&qexp,ilsp,wi,m);
 173
 174     pi=((pi*pi)>>16);
 175     qi=((qi*qi)>>16);
 176
 177     if(m&1){
 178       qexp= qexp*2-28*((m+1)>>1)+m;
 179       pi*=(1<<14)-((wi*wi)>>14);
 180       qi+=pi>>14;
 181     }else{
 182       qexp= qexp*2-13*m;
 183
 184       pi*=(1<<14)-wi;
 185       qi*=(1<<14)+wi;
 186
 187       qi=(qi+pi)>>14;
 188     }
 189
 190     if(qi&0xffff0000){ /* checks for 1.xxxxxxxxxxxxxxxx */
 191       qi>>=1; qexp++;
 192     }else
 193       lsp_norm_asm(&qi,&qexp);
 194
 195 #else
 196
 197     qi*=labs(ilsp[0]-wi);
 198     pi*=labs(ilsp[1]-wi);
 199
 200     for(j=3;j<m;j+=2){
 201       if(!(shift=MLOOP_1[(pi|qi)>>25]))
 202         if(!(shift=MLOOP_2[(pi|qi)>>19]))
 203           shift=MLOOP_3[(pi|qi)>>16];
 204       qi=(qi>>shift)*labs(ilsp[j-1]-wi);
 205       pi=(pi>>shift)*labs(ilsp[j]-wi);
 206       qexp+=shift;
 207     }
 208     if(!(shift=MLOOP_1[(pi|qi)>>25]))
 209       if(!(shift=MLOOP_2[(pi|qi)>>19]))
 210         shift=MLOOP_3[(pi|qi)>>16];
 211
 212     /* pi,qi normalized collectively, both tracked using qexp */
 213
 214     if(m&1){
 215       /* odd order filter; slightly assymetric */
 216       /* the last coefficient */
 217       qi=(qi>>shift)*labs(ilsp[j-1]-wi);
 218       pi=(pi>>shift)<<14;
 219       qexp+=shift;
 220
 221       if(!(shift=MLOOP_1[(pi|qi)>>25]))
 222         if(!(shift=MLOOP_2[(pi|qi)>>19]))
 223           shift=MLOOP_3[(pi|qi)>>16];
 224
 225       pi>>=shift;
 226       qi>>=shift;
 227       qexp+=shift-14*((m+1)>>1);
 228
 229       pi=((pi*pi)>>16);
 230       qi=((qi*qi)>>16);
 231       qexp=qexp*2+m;
 232
 233       pi*=(1<<14)-((wi*wi)>>14);
 234       qi+=pi>>14;
 235
 236     }else{
 237       /* even order filter; still symmetric */
 238
 239       /* p*=p(1-w), q*=q(1+w), let normalization drift because it isn't
 240          worth tracking step by step */
 241
 242       pi>>=shift;
 243       qi>>=shift;
 244       qexp+=shift-7*m;
 245
 246       pi=((pi*pi)>>16);
 247       qi=((qi*qi)>>16);
 248       qexp=qexp*2+m;
 249
 250       pi*=(1<<14)-wi;
 251       qi*=(1<<14)+wi;
 252       qi=(qi+pi)>>14;
 253
 254     }
 255
 256
 257     /* we've let the normalization drift because it wasn't important;
 258        however, for the lookup, things must be normalized again.  We
 259        need at most one right shift or a number of left shifts */
 260
 261     if(qi&0xffff0000){ /* checks for 1.xxxxxxxxxxxxxxxx */
 262       qi>>=1; qexp++;
 263     }else
 264       while(qi && !(qi&0x8000)){ /* checks for 0.0xxxxxxxxxxxxxxx or less*/
 265         qi<<=1; qexp--;
 266       }
 267
 268 #endif
 269
 270     amp=vorbis_fromdBlook_i(ampi*                     /*  n.4         */
 271                             vorbis_invsqlook_i(qi,qexp)-
 272                                                       /*  m.8, m+n<=8 */
 273                             ampoffseti);              /*  8.12[0]     */
 274
 275     curve[i]= MULT31_SHIFT15(curve[i],amp);
 276     while(map[++i]==k) curve[i]= MULT31_SHIFT15(curve[i],amp);
 277   }
 278 }
 279
 280 /*************** vorbis decode glue ************/
 281
 282 static void floor0_free_info(vorbis_info_floor *i){
 283   vorbis_info_floor0 *info=(vorbis_info_floor0 *)i;
 284   if(info){
 285     memset(info,0,sizeof(*info));
 286     _ogg_free(info);
 287   }
 288 }
 289
 290 static void floor0_free_look(vorbis_look_floor *i){
 291   vorbis_look_floor0 *look=(vorbis_look_floor0 *)i;
 292   if(look){
 293
 294     if(look->linearmap)_ogg_free(look->linearmap);
 295     if(look->lsp_look)_ogg_free(look->lsp_look);
 296     memset(look,0,sizeof(*look));
 297     _ogg_free(look);
 298   }
 299 }
 300
 301 static vorbis_info_floor *floor0_unpack (vorbis_info *vi,oggpack_buffer *opb){
 302   codec_setup_info     *ci=(codec_setup_info *)vi->codec_setup;
 303   int j;
 304
 305   vorbis_info_floor0 *info=(vorbis_info_floor0 *)_ogg_malloc(sizeof(*info));
 306   info->order=oggpack_read(opb,8);
 307   info->rate=oggpack_read(opb,16);
 308   info->barkmap=oggpack_read(opb,16);
 309   info->ampbits=oggpack_read(opb,6);
 310   info->ampdB=oggpack_read(opb,8);
 311   info->numbooks=oggpack_read(opb,4)+1;
 312
 313   if(info->order<1)goto err_out;
 314   if(info->rate<1)goto err_out;
 315   if(info->barkmap<1)goto err_out;
 316   if(info->numbooks<1)goto err_out;
 317
 318   for(j=0;j<info->numbooks;j++){
 319     info->books[j]=oggpack_read(opb,8);
 320     if(info->books[j]<0 || info->books[j]>=ci->books)goto err_out;
 321   }
 322   return(info);
 323
 324  err_out:
 325   floor0_free_info(info);
 326   return(NULL);
 327 }
 328
 329 /* initialize Bark scale and normalization lookups.  We could do this
 330    with static tables, but Vorbis allows a number of possible
 331    combinations, so it's best to do it computationally.
 332
 333    The below is authoritative in terms of defining scale mapping.
 334    Note that the scale depends on the sampling rate as well as the
 335    linear block and mapping sizes */
 336
 337 static vorbis_look_floor *floor0_look (vorbis_dsp_state *vd,vorbis_info_mode *mi,
 338                               vorbis_info_floor *i){
 339   int j;
 340   ogg_int32_t scale;
 341   vorbis_info        *vi=vd->vi;
 342   codec_setup_info   *ci=(codec_setup_info *)vi->codec_setup;
 343   vorbis_info_floor0 *info=(vorbis_info_floor0 *)i;
 344   vorbis_look_floor0 *look=(vorbis_look_floor0 *)_ogg_calloc(1,sizeof(*look));
 345   look->m=info->order;
 346   look->n=ci->blocksizes[mi->blockflag]/2;
 347   look->ln=info->barkmap;
 348   look->vi=info;
 349
 350   /* the mapping from a linear scale to a smaller bark scale is
 351      straightforward.  We do *not* make sure that the linear mapping
 352      does not skip bark-scale bins; the decoder simply skips them and
 353      the encoder may do what it wishes in filling them.  They're
 354      necessary in some mapping combinations to keep the scale spacing
 355      accurate */
 356   look->linearmap=(int *)_ogg_malloc((look->n+1)*sizeof(*look->linearmap));
 357   for(j=0;j<look->n;j++){
 358
 359     int val=(look->ln*
 360              ((toBARK(info->rate/2*j/look->n)<<11)/toBARK(info->rate/2)))>>11;
 361
 362     if(val>=look->ln)val=look->ln-1; /* guard against the approximation */
 363     look->linearmap[j]=val;
 364   }
 365   look->linearmap[j]=-1;
 366
 367   look->lsp_look=(ogg_int32_t *)_ogg_malloc(look->ln*sizeof(*look->lsp_look));
 368   for(j=0;j<look->ln;j++)
 369     look->lsp_look[j]=vorbis_coslook2_i(0x10000*j/look->ln);
 370
 371   return look;
 372 }
 373
 374 static void *floor0_inverse1(vorbis_block *vb,vorbis_look_floor *i){
 375   vorbis_look_floor0 *look=(vorbis_look_floor0 *)i;
 376   vorbis_info_floor0 *info=look->vi;
 377   int j,k;
 378
 379   int ampraw=oggpack_read(&vb->opb,info->ampbits);
 380   if(ampraw>0){ /* also handles the -1 out of data case */
 381     long maxval=(1<<info->ampbits)-1;
 382     int amp=((ampraw*info->ampdB)<<4)/maxval;
 383     int booknum=oggpack_read(&vb->opb,_ilog(info->numbooks));
 384
 385     if(booknum!=-1 && booknum<info->numbooks){ /* be paranoid */
 386       codec_setup_info  *ci=(codec_setup_info *)vb->vd->vi->codec_setup;
 387       codebook *b=ci->fullbooks+info->books[booknum];
 388       ogg_int32_t last=0;
 389       ogg_int32_t *lsp=(ogg_int32_t *)_vorbis_block_alloc(vb,sizeof(*lsp)*(look->m+1));
 390
 391       for(j=0;j<look->m;j+=b->dim)
 392         if(vorbis_book_decodev_set(b,lsp+j,&vb->opb,b->dim,-24)==-1)goto eop;
 393       for(j=0;j<look->m;){
 394         for(k=0;k<b->dim;k++,j++)lsp[j]+=last;
 395         last=lsp[j-1];
 396       }
 397
 398       lsp[look->m]=amp;
 399       return(lsp);
 400     }
 401   }
 402  eop:
 403   return(NULL);
 404 }
 405
 406 static int floor0_inverse2(vorbis_block *vb,vorbis_look_floor *i,
 407                            void *memo,ogg_int32_t *out){
 408   vorbis_look_floor0 *look=(vorbis_look_floor0 *)i;
 409   vorbis_info_floor0 *info=look->vi;
 410
 411   if(memo){
 412     ogg_int32_t *lsp=(ogg_int32_t *)memo;
 413     ogg_int32_t amp=lsp[look->m];
 414
 415     /* take the coefficients back to a spectral envelope curve */
 416     vorbis_lsp_to_curve(out,look->linearmap,look->n,look->ln,
 417                         lsp,look->m,amp,info->ampdB,look->lsp_look);
 418     return(1);
 419   }
 420   memset(out,0,sizeof(*out)*look->n);
 421   return(0);
 422 }
 423
 424 /* export hooks */
 425 vorbis_func_floor floor0_exportbundle={
 426   &floor0_unpack,&floor0_look,&floor0_free_info,
 427   &floor0_free_look,&floor0_inverse1,&floor0_inverse2
 428 };
 429
 430