function: PCM data vector blocking, windowing and dis/reassembly
author: Monty <xiphmont@mit.edu>
modifications by: Monty
- last modification date: Jul 29 1999
+ last modification date: Aug 05 1999
Handle windowing, overlap-add, etc of the PCM vectors. This is made
- more amusing by Vorbis' current two allowed block sizes (512 and 2048
- elements/channel).
+ more amusing by Vorbis' current two allowed block sizes.
Vorbis manipulates the dynamic range of the incoming PCM data
envelope to minimise time-domain energy leakage from percussive and
static int _vds_shared_init(vorbis_dsp_state *v,vorbis_info *vi){
memset(v,0,sizeof(vorbis_dsp_state));
+
+ memcpy(&v->vi,vi,sizeof(vorbis_info));
+ _ve_envelope_init(&v->ve,vi->envelopesa);
+
v->samples_per_envelope_step=vi->envelopesa;
v->block_size[0]=vi->smallblock;
v->block_size[1]=vi->largeblock;
if(vi->envelopech){
v->envelope_storage=v->pcm_storage/v->samples_per_envelope_step;
v->envelope_channels=vi->envelopech;
- v->multipliers=calloc(v->envelope_channels,sizeof(int *));
+ v->multipliers=calloc(v->envelope_channels,sizeof(double *));
{
int i;
for(i=0;i<v->envelope_channels;i++){
- v->multipliers[i]=calloc(v->envelope_storage,sizeof(int));
+ v->multipliers[i]=calloc(v->envelope_storage,sizeof(double));
}
}
}
vi->largeblock=2048;
vi->envelopesa=64;
vi->envelopech=vi->channels;
+ vi->preecho_thresh=10.;
+ vi->preecho_thresh=4.;
+ vi->envelopemap=calloc(2,sizeof(int));
+ vi->envelopemap[0]=0;
+ vi->envelopemap[1]=1;
_vds_shared_init(v,vi);
for(i=0;i<vb->pcm_channels;i++)
vb->pcm[i]=malloc(vb->pcm_storage*sizeof(double));
- vb->mult=malloc(vb->mult_channels*sizeof(int *));
+ vb->mult=malloc(vb->mult_channels*sizeof(double *));
for(i=0;i<vb->mult_channels;i++)
- vb->mult[i]=malloc(vb->mult_storage*sizeof(int));
+ vb->mult[i]=malloc(vb->mult_storage*sizeof(double));
return(0);
}
data, fill them up in before proceeding. */
if(v->pcm_current/v->samples_per_envelope_step>v->envelope_current){
+ /* This generates the multipliers, but does not sparsify the vector.
+ That's done by block before coding */
_ve_envelope_multipliers(v);
}
for(;i<v->envelope_current;i++){
for(j=0;j<v->envelope_channels;j++)
- if(v->multipliers[j][i])break;
+ if(v->multipliers[j][i-1]*v->vi.preecho_thresh<
+ v->multipliers[j][i])break;
if(j<v->envelope_channels)break;
}
memcpy(vb->pcm[i],v->pcm[i]+beginW,v->block_size[v->W]*sizeof(double));
for(i=0;i<v->envelope_channels;i++)
memcpy(vb->mult[i],v->multipliers[i]+beginM,v->block_size[v->W]/
- v->samples_per_envelope_step*sizeof(int));
+ v->samples_per_envelope_step*sizeof(double));
vb->frameno=v->frame;
int movementW=centerNext-new_centerNext;
int movementM=movementW/v->samples_per_envelope_step;
+ /* the multipliers and pcm stay synced up because the blocksizes
+ must be multiples of samples_per_envelope_step (minimum
+ multiple is 2) */
+
for(i=0;i<v->pcm_channels;i++)
memmove(v->pcm[i],v->pcm[i]+movementW,
(v->pcm_current-movementW)*sizeof(double));
for(i=0;i<v->envelope_channels;i++){
memmove(v->multipliers[i],v->multipliers[i]+movementM,
- (v->envelope_current-movementM)*sizeof(int));
+ (v->envelope_current-movementM)*sizeof(double));
}
v->pcm_current-=movementW;
}
int vorbis_synthesis_init(vorbis_dsp_state *v,vorbis_info *vi){
+ int temp=vi->envelopech;
+ vi->envelopech=0; /* we don't need multiplier buffering in syn */
_vds_shared_init(v,vi);
+ vi->envelopech=temp;
+
/* Adjust centerW to allow an easier mechanism for determining output */
v->pcm_returned=v->centerW;
v->centerW-= v->block_size[v->W]/4+v->block_size[v->lW]/4;
return(0);
}
-/* Unike in analysis, the window is only partially applied, and
- envelope *not* applied. */
+/* Unike in analysis, the window is only partially applied. Envelope
+ is previously applied (the whole envelope, if any, is shipped in
+ each block) */
int vorbis_synthesis_blockin(vorbis_dsp_state *v,vorbis_block *vb){
}
{
- int envperlong=v->block_size[1]/v->samples_per_envelope_step;
- int envperW=v->block_size[v->W]/v->samples_per_envelope_step;
int sizeW=v->block_size[vb->W];
int centerW=v->centerW+v->block_size[vb->lW]/4+sizeW/4;
int beginW=centerW-sizeW/2;
int endW=beginW+sizeW;
int beginSl;
int endSl;
- int i,j,k;
+ int i,j;
double *window;
/* Do we have enough PCM storage for the block? */
v->pcm[i]=realloc(v->pcm[i],v->pcm_storage*sizeof(double));
}
- /* multiplier storage works differently in decode than it does in
- encode; we only need to buffer the last 'half largeblock' and
- we don't need to keep it aligned with the PCM. */
-
- /* fill in the first half of the block's multipliers */
- i=v->envelope_current-1;
- j=envperW/2-1;
- for(;j>=0;j--)
- for(k=0;k<v->envelope_channels;k++)
- vb->mult[k][j]=v->multipliers[k][i];
-
- /* shift out unneeded buffered multipliers */
- {
- int needed=(envperlong-envperW)/2;
- if(needed){
- /* We need to keep some of the buffered ones */
- for(k=0;k<v->envelope_channels;k++)
- memmove(v->multipliers[k],
- v->multipliers[k]+v->envelope_current-needed,
- needed*sizeof(int));
- }
- v->envelope_current=needed;
- }
-
- /* add block's second half to the multiplier buffer */
- /* init makes certain we have enough storage; we only buffer a
- half longblock */
-
- for(i=envperW/2;i<envperW;i++){
- j=v->envelope_current;
- for(k=0;k<v->envelope_channels;k++)
- v->multipliers[k][j++]=vb->mult[k][i];
- }
- v->envelope_current+=envperW/2;
-
- /* manufacture/apply multiplier vector */
-
- _ve_envelope_apply(vb);
-
/* Overlap/add */
switch(vb->W){
case 0:
#include <stdio.h>
#include <math.h>
-void _ve_envelope_multipliers(vorbis_dsp_state *v){
- /* set 'random' deltas... */
- int new_current=v->pcm_current/v->samples_per_envelope_step;
- int i,j;
-
- for(i=v->envelope_current;i<new_current;i++){
- int flag=0;
- for(j=0;j<v->samples_per_envelope_step;j++)
- if(v->pcm[0][j+i*v->samples_per_envelope_step]>5)flag=1;
-
- for(j=0;j<v->envelope_channels;j++)
- v->multipliers[j][i]=flag;
- }
- v->envelope_current=i;
-}
-
-void _ve_envelope_apply(vorbis_block *vb){
-
-}
-
/* basic test of PCM blocking:
construct a PCM vector and block it using preset sizing in our fake
vi.rate=44100;
vi.version=0;
vi.mode=0;
-
vi.user_comments=temp;
vi.vendor="Xiphophorus";
double *window=encode.window[vb.W][vb.lW][vb.nW];
+ /* apply envelope */
+ _ve_envelope_sparsify(&vb);
+ _ve_envelope_apply(&vb,0);
+
for(i=0;i<vb.pcm_channels;i++)
MDCT(vb.pcm[i],vb.pcm[i],ml[vb.W],window);
out=fopen(path,"w");
for(i=0;i<avail;i++)
- fprintf(out,"%ld %g\n",i+beginW,vb.pcm[0][i]);
+ fprintf(out,"%d %g\n",i+beginW,vb.pcm[0][i]);
fprintf(out,"\n");
for(i=0;i<avail;i++)
- fprintf(out,"%ld %g\n",i+beginW,window[i]);
+ fprintf(out,"%d %g\n",i+beginW,window[i]);
fclose(out);
countermid+=encode.block_size[vb.W]/4+encode.block_size[vb.nW]/4;
}
+ _ve_envelope_apply(&vb,1);
vorbis_synthesis_blockin(&decode,&vb);
function: PCM data envelope analysis and manipulation
author: Monty <xiphmont@mit.edu>
modifications by: Monty
- last modification date: Jun 17 1999
+ last modification date: Aug 05 1999
Vorbis manipulates the dynamic range of the incoming PCM data
envelope to minimise time-domain energy leakage from percussive and
********************************************************************/
#include <stdlib.h>
+#include <string.h>
#include <stdio.h>
#include <math.h>
-static typedef struct {
- int divisor;
- double *window;
-} envelope_lookup;
+#include "codec.h"
+#include "envelope.h"
-static double oPI = 3.14159265358979323846;
-
-envelope_lookup *init_envelope(int length,int divleng){
- envelope_lookup *ret=malloc(sizeof(envelope_lookup));
+void _ve_envelope_init(envelope_lookup *e,int samples_per){
int i;
- ret->length=divleng;
- ret->window=malloc(divleng*sizeof(double)*2);
+ e->winlen=samples_per*2;
+ e->window=malloc(e->winlen*sizeof(double));
/* We just use a straight sin^2(x) window for this */
- for(i=0;i<divleng*2;i++){
- double temp=sin((i+.5)/divleng*oPI);
- ret->window[i]=temp*temp;
+ for(i=0;i<e->winlen;i++){
+ double temp=sin((i+.5)/e->winlen*M_PI);
+ e->window[i]=temp*temp;
}
}
-/* right now, we do things simple and dirty. Should this prove
- inadequate, then we'll think of something different. The details
- of the encoding format do not depend on the exact behavior, only
- the format of the bits that come out.
-
- Using residual from an LPC whitening filter to judge envelope
- energy would probably yield cleaner results, but that's slow.
- Let's see if simple delta analysis gives us acceptible results. */
+/* initial and final blocks are special cases. Eg:
+ ______
+ `--_
+ |_______|_`-.___|_______|_______|
-int analyze_envelope0(double *vector, envelope_lookup *init, int n,
- double *deltas){
+ ___
+ _--' `--_
+ |___.-'_|_______|_`-.___|_______|
- int divisor=init->length;
- int divs=n/divisor-1;
- double *win=init->window;
- int i,j,count=0;
+ ___
+ _--' `--_
+ |_______|___.-'_|_______|_`-.___|
- double max,spanlo,spanhi;
-
- /* initial and final blocks are special cases. Eg:
- ______________
- \
- |_______|______\|_______|_______|
-
- ___________
- / \
- |_______|/______|______\|_______|
-
- _____________
- /
- |_______|_______|/______|_______|
+ _____
+ _--'
+ |_______|_______|____.-'|_______|
- as we go block by block, we watch the collective metrics span. If we
- span the threshhold (assuming the threshhold is active), we use an
- abbreviated vector */
+ as we go block by block, we watch the collective metrics span. If we
+ span the threshhold (assuming the threshhold is active), we use an
+ abbreviated vector */
+
+static void _ve_envelope_generate(double *mult,double *env,double *look,
+ int n,int step){
+ int i,j,p;
+ double m;
+ n*=step;
+
+ /* first multiplier special case */
+ m=ldexp(2,mult[0]-1);
+ for(i=0;i<step/2;i++)env[i]=m;
+ p=i;
+ for(i=step;i<step*2;i++,p++)env[p]=m*look[i];
- /* initial frame */
- max=0;
- for(i=1;i<divisor;i++){
- double temp=abs(vector[i-1]-vector[i]);
- if(max<temp)max=temp;
+ /* mid multipliers normal case */
+ for(j=1;p<n-step/2;j++){
+ p-=step;
+ m=ldexp(2,mult[j]-1);
+ for(i=0;i<step;i++,p++)env[p]+=m*look[i];
+ for(;i<step*2;i++,p++)env[p]=m*look[i];
}
- for(;i<divisor*2;i++){
- double temp=abs(win[i-1]*vector[i-1]-win[i]*vector[i]);
- if(max<temp)max=temp;
+
+ /* last multiplier special case */
+ p-=step;
+ m=ldexp(2,mult[j]-1);
+ for(i=0;i<step;i++,p++)env[p]+=m*look[i];
+ for(;p<n;p++)env[p]=m;
+
+ {
+ static int frameno=0;
+ FILE *out;
+ char path[80];
+ int i;
+
+ sprintf(path,"env%d",frameno);
+ out=fopen(path,"w");
+ for(i=0;i<n;i++)
+ fprintf(out,"%g\n",env[i]);
+ fclose(out);
+
+ frameno++;
}
- spanlo=spanhi=deltas[count++]=max;
-
- /* mid frames */
- for(j=divisor;j<n-divisor*2;j+=divisor){
- max=0;
- for(i=1;i<divisor*2;i++){
- double temp=abs(win[i-1]*vector[j+i-1]-win[i]*vector[j+i]);
- if(max<temp)max=temp;
+}
+
+/* right now, we do things simple and dirty (read: our current preecho
+ is a joke). Should this prove inadequate, then we'll think of
+ something different. The details of the encoding format do not
+ depend on the exact behavior, only the format of the bits that come
+ out.
+
+ Mark Taylor probably has much witter ways of doing this... Let's
+ see if simple delta analysis gives us acceptible results for now. */
+
+static void _ve_deltas(double *deltas,double *pcm,int n,double *win,
+ int winsize){
+ int i,j,p=winsize/2;
+ for(j=0;j<n;j++){
+ p-=winsize/2;
+ for(i=0;i<winsize-1;i++,p++){
+ double temp=fabs(win[i]*pcm[p]-win[i+1]*pcm[p+1]);
+ if(deltas[j]<temp)deltas[j]=temp;
}
- deltas[count++]=max;
- if(max<spanlo)spanlo=max;
- if(max>spanhi)spanhi=max;
- if(threshhold>1 && spanlo*threshhold<spanhi)
- abbrevflag=1;
- if(abbrevflag && j>n0-divisor/2)break;
+ p++;
}
+}
- /* last frame */
- if(!abbrevflag){
- max=0;
- for(i=1;i<divisor;i++){
- double temp=abs(win[i-1]*vector[j+i-1]-win[i]*vector[j+i]);
- if(max<temp)max=temp;
- }
- for(;i<divisor*2;i++){
- double temp=abs(vector[j+i-1]-vector[j+i]);
- if(max<temp)max=temp;
+void _ve_envelope_multipliers(vorbis_dsp_state *v){
+ int step=v->samples_per_envelope_step;
+ static int frame=0;
+
+ /* we need a 1-1/4 envelope window overlap begin and 1/4 end */
+ int dtotal=(v->pcm_current-step/2)/v->samples_per_envelope_step;
+ int dcurr=v->envelope_current;
+ double *window=v->ve.window;
+ int winlen=v->ve.winlen;
+ int pch,ech;
+ vorbis_info *vi=&v->vi;
+
+ if(dtotal>dcurr){
+ for(ech=0;ech<vi->envelopech;ech++){
+ double *mult=v->multipliers[ech]+dcurr;
+ memset(mult,0,sizeof(double)*(dtotal-dcurr));
+
+ for(pch=0;pch<vi->channels;pch++){
+
+ /* does this channel contribute to the envelope analysis */
+ if(vi->envelopemap[pch]==ech){
+
+ /* we need a 1/4 envelope window overlap front and back */
+ double *pcm=v->pcm[pch]+dcurr*step-step/2;
+ _ve_deltas(mult,pcm,dtotal-dcurr,window,winlen);
+
+ }
+ }
}
- deltas[count++]=max;
- if(max<spanlo)spanlo=max;
- if(max>spanhi)spanhi=max;
- if(threshhold>1 && spanlo*threshhold<spanhi)
- abbrevflag=1;
+ v->envelope_current=dtotal;
+ frame++;
}
+}
- if(abbrevflag)return(n0);
- return(n);
-}
-
-/* also decide if we're going with a full sized or abbreviated
- vector. Some encoding tactics might want to use envelope massaging
- fully and discard abbreviated vectors entriely. We make that
- decision here */
-
-int analyze_envelope1(envelope_lookup *init,int n,
- double triggerthresh,double spanthresh,
- double *deltas){
-
- /* Look at the delta values; decide if we need to do any envelope
- manipulation at all on this vector; if so, choose the
- multipliers and placeholders.
-
- '0' is a placeholder. Other values specify a
- multiplier/divisor. Multipliers are used by the decoder, divisors
- in the encoder. The mapped m/d value for each segment is
- 2^(n-1). Placeholders (zeros) take on the value of the last
- non-zero multiplier/divisor. When the placeholder is not
- preceeded by a non-placeholder value in the current vector, it
- assumes the value of the *next* non-zero value. In this way, the
- vector manipulation is local to the current vector and does not
- rely on preceeding vectors.
-
- */
-
- /* scan forward with sliding windows; we start manipulating envelopes
- when the collective deltas span over a threshhold. If in fact we
- begin manipulating, we can manage on a finer scale than the
- original threshhold. first look for the larger threshhold and if
- we span it, manipulate the vector to hold within the smaller span
- threshhold. */
-
- /* scan for the trigger */
-
- int divisor=init->length;
- int divs=n/divisor-1;
- int i,triggerflag=0;
- double spanlo,spanhi;
-
- spanlo=spanhi=deltas[0];
-
- for(i=1;i<divs;i++){
- double max=deltas[i];
- if(max<spanlo)spanlo=max;
- if(max>spanhi)spanhi=max;
- if(spanlo*triggerthresh<spanhi){
- triggerflag=1;
- break;
+/* This readies the multiplier vector for use/coding. Clamp/adjust
+ the multipliers to the allowed range and eliminate unneeded
+ coefficients */
+
+void _ve_envelope_sparsify(vorbis_block *vb){
+ int ch;
+ for(ch=0;ch<vb->vd->vi.envelopech;ch++){
+ int flag=0;
+ double *mult=vb->mult[ch];
+ int n=vb->multend;
+ double first=mult[0];
+ double last=first;
+ double clamp;
+ int i;
+
+ /* are we going to multiply anything? */
+
+ for(i=1;i<n;i++){
+ if(mult[i]>=last*vb->vd->vi.preecho_thresh){
+ flag=1;
+ break;
+ }
+ if(i<n-1 && mult[i+1]>=last*vb->vd->vi.preecho_thresh){
+ flag=1;
+ break;
+ }
+ last=mult[i];
}
+
+ if(flag){
+ /* we need to adjust, so we might as well go nuts */
+
+ int begin=i;
+ clamp=last?last:1;
+
+ for(i=0;i<begin;i++)mult[i]=0;
+
+ last=1;
+ for(;i<n;i++){
+ if(mult[i]/last>clamp*vb->vd->vi.preecho_thresh){
+ last=mult[i]/vb->vd->vi.preecho_clamp;
+
+ mult[i]=floor(log(mult[i]/clamp/vb->vd->vi.preecho_clamp)/log(2))-1;
+ if(mult[i]>15)mult[i]=15;
+ }else{
+ mult[i]=0;
+ }
+ }
+ }else
+ memset(mult,0,sizeof(double)*n);
}
+}
- if(triggerflag){
- /* choose divisors/multipliers to fit the vector into the
- specified span. In the decoder, these values are *multipliers*, so */
-
-
-
-
-
-
-
-
-
+void _ve_envelope_apply(vorbis_block *vb,int multp){
+ vorbis_info *vi=&vb->vd->vi;
+ double env[vb->multend*vi->envelopesa];
+ envelope_lookup *look=&vb->vd->ve;
+ int i,j,k;
+
+ for(i=0;i<vi->envelopech;i++){
+ double *mult=vb->mult[i];
+ double last=1.;
+
+ /* fill in the multiplier placeholders */
+
+ for(j=0;j<vb->multend;j++){
+ if(mult[j]){
+ last=mult[j];
+ }else
+ mult[j]=last;
+ }
+ /* compute the envelope curve */
+ _ve_envelope_generate(mult,env,look->window,vb->multend,vi->envelopesa);
+ /* apply the envelope curve */
+ for(j=0;j<vi->channels;j++){
+ /* check to see if the generated envelope applies to this channel */
+ if(vi->envelopemap[j]==i){
+
+ if(multp)
+ for(k=0;k<vb->multend*vi->envelopesa;k++)
+ vb->pcm[j][k]*=env[k];
+ else
+ for(k=0;k<vb->multend*vi->envelopesa;k++)
+ vb->pcm[j][k]/=env[k];
+ }
+ }
}
- return(triggerflag);
}
projects / platform / upstream / libvorbis.git / commitdiff