summaryrefslogtreecommitdiff
diff options
context:
space:
mode:
authorMonty <xiphmont@xiph.org>2001-04-29 22:21:06 +0000
committerMonty <xiphmont@xiph.org>2001-04-29 22:21:06 +0000
commit2b1605de070613fa21f435844d9c5ea8a3688f19 (patch)
tree26c9fdaad52bcd1c3e223128948753f711761f7d
parent78540c3db714ab4dda28e11659f2a69e9961bf69 (diff)
downloadlibvorbis-git-2b1605de070613fa21f435844d9c5ea8a3688f19.tar.gz
...such that I don;t lose anything (floor 1 not yet tested/functional)
Monty svn path=/branches/monty-branch-20010404/vorbis/; revision=1435
-rw-r--r--lib/backends.h32
-rw-r--r--lib/barkmel.c64
-rw-r--r--lib/floor0.c10
-rw-r--r--lib/floor1.c949
-rw-r--r--lib/mapping0.c401
-rw-r--r--lib/modes/mode_A.h316
-rw-r--r--lib/psy.c727
-rw-r--r--lib/psy.h101
8 files changed, 2477 insertions, 123 deletions
diff --git a/lib/backends.h b/lib/backends.h
index 72097000..8c64ca41 100644
--- a/lib/backends.h
+++ b/lib/backends.h
@@ -12,7 +12,7 @@
function: libvorbis backend and mapping structures; needed for
static mode headers
- last mod: $Id: backends.h,v 1.6.4.1 2001/04/05 00:22:48 xiphmont Exp $
+ last mod: $Id: backends.h,v 1.6.4.2 2001/04/29 22:21:04 xiphmont Exp $
********************************************************************/
@@ -61,7 +61,7 @@ typedef struct{
void (*free_info) (vorbis_info_floor *);
void (*free_look) (vorbis_look_floor *);
int (*forward) (struct vorbis_block *,vorbis_look_floor *,
- float *, float *);
+ float *);
int (*inverse) (struct vorbis_block *,vorbis_look_floor *,
float *);
} vorbis_func_floor;
@@ -82,21 +82,27 @@ typedef struct{
} vorbis_info_floor0;
+#define VIF_POSIT 63
+#define VIF_CLASS 16
+#define VIF_PARTS 31
typedef struct{
- long maxrange;
- int rangeres; /* 6, 7 or 8 */
+ int partitions; /* 0 to 31 */
+ int partitionclass[VIF_PARTS]; /* 0 to 15 */
- int ranges; /* up to 8 */
- int rangesizes[8]; /* 1 to 8 */
- int rangelist[64];
-
- int positionbook[8];
- int ampbook[8];
+ int class_dim[VIF_CLASS]; /* 1 to 8 */
+ int class_subs[VIF_CLASS]; /* 0,1,2,3 (bits: 1<<n poss) */
+ int class_book[VIF_CLASS]; /* subs ^ dim entries */
+ int class_subbook[VIF_CLASS][8]; /* [VIF_CLASS][subs] */
+
+
+ int mult; /* 1 2 3 or 4 */
+ int postlist[VIF_POSIT+2]; /* first two implicit */
- float maxover; /* encode side analysis parameter */
- float maxunder; /* encode side analysis parameter */
- float maxerr; /* encode side analysis parameter */
+ float maxover; /* encode side analysis parameter */
+ float maxunder; /* encode side analysis parameter */
+ float maxerr; /* encode side analysis parameter */
+ int searchstart; /* encode side analysis parameter */
} vorbis_info_floor1;
/* Residue backend generic *****************************************/
diff --git a/lib/barkmel.c b/lib/barkmel.c
new file mode 100644
index 00000000..1f3ba38c
--- /dev/null
+++ b/lib/barkmel.c
@@ -0,0 +1,64 @@
+/********************************************************************
+ * *
+ * THIS FILE IS PART OF THE OggVorbis SOFTWARE CODEC SOURCE CODE. *
+ * USE, DISTRIBUTION AND REPRODUCTION OF THIS LIBRARY SOURCE IS *
+ * GOVERNED BY A BSD-STYLE SOURCE LICENSE INCLUDED WITH THIS SOURCE *
+ * IN 'COPYING'. PLEASE READ THESE TERMS BEFORE DISTRIBUTING. *
+ * *
+ * THE OggVorbis SOURCE CODE IS (C) COPYRIGHT 1994-2001 *
+ * by the XIPHOPHORUS Company http://www.xiph.org/ *
+
+ ********************************************************************
+
+ function: bark scale utility
+ last mod: $Id: barkmel.c,v 1.5.4.1 2001/04/29 22:21:04 xiphmont Exp $
+
+ ********************************************************************/
+
+#include <stdio.h>
+#include "scales.h"
+int main(){
+ int i;
+ double rate;
+ for(i=64;i<32000;i*=2){
+ rate=48000.f;
+ fprintf(stderr,"rate=%gHz, block=%d, f(1)=%.2gHz bark(1)=%.2g (of %.2g)\n",
+ rate,i,rate/2 / (i/2),toBARK(rate/2 /(i/2)),toBARK(rate/2));
+
+ rate=44100.f;
+ fprintf(stderr,"rate=%gHz, block=%d, f(1)=%.2gHz bark(1)=%.2g (of %.2g)\n",
+ rate,i,rate/2 / (i/2),toBARK(rate/2 /(i/2)),toBARK(rate/2));
+
+ rate=32000.f;
+ fprintf(stderr,"rate=%gHz, block=%d, f(1)=%.2gHz bark(1)=%.2g (of %.2g)\n",
+ rate,i,rate/2 / (i/2),toBARK(rate/2 /(i/2)),toBARK(rate/2));
+
+ rate=22050.f;
+ fprintf(stderr,"rate=%gHz, block=%d, f(1)=%.2gHz bark(1)=%.2g (of %.2g)\n",
+ rate,i,rate/2 / (i/2),toBARK(rate/2 /(i/2)),toBARK(rate/2));
+
+ rate=16000.f;
+ fprintf(stderr,"rate=%gHz, block=%d, f(1)=%.2gHz bark(1)=%.2g (of %.2g)\n",
+ rate,i,rate/2 / (i/2),toBARK(rate/2 /(i/2)),toBARK(rate/2));
+
+ rate=11025.f;
+ fprintf(stderr,"rate=%gHz, block=%d, f(1)=%.2gHz bark(1)=%.2g (of %.2g)\n",
+ rate,i,rate/2 / (i/2),toBARK(rate/2 /(i/2)),toBARK(rate/2));
+
+ rate=8000.f;
+ fprintf(stderr,"rate=%gHz, block=%d, f(1)=%.2gHz bark(1)=%.2g (of %.2g)\n\n",
+ rate,i,rate/2 / (i/2),toBARK(rate/2 /(i/2)),toBARK(rate/2));
+
+
+ }
+ {
+ float i;
+ int j;
+ for(i=0.,j=0;i<28;i+=.816,j++){
+ fprintf(stderr,"(%d) bark=%f %gHz (%d of 1024)\n",
+ j,i,fromBARK(i),(int)(fromBARK(i)/22050.*1024.));
+ }
+ }
+ return(0);
+}
+
diff --git a/lib/floor0.c b/lib/floor0.c
index 6be7fa90..e2f443f1 100644
--- a/lib/floor0.c
+++ b/lib/floor0.c
@@ -11,7 +11,7 @@
********************************************************************
function: floor backend 0 implementation
- last mod: $Id: floor0.c,v 1.39.4.1 2001/04/05 00:22:48 xiphmont Exp $
+ last mod: $Id: floor0.c,v 1.39.4.2 2001/04/29 22:21:04 xiphmont Exp $
********************************************************************/
@@ -240,8 +240,8 @@ float _curve_to_lpc(float *curve,float *lpc,
for(i=0;i<l->n;i++)
curve[i]-=150;
- _analysis_output_always("barkfloor",seq,work,bark,0,0);
- _analysis_output_always("barkcurve",seq++,curve,l->n,1,0);
+ _analysis_output("barkfloor",seq,work,bark,0,0);
+ _analysis_output("barkcurve",seq++,curve,l->n,1,0);
for(i=0;i<l->n;i++)
curve[i]+=150;
@@ -255,7 +255,7 @@ float _curve_to_lpc(float *curve,float *lpc,
/* didn't need in->out seperation, modifies the flr[] vector; takes in
a dB scale floor, puts out linear */
static int floor0_forward(vorbis_block *vb,vorbis_look_floor *i,
- float *flr,float *dummy){
+ float *flr){
long j;
vorbis_look_floor0 *look=(vorbis_look_floor0 *)i;
vorbis_info_floor0 *info=look->vi;
@@ -385,7 +385,7 @@ static int floor0_forward(vorbis_block *vb,vorbis_look_floor *i,
vorbis_lsp_to_curve(flr,look->linearmap,look->n,look->ln,
lspwork,look->m,amp,info->ampdB);
- _analysis_output_always("barklsp",seq-1,flr,look->n,1,1);
+ _analysis_output("barklsp",seq-1,flr,look->n,1,1);
_analysis_output("lsp3",seq-1,flr,look->n,0,1);
return(val);
}
diff --git a/lib/floor1.c b/lib/floor1.c
index 6457f506..84128f80 100644
--- a/lib/floor1.c
+++ b/lib/floor1.c
@@ -11,7 +11,7 @@
********************************************************************
function: floor backend 1 implementation
- last mod: $Id: floor1.c,v 1.1.2.1 2001/04/05 00:22:48 xiphmont Exp $
+ last mod: $Id: floor1.c,v 1.1.2.2 2001/04/29 22:21:04 xiphmont Exp $
********************************************************************/
@@ -24,21 +24,43 @@
#include "registry.h"
#include "codebook.h"
#include "misc.h"
+#include "scales.h"
#include <stdio.h>
+#define VORBIS_IEEE_FLOAT32
+
#define floor1_rangedB 140 /* floor 0 fixed at -140dB to 0dB range */
typedef struct {
- int lo_neighbor[66];
- int hi_neighbor[66];
+ int sorted_index[VIF_POSIT+2]; /* en/de */
+ int forward_index[VIF_POSIT+2]; /* en/de */
+ int reverse_index[VIF_POSIT+2]; /* en */
+
+ int hineighbor[VIF_POSIT]; /* de */
+ int loneighbor[VIF_POSIT]; /* de */
+ int posts;
+
int n;
+ int quant_q;
vorbis_info_floor1 *vi;
} vorbis_look_floor1;
-/***********************************************/
+typedef struct lsfit_acc{
+ long x0;
+ long x1;
+ long xa;
+ long ya;
+ long x2a;
+ long y2a;
+ long xya;
+ long n;
+} lsfit_acc;
+
+/***********************************************/
+
static vorbis_info_floor *floor1_copy_info (vorbis_info_floor *i){
vorbis_info_floor1 *info=(vorbis_info_floor1 *)i;
vorbis_info_floor1 *ret=_ogg_malloc(sizeof(vorbis_info_floor1));
@@ -61,6 +83,15 @@ static void floor1_free_look(vorbis_look_floor *i){
}
}
+static int ilog(unsigned int v){
+ int ret=0;
+ while(v){
+ ret++;
+ v>>=1;
+ }
+ return(ret);
+}
+
static int ilog2(unsigned int v){
int ret=0;
while(v>1){
@@ -75,73 +106,103 @@ static void floor1_pack (vorbis_info_floor *i,oggpack_buffer *opb){
int j,k;
int count=0;
int rangebits;
+ int maxposit=info->postlist[1];
+ int maxclass=-1;
- oggpack_write(opb,info->maxrange,24);
- rangebits=ilog2(info->maxrange);
-
- oggpack_write(opb,info->ranges,4); /* only 0 to 8 legal now */
- oggpack_write(opb,info->rangeres-1,3); /* only 6,7 or 8 legal now */
- for(j=0,k=0;j<info->ranges;j++){
- oggpack_write(opb,info->rangesizes[j]-1,3);
- oggpack_write(opb,info->positionbook[j],8);
- oggpack_write(opb,info->ampbook[j],8);
- count+=info->rangesizes[j];
+ /* save out partitions */
+ oggpack_write(opb,info->partitions,5); /* only 0 to 31 legal */
+ for(j=0;j<info->partitions;j++){
+ oggpack_write(opb,info->partitionclass[j],4); /* only 0 to 15 legal */
+ if(maxclass<info->partitionclass[j])maxclass=info->partitionclass[j];
+ }
+
+ /* save out partition classes */
+ for(j=0;j<maxclass+1;j++){
+ oggpack_write(opb,info->class_dim[j]-1,3); /* 1 to 8 */
+ oggpack_write(opb,info->class_subs[j],2); /* 0 to 3 */
+ if(info->class_subs[j])oggpack_write(opb,info->class_book[j],8);
+ for(k=0;k<info->class_subs[j];k++)
+ oggpack_write(opb,info->class_subbook[j][k]+1,8);
+ }
+
+ /* save out the post list */
+ oggpack_write(opb,info->mult-1,2); /* only 1,2,3,4 legal now */
+ oggpack_write(opb,ilog2(maxposit),4);
+ rangebits=ilog2(maxposit);
+
+ for(j=0,k=0;j<info->partitions;j++){
+ count+=info->class_dim[info->partitionclass[j]];
for(;k<count;k++)
- oggpack_write(opb,info->rangelist[k],rangebits);
+ oggpack_write(opb,info->postlist[k+2],rangebits);
}
}
+
static vorbis_info_floor *floor1_unpack (vorbis_info *vi,oggpack_buffer *opb){
codec_setup_info *ci=vi->codec_setup;
- int j,k,count;
-
- vorbis_info_floor1 *info=_ogg_malloc(sizeof(vorbis_info_floor1));
- int rangebits;
-
- info->maxrange=oggpack_read(opb,24);
- rangebits=ilog2(info->maxrange);
-
- info->ranges=oggpack_read(opb,4);
- if(info->ranges>8)goto err_out;
+ int j,k,count,maxclass=-1,rangebits;
- info->rangeres=oggpack_read(opb,3)+1;
- if(info->rangeres<6)goto err_out;
- if(info->rangeres>8)goto err_out;
+ vorbis_info_floor1 *info=_ogg_calloc(1,sizeof(vorbis_info_floor1));
+ /* read partitions */
+ info->partitions=oggpack_read(opb,5); /* only 0 to 31 legal */
+ for(j=0;j<info->partitions;j++){
+ info->partitionclass[j]=oggpack_read(opb,4); /* only 0 to 15 legal */
+ if(maxclass<info->partitionclass[j])maxclass=info->partitionclass[j];
+ }
- for(j=0,k=0;j<info->ranges;j++){
- info->rangesizes[j]=oggpack_read(opb,3)+1;
- info->positionbook[j]=oggpack_read(opb,8);
- if(info->positionbook[j]<0 || info->positionbook[j]>=ci->books)
+ /* read partition classes */
+ for(j=0;j<maxclass+1;j++){
+ info->class_dim[j]=oggpack_read(opb,3)+1; /* 1 to 8 */
+ info->class_subs[j]=oggpack_read(opb,2); /* 0,1,2,3 bits */
+ if(info->class_subs[j]<0)
goto err_out;
-
- info->ampbook[j]=oggpack_read(opb,8);
- if(info->ampbook[j]<0 || info->ampbook[j]>=ci->books)
+ if(info->class_subs[j])info->class_book[j]=oggpack_read(opb,8);
+ if(info->class_book[j]<0 || info->class_book[j]>=ci->books)
goto err_out;
- count+=info->rangesizes[j];
+ for(k=0;k<info->class_subs[j];k++){
+ info->class_subbook[j][k]=oggpack_read(opb,8)-1;
+ if(info->class_subbook[j][k]<-1 || info->class_subbook[j][k]>=ci->books)
+ goto err_out;
+ }
+ }
+
+ /* read the post list */
+ info->mult=oggpack_read(opb,2)+1; /* only 1,2,3,4 legal now */
+ rangebits=oggpack_read(opb,4);
+
+ for(j=0,k=0;j<info->partitions;j++){
+ count+=info->class_dim[info->partitionclass[j]];
for(;k<count;k++){
- info->rangelist[k]=oggpack_read(opb,rangebits);
- if(info->rangelist[k]<0 || info->rangelist[k]>=maxrange)
+ int t=info->postlist[k+2]=oggpack_read(opb,rangebits);
+ if(t<0 || t>=(1<<rangebits))
goto err_out;
}
}
+ info->postlist[0]=0;
+ info->postlist[1]=1<<rangebits;
+
+ return(info);
err_out:
floor1_free_info(info);
return(NULL);
}
+static int icomp(const void *a,const void *b){
+ return(**(int **)a-**(int **)b);
+}
+
static vorbis_look_floor *floor1_look(vorbis_dsp_state *vd,vorbis_info_mode *mi,
- vorbis_info_floor *i){
+ vorbis_info_floor *in){
- vorbis_info *vi=vd->vi;
- codec_setup_info *ci=vi->codec_setup;
- vorbis_info_floor1 *info=(vorbis_info_floor1 *)i;
+ int *sortpointer[VIF_POSIT+2];
+ vorbis_info_floor1 *info=(vorbis_info_floor1 *)in;
vorbis_look_floor1 *look=_ogg_calloc(1,sizeof(vorbis_look_floor1));
int i,j,n=0;
- look->vi=vi;
- look->n=ci->blocksizes[mi->blockflag]/2;
-
+ look->vi=info;
+ look->n=info->postlist[1];
+
/* we drop each position value in-between already decoded values,
and use linear interpolation to predict each new value past the
edges. The positions are read in the order of the position
@@ -149,37 +210,232 @@ static vorbis_look_floor *floor1_look(vorbis_dsp_state *vd,vorbis_info_mode *mi,
course, the neighbors can change (if a position is declined), but
this is an initial mapping */
- for(i=0;i<info->ranges;i++)n+=info->rangesizes[j];
- for(i=0;i<info->ranges;i++){
- int hi=-1;
- int lo=-1;
- for(j=0;j<i;j++){
- if(info->rangelist[j]<info->rangelist[i] &&
- (lo==-1 || info->rangelist[j]>info->rangelist[lo]))
+ for(i=0;i<info->partitions;i++)n+=info->class_dim[info->partitionclass[i]];
+ n+=2;
+ look->posts=n;
+
+ /* also store a sorted position index */
+ for(i=0;i<n;i++)sortpointer[i]=info->postlist+i;
+ qsort(sortpointer,n,sizeof(int),icomp);
+
+ /* points from sort order back to range number */
+ for(i=0;i<n;i++)look->forward_index[i]=sortpointer[i]-info->postlist;
+ /* points from range order to sorted position */
+ for(i=0;i<n;i++)look->reverse_index[look->forward_index[i]]=i;
+ /* we actually need the post values too */
+ for(i=0;i<n;i++)look->sorted_index[i]=info->postlist[look->forward_index[i]];
+
+ /* quantize values to multiplier spec */
+ switch(info->mult){
+ case 1: /* 1024 -> 256 */
+ look->quant_q=256;
+ break;
+ case 2: /* 1024 -> 128 */
+ look->quant_q=128;
+ break;
+ case 3: /* 1024 -> 86 */
+ look->quant_q=86;
+ break;
+ case 4: /* 1024 -> 64 */
+ look->quant_q=64;
+ break;
+ }
+
+ /* discover our neighbors for decode where we don't use fit flags
+ (that would push the neighbors outward) */
+ for(i=0;i<n-2;i++){
+ int lo=0;
+ int hi=1;
+ int lx=0;
+ int hx=look->n;
+ int currentx=info->postlist[i+2];
+ for(j=0;j<i+2;j++){
+ int x=info->postlist[j];
+ if(x>lx && x<currentx){
lo=j;
- if(info->rangelist[j]>info->rangelist[i] &&
- (hi==-1 || info->rangelist[j]<info->rangelist[hi]))
+ lx=x;
+ }
+ if(x<hx && x>currentx){
hi=j;
+ hx=x;
+ }
}
- look->lo_neighbor[i]=lo;
- look->hi_neighbor[i]=hi;
+ look->loneighbor[i]=lo;
+ look->hineighbor[i]=hi;
}
return(look);
}
-static int render_point_dB(int x0,int x1,int y0,int y1,int x){
+static int render_point(int x0,int x1,int y0,int y1,int x){
int dy=y1-y0;
int adx=x1-x0;
int ady=abs(dy);
- int err=(ady>>1)+ady*(x-x0);
+ int err=ady*(x-x0);
- int off=err%adx;
+ int off=err/adx;
if(dy<0)return(y0-off);
return(y0+off);
}
-static int render_line_dB(int x0,int x1,int y0,int y1,int *d){
+static int FLOOR_todB_LOOKUP[560]={
+ 1023, 1021, 1019, 1018, 1016, 1014, 1012, 1010,
+ 1008, 1007, 1005, 1003, 1001, 999, 997, 996,
+ 994, 992, 990, 988, 986, 985, 983, 981,
+ 979, 977, 975, 974, 972, 970, 968, 966,
+ 964, 963, 961, 959, 957, 955, 954, 952,
+ 950, 948, 946, 944, 943, 941, 939, 937,
+ 935, 933, 932, 930, 928, 926, 924, 922,
+ 921, 919, 917, 915, 913, 911, 910, 908,
+ 906, 904, 902, 900, 899, 897, 895, 893,
+ 891, 890, 888, 886, 884, 882, 880, 879,
+ 877, 875, 873, 871, 869, 868, 866, 864,
+ 862, 860, 858, 857, 855, 853, 851, 849,
+ 847, 846, 844, 842, 840, 838, 836, 835,
+ 833, 831, 829, 827, 826, 824, 822, 820,
+ 818, 816, 815, 813, 811, 809, 807, 805,
+ 804, 802, 800, 798, 796, 794, 793, 791,
+ 789, 787, 785, 783, 782, 780, 778, 776,
+ 774, 772, 771, 769, 767, 765, 763, 762,
+ 760, 758, 756, 754, 752, 751, 749, 747,
+ 745, 743, 741, 740, 738, 736, 734, 732,
+ 730, 729, 727, 725, 723, 721, 719, 718,
+ 716, 714, 712, 710, 708, 707, 705, 703,
+ 701, 699, 698, 696, 694, 692, 690, 688,
+ 687, 685, 683, 681, 679, 677, 676, 674,
+ 672, 670, 668, 666, 665, 663, 661, 659,
+ 657, 655, 654, 652, 650, 648, 646, 644,
+ 643, 641, 639, 637, 635, 634, 632, 630,
+ 628, 626, 624, 623, 621, 619, 617, 615,
+ 613, 612, 610, 608, 606, 604, 602, 601,
+ 599, 597, 595, 593, 591, 590, 588, 586,
+ 584, 582, 580, 579, 577, 575, 573, 571,
+ 570, 568, 566, 564, 562, 560, 559, 557,
+ 555, 553, 551, 549, 548, 546, 544, 542,
+ 540, 538, 537, 535, 533, 531, 529, 527,
+ 526, 524, 522, 520, 518, 516, 515, 513,
+ 511, 509, 507, 506, 504, 502, 500, 498,
+ 496, 495, 493, 491, 489, 487, 485, 484,
+ 482, 480, 478, 476, 474, 473, 471, 469,
+ 467, 465, 463, 462, 460, 458, 456, 454,
+ 452, 451, 449, 447, 445, 443, 442, 440,
+ 438, 436, 434, 432, 431, 429, 427, 425,
+ 423, 421, 420, 418, 416, 414, 412, 410,
+ 409, 407, 405, 403, 401, 399, 398, 396,
+ 394, 392, 390, 388, 387, 385, 383, 381,
+ 379, 378, 376, 374, 372, 370, 368, 367,
+ 365, 363, 361, 359, 357, 356, 354, 352,
+ 350, 348, 346, 345, 343, 341, 339, 337,
+ 335, 334, 332, 330, 328, 326, 324, 323,
+ 321, 319, 317, 315, 314, 312, 310, 308,
+ 306, 304, 303, 301, 299, 297, 295, 293,
+ 292, 290, 288, 286, 284, 282, 281, 279,
+ 277, 275, 273, 271, 270, 268, 266, 264,
+ 262, 260, 259, 257, 255, 253, 251, 250,
+ 248, 246, 244, 242, 240, 239, 237, 235,
+ 233, 231, 229, 228, 226, 224, 222, 220,
+ 218, 217, 215, 213, 211, 209, 207, 206,
+ 204, 202, 200, 198, 196, 195, 193, 191,
+ 189, 187, 186, 184, 182, 180, 178, 176,
+ 175, 173, 171, 169, 167, 165, 164, 162,
+ 160, 158, 156, 154, 153, 151, 149, 147,
+ 145, 143, 142, 140, 138, 136, 134, 132,
+ 131, 129, 127, 125, 123, 122, 120, 118,
+ 116, 114, 112, 111, 109, 107, 105, 103,
+ 101, 100, 98, 96, 94, 92, 90, 89,
+ 87, 85, 83, 81, 79, 78, 76, 74,
+ 72, 70, 68, 67, 65, 63, 61, 59,
+ 58, 56, 54, 52, 50, 48, 47, 45,
+ 43, 41, 39, 37, 36, 34, 32, 30,
+ 28, 26, 25, 23, 21, 19, 17, 15,
+ 14, 12, 10, 8, 6, 4, 3, 1,
+};
+
+static float FLOOR_fromdB_LOOKUP[256]={
+ 1.0649863e-07F, 1.1341951e-07F, 1.2079015e-07F, 1.2863978e-07F,
+ 1.3699951e-07F, 1.4590251e-07F, 1.5538408e-07F, 1.6548181e-07F,
+ 1.7623575e-07F, 1.8768855e-07F, 1.9988561e-07F, 2.128753e-07F,
+ 2.2670913e-07F, 2.4144197e-07F, 2.5713223e-07F, 2.7384213e-07F,
+ 2.9163793e-07F, 3.1059021e-07F, 3.3077411e-07F, 3.5226968e-07F,
+ 3.7516214e-07F, 3.9954229e-07F, 4.2550680e-07F, 4.5315863e-07F,
+ 4.8260743e-07F, 5.1396998e-07F, 5.4737065e-07F, 5.8294187e-07F,
+ 6.2082472e-07F, 6.6116941e-07F, 7.0413592e-07F, 7.4989464e-07F,
+ 7.9862701e-07F, 8.5052630e-07F, 9.0579828e-07F, 9.6466216e-07F,
+ 1.0273513e-06F, 1.0941144e-06F, 1.1652161e-06F, 1.2409384e-06F,
+ 1.3215816e-06F, 1.4074654e-06F, 1.4989305e-06F, 1.5963394e-06F,
+ 1.7000785e-06F, 1.8105592e-06F, 1.9282195e-06F, 2.0535261e-06F,
+ 2.1869758e-06F, 2.3290978e-06F, 2.4804557e-06F, 2.6416497e-06F,
+ 2.8133190e-06F, 2.9961443e-06F, 3.1908506e-06F, 3.3982101e-06F,
+ 3.6190449e-06F, 3.8542308e-06F, 4.1047004e-06F, 4.3714470e-06F,
+ 4.6555282e-06F, 4.9580707e-06F, 5.2802740e-06F, 5.6234160e-06F,
+ 5.9888572e-06F, 6.3780469e-06F, 6.7925283e-06F, 7.2339451e-06F,
+ 7.7040476e-06F, 8.2047000e-06F, 8.7378876e-06F, 9.3057248e-06F,
+ 9.9104632e-06F, 1.0554501e-05F, 1.1240392e-05F, 1.1970856e-05F,
+ 1.2748789e-05F, 1.3577278e-05F, 1.4459606e-05F, 1.5399272e-05F,
+ 1.6400004e-05F, 1.7465768e-05F, 1.8600792e-05F, 1.9809576e-05F,
+ 2.1096914e-05F, 2.2467911e-05F, 2.3928002e-05F, 2.5482978e-05F,
+ 2.7139006e-05F, 2.8902651e-05F, 3.0780908e-05F, 3.2781225e-05F,
+ 3.4911534e-05F, 3.7180282e-05F, 3.9596466e-05F, 4.2169667e-05F,
+ 4.4910090e-05F, 4.7828601e-05F, 5.0936773e-05F, 5.4246931e-05F,
+ 5.7772202e-05F, 6.1526565e-05F, 6.5524908e-05F, 6.9783085e-05F,
+ 7.4317983e-05F, 7.9147585e-05F, 8.4291040e-05F, 8.9768747e-05F,
+ 9.5602426e-05F, 0.00010181521F, 0.00010843174F, 0.00011547824F,
+ 0.00012298267F, 0.00013097477F, 0.00013948625F, 0.00014855085F,
+ 0.00015820453F, 0.00016848555F, 0.00017943469F, 0.00019109536F,
+ 0.00020351382F, 0.00021673929F, 0.00023082423F, 0.00024582449F,
+ 0.00026179955F, 0.00027881276F, 0.00029693158F, 0.00031622787F,
+ 0.00033677814F, 0.00035866388F, 0.00038197188F, 0.00040679456F,
+ 0.00043323036F, 0.00046138411F, 0.00049136745F, 0.00052329927F,
+ 0.00055730621F, 0.00059352311F, 0.00063209358F, 0.00067317058F,
+ 0.00071691700F, 0.00076350630F, 0.00081312324F, 0.00086596457F,
+ 0.00092223983F, 0.00098217216F, 0.0010459992F, 0.0011139742F,
+ 0.0011863665F, 0.0012634633F, 0.0013455702F, 0.0014330129F,
+ 0.0015261382F, 0.0016253153F, 0.0017309374F, 0.0018434235F,
+ 0.0019632195F, 0.0020908006F, 0.0022266726F, 0.0023713743F,
+ 0.0025254795F, 0.0026895994F, 0.0028643847F, 0.0030505286F,
+ 0.0032487691F, 0.0034598925F, 0.0036847358F, 0.0039241906F,
+ 0.0041792066F, 0.0044507950F, 0.0047400328F, 0.0050480668F,
+ 0.0053761186F, 0.0057254891F, 0.0060975636F, 0.0064938176F,
+ 0.0069158225F, 0.0073652516F, 0.0078438871F, 0.0083536271F,
+ 0.0088964928F, 0.009474637F, 0.010090352F, 0.010746080F,
+ 0.011444421F, 0.012188144F, 0.012980198F, 0.013823725F,
+ 0.014722068F, 0.015678791F, 0.016697687F, 0.017782797F,
+ 0.018938423F, 0.020169149F, 0.021479854F, 0.022875735F,
+ 0.024362330F, 0.025945531F, 0.027631618F, 0.029427276F,
+ 0.031339626F, 0.033376252F, 0.035545228F, 0.037855157F,
+ 0.040315199F, 0.042935108F, 0.045725273F, 0.048696758F,
+ 0.051861348F, 0.055231591F, 0.058820850F, 0.062643361F,
+ 0.066714279F, 0.071049749F, 0.075666962F, 0.080584227F,
+ 0.085821044F, 0.091398179F, 0.097337747F, 0.10366330F,
+ 0.11039993F, 0.11757434F, 0.12521498F, 0.13335215F,
+ 0.14201813F, 0.15124727F, 0.16107617F, 0.17154380F,
+ 0.18269168F, 0.19456402F, 0.20720788F, 0.22067342F,
+ 0.23501402F, 0.25028656F, 0.26655159F, 0.28387361F,
+ 0.30232132F, 0.32196786F, 0.34289114F, 0.36517414F,
+ 0.38890521F, 0.41417847F, 0.44109412F, 0.46975890F,
+ 0.50028648F, 0.53279791F, 0.56742212F, 0.60429640F,
+ 0.64356699F, 0.68538959F, 0.72993007F, 0.77736504F,
+ 0.82788260F, 0.88168307F, 0.9389798F, 1.F,
+};
+
+#ifdef VORBIS_IEEE_FLOAT32
+static int vorbis_floor1_dBquant(float *x){
+ float temp=*x-256.f;
+ ogg_uint32_t *i=(ogg_uint32_t *)(&temp);
+ if(*i<(ogg_uint32_t)0xc3800000)return(1023);
+ if(*i>(ogg_uint32_t)0xc3c5e000)return(0);
+ return FLOOR_quantdB_LOOKUP[(*i-0xc3800000)>>13];
+}
+#else
+static int vorbis_floor1_dBquant(float *x){
+ int i= ((*x)+140.)/140.*1024.+.5;
+ if(i>1023)return(1023);
+ if(i<0)return(0);
+ return i;
+}
+#endif
+
+static void render_line(int x0,int x1,int y0,int y1,float *d){
int dy=y1-y0;
int adx=x1-x0;
int ady=abs(dy);
@@ -187,95 +443,579 @@ static int render_line_dB(int x0,int x1,int y0,int y1,int *d){
int sy=(dy<0?base-1:base+1);
int x=x0;
int y=y0;
- int err;
+ int err=0;
- ady-=base*adx;
- err=ady>>1;
+ ady-=abs(base*adx);
- d[x]=y;
+ d[x]=FLOOR_fromdB_LOOKUP[y];
while(++x<x1){
err=err+ady;
- if(err>adx){
+ if(err>=adx){
err-=adx;
y+=sy;
}else{
y+=base;
}
- d[x]=y;
+ d[x]=FLOOR_fromdB_LOOKUP[y];
}
}
/* the floor has already been filtered to only include relevant sections */
-static int fit_line(int *floor,int x0, int x1,int *y0, int *y1){
- long i,n=0;
- long x=0,y=0,x2=0,y2=0,xy=0;
-
- if(*y0>-900){ /* hint used to break degenerate cases */
- x+=x0;
- y+=*y0;
- x2+=(x0*x0);
+static int accumulate_fit(float *floor,int x0, int x1,lsfit_acc *a){
+ long i;
+
+ memset(a,0,sizeof(lsfit_acc));
+ a->x0=x0;
+ a->x1=x1;
+
+ for(i=x0;i<x1;i++){
+ int quantized=vorbis_floor1_dBquant(floor+i);
+ if(quantized){
+ a->xa += i;
+ a->ya += quantized;
+ a->x2a += (i*i);
+ a->y2a += quantized*quantized;
+ a->xya += i*quantized;
+ a->n++;
+ }
+ }
+ return(a->n);
+}
+
+/* returns < 0 on too few points to fit, >=0 (meansq error) on success */
+static int fit_line(lsfit_acc *a,int fits,int *y0,int *y1){
+ long x=0,y=0,x2=0,y2=0,xy=0,n=0,i;
+ long x0=a[0].x0;
+ long x1=a[fits-1].x1;
+
+ if(*y0>=0){ /* hint used to break degenerate cases */
+ x+= x0;
+ y+= *y0;
+ x2+= x0 * x0;
y2+= *y0 * *y0;
- xy= *y0*x0;
+ xy+= *y0 * x0;
n++;
}
- if(*y1>-900){ /* hint used to break degenerate cases */
- x+=x1;
- y+=*y1;
- x2+=(x1*x1);
+ if(*y1>=0){ /* hint used to break degenerate cases */
+ x+= x1;
+ y+= *y1;
+ x2+= x1 * x1;
y2+= *y1 * *y1;
- xy= *y1*x1;
+ xy+= *y1 * x1;
n++;
}
- for(i=x0;i<x1;i++)
- if(floor[i]>-900){
- x+=i;
- y+=floor[i];
- x2+=(i*i);
- y2+=floor[i]*floor[i];
- xy=i*floor[i];
- n++;
- }
+ for(i=0;i<fits;i++){
+ x+=a[i].xa;
+ y+=a[i].ya;
+ x2+=a[i].x2a;
+ y2+=a[i].y2a;
+ xy+=a[i].xya;
+ n+=a[i].n;
+ }
if(n<2)return(-1);
{
- float a=(float)(y*x2-xy*x)/(n*x2-x*x);
- float b=(float)(n*xy-x*y)/(n*x2-x*x);
+ /* need 64 bit multiplies, which C doesn't give portably as int */
+ double fx=x;
+ double fy=y;
+ double fx2=x2;
+ double fy2=y2;
+ double fxy=xy;
+ double a=(fy*fx2-fxy*fx)/(n*fx2-fx*fx);
+ double b=(n*fxy-fx*fy)/(n*fx2-fx*fx);
+ int s=rint((a*a*n + a*b*(2*fx) - a*(2*fy) + b*b*fx2 - b*(2*fxy) + fy2)/(n*n));
+ if(s<0)s=0;
*y0=rint(a+b*x0);
*y1=rint(a+b*x1);
+ return(s);
}
+}
+
+static int inspect_error(int x0,int x1,int y0,int y1,float *flr,
+ vorbis_info_floor1 *info){
+ int dy=y1-y0;
+ int adx=x1-x0;
+ int ady=abs(dy);
+ int base=dy/adx;
+ int sy=(dy<0?base-1:base+1);
+ int x=x0;
+ int y=y0;
+ int err=0;
+ int mse=0;
+ int val=vorbis_floor1_dBquant(flr+x);
+ int n=0;
+
+ ady-=abs(base*adx);
+
+ if(val){
+ if(y-info->maxover>val)return(1);
+ if(y+info->maxunder<val)return(1);
+ mse=(y-val);
+ mse*=mse;
+ n++;
+ }
+
+ while(++x<x1){
+ err=err+ady;
+ if(err>=adx){
+ err-=adx;
+ y+=sy;
+ }else{
+ y+=base;
+ }
+
+ val=vorbis_floor1_dBquant(flr+x);
+ if(val){
+ if(y-info->maxover>val)return(1);
+ if(y+info->maxunder<val)return(1);
+ mse+=((y-val)*(y-val));
+ n++;
+ }
+ }
+
+ if(mse/n>info->maxerr)return(1);
return(0);
}
+int post_Y(int *A,int *B,int pos){
+ if(A[pos]==-1)
+ return B[pos];
+ if(B[pos]==-1)
+ return A[pos];
+ return (A[pos]+B[pos])>>1;
+}
/* didn't need in->out seperation, modifies the flr[] vector; takes in
a dB scale floor, puts out linear */
-static int floor1_forward(vorbis_block *vb,vorbis_look_floor *i,
- float *flr,float *mdct){
- long i,j;
- vorbis_look_floor1 *look=(vorbis_look_floor1 *)i;
+static int floor1_forward(vorbis_block *vb,vorbis_look_floor *in,
+ float *flr){
+ long i,j,k,l;
+ vorbis_look_floor1 *look=(vorbis_look_floor1 *)in;
vorbis_info_floor1 *info=look->vi;
long n=look->n;
- /* first, quantize the floor into the units we're using,
- silmultaneously weeding out values for x positions where the
- residue will likely quantize to zero */
- int workfloor=alloca(sizeof(int)*n);
- memset(workfloor,0,sizeof(int)*n);
- for(i=0;i<n;i++){
- if
+ long posts=look->posts;
+ long nonzero=0;
+ lsfit_acc fits[VIF_POSIT+1];
+ int fit_valueA[VIF_POSIT+2]; /* index by range list position */
+ int fit_valueB[VIF_POSIT+2]; /* index by range list position */
+ int fit_flag[VIF_POSIT+2];
+
+ int loneighbor[VIF_POSIT+2]; /* sorted index of range list position (+2) */
+ int hineighbor[VIF_POSIT+2];
+ int memo[VIF_POSIT+2];
+ static_codebook *sbooks=vb->vd->vi->codec_setup->book_param;
+ codebook *books=vb->vd->backend_state->fullbooks;
+
+ memset(fit_flag,0,sizeof(fit_flag));
+ for(i=0;i<posts;i++)loneighbor[i]=0; /* 0 for the implicit 0 post */
+ for(i=0;i<posts;i++)hineighbor[i]=1; /* 1 for the implicit post at n */
+ for(i=0;i<posts;i++)memo[i]=-1; /* no neighbor yet */
+
+ /* quantize the relevant floor points and collect them into line fit
+ structures (one per minimal division) at the same time */
+ if(posts==0){
+ nonzero+=accumulate_fit(flr,0,n,fits);
+ }else{
+ for(i=0;i<posts-1;i++)
+ nonzero+=accumulate_fit(flr,look->sorted_index[i],
+ look->sorted_index[i+1],fits+i);
+ }
+
+ if(nonzero){
+ /* start by fitting the implicit base case.... */
+ int y0=-999;
+ int y1=-999;
+ int mse=fit_line(fits,posts-1,&y0,&y1);
+ if(mse<0){
+ /* Only a single nonzero point */
+ y0=-999;
+ y1=0;
+ fit_line(fits,posts-1,&y0,&y1);
+ }
+
+ fit_flag[0]=1;
+ fit_flag[1]=1;
+ fit_valueA[0]=y0;
+ fit_valueB[0]=y0;
+ fit_valueB[1]=y1;
+ fit_valueA[1]=y1;
+ if(mse>=0){
+ /* Non degenerate case */
+ /* start progressive splitting. This is a greedy, non-optimal
+ algorithm, but simple and close enough to the best
+ answer. */
+ for(i=2;i<posts;i++){
+ int sortpos=look->reverse_index[i];
+ int ln=loneighbor[sortpos];
+ int hn=hineighbor[sortpos];
- return(val);
+ /* eliminate repeat searches of a particular range with a memo */
+ if(memo[ln]!=hn){
+ /* haven't performed this error search yet */
+ int lsortpos=look->reverse_index[ln];
+ int hsortpos=look->reverse_index[hn];
+
+ /* if this is an empty segment, its endpoints don't matter.
+ Mark as such */
+ for(j=lsortpos;j<hsortpos;j++)
+ if(fits[j].n)break;
+ if(j==hsortpos){
+ /* empty segment; important to note that this does not
+ break 0/n post case */
+ fit_valueB[ln]=-1;
+ fit_valueA[hn]=-1;
+ }else{
+ /* A note: we want to bound/minimize *local*, not global, error */
+ int lx=info->postlist[ln];
+ int hx=info->postlist[hn];
+ int ly=post_Y(fit_valueA,fit_valueB,ln);
+ int hy=post_Y(fit_valueA,fit_valueB,hn);
+ memo[ln]=hn;
+
+ if(i<info->searchstart ||
+ inspect_error(lx,hx,ly,hy,flr,info)){
+ /* outside error bounds/begin search area. Split it. */
+ int ly0=-999;
+ int ly1=-999;
+ int hy0=-999;
+ int hy1=-999;
+ int lmse=fit_line(fits+lsortpos,sortpos-lsortpos,&ly0,&ly1);
+ int hmse=fit_line(fits+sortpos,hsortpos-sortpos,&hy0,&hy1);
+
+ /* Handle degeneracy */
+ if(lmse<0 && hmse<0){
+ ly0=fit_valueA[ln];
+ hy1=fit_valueB[hn];
+ lmse=fit_line(fits+lsortpos,sortpos-lsortpos,&ly0,&ly1);
+ hmse=fit_line(fits+sortpos,hsortpos-sortpos,&hy0,&hy1);
+ }
+
+ if(lmse<0 && hmse<0) continue;
+
+ if(lmse<0){
+ ly0=fit_valueA[ln];
+ ly1=hy0;
+ lmse=fit_line(fits+lsortpos,sortpos-lsortpos,&ly0,&ly1);
+ }
+ if(hmse<0){
+ hy1=fit_valueB[hn];
+ hy0=ly1;
+ hmse=fit_line(fits+sortpos,hsortpos-sortpos,&hy0,&hy1);
+ }
+
+ /* store new edge values */
+ fit_valueB[ln]=ly0;
+ if(ln==0)fit_valueA[ln]=ly0;
+ fit_valueA[i]=ly1;
+ fit_valueB[i]=hy0;
+ fit_valueA[hn]=hy1;
+ if(hn==1)fit_valueB[hn]=hy1;
+
+ /* store new neighbor values */
+ for(j=sortpos-1;j>=0;j--)
+ if(hineighbor[j]==hn)
+ hineighbor[j]=i;
+ else
+ break;
+ for(j=sortpos+1;j<posts;j++)
+ if(loneighbor[j]==ln)
+ loneighbor[j]=i;
+ else
+ break;
+
+ /* store flag (set) */
+ fit_flag[i]=1;
+
+ }
+ }
+ }
+ }
+ }
+
+ /* quantize values to multiplier spec */
+ switch(info->mult){
+ case 1: /* 1024 -> 256 */
+ for(i=0;i<posts;i++)
+ if(fit_flag[i])
+ fit_valueA[i]=(post_Y(fit_valueA,fit_valueB,i)+2)>>2;
+ break;
+ case 2: /* 1024 -> 128 */
+ for(i=0;i<posts;i++)
+ if(fit_flag[i])
+ fit_valueA[i]=(post_Y(fit_valueA,fit_valueB,i)+4)>>3;
+ break;
+ case 3: /* 1024 -> 86 */
+ for(i=0;i<posts;i++)
+ if(fit_flag[i])
+ fit_valueA[i]=(post_Y(fit_valueA,fit_valueB,i)+6)/12;
+ break;
+ case 4: /* 1024 -> 64 */
+ for(i=0;i<posts;i++)
+ if(fit_flag[i])
+ fit_valueA[i]=(post_Y(fit_valueA,fit_valueB,i)+8)>>4;
+ break;
+ }
+
+ /* find prediction values for each post and subtract them */
+ /* work backward to avoid a copy; unwind the neighbor arrays */
+ for(i=posts-1;i>1;i--){
+ int sp=look->reverse_index[i];
+ int ln=loneighbor[i];
+ int hn=hineighbor[i];
+
+ int x0=info->postlist[ln];
+ int x1=info->postlist[hn];
+ int y0=fit_valueA[ln];
+ int y1=fit_valueA[hn];
+
+ int predicted=render_point(x0,x1,y0,y1,info->postlist[i]);
+
+ if(fit_flag[i]){
+ int headroom=(look->quant_q-predicted<predicted?
+ look->quant_q-predicted:predicted);
+
+ int val=fit_valueA[i]-predicted;
+
+ /* at this point the 'deviation' value is in the range +/- max
+ range, but the real, unique range can be mapped to only
+ [0-maxrange). So we want to wrap the deviation into this
+ limited range, but do it in the way that least screws an
+ essentially gaussian probability distribution. */
+
+ if(val<0)
+ if(val<-headroom)
+ val=headroom-val-1;
+ else
+ val=1-(val<<1);
+ else
+ if(val>=headroom)
+ val= val+headroom;
+ else
+ val>>=1;
+
+ fit_valueB[i]=val;
+
+ /* unroll the neighbor arrays */
+ for(j=sp+1;j<posts;j++)
+ if(loneighbor[j]==i)
+ loneighbor[j]=loneighbor[sp];
+ else
+ break;
+ for(j=sp-1;j>=0;j--)
+ if(hineighbor[j]==i)
+ hineighbor[j]=hineighbor[sp];
+ else
+ break;
+
+ }else{
+ fit_valueA[i]=predicted;
+ fit_valueB[i]=0;
+ }
+ }
+
+ /* we have everything we need. pack it out */
+ /* mark nontrivial floor */
+ oggpack_write(&vb->opb,1,1);
+
+ /* beginning/end post */
+ oggpack_write(&vb->opb,fit_valueA[0],ilog(look->quant_q-1));
+ oggpack_write(&vb->opb,fit_valueA[1],ilog(look->quant_q-1));
+
+ /* partition by partition */
+ for(i=0,j=2;i<info->partitions;i++){
+ int class=info->partitionclass[i];
+ int cdim=info->class_dim[class];
+ int csubbits=info->class_subs[class];
+ int csub=1<<csubbits;
+ int bookas[8]={0,0,0,0,0,0,0,0};
+ int cval=0;
+ int cshift=0;
+
+ /* generate the partition's first stage cascade value */
+ if(csubbits){
+ int maxval[8];
+ for(k=0;k<csub;k++){
+ int booknum=info->class_subbook[class][k];
+ if(booknum<0){
+ maxval[k]=0;
+ }else{
+ maxval[k]=sbooks[info->class_subbook[class][k]].entries;
+ }
+ }
+ for(k=0;k<cdim;k++){
+ for(l=0;l<csub;l++){
+ val=fit_valueB[j+k];
+ if(val<maxval[k]){
+ bookas[k]=l;
+ break;
+ }
+ }
+ cval|= bookas[k]<<cshift;
+ cshift+=csubbits;
+ }
+ /* write it */
+ vorbis_book_encode(books+info->class_book[class],cval,&vb->opb);
+
+#ifdef TRAIN_FLOOR1
+ {
+ FILE *of;
+ char buffer[80];
+ sprintf(buffer,"floor1_class%d.vqd",class);
+ of=fopen(buffer,"a");
+ fprintf(of,"%d\n",cval);
+ fclose(of);
+ }
+#endif
+ }
+
+ /* write post values */
+ for(k=0;k<cdim;k++){
+ int book=info->class_subbook[class][bookas[k]];
+ if(book>=0){
+ vorbis_book_encode(books+book,
+ fit_valueB[j+k],&vb->opb);
+#ifdef TRAIN_FLOOR1
+ {
+ FILE *of;
+ char buffer[80];
+ sprintf(buffer,"floor1_class%dsub%d.vqd",class,bookas[k]);
+ of=fopen(buffer,"a");
+ fprintf(of,"%d\n",fit_valueB[j+k]);
+ fclose(of);
+ }
+#endif
+ }
+ }
+ j+=cdim;
+ }
+
+ /* generate quantized floor equivalent to what we'd unpack in decode */
+ {
+ int hx;
+ int lx=0;
+ int ly=fit_valueA[0]*info->mult;
+ for(j=1;j<posts;j++){
+ int current=look->forward_index[j];
+ int hy=fit_valueA[current]*info->mult;
+ hx=info->postlist[current];
+
+ render_line(lx,hx,ly,hy,out);
+
+ lx=hx;
+ ly=hy;
+ }
+ for(j=hx;j<n;j++)out[j]=0.f; /* be certain */
+ }
+
+ }else{
+ oggpack_write(&vb->opb,0,1);
+ memset(flr,0,n*sizeof(float));
+ }
+ return(nonzero);
}
static int floor1_inverse(vorbis_block *vb,vorbis_look_floor *i,float *out){
vorbis_look_floor1 *look=(vorbis_look_floor1 *)i;
vorbis_info_floor1 *info=look->vi;
+ codec_setup_info *ci=vb->vd->vi->codec_setup;
+ int n=ci->blocksizes[vb->mode]/2;
+ int i,j;
+ codebook *books=vb->vd->backend_state->fullbooks;
+
+ /* unpack wrapped/predicted values from stream */
+ if(oggpack_read(&vb->opb,1)==1){
+ int fit_value[VIF_POSIT+2];
+
+ fit_value[0]=oggpack_read(&vb->opb,ilog(look->quant_q-1));
+ fit_value[1]=oggpack_read(&vb->opb,ilog(look->quant_q-1));
+
+ /* partition by partition */
+ /* partition by partition */
+ for(i=0,j=2;i<info->partitions;i++){
+ int class=info->partitionclass[i];
+ int cdim=info->class_dim[class];
+ int csubbits=info->class_subs[class];
+ int csub=1<<csubbits;
+ int bookas[8]={0,0,0,0,0,0,0,0};
+
+ /* decode the partition's first stage cascade value */
+ if(csubbits){
+ int cval=vorbis_book_decode(books+info->class_book[class],&vb->opb);
+ if(cval==-1)goto eop;
+
+ for(k=0;k<cdim;k++){
+ int book=info->class_subbook[class][val&(csub-1)];
+ cval>>=csubbits;
+ if(book>=0){
+ if((fit_value[j+k]=vorbis_book_decode(books+book,&vb->opb))==-1)
+ goto eop;
+ }else{
+ fit_value[j+k]=0;
+ }
+ }
+ }else{
+ for(k=0;k<cdim;k++)
+ fit_value[j+k]=0;
+ }
+ j+=cdim;
+ }
+ /* unwrap positive values and reconsitute via linear interpolation */
+ for(i=2;i<posts;i++){
+ int predicted=render_point(info->postlist[look->loneighbor[i-2]],
+ info->postlist[look->hineighbor[i-2]],
+ fit_value[look->loneighbor[i-2]],
+ fit_value[look->hineighbor[i-2]],
+ info->postlist[i]);
+ int hiroom=look->quant_q-predicted;
+ int loroom=predicted;
+ int room=(hiroom<loroom?hiroom:loroom)<<1;
+ int val=fit_value[i];
+
+ if(val>=room){
+ if(hiroom>loroom){
+ val = val-loroom;
+ }else{
+ val = -1-(val-hiroom);
+ }
+ }else{
+ if(val&1){
+ val= -1 -(val>>1);
+ }else{
+ val>>=1;
+ }
+ }
+
+ fit_value[i]=val+predicted;
+ }
+
+ /* render the lines */
+ {
+ int hx;
+ int lx=0;
+ int ly=fit_value[0]*info->mult;
+ for(j=1;j<posts;j++){
+ int current=look->forward_index[j];
+ int hy=fit_value[current]*info->mult;
+ hx=info->postlist[current];
+
+ render_line(lx,hx,ly,hy,out);
+
+ lx=hx;
+ ly=hy;
+ }
+ for(j=hx;j<n;j++)out[j]=0.f; /* be certain */
+ }
+ return(1);
+ }
+
+ /* fall through */
eop:
- memset(out,0,sizeof(float)*look->n);
+ memset(out,0,sizeof(float)*n);
return(0);
}
@@ -285,4 +1025,3 @@ vorbis_func_floor floor1_exportbundle={
&floor1_free_look,&floor1_forward,&floor1_inverse
};
-
diff --git a/lib/mapping0.c b/lib/mapping0.c
new file mode 100644
index 00000000..0e02d9f1
--- /dev/null
+++ b/lib/mapping0.c
@@ -0,0 +1,401 @@
+/********************************************************************
+ * *
+ * THIS FILE IS PART OF THE OggVorbis SOFTWARE CODEC SOURCE CODE. *
+ * USE, DISTRIBUTION AND REPRODUCTION OF THIS LIBRARY SOURCE IS *
+ * GOVERNED BY A BSD-STYLE SOURCE LICENSE INCLUDED WITH THIS SOURCE *
+ * IN 'COPYING'. PLEASE READ THESE TERMS BEFORE DISTRIBUTING. *
+ * *
+ * THE OggVorbis SOURCE CODE IS (C) COPYRIGHT 1994-2001 *
+ * by the XIPHOPHORUS Company http://www.xiph.org/ *
+
+ ********************************************************************
+
+ function: channel mapping 0 implementation
+ last mod: $Id: mapping0.c,v 1.27.4.1 2001/04/29 22:21:04 xiphmont Exp $
+
+ ********************************************************************/
+
+#include <stdlib.h>
+#include <stdio.h>
+#include <string.h>
+#include <math.h>
+#include <ogg/ogg.h>
+#include "vorbis/codec.h"
+#include "codec_internal.h"
+#include "codebook.h"
+#include "bitbuffer.h"
+#include "registry.h"
+#include "psy.h"
+#include "misc.h"
+
+/* simplistic, wasteful way of doing this (unique lookup for each
+ mode/submapping); there should be a central repository for
+ identical lookups. That will require minor work, so I'm putting it
+ off as low priority.
+
+ Why a lookup for each backend in a given mode? Because the
+ blocksize is set by the mode, and low backend lookups may require
+ parameters from other areas of the mode/mapping */
+
+typedef struct {
+ drft_lookup fft_look;
+ vorbis_info_mode *mode;
+ vorbis_info_mapping0 *map;
+
+ vorbis_look_time **time_look;
+ vorbis_look_floor **floor_look;
+
+ vorbis_look_residue **residue_look;
+ vorbis_look_psy *psy_look;
+
+ vorbis_func_time **time_func;
+ vorbis_func_floor **floor_func;
+ vorbis_func_residue **residue_func;
+
+ int ch;
+ long lastframe; /* if a different mode is called, we need to
+ invalidate decay */
+} vorbis_look_mapping0;
+
+static vorbis_info_mapping *mapping0_copy_info(vorbis_info_mapping *vm){
+ vorbis_info_mapping0 *info=(vorbis_info_mapping0 *)vm;
+ vorbis_info_mapping0 *ret=_ogg_malloc(sizeof(vorbis_info_mapping0));
+ memcpy(ret,info,sizeof(vorbis_info_mapping0));
+ return(ret);
+}
+
+static void mapping0_free_info(vorbis_info_mapping *i){
+ if(i){
+ memset(i,0,sizeof(vorbis_info_mapping0));
+ _ogg_free(i);
+ }
+}
+
+static void mapping0_free_look(vorbis_look_mapping *look){
+ int i;
+ vorbis_look_mapping0 *l=(vorbis_look_mapping0 *)look;
+ if(l){
+ drft_clear(&l->fft_look);
+
+ for(i=0;i<l->map->submaps;i++){
+ l->time_func[i]->free_look(l->time_look[i]);
+ l->floor_func[i]->free_look(l->floor_look[i]);
+ l->residue_func[i]->free_look(l->residue_look[i]);
+ if(l->psy_look)_vp_psy_clear(l->psy_look+i);
+ }
+
+ _ogg_free(l->time_func);
+ _ogg_free(l->floor_func);
+ _ogg_free(l->residue_func);
+ _ogg_free(l->time_look);
+ _ogg_free(l->floor_look);
+ _ogg_free(l->residue_look);
+ if(l->psy_look)_ogg_free(l->psy_look);
+ memset(l,0,sizeof(vorbis_look_mapping0));
+ _ogg_free(l);
+ }
+}
+
+static vorbis_look_mapping *mapping0_look(vorbis_dsp_state *vd,vorbis_info_mode *vm,
+ vorbis_info_mapping *m){
+ int i;
+ vorbis_info *vi=vd->vi;
+ codec_setup_info *ci=vi->codec_setup;
+ vorbis_look_mapping0 *look=_ogg_calloc(1,sizeof(vorbis_look_mapping0));
+ vorbis_info_mapping0 *info=look->map=(vorbis_info_mapping0 *)m;
+ look->mode=vm;
+
+ look->time_look=_ogg_calloc(info->submaps,sizeof(vorbis_look_time *));
+ look->floor_look=_ogg_calloc(info->submaps,sizeof(vorbis_look_floor *));
+
+ look->residue_look=_ogg_calloc(info->submaps,sizeof(vorbis_look_residue *));
+ if(ci->psys)look->psy_look=_ogg_calloc(info->submaps,sizeof(vorbis_look_psy));
+
+ look->time_func=_ogg_calloc(info->submaps,sizeof(vorbis_func_time *));
+ look->floor_func=_ogg_calloc(info->submaps,sizeof(vorbis_func_floor *));
+ look->residue_func=_ogg_calloc(info->submaps,sizeof(vorbis_func_residue *));
+
+ for(i=0;i<info->submaps;i++){
+ int timenum=info->timesubmap[i];
+ int floornum=info->floorsubmap[i];
+ int resnum=info->residuesubmap[i];
+
+ look->time_func[i]=_time_P[ci->time_type[timenum]];
+ look->time_look[i]=look->time_func[i]->
+ look(vd,vm,ci->time_param[timenum]);
+ look->floor_func[i]=_floor_P[ci->floor_type[floornum]];
+ look->floor_look[i]=look->floor_func[i]->
+ look(vd,vm,ci->floor_param[floornum]);
+ look->residue_func[i]=_residue_P[ci->residue_type[resnum]];
+ look->residue_look[i]=look->residue_func[i]->
+ look(vd,vm,ci->residue_param[resnum]);
+
+ if(ci->psys && vd->analysisp){
+ int psynum=info->psysubmap[i];
+ _vp_psy_init(look->psy_look+i,ci->psy_param[psynum],
+ ci->blocksizes[vm->blockflag]/2,vi->rate);
+ }
+ }
+
+ look->ch=vi->channels;
+
+ if(vd->analysisp)drft_init(&look->fft_look,ci->blocksizes[vm->blockflag]);
+ return(look);
+}
+
+static void mapping0_pack(vorbis_info *vi,vorbis_info_mapping *vm,oggpack_buffer *opb){
+ int i;
+ vorbis_info_mapping0 *info=(vorbis_info_mapping0 *)vm;
+
+ oggpack_write(opb,info->submaps-1,4);
+ /* we don't write the channel submappings if we only have one... */
+ if(info->submaps>1){
+ for(i=0;i<vi->channels;i++)
+ oggpack_write(opb,info->chmuxlist[i],4);
+ }
+ for(i=0;i<info->submaps;i++){
+ oggpack_write(opb,info->timesubmap[i],8);
+ oggpack_write(opb,info->floorsubmap[i],8);
+ oggpack_write(opb,info->residuesubmap[i],8);
+ }
+}
+
+/* also responsible for range checking */
+static vorbis_info_mapping *mapping0_unpack(vorbis_info *vi,oggpack_buffer *opb){
+ int i;
+ vorbis_info_mapping0 *info=_ogg_calloc(1,sizeof(vorbis_info_mapping0));
+ codec_setup_info *ci=vi->codec_setup;
+ memset(info,0,sizeof(vorbis_info_mapping0));
+
+ info->submaps=oggpack_read(opb,4)+1;
+
+ if(info->submaps>1){
+ for(i=0;i<vi->channels;i++){
+ info->chmuxlist[i]=oggpack_read(opb,4);
+ if(info->chmuxlist[i]>=info->submaps)goto err_out;
+ }
+ }
+ for(i=0;i<info->submaps;i++){
+ info->timesubmap[i]=oggpack_read(opb,8);
+ if(info->timesubmap[i]>=ci->times)goto err_out;
+ info->floorsubmap[i]=oggpack_read(opb,8);
+ if(info->floorsubmap[i]>=ci->floors)goto err_out;
+ info->residuesubmap[i]=oggpack_read(opb,8);
+ if(info->residuesubmap[i]>=ci->residues)goto err_out;
+ }
+
+ return info;
+
+ err_out:
+ mapping0_free_info(info);
+ return(NULL);
+}
+
+#include "os.h"
+#include "lpc.h"
+#include "lsp.h"
+#include "envelope.h"
+#include "mdct.h"
+#include "psy.h"
+#include "scales.h"
+
+/* no time mapping implementation for now */
+static long seq=0;
+static int mapping0_forward(vorbis_block *vb,vorbis_look_mapping *l){
+ vorbis_dsp_state *vd=vb->vd;
+ vorbis_info *vi=vd->vi;
+ backend_lookup_state *b=vb->vd->backend_state;
+ vorbis_look_mapping0 *look=(vorbis_look_mapping0 *)l;
+ vorbis_info_mapping0 *info=look->map;
+ vorbis_info_mode *mode=look->mode;
+ vorbis_block_internal *vbi=(vorbis_block_internal *)vb->internal;
+ int n=vb->pcmend;
+ int i,j;
+ float *window=b->window[vb->W][vb->lW][vb->nW][mode->windowtype];
+
+ float **pcmbundle=alloca(sizeof(float *)*vi->channels);
+
+ int *nonzero=alloca(sizeof(int)*vi->channels);
+
+ float **floor=_vorbis_block_alloc(vb,vi->channels*sizeof(float *));
+ float *additional=_vorbis_block_alloc(vb,n*sizeof(float));
+ float newmax=vbi->ampmax;
+
+ for(i=0;i<vi->channels;i++){
+ float *pcm=vb->pcm[i];
+ float scale=4.f/n;
+ int submap=info->chmuxlist[i];
+ float ret;
+
+ _analysis_output("pcm",seq,pcm,n,0,0);
+
+ /* window the PCM data */
+ for(j=0;j<n;j++)
+ additional[j]=pcm[j]*=window[j];
+
+ _analysis_output("windowed",seq,pcm,n,0,0);
+
+ /* transform the PCM data */
+ /* only MDCT right now.... */
+ mdct_forward(b->transform[vb->W][0],pcm,pcm);
+
+ /* FFT yields more accurate tonal estimation (not phase sensitive) */
+ drft_forward(&look->fft_look,additional);
+
+ additional[0]=fabs(additional[0]*scale);
+ for(j=1;j<n-1;j+=2)
+ additional[(j+1)>>1]=scale*FAST_HYPOT(additional[j],additional[j+1]);
+
+ /* set up our masking data working vector, and stash unquantized
+ data for later */
+ /*memcpy(pcm+n/2,pcm,n*sizeof(float)/2);*/
+ memcpy(additional+n/2,pcm,n*sizeof(float)/2);
+
+ /* begin masking work */
+ floor[i]=_vorbis_block_alloc(vb,n*sizeof(float)/2);
+
+ _analysis_output("fft",seq,additional,n/2,0,1);
+ _analysis_output("mdct",seq,additional+n/2,n/2,0,1);
+ _analysis_output("lfft",seq,additional,n/2,0,0);
+ _analysis_output("lmdct",seq,additional+n/2,n/2,0,0);
+
+ /* perform psychoacoustics; do masking */
+ ret=_vp_compute_mask(look->psy_look+submap,additional,additional+n/2,
+ floor[i],NULL,vbi->ampmax);
+ if(ret>newmax)newmax=ret;
+
+ _analysis_output_always("mask",seq,floor[i],n/2,0,0);
+
+ /* perform floor encoding */
+ nonzero[i]=look->floor_func[submap]->
+ forward(vb,look->floor_look[submap],floor[i]);
+
+ _analysis_output("floor",seq,floor[i],n/2,0,1);
+
+ /* apply the floor, do optional noise levelling */
+ _vp_apply_floor(look->psy_look+submap,pcm,floor[i]);
+
+ _analysis_output("res",seq++,pcm,n/2,0,0);
+
+#ifdef TRAIN_RES
+ if(nonzero[i]){
+ FILE *of;
+ char buffer[80];
+ int i;
+
+ sprintf(buffer,"residue_%d.vqd",vb->mode);
+ of=fopen(buffer,"a");
+ for(i=0;i<n/2;i++){
+ fprintf(of,"%.2f, ",pcm[i]);
+ if(fabs(pcm[i])>1000){
+ fprintf(stderr," %d ",seq-1);
+ }
+ }
+ fprintf(of,"\n");
+ fclose(of);
+ }
+#endif
+ }
+
+ seq++;
+ vbi->ampmax=newmax;
+
+ /* perform residue encoding with residue mapping; this is
+ multiplexed. All the channels belonging to one submap are
+ encoded (values interleaved), then the next submap, etc */
+
+ for(i=0;i<info->submaps;i++){
+ int ch_in_bundle=0;
+ for(j=0;j<vi->channels;j++){
+ if(info->chmuxlist[j]==i && nonzero[j]>0){
+ pcmbundle[ch_in_bundle++]=vb->pcm[j];
+ }
+ }
+
+ look->residue_func[i]->forward(vb,look->residue_look[i],
+ pcmbundle,ch_in_bundle);
+ }
+
+ look->lastframe=vb->sequence;
+ return(0);
+}
+
+static int mapping0_inverse(vorbis_block *vb,vorbis_look_mapping *l){
+ vorbis_dsp_state *vd=vb->vd;
+ vorbis_info *vi=vd->vi;
+ codec_setup_info *ci=vi->codec_setup;
+ backend_lookup_state *b=vd->backend_state;
+ vorbis_look_mapping0 *look=(vorbis_look_mapping0 *)l;
+ vorbis_info_mapping0 *info=look->map;
+ vorbis_info_mode *mode=look->mode;
+ int i,j;
+ long n=vb->pcmend=ci->blocksizes[vb->W];
+
+ float *window=b->window[vb->W][vb->lW][vb->nW][mode->windowtype];
+ float **pcmbundle=alloca(sizeof(float *)*vi->channels);
+ int *nonzero=alloca(sizeof(int)*vi->channels);
+
+ /* time domain information decode (note that applying the
+ information would have to happen later; we'll probably add a
+ function entry to the harness for that later */
+ /* NOT IMPLEMENTED */
+
+ /* recover the spectral envelope; store it in the PCM vector for now */
+ for(i=0;i<vi->channels;i++){
+ float *pcm=vb->pcm[i];
+ int submap=info->chmuxlist[i];
+ nonzero[i]=look->floor_func[submap]->
+ inverse(vb,look->floor_look[submap],pcm);
+ _analysis_output("ifloor",seq+i,pcm,n/2,0,1);
+ }
+
+ /* recover the residue, apply directly to the spectral envelope */
+
+ for(i=0;i<info->submaps;i++){
+ int ch_in_bundle=0;
+ for(j=0;j<vi->channels;j++){
+ if(info->chmuxlist[j]==i && nonzero[j])
+ pcmbundle[ch_in_bundle++]=vb->pcm[j];
+ }
+
+ look->residue_func[i]->inverse(vb,look->residue_look[i],pcmbundle,ch_in_bundle);
+ }
+
+ /* transform the PCM data; takes PCM vector, vb; modifies PCM vector */
+ /* only MDCT right now.... */
+ for(i=0;i<vi->channels;i++){
+ float *pcm=vb->pcm[i];
+ _analysis_output("out",seq+i,pcm,n/2,0,1);
+ mdct_backward(b->transform[vb->W][0],pcm,pcm);
+ }
+
+ /* now apply the decoded pre-window time information */
+ /* NOT IMPLEMENTED */
+
+ /* window the data */
+ for(i=0;i<vi->channels;i++){
+ float *pcm=vb->pcm[i];
+ if(nonzero[i])
+ for(j=0;j<n;j++)
+ pcm[j]*=window[j];
+ else
+ for(j=0;j<n;j++)
+ pcm[j]=0.f;
+ _analysis_output("final",seq++,pcm,n,0,0);
+ }
+
+ /* now apply the decoded post-window time information */
+ /* NOT IMPLEMENTED */
+
+ /* all done! */
+ return(0);
+}
+
+/* export hooks */
+vorbis_func_mapping mapping0_exportbundle={
+ &mapping0_pack,&mapping0_unpack,&mapping0_look,&mapping0_copy_info,
+ &mapping0_free_info,&mapping0_free_look,&mapping0_forward,&mapping0_inverse
+};
+
+
+
diff --git a/lib/modes/mode_A.h b/lib/modes/mode_A.h
new file mode 100644
index 00000000..09694c04
--- /dev/null
+++ b/lib/modes/mode_A.h
@@ -0,0 +1,316 @@
+/********************************************************************
+ * *
+ * THIS FILE IS PART OF THE OggVorbis SOFTWARE CODEC SOURCE CODE. *
+ * USE, DISTRIBUTION AND REPRODUCTION OF THIS LIBRARY SOURCE IS *
+ * GOVERNED BY A BSD-STYLE SOURCE LICENSE INCLUDED WITH THIS SOURCE *
+ * IN 'COPYING'. PLEASE READ THESE TERMS BEFORE DISTRIBUTING. *
+ * *
+ * THE OggVorbis SOURCE CODE IS (C) COPYRIGHT 1994-2001 *
+ * by the XIPHOPHORUS Company http://www.xiph.org/ *
+
+ ********************************************************************
+
+ function: predefined encoding modes
+ last mod: $Id: mode_A.h,v 1.14.4.1 2001/04/29 22:21:06 xiphmont Exp $
+
+ ********************************************************************/
+
+#ifndef _V_MODES_A_H_
+#define _V_MODES_A_H_
+
+#include <stdio.h>
+#include "vorbis/codec.h"
+#include "backends.h"
+
+#include "books/lsp12_0.vqh"
+#include "books/lsp30_0.vqh"
+#include "books/lsp12_1.vqh"
+#include "books/lsp30_1.vqh"
+
+#include "books/res0_128_128aux.vqh"
+#include "books/res0_128_1024aux.vqh"
+
+#include "books/res0_128_128_1.vqh"
+#include "books/res0_128_128_2.vqh"
+#include "books/res0_128_128_3.vqh"
+#include "books/res0_128_128_4.vqh"
+#include "books/res0_128_128_5.vqh"
+
+#include "books/res0_128_1024_1.vqh"
+#include "books/res0_128_1024_2.vqh"
+#include "books/res0_128_1024_3.vqh"
+#include "books/res0_128_1024_4.vqh"
+#include "books/res0_128_1024_5.vqh"
+#include "books/res0_128_1024_6.vqh"
+#include "books/res0_128_1024_7.vqh"
+#include "books/res0_128_1024_8.vqh"
+#include "books/res0_128_1024_9.vqh"
+
+
+static vorbis_info_psy _psy_set_A0={
+ 1,/*athp*/
+ 1,/*decayp*/
+
+ -100.,
+ -140.,
+
+ 8,
+
+ /* 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 */
+ /* x: 63 88 125 175 250 350 500 700 1k 1.4k 2k 2.8k 4k 5.6k 8k 11.5k 16k Hz */
+ /* y: 0 10 20 30 40 50 60 70 80 90 100 dB */
+
+ 1,/* tonemaskp */
+ /* 0 10 20 30 40 50 60 70 80 90 100 */
+ {
+ {-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.}, /*63*/
+ {-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.}, /*88*/
+ {-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.}, /*125*/
+ {-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.}, /*175*/
+ {-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.}, /*250*/
+ {-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.}, /*350*/
+ {-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.}, /*500*/
+ {-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.,-999.}, /*700*/
+
+ {-30.,-35.,-35.,-40.,-40.,-50.,-60.,-70.,-80.,-90.,-100.}, /*1000*/
+ {-30.,-35.,-35.,-40.,-40.,-50.,-60.,-70.,-80.,-90.,-100.}, /*1400*/
+ {-40.,-40.,-40.,-40.,-40.,-50.,-60.,-70.,-80.,-90.,-100.}, /*2000*/
+ {-40.,-40.,-40.,-40.,-40.,-50.,-60.,-70.,-80.,-90.,-100.}, /*2800*/
+ {-30.,-35.,-35.,-40.,-40.,-50.,-60.,-70.,-80.,-90.,-100.}, /*4000*/
+ {-30.,-35.,-35.,-40.,-40.,-50.,-60.,-70.,-80.,-90.,-100.}, /*5600*/
+ {-30.,-35.,-35.,-40.,-40.,-50.,-60.,-70.,-80.,-90.,-100.}, /*8000*/
+ {-30.,-35.,-35.,-40.,-40.,-50.,-60.,-70.,-80.,-90.,-100.}, /*11500*/
+ {-30.,-35.,-35.,-40.,-40.,-50.,-60.,-70.,-80.,-90.,-100.}, /*16000*/
+ },
+
+ 1,/* peakattp */
+ {{-14.f,-16.f,-18.f,-19.f,-20.f,-21.f,-22.f,-22.f,-22.f,-22.f,-22.f},/*63*/
+ {-14.f,-16.f,-18.f,-19.f,-20.f,-21.f,-22.f,-22.f,-22.f,-22.f,-22.f},/*88*/
+ {-14.f,-16.f,-18.f,-19.f,-20.f,-21.f,-22.f,-22.f,-22.f,-22.f,-22.f},/*125*/
+ {-10.f,-10.f,-10.f,-10.f,-16.f,-16.f,-18.f,-20.f,-22.f,-22.f,-22.f},/*175*/
+ {-10.f,-10.f,-10.f,-10.f,-16.f,-16.f,-18.f,-20.f,-22.f,-22.f,-22.f},/*250*/
+ {-10.f,-10.f,-10.f,-10.f,-16.f,-16.f,-18.f,-20.f,-22.f,-22.f,-22.f},/*350*/
+ {-10.f,-10.f,-10.f,-10.f,-16.f,-16.f,-18.f,-20.f,-22.f,-22.f,-22.f},/*500*/
+ {-10.f,-10.f,-10.f,-10.f,-14.f,-14.f,-14.f,-20.f,-22.f,-22.f,-22.f},/*700*/
+ {-10.f,-10.f,-10.f,-10.f,-14.f,-14.f,-14.f,-20.f,-22.f,-22.f,-22.f},/*1000*/
+ {-10.f,-10.f,-10.f,-10.f,-14.f,-14.f,-16.f,-20.f,-22.f,-22.f,-22.f},/*1400*/
+ {-10.f,-10.f,-10.f,-10.f,-14.f,-14.f,-16.f,-20.f,-22.f,-22.f,-22.f},/*2000*/
+ {-10.f,-10.f,-10.f,-10.f,-14.f,-14.f,-16.f,-20.f,-22.f,-22.f,-22.f},/*2400*/
+ {-10.f,-10.f,-10.f,-10.f,-14.f,-14.f,-16.f,-20.f,-22.f,-22.f,-22.f},/*4000*/
+ {-10.f,-10.f,-10.f,-10.f,-10.f,-11.f,-12.f,-13.f,-22.f,-22.f,-22.f},/*5600*/
+ {-10.f,-10.f,-10.f,-10.f,-10.f,-11.f,-12.f,-13.f,-22.f,-22.f,-22.f},/*8000*/
+ {-10.f,-10.f,-10.f,-10.f,-10.f,-10.f,-10.f,-11.f,-22.f,-22.f,-22.f},/*11500*/
+ {-10.f,-10.f,-10.f,-10.f,-10.f,-10.f,-10.f,-10.f,-20.f,-21.f,-22.f},/*16000*/
+ },
+
+ 1,/*noisemaskp */
+ -10.f, /* suppress any noise curve over maxspec+n */
+ .5f, /* low window */
+ .5f, /* high window */
+ 5,
+ 5,
+ {.000f, 0.f,/*63*/
+ .000f, 0.f,/*88*/
+ .000f, 0.f,/*125*/
+ .000f, 0.f,/*175*/
+ .000f, 0.f,/*250*/
+ .000f, 0.f,/*350*/
+ .000f, 0.f,/*500*/
+ .300f, 0.f,/*700*/
+ .500f, 0.f,/*1000*/
+ .500f, 0.f,/*1400*/
+ .500f, 0.f,/*2000*/
+ .500f, 0.f,/*2800*/
+ .600f, 0.f,/*4000*/
+ .700f, 0.f,/*5600*/
+ .850f, 0.f,/*8000*/
+ .900f, 0.f,/*11500*/
+ .900f, 0.f,/*16000*/
+ },
+
+ 95.f, /* even decade + 5 is important; saves an rint() later in a
+ tight loop) */
+ -22.,
+
+ 1,
+};
+
+static vorbis_info_psy _psy_set_A={
+ 1,/*athp*/
+ 1,/*decayp*/
+
+ -100.f,
+ -140.f,
+
+ 8,
+
+ /* 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 */
+ /* x: 63 88 125 175 250 350 500 700 1k 1.4k 2k 2.8k 4k 5.6k 8k 11.5k 16k Hz */
+ /* y: 0 10 20 30 40 50 60 70 80 90 100 dB */
+ 1,/* tonemaskp */
+ /* 0 10 20 30 40 50 60 70 80 90 100 */
+ {
+ {-999.f,-999.f,-999.f,-999.f,-999.f,-999.f,-999.f,-999.f,-999.f,-999.f,-999.f}, /*63*/
+ {-999.f,-999.f,-999.f,-999.f,-999.f,-999.f,-999.f,-999.f,-999.f,-999.f,-999.f}, /*88*/
+ {-999.f,-999.f,-999.f,-999.f,-999.f,-999.f,-999.f,-999.f,-999.f,-999.f,-999.f}, /*125*/
+
+ {-35.f,-35.f,-35.f,-40.f,-40.f,-50.f,-60.f,-70.f,-80.f,-90.f,-100.f}, /*175*/
+ {-35.f,-35.f,-35.f,-40.f,-40.f,-50.f,-60.f,-70.f,-80.f,-90.f,-100.f}, /*250*/
+ {-35.f,-35.f,-35.f,-40.f,-40.f,-50.f,-60.f,-70.f,-80.f,-90.f,-100.f}, /*350*/
+ {-35.f,-35.f,-35.f,-40.f,-40.f,-50.f,-60.f,-70.f,-80.f,-90.f,-100.f}, /*500*/
+ {-35.f,-35.f,-35.f,-40.f,-40.f,-50.f,-60.f,-70.f,-80.f,-90.f,-100.f}, /*700*/
+
+ {-35.f,-35.f,-35.f,-40.f,-40.f,-50.f,-60.f,-70.f,-80.f,-90.f,-100.f}, /*1000*/
+ {-35.f,-35.f,-35.f,-40.f,-40.f,-50.f,-60.f,-70.f,-80.f,-90.f,-100.f}, /*1400*/
+ {-40.f,-40.f,-40.f,-40.f,-40.f,-50.f,-60.f,-70.f,-80.f,-90.f,-100.f}, /*2000*/
+ {-40.f,-40.f,-40.f,-40.f,-40.f,-50.f,-60.f,-70.f,-80.f,-90.f,-100.f}, /*2800*/
+ {-35.f,-35.f,-35.f,-40.f,-40.f,-50.f,-60.f,-70.f,-80.f,-90.f,-100.f}, /*4000*/
+ {-35.f,-35.f,-35.f,-40.f,-40.f,-50.f,-60.f,-70.f,-80.f,-90.f,-100.f}, /*5600*/
+
+ {-30.f,-30.f,-33.f,-35.f,-40.f,-50.f,-60.f,-70.f,-80.f,-90.f,-100.f}, /*8000*/
+ {-30.f,-30.f,-33.f,-35.f,-35.f,-45.f,-50.f,-60.f,-70.f,-90.f,-100.f}, /*11500*/
+ {-24.f,-24.f,-26.f,-32.f,-32.f,-42.f,-50.f,-60.f,-70.f,-90.f,-100.f}, /*16000*/
+
+ },
+
+ 1,/* peakattp */
+ {{-14.f,-16.f,-18.f,-19.f,-20.f,-21.f,-22.f,-22.f,-24.f,-24.f,-24.f},/*63*/
+ {-14.f,-16.f,-18.f,-19.f,-20.f,-21.f,-22.f,-22.f,-24.f,-24.f,-24.f},/*88*/
+ {-14.f,-16.f,-18.f,-19.f,-20.f,-21.f,-22.f,-22.f,-24.f,-24.f,-24.f},/*125*/
+ {-10.f,-10.f,-10.f,-10.f,-16.f,-16.f,-18.f,-20.f,-24.f,-24.f,-24.f},/*175*/
+ {-10.f,-10.f,-10.f,-10.f,-16.f,-16.f,-18.f,-20.f,-24.f,-24.f,-24.f},/*250*/
+ {-10.f,-10.f,-10.f,-10.f,-16.f,-16.f,-18.f,-20.f,-22.f,-24.f,-24.f},/*350*/
+ {-10.f,-10.f,-10.f,-10.f,-16.f,-16.f,-18.f,-20.f,-22.f,-24.f,-24.f},/*500*/
+ {-10.f,-10.f,-10.f,-10.f,-14.f,-14.f,-16.f,-20.f,-22.f,-24.f,-24.f},/*700*/
+ {-10.f,-10.f,-10.f,-10.f,-14.f,-14.f,-16.f,-20.f,-22.f,-24.f,-24.f},/*1000*/
+ {-10.f,-10.f,-10.f,-10.f,-14.f,-14.f,-16.f,-20.f,-22.f,-24.f,-24.f},/*1400*/
+ {-10.f,-10.f,-10.f,-10.f,-14.f,-14.f,-16.f,-20.f,-22.f,-24.f,-24.f},/*2000*/
+ {-10.f,-10.f,-10.f,-12.f,-16.f,-16.f,-16.f,-20.f,-22.f,-24.f,-24.f},/*2400*/
+ {-10.f,-10.f,-10.f,-12.f,-16.f,-16.f,-16.f,-20.f,-22.f,-24.f,-24.f},/*4000*/
+ {-10.f,-10.f,-10.f,-12.f,-12.f,-14.f,-16.f,-18.f,-22.f,-24.f,-24.f},/*5600*/
+ {-10.f,-10.f,-10.f,-10.f,-10.f,-14.f,-16.f,-18.f,-22.f,-24.f,-24.f},/*8000*/
+ {-10.f,-10.f,-10.f,-10.f,-10.f,-14.f,-16.f,-18.f,-22.f,-24.f,-24.f},/*11500*/
+ {-10.f,-10.f,-10.f,-10.f,-10.f,-12.f,-16.f,-18.f,-22.f,-24.f,-24.f},/*16000*/
+ },
+
+ 1,/*noisemaskp */
+ -24.f, /* suppress any noise curve over maxspec+n */
+ .5f, /* low window */
+ .5f, /* high window */
+ 25,
+ 25,
+ {.000f, 0.f, /*63*/
+ .000f, 0.f, /*88*/
+ .000f, 0.f, /*125*/
+ .000f, 0.f, /*175*/
+ .000f, 0.f, /*250*/
+ .000f, 0.f, /*350*/
+ .000f, 0.f, /*500*/
+ .200f, 0.f, /*700*/
+ .300f, 0.f, /*1000*/
+ .400f, 0.f, /*1400*/
+ .400f, 0.f, /*2000*/
+ .400f, 0.f, /*2800*/
+ .700f, 0.f, /*4000*/
+ .850f, 0.f, /*5600*/
+ .900f, 0.f, /*8000*/
+ .900f, 0.f, /*11500*/
+ .900f, 1.f, /*16000*/
+ },
+
+ 95.f, /* even decade + 5 is important; saves an rint() later in a
+ tight loop) */
+ -28.,
+
+ 1,
+
+};
+
+/* with GNUisms, this could be short and readable. Oh well */
+static vorbis_info_time0 _time_set0A={0};
+static vorbis_info_floor0 _floor_set0A={12, 44100, 64, 10,130, 2, {0,1},
+ 0.199f, .285f};
+static vorbis_info_floor0 _floor_set1A={30, 44100, 256, 12,150, 2, {2,3},
+ .082f, .126f};
+static vorbis_info_residue0 _residue_set0A={0,96,16,6,4,
+ {0,1,1,1,1,1},
+ {6,7,8,9,10},
+
+ {0,99999,9999,9999,9999},
+ {999.f,1.5f,3.5f,15.5f,26.5f},
+ {4,4,4,4,4},
+ {99,99,99,99,99}};
+
+static vorbis_info_residue0 _residue_set1A={0, 960, 32,10,5,
+ {0,1,1,1,1,1,1,1,1,1},
+ {11,12,13,14,15,16,17,18,19},
+
+ {0,8,9999,16,9999,
+ 24,9999,9999,9999},
+ {999.f,1.5f,1.5f,2.5f,2.5f,
+ 6.5f,6.5f,12.5f,22.5f},
+ {5,5,5,5,5,5,5,5,5},
+ {99,99,99,99,99,99,99,99,99}};
+
+static vorbis_info_mapping0 _mapping_set0A={1, {0,0}, {0}, {0}, {0}, {0}};
+static vorbis_info_mapping0 _mapping_set1A={1, {0,0}, {0}, {1}, {1}, {1}};
+static vorbis_info_mode _mode_set0A={0,0,0,0};
+static vorbis_info_mode _mode_set1A={1,0,0,1};
+
+/* CD quality stereo, no channel coupling */
+codec_setup_info info_A={
+
+ /* smallblock, largeblock */
+ {256, 2048},
+ /* modes,maps,times,floors,residues,books,psys */
+ 2, 2, 1, 2, 2, 20, 2,
+ /* modes */
+ {&_mode_set0A,&_mode_set1A},
+ /* maps */
+ {0,0},{&_mapping_set0A,&_mapping_set1A},
+ /* times */
+ {0,0},{&_time_set0A},
+ /* floors */
+ {0,0},{&_floor_set0A,&_floor_set1A},
+ /* residue */
+ {0,0},{&_residue_set0A,&_residue_set1A},
+ /* books */
+ {&_vq_book_lsp12_0, /* 0 */
+ &_vq_book_lsp12_1, /* 1 */
+ &_vq_book_lsp30_0, /* 2 */
+ &_vq_book_lsp30_1, /* 3 */
+
+ &_huff_book_res0_128_128aux,
+ &_huff_book_res0_128_1024aux,
+
+ &_vq_book_res0_128_128_1,
+ &_vq_book_res0_128_128_2,
+ &_vq_book_res0_128_128_3,
+ &_vq_book_res0_128_128_4,
+ &_vq_book_res0_128_128_5,
+
+ &_vq_book_res0_128_1024_1,
+ &_vq_book_res0_128_1024_2,
+ &_vq_book_res0_128_1024_3,
+ &_vq_book_res0_128_1024_4,
+ &_vq_book_res0_128_1024_5,
+ &_vq_book_res0_128_1024_6,
+ &_vq_book_res0_128_1024_7,
+ &_vq_book_res0_128_1024_8,
+ &_vq_book_res0_128_1024_9,
+
+ },
+ /* psy */
+ {&_psy_set_A0,&_psy_set_A},
+
+ /* thresh sample period, preecho clamp trigger threshhold, range, minenergy */
+ 256, {26.f,26.f,26.f,30.f}, {-90.f,-90.f,-90.f,-90.f}, -90.f,
+
+ -10.,
+
+ 0,
+};
+
+#define PREDEF_INFO_MAX 0
+
+#endif
diff --git a/lib/psy.c b/lib/psy.c
new file mode 100644
index 00000000..8cd4272e
--- /dev/null
+++ b/lib/psy.c
@@ -0,0 +1,727 @@
+/********************************************************************
+ * *
+ * THIS FILE IS PART OF THE OggVorbis SOFTWARE CODEC SOURCE CODE. *
+ * USE, DISTRIBUTION AND REPRODUCTION OF THIS LIBRARY SOURCE IS *
+ * GOVERNED BY A BSD-STYLE SOURCE LICENSE INCLUDED WITH THIS SOURCE *
+ * IN 'COPYING'. PLEASE READ THESE TERMS BEFORE DISTRIBUTING. *
+ * *
+ * THE OggVorbis SOURCE CODE IS (C) COPYRIGHT 1994-2001 *
+ * by the XIPHOPHORUS Company http://www.xiph.org/ *
+
+ ********************************************************************
+
+ function: psychoacoustics not including preecho
+ last mod: $Id: psy.c,v 1.44.4.1 2001/04/29 22:21:04 xiphmont Exp $
+
+ ********************************************************************/
+
+#include <stdlib.h>
+#include <math.h>
+#include <string.h>
+#include "vorbis/codec.h"
+#include "codec_internal.h"
+
+#include "masking.h"
+#include "psy.h"
+#include "os.h"
+#include "lpc.h"
+#include "smallft.h"
+#include "scales.h"
+#include "misc.h"
+
+#define NEGINF -9999.f
+
+/* Why Bark scale for encoding but not masking computation? Because
+ masking has a strong harmonic dependancy */
+
+/* the beginnings of real psychoacoustic infrastructure. This is
+ still not tightly tuned */
+void _vi_psy_free(vorbis_info_psy *i){
+ if(i){
+ memset(i,0,sizeof(vorbis_info_psy));
+ _ogg_free(i);
+ }
+}
+
+vorbis_info_psy *_vi_psy_copy(vorbis_info_psy *i){
+ vorbis_info_psy *ret=_ogg_malloc(sizeof(vorbis_info_psy));
+ memcpy(ret,i,sizeof(vorbis_info_psy));
+ return(ret);
+}
+
+/* Set up decibel threshold slopes on a Bark frequency scale */
+/* ATH is the only bit left on a Bark scale. No reason to change it
+ right now */
+static void set_curve(float *ref,float *c,int n, float crate){
+ int i,j=0;
+
+ for(i=0;i<MAX_BARK-1;i++){
+ int endpos=rint(fromBARK(i+1)*2*n/crate);
+ float base=ref[i];
+ if(j<endpos){
+ float delta=(ref[i+1]-base)/(endpos-j);
+ for(;j<endpos && j<n;j++){
+ c[j]=base;
+ base+=delta;
+ }
+ }
+ }
+}
+
+static void min_curve(float *c,
+ float *c2){
+ int i;
+ for(i=0;i<EHMER_MAX;i++)if(c2[i]<c[i])c[i]=c2[i];
+}
+static void max_curve(float *c,
+ float *c2){
+ int i;
+ for(i=0;i<EHMER_MAX;i++)if(c2[i]>c[i])c[i]=c2[i];
+}
+
+static void attenuate_curve(float *c,float att){
+ int i;
+ for(i=0;i<EHMER_MAX;i++)
+ c[i]+=att;
+}
+
+static void interp_curve(float *c,float *c1,float *c2,float del){
+ int i;
+ for(i=0;i<EHMER_MAX;i++)
+ c[i]=c2[i]*del+c1[i]*(1.f-del);
+}
+
+static void setup_curve(float **c,
+ int band,
+ float *curveatt_dB){
+ int i,j;
+ float ath[EHMER_MAX];
+ float tempc[P_LEVELS][EHMER_MAX];
+
+ memcpy(c[0]+2,c[4]+2,sizeof(float)*EHMER_MAX);
+ memcpy(c[2]+2,c[4]+2,sizeof(float)*EHMER_MAX);
+
+ /* we add back in the ATH to avoid low level curves falling off to
+ -infinity and unneccessarily cutting off high level curves in the
+ curve limiting (last step). But again, remember... a half-band's
+ settings must be valid over the whole band, and it's better to
+ mask too little than too much, so be pessimal. */
+
+ for(i=0;i<EHMER_MAX;i++){
+ float oc_min=band*.5+(i-EHMER_OFFSET)*.125;
+ float oc_max=band*.5+(i-EHMER_OFFSET+1)*.125;
+ float bark=toBARK(fromOC(oc_min));
+ int ibark=floor(bark);
+ float del=bark-ibark;
+ float ath_min,ath_max;
+
+ if(ibark<26)
+ ath_min=ATH_Bark_dB[ibark]*(1.f-del)+ATH_Bark_dB[ibark+1]*del;
+ else
+ ath_min=ATH_Bark_dB[25];
+
+ bark=toBARK(fromOC(oc_max));
+ ibark=floor(bark);
+ del=bark-ibark;
+
+ if(ibark<26)
+ ath_max=ATH_Bark_dB[ibark]*(1.f-del)+ATH_Bark_dB[ibark+1]*del;
+ else
+ ath_max=ATH_Bark_dB[25];
+
+ ath[i]=min(ath_min,ath_max);
+ }
+
+ /* The c array is comes in as dB curves at 20 40 60 80 100 dB.
+ interpolate intermediate dB curves */
+ for(i=1;i<P_LEVELS;i+=2){
+ interp_curve(c[i]+2,c[i-1]+2,c[i+1]+2,.5);
+ }
+
+ /* normalize curves so the driving amplitude is 0dB */
+ /* make temp curves with the ATH overlayed */
+ for(i=0;i<P_LEVELS;i++){
+ attenuate_curve(c[i]+2,curveatt_dB[i]);
+ memcpy(tempc[i],ath,EHMER_MAX*sizeof(float));
+ attenuate_curve(tempc[i],-i*10.f);
+ max_curve(tempc[i],c[i]+2);
+ }
+
+ /* Now limit the louder curves.
+
+ the idea is this: We don't know what the playback attenuation
+ will be; 0dB SL moves every time the user twiddles the volume
+ knob. So that means we have to use a single 'most pessimal' curve
+ for all masking amplitudes, right? Wrong. The *loudest* sound
+ can be in (we assume) a range of ...+100dB] SL. However, sounds
+ 20dB down will be in a range ...+80], 40dB down is from ...+60],
+ etc... */
+
+ for(j=1;j<P_LEVELS;j++){
+ min_curve(tempc[j],tempc[j-1]);
+ min_curve(c[j]+2,tempc[j]);
+ }
+
+ /* add fenceposts */
+ for(j=0;j<P_LEVELS;j++){
+
+ for(i=0;i<EHMER_MAX;i++)
+ if(c[j][i+2]>-200.f)break;
+ c[j][0]=i;
+
+ for(i=EHMER_MAX-1;i>=0;i--)
+ if(c[j][i+2]>-200.f)
+ break;
+ c[j][1]=i;
+
+ }
+}
+
+void _vp_psy_init(vorbis_look_psy *p,vorbis_info_psy *vi,int n,long rate){
+ long i,j;
+ long maxoc;
+ memset(p,0,sizeof(vorbis_look_psy));
+
+
+ p->eighth_octave_lines=vi->eighth_octave_lines;
+ p->shiftoc=rint(log(vi->eighth_octave_lines*8)/log(2))-1;
+
+ p->firstoc=toOC(.25f*rate/n)*(1<<(p->shiftoc+1))-vi->eighth_octave_lines;
+ maxoc=toOC((n*.5f-.25f)*rate/n)*(1<<(p->shiftoc+1))+.5f;
+ p->total_octave_lines=maxoc-p->firstoc+1;
+
+ p->ath=_ogg_malloc(n*sizeof(float));
+ p->octave=_ogg_malloc(n*sizeof(long));
+ p->bark=_ogg_malloc(n*sizeof(float));
+ p->vi=vi;
+ p->n=n;
+
+ /* set up the lookups for a given blocksize and sample rate */
+ /* Vorbis max sample rate is currently limited by 26 Bark (54kHz) */
+ set_curve(ATH_Bark_dB, p->ath,n,rate);
+ for(i=0;i<n;i++)
+ p->bark[i]=toBARK(rate/(2*n)*i);
+
+ for(i=0;i<n;i++)
+ p->octave[i]=toOC((i*.5f+.25f)*rate/n)*(1<<(p->shiftoc+1))+.5f;
+
+ p->tonecurves=_ogg_malloc(P_BANDS*sizeof(float **));
+ p->noisemedian=_ogg_malloc(n*sizeof(float));
+ p->noiseoffset=_ogg_malloc(n*sizeof(float));
+ p->peakatt=_ogg_malloc(P_BANDS*sizeof(float *));
+ for(i=0;i<P_BANDS;i++){
+ p->tonecurves[i]=_ogg_malloc(P_LEVELS*sizeof(float *));
+ p->peakatt[i]=_ogg_malloc(P_LEVELS*sizeof(float));
+ }
+
+ for(i=0;i<P_BANDS;i++)
+ for(j=0;j<P_LEVELS;j++){
+ p->tonecurves[i][j]=_ogg_malloc((EHMER_MAX+2)*sizeof(float));
+ }
+
+ /* OK, yeah, this was a silly way to do it */
+ memcpy(p->tonecurves[0][4]+2,tone_125_40dB_SL,sizeof(float)*EHMER_MAX);
+ memcpy(p->tonecurves[0][6]+2,tone_125_60dB_SL,sizeof(float)*EHMER_MAX);
+ memcpy(p->tonecurves[0][8]+2,tone_125_80dB_SL,sizeof(float)*EHMER_MAX);
+ memcpy(p->tonecurves[0][10]+2,tone_125_100dB_SL,sizeof(float)*EHMER_MAX);
+
+ memcpy(p->tonecurves[2][4]+2,tone_125_40dB_SL,sizeof(float)*EHMER_MAX);
+ memcpy(p->tonecurves[2][6]+2,tone_125_60dB_SL,sizeof(float)*EHMER_MAX);
+ memcpy(p->tonecurves[2][8]+2,tone_125_80dB_SL,sizeof(float)*EHMER_MAX);
+ memcpy(p->tonecurves[2][10]+2,tone_125_100dB_SL,sizeof(float)*EHMER_MAX);
+
+ memcpy(p->tonecurves[4][4]+2,tone_250_40dB_SL,sizeof(float)*EHMER_MAX);
+ memcpy(p->tonecurves[4][6]+2,tone_250_60dB_SL,sizeof(float)*EHMER_MAX);
+ memcpy(p->tonecurves[4][8]+2,tone_250_80dB_SL,sizeof(float)*EHMER_MAX);
+ memcpy(p->tonecurves[4][10]+2,tone_250_100dB_SL,sizeof(float)*EHMER_MAX);
+
+ memcpy(p->tonecurves[6][4]+2,tone_500_40dB_SL,sizeof(float)*EHMER_MAX);
+ memcpy(p->tonecurves[6][6]+2,tone_500_60dB_SL,sizeof(float)*EHMER_MAX);
+ memcpy(p->tonecurves[6][8]+2,tone_500_80dB_SL,sizeof(float)*EHMER_MAX);
+ memcpy(p->tonecurves[6][10]+2,tone_500_100dB_SL,sizeof(float)*EHMER_MAX);
+
+ memcpy(p->tonecurves[8][4]+2,tone_1000_40dB_SL,sizeof(float)*EHMER_MAX);
+ memcpy(p->tonecurves[8][6]+2,tone_1000_60dB_SL,sizeof(float)*EHMER_MAX);
+ memcpy(p->tonecurves[8][8]+2,tone_1000_80dB_SL,sizeof(float)*EHMER_MAX);
+ memcpy(p->tonecurves[8][10]+2,tone_1000_100dB_SL,sizeof(float)*EHMER_MAX);
+
+ memcpy(p->tonecurves[10][4]+2,tone_2000_40dB_SL,sizeof(float)*EHMER_MAX);
+ memcpy(p->tonecurves[10][6]+2,tone_2000_60dB_SL,sizeof(float)*EHMER_MAX);
+ memcpy(p->tonecurves[10][8]+2,tone_2000_80dB_SL,sizeof(float)*EHMER_MAX);
+ memcpy(p->tonecurves[10][10]+2,tone_2000_100dB_SL,sizeof(float)*EHMER_MAX);
+
+ memcpy(p->tonecurves[12][4]+2,tone_4000_40dB_SL,sizeof(float)*EHMER_MAX);
+ memcpy(p->tonecurves[12][6]+2,tone_4000_60dB_SL,sizeof(float)*EHMER_MAX);
+ memcpy(p->tonecurves[12][8]+2,tone_4000_80dB_SL,sizeof(float)*EHMER_MAX);
+ memcpy(p->tonecurves[12][10]+2,tone_4000_100dB_SL,sizeof(float)*EHMER_MAX);
+
+ memcpy(p->tonecurves[14][4]+2,tone_8000_40dB_SL,sizeof(float)*EHMER_MAX);
+ memcpy(p->tonecurves[14][6]+2,tone_8000_60dB_SL,sizeof(float)*EHMER_MAX);
+ memcpy(p->tonecurves[14][8]+2,tone_8000_80dB_SL,sizeof(float)*EHMER_MAX);
+ memcpy(p->tonecurves[14][10]+2,tone_8000_100dB_SL,sizeof(float)*EHMER_MAX);
+
+ memcpy(p->tonecurves[16][4]+2,tone_8000_40dB_SL,sizeof(float)*EHMER_MAX);
+ memcpy(p->tonecurves[16][6]+2,tone_8000_60dB_SL,sizeof(float)*EHMER_MAX);
+ memcpy(p->tonecurves[16][8]+2,tone_8000_80dB_SL,sizeof(float)*EHMER_MAX);
+ memcpy(p->tonecurves[16][10]+2,tone_8000_100dB_SL,sizeof(float)*EHMER_MAX);
+
+ /* interpolate curves between */
+ for(i=1;i<P_BANDS;i+=2)
+ for(j=4;j<P_LEVELS;j+=2){
+ memcpy(p->tonecurves[i][j]+2,p->tonecurves[i-1][j]+2,EHMER_MAX*sizeof(float));
+ /*interp_curve(p->tonecurves[i][j],
+ p->tonecurves[i-1][j],
+ p->tonecurves[i+1][j],.5);*/
+ min_curve(p->tonecurves[i][j]+2,p->tonecurves[i+1][j]+2);
+ }
+
+ /* set up the final curves */
+ for(i=0;i<P_BANDS;i++)
+ setup_curve(p->tonecurves[i],i,vi->toneatt[i]);
+
+ /* set up attenuation levels */
+ for(i=0;i<P_BANDS;i++)
+ for(j=0;j<P_LEVELS;j++){
+ p->peakatt[i][j]=p->vi->peakatt[i][j];
+ }
+
+ /* set up rolling noise median */
+ for(i=0;i<n;i++){
+ float halfoc=toOC((i+.5)*rate/(2.*n))*2.+2.;
+ int inthalfoc;
+ float del;
+
+ if(halfoc<0)halfoc=0;
+ if(halfoc>=P_BANDS-1)halfoc=P_BANDS-1;
+ inthalfoc=(int)halfoc;
+ del=halfoc-inthalfoc;
+
+ p->noisemedian[i]=
+ p->vi->noisemedian[inthalfoc*2]*(1.-del) +
+ p->vi->noisemedian[inthalfoc*2+2]*del;
+ p->noiseoffset[i]=
+ p->vi->noisemedian[inthalfoc*2+1]*(1.-del) +
+ p->vi->noisemedian[inthalfoc*2+3]*del;
+ }
+ /*_analysis_output("mediancurve",0,p->noisemedian,n,0,0);*/
+}
+
+void _vp_psy_clear(vorbis_look_psy *p){
+ int i,j;
+ if(p){
+ if(p->ath)_ogg_free(p->ath);
+ if(p->octave)_ogg_free(p->octave);
+ if(p->bark)_ogg_free(p->bark);
+ if(p->tonecurves){
+ for(i=0;i<P_BANDS;i++){
+ for(j=0;j<P_LEVELS;j++){
+ _ogg_free(p->tonecurves[i][j]);
+ }
+ _ogg_free(p->tonecurves[i]);
+ _ogg_free(p->peakatt[i]);
+ }
+ _ogg_free(p->tonecurves);
+ _ogg_free(p->noisemedian);
+ _ogg_free(p->noiseoffset);
+ _ogg_free(p->peakatt);
+ }
+ memset(p,0,sizeof(vorbis_look_psy));
+ }
+}
+
+/* octave/(8*eighth_octave_lines) x scale and dB y scale */
+static void seed_curve(float *seed,
+ float **curves,
+ float amp,
+ int oc,int n,int linesper,float dBoffset){
+ int i;
+ long seedptr;
+ float *posts,*curve;
+
+ int choice=(int)((amp+dBoffset)*.1f);
+ choice=max(choice,0);
+ choice=min(choice,P_LEVELS-1);
+ posts=curves[choice];
+ curve=posts+2;
+ seedptr=oc+(posts[0]-16)*linesper-(linesper>>1);
+
+ for(i=posts[0];i<posts[1];i++){
+ if(seedptr>0){
+ float lin=amp+curve[i];
+ if(seed[seedptr]<lin)seed[seedptr]=lin;
+ }
+ seedptr+=linesper;
+ if(seedptr>=n)break;
+ }
+}
+
+static void seed_peak(float *seed,
+ float *att,
+ float amp,
+ int oc,
+ int linesper,
+ float dBoffset){
+ long seedptr;
+
+ int choice=(int)((amp+dBoffset)*.1f);
+ choice=max(choice,0);
+ choice=min(choice,P_LEVELS-1);
+ seedptr=oc-(linesper>>1);
+
+ amp+=att[choice];
+ if(seed[seedptr]<amp)seed[seedptr]=amp;
+
+}
+
+static void seed_loop(vorbis_look_psy *p,
+ float ***curves,
+ float **att,
+ float *f,
+ float *flr,
+ float *minseed,
+ float *maxseed,
+ float specmax){
+ vorbis_info_psy *vi=p->vi;
+ long n=p->n,i;
+ float dBoffset=vi->max_curve_dB-specmax;
+
+ /* prime the working vector with peak values */
+
+ for(i=0;i<n;i++){
+ float max=f[i];
+ long oc=p->octave[i];
+ while(i+1<n && p->octave[i+1]==oc){
+ i++;
+ if(f[i]>max)max=f[i];
+ }
+
+ if(max>flr[i]){
+ oc=oc>>p->shiftoc;
+ if(oc>=P_BANDS)oc=P_BANDS-1;
+ if(oc<0)oc=0;
+ if(vi->tonemaskp)
+ seed_curve(minseed,
+ curves[oc],
+ max,
+ p->octave[i]-p->firstoc,
+ p->total_octave_lines,
+ p->eighth_octave_lines,
+ dBoffset);
+ if(vi->peakattp)
+ seed_peak(maxseed,
+ att[oc],
+ max,
+ p->octave[i]-p->firstoc,
+ p->eighth_octave_lines,
+ dBoffset);
+ }
+ }
+}
+
+static void bound_loop(vorbis_look_psy *p,
+ float *f,
+ float *seeds,
+ float *flr,
+ float att){
+ long n=p->n,i;
+
+ long off=(p->eighth_octave_lines>>1)+p->firstoc;
+ long *ocp=p->octave;
+
+ for(i=0;i<n;i++){
+ long oc=ocp[i]-off;
+ float v=f[i]+att;
+ if(seeds[oc]<v)seeds[oc]=v;
+ }
+}
+
+static void seed_chase(float *seeds, int linesper, long n){
+ long *posstack=alloca(n*sizeof(long));
+ float *ampstack=alloca(n*sizeof(float));
+ long stack=0;
+ long pos=0;
+ long i;
+
+ for(i=0;i<n;i++){
+ if(stack<2){
+ posstack[stack]=i;
+ ampstack[stack++]=seeds[i];
+ }else{
+ while(1){
+ if(seeds[i]<ampstack[stack-1]){
+ posstack[stack]=i;
+ ampstack[stack++]=seeds[i];
+ break;
+ }else{
+ if(i<posstack[stack-1]+linesper){
+ if(stack>1 && ampstack[stack-1]<=ampstack[stack-2] &&
+ i<posstack[stack-2]+linesper){
+ /* we completely overlap, making stack-1 irrelevant. pop it */
+ stack--;
+ continue;
+ }
+ }
+ posstack[stack]=i;
+ ampstack[stack++]=seeds[i];
+ break;
+
+ }
+ }
+ }
+ }
+
+ /* the stack now contains only the positions that are relevant. Scan
+ 'em straight through */
+
+ for(i=0;i<stack;i++){
+ long endpos;
+ if(i<stack-1 && ampstack[i+1]>ampstack[i]){
+ endpos=posstack[i+1];
+ }else{
+ endpos=posstack[i]+linesper+1; /* +1 is important, else bin 0 is
+ discarded in short frames */
+ }
+ if(endpos>n)endpos=n;
+ for(;pos<endpos;pos++)
+ seeds[pos]=ampstack[i];
+ }
+
+ /* there. Linear time. I now remember this was on a problem set I
+ had in Grad Skool... I didn't solve it at the time ;-) */
+
+}
+
+/* bleaugh, this is more complicated than it needs to be */
+static void max_seeds(vorbis_look_psy *p,float *minseed,float *maxseed,
+ float *flr){
+ long n=p->total_octave_lines;
+ int linesper=p->eighth_octave_lines;
+ long linpos=0;
+ long pos;
+
+ seed_chase(minseed,linesper,n); /* for masking */
+ seed_chase(maxseed,linesper,n); /* for peak att */
+
+ pos=p->octave[0]-p->firstoc-(linesper>>1);
+ while(linpos+1<p->n){
+ float min=minseed[pos];
+ float max=maxseed[pos];
+ long end=((p->octave[linpos]+p->octave[linpos+1])>>1)-p->firstoc;
+ while(pos+1<=end){
+ pos++;
+ if((minseed[pos]>NEGINF && minseed[pos]<min) || min==NEGINF)
+ min=minseed[pos];
+ if(maxseed[pos]>max)max=maxseed[pos];
+ }
+ if(max<min)max=min;
+
+ /* seed scale is log. Floor is linear. Map back to it */
+ end=pos+p->firstoc;
+ for(;linpos<p->n && p->octave[linpos]<=end;linpos++)
+ if(flr[linpos]<max)flr[linpos]=max;
+ }
+
+ {
+ float min=minseed[p->total_octave_lines-1];
+ float max=maxseed[p->total_octave_lines-1];
+ if(max<min)max=min;
+ for(;linpos<p->n;linpos++)
+ if(flr[linpos]<max)flr[linpos]=max;
+ }
+
+}
+
+/* quarter-dB bins */
+#define BIN(x) ((int)((x)*negFour))
+#define BINdB(x) ((x)*negQuarter)
+#define BINCOUNT (200*4)
+#define LASTBIN (BINCOUNT-1)
+
+static void bark_noise_median(long n,float *b,float *f,float *noise,
+ float lowidth,float hiwidth,
+ int lomin,int himin,
+ float *thresh,float *off){
+ long i=0,lo=0,hi=0;
+ float bi,threshi;
+ long median=LASTBIN;
+ float negFour = -4.0f;
+ float negQuarter = -0.25f;
+
+ /* these are really integral values, but we store them in floats to
+ avoid excessive float/int conversions, which GCC and MSVC are
+ farily poor at optimizing. */
+
+ float radix[BINCOUNT];
+ float countabove=0;
+ float countbelow=0;
+
+ memset(radix,0,sizeof(radix));
+
+ for(i=0;i<n;i++){
+ /* find new lo/hi */
+ bi=b[i]+hiwidth;
+ for(;hi<n && (hi<i+himin || b[hi]<=bi);hi++){
+ int bin=BIN(f[hi]);
+ if(bin>LASTBIN)bin=LASTBIN;
+ if(bin<0)bin=0;
+ radix[bin]++;
+ if(bin<median)
+ countabove++;
+ else
+ countbelow++;
+ }
+ bi=b[i]-lowidth;
+ for(;lo<i && lo+lomin<i && b[lo]<=bi;lo++){
+ int bin=BIN(f[lo]);
+ if(bin>LASTBIN)bin=LASTBIN;
+ if(bin<0)bin=0;
+ radix[bin]--;
+ if(bin<median)
+ countabove--;
+ else
+ countbelow--;
+ }
+
+ /* move the median if needed */
+ if(countabove+countbelow){
+ threshi = thresh[i]*(countabove+countbelow);
+
+ while(threshi>countbelow && median>0){
+ median--;
+ countabove-=radix[median];
+ countbelow+=radix[median];
+ }
+
+ while(threshi<(countbelow-radix[median]) &&
+ median<LASTBIN){
+ countabove+=radix[median];
+ countbelow-=radix[median];
+ median++;
+ }
+ }
+ noise[i]=BINdB(median)+off[i];
+ }
+
+}
+
+float _vp_compute_mask(vorbis_look_psy *p,
+ float *fft,
+ float *mdct,
+ float *flr,
+ float *decay,
+ float specmax){
+ int i,n=p->n;
+ float localmax=NEGINF;
+ static int seq=0;
+
+ float *minseed=alloca(sizeof(float)*p->total_octave_lines);
+ float *maxseed=alloca(sizeof(float)*p->total_octave_lines);
+ for(i=0;i<p->total_octave_lines;i++)minseed[i]=maxseed[i]=NEGINF;
+
+ /* go to dB scale. Also find the highest peak so we know the limits */
+ for(i=0;i<n;i++){
+ fft[i]=todB_nn(fft[i]);
+ if(fft[i]>localmax)localmax=fft[i];
+ }
+ if(specmax<localmax)specmax=localmax;
+
+
+ for(i=0;i<n;i++){
+ mdct[i]=todB(mdct[i]);
+ }
+
+ _analysis_output("mdct",seq,mdct,n,0,0);
+ _analysis_output("fft",seq,fft,n,0,0);
+
+ /* noise masking */
+ if(p->vi->noisemaskp){
+ bark_noise_median(n,p->bark,mdct,flr,
+ p->vi->noisewindowlo,
+ p->vi->noisewindowhi,
+ p->vi->noisewindowlomin,
+ p->vi->noisewindowhimin,
+ p->noisemedian,
+ p->noiseoffset);
+ /* suppress any noise curve > specmax+p->vi->noisemaxsupp */
+ for(i=0;i<n;i++)
+ if(flr[i]>specmax+p->vi->noisemaxsupp)
+ flr[i]=specmax+p->vi->noisemaxsupp;
+ _analysis_output("noise",seq,flr,n,0,0);
+ }else{
+ for(i=0;i<n;i++)flr[i]=NEGINF;
+ }
+
+ /* set the ATH (floating below localmax, not global max by a
+ specified att) */
+ if(p->vi->athp){
+ float att=localmax+p->vi->ath_adjatt;
+ if(att<p->vi->ath_maxatt)att=p->vi->ath_maxatt;
+
+ for(i=0;i<n;i++){
+ float av=p->ath[i]+att;
+ if(av>flr[i])flr[i]=av;
+ }
+ }
+
+ _analysis_output("ath",seq,flr,n,0,0);
+
+ /* tone/peak masking */
+
+ /* XXX apply decay to the fft here */
+
+ seed_loop(p,p->tonecurves,p->peakatt,fft,flr,minseed,maxseed,specmax);
+ bound_loop(p,mdct,maxseed,flr,p->vi->bound_att_dB);
+ _analysis_output("minseed",seq,minseed,p->total_octave_lines,0,0);
+ _analysis_output("maxseed",seq,maxseed,p->total_octave_lines,0,0);
+ max_seeds(p,minseed,maxseed,flr);
+ _analysis_output("final",seq,flr,n,0,0);
+
+ /* doing this here is clean, but we need to find a faster way to do
+ it than to just tack it on */
+
+ if(p->vi->floor_cullp){
+ for(i=0;i<n;i++)
+ if(mdct[i]+6.<flr[i])
+ flr[i]=NEGINF;
+ }else{
+ for(i=0;i<n;i++)if(mdct[i]>=flr[i])break;
+ if(i==n)for(i=0;i<n;i++)flr[i]=NEGINF;
+ }
+
+ seq++;
+
+ return(specmax);
+}
+
+
+/* this applies the floor and (optionally) tries to preserve noise
+ energy in low resolution portions of the spectrum */
+/* f and flr are *linear* scale, not dB */
+void _vp_apply_floor(vorbis_look_psy *p,float *f, float *flr){
+ float *work=alloca(p->n*sizeof(float));
+ int j;
+
+ /* subtract the floor */
+ for(j=0;j<p->n;j++){
+ if(flr[j]<=0)
+ work[j]=0.f;
+ else
+ work[j]=f[j]/flr[j];
+ }
+
+ memcpy(f,work,p->n*sizeof(float));
+}
+
+float _vp_ampmax_decay(float amp,vorbis_dsp_state *vd){
+ vorbis_info *vi=vd->vi;
+ codec_setup_info *ci=vi->codec_setup;
+ int n=ci->blocksizes[vd->W]/2;
+ float secs=(float)n/vi->rate;
+
+ amp+=secs*ci->ampmax_att_per_sec;
+ if(amp<-9999)amp=-9999;
+ return(amp);
+}
+
+
+
diff --git a/lib/psy.h b/lib/psy.h
new file mode 100644
index 00000000..1e0afccd
--- /dev/null
+++ b/lib/psy.h
@@ -0,0 +1,101 @@
+/********************************************************************
+ * *
+ * THIS FILE IS PART OF THE OggVorbis SOFTWARE CODEC SOURCE CODE. *
+ * USE, DISTRIBUTION AND REPRODUCTION OF THIS LIBRARY SOURCE IS *
+ * GOVERNED BY A BSD-STYLE SOURCE LICENSE INCLUDED WITH THIS SOURCE *
+ * IN 'COPYING'. PLEASE READ THESE TERMS BEFORE DISTRIBUTING. *
+ * *
+ * THE OggVorbis SOURCE CODE IS (C) COPYRIGHT 1994-2001 *
+ * by the XIPHOPHORUS Company http://www.xiph.org/ *
+
+ ********************************************************************
+
+ function: random psychoacoustics (not including preecho)
+ last mod: $Id: psy.h,v 1.19.4.1 2001/04/29 22:21:05 xiphmont Exp $
+
+ ********************************************************************/
+
+#ifndef _V_PSY_H_
+#define _V_PSY_H_
+#include "smallft.h"
+
+#ifndef EHMER_MAX
+#define EHMER_MAX 56
+#endif
+
+/* psychoacoustic setup ********************************************/
+#define MAX_BARK 27
+#define P_BANDS 17
+#define P_LEVELS 11
+typedef struct vorbis_info_psy{
+ int athp;
+ int decayp;
+
+ float ath_adjatt;
+ float ath_maxatt;
+
+ int eighth_octave_lines;
+
+ /* 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 */
+ /* x: 63 88 125 175 250 350 500 700 1k 1.4k 2k 2.8k 4k 5.6k 8k 11.5k 16k Hz */
+ /* y: 0 10 20 30 40 50 60 70 80 90 100 dB */
+
+ int tonemaskp;
+ float toneatt[P_BANDS][P_LEVELS];
+
+ int peakattp;
+ float peakatt[P_BANDS][P_LEVELS];
+
+ int noisemaskp;
+ float noisemaxsupp;
+ float noisewindowlo;
+ float noisewindowhi;
+ int noisewindowlomin;
+ int noisewindowhimin;
+ float noisemedian[P_BANDS*2];
+
+ float max_curve_dB;
+ float bound_att_dB;
+
+ int floor_cullp;
+} vorbis_info_psy;
+
+typedef struct {
+ int n;
+ struct vorbis_info_psy *vi;
+
+ float ***tonecurves;
+ float **peakatt;
+ float *noisemedian;
+ float *noiseoffset;
+
+ float *ath;
+ long *octave; /* in n.ocshift format */
+ float *bark;
+
+ long firstoc;
+ long shiftoc;
+ int eighth_octave_lines; /* power of two, please */
+ int total_octave_lines;
+
+} vorbis_look_psy;
+
+extern void _vp_psy_init(vorbis_look_psy *p,vorbis_info_psy *vi,int n,long rate);
+extern void _vp_psy_clear(vorbis_look_psy *p);
+extern void *_vi_psy_dup(void *source);
+
+extern void _vi_psy_free(vorbis_info_psy *i);
+extern vorbis_info_psy *_vi_psy_copy(vorbis_info_psy *i);
+
+extern float _vp_compute_mask(vorbis_look_psy *p,
+ float *fft,
+ float *mdct,
+ float *floor,
+ float *decay,
+ float prev_maxamp);
+extern void _vp_apply_floor(vorbis_look_psy *p,float *f,float *flr);
+extern float _vp_ampmax_decay(float amp,vorbis_dsp_state *vd);
+
+#endif
+
+