lib/envelope.c

   1 /********************************************************************
   2  *                                                                  *
   3  * THIS FILE IS PART OF THE Ogg Vorbis SOFTWARE CODEC SOURCE CODE.  *
   4  * USE, DISTRIBUTION AND REPRODUCTION OF THIS SOURCE IS GOVERNED BY *
   5  * THE GNU PUBLIC LICENSE 2, WHICH IS INCLUDED WITH THIS SOURCE.    *
   6  * PLEASE READ THESE TERMS DISTRIBUTING.                            *
   7  *                                                                  *
   8  * THE OggSQUISH SOURCE CODE IS (C) COPYRIGHT 1994-1999             *
   9  * by 1999 Monty <monty@xiph.org> and The XIPHOPHORUS Company       *
  10  * http://www.xiph.org/                                             *
  11  *                                                                  *
  12  ********************************************************************
  13
  14  function: PCM data envelope analysis and manipulation
  15  author: Monty <xiphmont@mit.edu>
  16  modifications by: Monty
  17  last modification date: Aug 05 1999
  18
  19  Vorbis manipulates the dynamic range of the incoming PCM data
  20  envelope to minimise time-domain energy leakage from percussive and
  21  plosive waveforms being quantized in the MDCT domain.
  22
  23  ********************************************************************/
  24
  25 #include <stdlib.h>
  26 #include <string.h>
  27 #include <stdio.h>
  28 #include <math.h>
  29
  30 #include "codec.h"
  31 #include "envelope.h"
  32
  33 void _ve_envelope_init(envelope_lookup *e,int samples_per){
  34   int i;
  35
  36   e->winlen=samples_per*2;
  37   e->window=malloc(e->winlen*sizeof(double));
  38
  39   /* We just use a straight sin^2(x) window for this */
  40   for(i=0;i<e->winlen;i++){
  41     double temp=sin((i+.5)/e->winlen*M_PI);
  42     e->window[i]=temp*temp;
  43   }
  44 }
  45
  46 /* initial and final blocks are special cases. Eg:
  47    ______
  48          `--_
  49    |_______|_`-.___|_______|_______|
  50
  51               ___
  52           _--'   `--_
  53    |___.-'_|_______|_`-.___|_______|
  54
  55                       ___
  56                   _--'   `--_
  57    |_______|___.-'_|_______|_`-.___|
  58
  59                                _____
  60                            _--'
  61    |_______|_______|____.-'|_______|
  62
  63    as we go block by block, we watch the collective metrics span. If we
  64    span the threshhold (assuming the threshhold is active), we use an
  65    abbreviated vector */
  66
  67 static void _ve_envelope_generate(double *mult,double *env,double *look,
  68                                   int n,int step){
  69   int i,j,p;
  70   double m;
  71   n*=step;
  72
  73   /* first multiplier special case */
  74   m=ldexp(2,mult[0]-1);
  75   for(i=0;i<step/2;i++)env[i]=m;
  76   p=i;
  77   for(i=step;i<step*2;i++,p++)env[p]=m*look[i];
  78
  79   /* mid multipliers normal case */
  80   for(j=1;p<n-step/2;j++){
  81     p-=step;
  82     m=ldexp(2,mult[j]-1);
  83     for(i=0;i<step;i++,p++)env[p]+=m*look[i];
  84     for(;i<step*2;i++,p++)env[p]=m*look[i];
  85   }
  86
  87   /* last multiplier special case */
  88   p-=step;
  89   m=ldexp(2,mult[j]-1);
  90   for(i=0;i<step;i++,p++)env[p]+=m*look[i];
  91   for(;p<n;p++)env[p]=m;
  92
  93   {
  94     static int frameno=0;
  95     FILE *out;
  96     char path[80];
  97     int i;
  98
  99     sprintf(path,"env%d",frameno);
 100     out=fopen(path,"w");
 101     for(i=0;i<n;i++)
 102       fprintf(out,"%g\n",env[i]);
 103     fclose(out);
 104
 105     frameno++;
 106   }
 107 }
 108
 109 /* right now, we do things simple and dirty (read: our current preecho
 110    is a joke).  Should this prove inadequate, then we'll think of
 111    something different.  The details of the encoding format do not
 112    depend on the exact behavior, only the format of the bits that come
 113    out.
 114
 115    Mark Taylor probably has much witter ways of doing this...  Let's
 116    see if simple delta analysis gives us acceptible results for now.  */
 117
 118 static void _ve_deltas(double *deltas,double *pcm,int n,double *win,
 119                        int winsize){
 120   int i,j,p=winsize/2;
 121   for(j=0;j<n;j++){
 122     p-=winsize/2;
 123     for(i=0;i<winsize-1;i++,p++){
 124       double temp=fabs(win[i]*pcm[p]-win[i+1]*pcm[p+1]);
 125       if(deltas[j]<temp)deltas[j]=temp;
 126     }
 127     p++;
 128   }
 129 }
 130
 131 void _ve_envelope_multipliers(vorbis_dsp_state *v){
 132   int step=v->samples_per_envelope_step;
 133   static int frame=0;
 134
 135   /* we need a 1-1/4 envelope window overlap begin and 1/4 end */
 136   int dtotal=(v->pcm_current-step/2)/v->samples_per_envelope_step;
 137   int dcurr=v->envelope_current;
 138   double *window=v->ve.window;
 139   int winlen=v->ve.winlen;
 140   int pch,ech;
 141   vorbis_info *vi=&v->vi;
 142
 143   if(dtotal>dcurr){
 144     for(ech=0;ech<vi->envelopech;ech++){
 145       double *mult=v->multipliers[ech]+dcurr;
 146       memset(mult,0,sizeof(double)*(dtotal-dcurr));
 147
 148       for(pch=0;pch<vi->channels;pch++){
 149
 150         /* does this channel contribute to the envelope analysis */
 151         if(vi->envelopemap[pch]==ech){
 152
 153           /* we need a 1/4 envelope window overlap front and back */
 154           double *pcm=v->pcm[pch]+dcurr*step-step/2;
 155           _ve_deltas(mult,pcm,dtotal-dcurr,window,winlen);
 156
 157         }
 158       }
 159     }
 160     v->envelope_current=dtotal;
 161     frame++;
 162   }
 163 }
 164
 165 /* This readies the multiplier vector for use/coding.  Clamp/adjust
 166    the multipliers to the allowed range and eliminate unneeded
 167    coefficients */
 168
 169 void _ve_envelope_sparsify(vorbis_block *vb){
 170   int ch;
 171   for(ch=0;ch<vb->vd->vi.envelopech;ch++){
 172     int flag=0;
 173     double *mult=vb->mult[ch];
 174     int n=vb->multend;
 175     double first=mult[0];
 176     double last=first;
 177     double clamp;
 178     int i;
 179
 180     /* are we going to multiply anything? */
 181
 182     for(i=1;i<n;i++){
 183       if(mult[i]>=last*vb->vd->vi.preecho_thresh){
 184         flag=1;
 185         break;
 186       }
 187       if(i<n-1 && mult[i+1]>=last*vb->vd->vi.preecho_thresh){
 188         flag=1;
 189         break;
 190       }
 191       last=mult[i];
 192     }
 193
 194     if(flag){
 195       /* we need to adjust, so we might as well go nuts */
 196
 197       int begin=i;
 198       clamp=last?last:1;
 199
 200       for(i=0;i<begin;i++)mult[i]=0;
 201
 202       last=1;
 203       for(;i<n;i++){
 204         if(mult[i]/last>clamp*vb->vd->vi.preecho_thresh){
 205           last=mult[i]/vb->vd->vi.preecho_clamp;
 206
 207           mult[i]=floor(log(mult[i]/clamp/vb->vd->vi.preecho_clamp)/log(2))-1;
 208           if(mult[i]>15)mult[i]=15;
 209         }else{
 210           mult[i]=0;
 211         }
 212       }
 213     }else
 214       memset(mult,0,sizeof(double)*n);
 215   }
 216 }
 217
 218 void _ve_envelope_apply(vorbis_block *vb,int multp){
 219   vorbis_info *vi=&vb->vd->vi;
 220   double env[vb->multend*vi->envelopesa];
 221   envelope_lookup *look=&vb->vd->ve;
 222   int i,j,k;
 223
 224   for(i=0;i<vi->envelopech;i++){
 225     double *mult=vb->mult[i];
 226     double last=1.;
 227
 228     /* fill in the multiplier placeholders */
 229
 230     for(j=0;j<vb->multend;j++){
 231       if(mult[j]){
 232         last=mult[j];
 233       }else
 234         mult[j]=last;
 235     }
 236
 237     /* compute the envelope curve */
 238     _ve_envelope_generate(mult,env,look->window,vb->multend,vi->envelopesa);
 239
 240     /* apply the envelope curve */
 241     for(j=0;j<vi->channels;j++){
 242
 243       /* check to see if the generated envelope applies to this channel */
 244       if(vi->envelopemap[j]==i){
 245
 246         if(multp)
 247           for(k=0;k<vb->multend*vi->envelopesa;k++)
 248             vb->pcm[j][k]*=env[k];
 249         else
 250           for(k=0;k<vb->multend*vi->envelopesa;k++)
 251             vb->pcm[j][k]/=env[k];
 252
 253       }
 254     }
 255   }
 256 }
 257